texts: Splitting polyline into list of subpolylines by distance using PyQGIS?

texts:
I'm actually working on a way to print parts of a Polyline according to both specified paper size and scale. The paper should be rotated in order to maximize the visualisation comfort, and this must be automated.

Note : I'm using QGIS 2.18.x, and can't migrate to 3.x due to major internal extensions that haven't been switched from Python 2.x to 3.x yet.

Here are the steps I followed so far:

Get the line length, and divide it by the page width in order to know the number of pages to print.
Divide line length by this number to get the average distance of cable to print on each page (in order to avoid the last page to be partially used)
Subdivide the polyline using the average distance. This is where I'm actually stuck

I found a way to get the points where I should cut, but I didn't found any function to subdivide a line using a distance. (Apart from splitGeometry(), which doesn't works with a list of points. It doesn't make them collide the line. Also giving a polyline version of this list cuts too much, as it splits at each interesection of the lines, not only the specified points)

I only succeed to get the points that should cut the line with this code :

pointsList = []
for i in range (1, pagesNb):
    point = line.interpolate(dist * i)
    pointsList.append(point.asPoint())

Make a Convex Hull of each subline and search for the longuest included line
Make the orientation of the page parallel to this line and adjust the page using the centroid and the page dimensions

Is there a function to cut a line by a ratio, a number of segments or by a list of points?

I'm also interested in any hints/improvements on this method.

I am aware that a minimumOrientedBoundingBox() function exists, the problem is that it's not present in the QGIS version I use.

Nethertheless, if there is a way to catch the sourcecode of it and use it, I will also take that as a solution (not the one present on the toolbox and that can be called with processing, only the PyQGIS function)


texts: pgrouting 2: missing pgr_findnearestnodedwithin function

texts:
With the PgRouting 1.5 I was using the find_nearest_node_within_distance function and it works, but now I'm moving the project to PgRouting 2 + Postgis 2 + PostgreSQL 9.1, and now that function, I think, was renamed to pgr_findnearestnodedwithin but when I try to execute the query I get the following error

ERROR:  function pgr_findnearestnodedwithin(unknown, numeric, unknown) does not exist

I was looking if there is a new way to import this function but I didn't have success. please help. there are something I'm missing?

Thanks in advance


texts: how to not print scalebar and how to print the zoom bar that shows on map

texts:
I am using javascript api for ArcGIS and need to print the zoom bar on map, but don't want the scale information to be printed.

How can I achieve this?


texts:
I am using "date '2015-8-10'" in ArcGIS Server REST API to retreive records from map services. But it returns the error "Fail to Execute". My date is stored as "10-AUG-15".

When it is in the right date format, as it is stored in map service, it returns an accurate result. It is not matching the Esri date format with date format in the attribute table of Map service. Column type is also specified as date in the Map service attribute table.

Do I need to make some changes in ArcGIS for Server to convert the Esri date format to my stored format?

I am using Oracle Spatial feature class in ArcMap and then creating map service on ArcGIS Server. Type of column that contains dates is also date in the Oracle table.


texts: How to put a roads layer (OpenStreetMap) on top of KML polygons?

texts: I need to develop a map where KMLs are overlayed with road layers. So, KMLs are polygons, and I would need to use a map tiles like OpenStreetMaps. But the problem is how to put roads on top of KML areas? What do you suggest me? I would like to use some open source software.


texts:
I'm digitising a large number of maps from hand drawn mud maps and digital images, and although I have used the snapping tools ("snap to vertex and segment", Tolerance 10 pixels, avoid intersections) I consistently get ~30 geometry errors. My workflow is as follows:

digitise (polygons) using snapping settings as above
Use GRASS tools to remove most bad geometry (v.in.ogr.qgis, snapping threshold 0.02 min size 0.0001)
Dissolve polygons with same attribute
Check geometry validity

I consistently get around 30 errors that persist, but if I tweak the GRASS setting anymore it changes the map into an unrecognisable blob. Tutorials like this one tell you how to find the errors and then it say "fix them", but I guess I haven't figured out a good way to fix them? The errors are always of this type, a green cross where polygons meet and their handles overlap:

It seems the only way I can figure out to fix them is to nudge one polygon slightly to the side, so it's no longer overlapping that one point anymore, but I worry that will mess with any spatial analyses I run?

My questions are:

How do I avoid getting all these errors?
How do I fix them if they do occur?
Are they important?

I'm not a QGIS power user.

This is what the problem looks like when I nudge the polygons apart, the top three polygons all join in the same node, while the fourth (bottom right) is the one that is overlapping and causing the problem.

zoomed out (area ~1x1m), error highlighted with red cross


texts: How to show hierarchical clusters in google maps or other tool?

texts:
Google maps marker can be clustered but those clusters expand on different levels.

I want to stay on a specific zoom level and cluster a set of markers. Then I want to show sub clusters inside that cluster. I can not find any way to do that in Google Maps. Is there any other tool which I can use to do that?

I want to do something like these images


texts: In Mapbasic, How to populate my table row data in a dialog box with RadioGroup and EditText

texts:
Good day!

Im trying to pupulate my table row data in a dialog box with Radiogroup and editText control.

I successfully get the column value list from popupMenu control after selecting the table and column popupmenu.

But i wanted to populate from the column value list i selected, the row data from my selected table and update the Radiogroup and EditText every time i select a column value list.

br,

Ivan


texts:
I am creating a plugin that digitizes a line by "add line feature " in the "Toggle Editing" toolbar, I want the digitized line to be between two points of a point layer which is a layer different than the line layer to avoid adding any random line.

Is controlling the line digitizing possible in PyQGIS?

The digitized data are in PostgreSQL tables


texts:
I'm using a webmap on my page. The webmap contains layers with popups. These popups contain a link with url GET parameters (ex: www.example.com?id=123)

Is there a way to force the link to open in the same window? How do I capture the GET parameter on this page?


texts: We have some geometry records - water bodies in a DB which does not have name. I have tried to use Google Geocoding to extract them but it seems that it is more suitable for addresses. Is there any API to get water body names by coordinates?


texts: Not able to overwrite/delete an existing Service in ArcGIS Server?

texts:
I am using ArcGIS 10.3 for Server.

I have an existing service called xyz in a folder called DATA published using xyz.mxd. Now when I am trying to overwrite the existing service, it is giving me an error "ERROR 001369: Failed to create the service". So I thought to delete the service and create a new one, but it doesn't allow me to delete the service either and gives me an error "Access Denied" even though I have admin rights.

I checked the logs on ArcGIS Server Manager while trying to publish the service and the following error came up:

java.lang.RuntimeException: Unable to rename old path 'E:\arcgisserver\directories\arcgissystem\arcgisinput\Data/OS.MapServer' to new path 'E:\arcgisserver\directories\arcgissystem\arcgisinput\Data/OS_201604151509.MapServer'.
at com.esri.arcgis.discovery.util.FileUtil.renameTo(FileUtil.java:81)
at com.esri.arcgis.discovery.admin.services.ServiceManager.renameService(ServiceManager.java:856)
at com.esri.arcgis.discovery.admin.rest.services.FolderServlet.b(FolderServlet.java:444) at

I am not sure why it is trying to rename the file!!! The issue is with this particular folder because when I create a new service in this folder it again give an error.

On deletion, the following error comes:

Could not connect to the identity store as one or more of the connection parameters is incorrect. Verify that you can connect to the identity store outside of ArcGIS Server using the same parameters. [LDAP: error code 49 - 80090308: LdapErr: DSID-0C0903A9, comment: AcceptSecurityContext error, data 52e, v1db1]


texts: ogr2ogr polygon to point with point inside the polygon

texts: Can ogr2ogr be used to convert a polygon to a point but the point must remain inside the polygon? It does not matter where the point is as long as it is inside the object polygon.


texts: How draw coordinate points automatically by Python in QGIS

texts:
I am extracting the coordinates of the vertex lines, but I could not have created them, I tried with itself code:

# create a memory layer with two points
layer =  QgsVectorLayer('Point', 'points' , "memory")
pr = layer.dataProvider()
# add the first point
pt = QgsFeature()
for i in list_glob:
point1 = QgsPoint(i[0],i[1])
pt.setGeometry(QgsGeometry.fromPoint(point1))
pr.addFeatures([pt])
# update extent of the layer
layer.updateExtents()
# update extent
layer.updateExtents()
# add the layer to the canvas
QgsMapLayerRegistry.instance().addMapLayers([layer])

but it does not work


texts: Estimating the Distance Between Two Points

texts:
I am wondering how to estimate the distance between two points on a map. For example, consider the distance from A to B on the map below.

If the latitude and longitude lines shown are separated by 10°, how is it possible to estimate the distance AB.


texts:
I followed the answer given in this topic for searching the distance to the nearest coastline given a lat/lon coordinate. The distance output obtained is in degree (I think) but I need the distance expressed in meters.

By citing the answer given:

  There will be a units problem if you create the distance raster in
  geographic coordinates. In order to get kilimetres you will need to
  rasterize the sea in a suitable projected coordinate system and then
  project the distance raster to geograpic coordinates to do the lookup.
  This will give you the fastest response time as you're not projecting the
  query coordinate at the time of lookup.

At the moment my python code is the following:

from osgeo import gdal
ds = gdal.Open("/home/user/distance.tif")

gt = ds.GetGeoTransform()

topLeftX = gt[0]
dimX = gt[1]
topLeftY = gt[3]
dimY = gt[5]

band = ds.GetRasterBand(1)
#Numpy array of the distances in the whole area
data = band.ReadAsArray()
#Indexing by longitude/latitude position

def getDistance(coord):
    x = int((coord[0] - topLeftX)/dimX)
    y = int((coord[1] - topLeftY)/dimY)
    distance = data[y,x]
    return distance

The coordinate system is WGS84. I think I have to rasterize the shape file in a Projected Coordinate System but then I don't know what to do. Any suggestions?


texts: QGIS crashing when opening Layer Properties in v2.14 and after upgrading to 2.18

texts:
I have a similar problem to this post I believe the issue began after I re-arranged my directory structure (which included some .shp files).

Now I am unable to open any layer properties in any .qgs document, even after upgrading to 2.18 as per these instructions as well as Kyngchaos.

I am running OSX 10.11.5


texts:
I would like to plot a tilted map, similar to the image below. I would prefer to do this using ggplot2 but I'm glad if there is a solution using base R. Any ideas on how to do this?

Here is a dataset for a reproducible example:

library(rgeos)
library(UScensus2000tract)
library(ggplot2)

# load data
  data("oregon.tract")

# plot Census Tract map
  plot(oregon.tract)


texts: Geopandas and significant digits of coords


texts: QGIS GDAL Build VRT results do not reflect input raster values

texts:
I am having trouble successfully creating a VRT.

My Build VRT gdal command results in values of ±1.79769e+308. The values should be 1,2,3, or 4.

I've verified projection info is the same across the input rasters.

Input rasters (.tifs) are:

FLT4S (also tried INT1U)
LZW compressed
1 m pixel size

There are large nodata spaces in between the rasters, but setting no data to -9999, does not change results.

±1.79769e+308 is strangely the range of values for 32 bit floating data types, (FLT8S).

Is this a clue?

I am using QGIS 3.4.14 on Mac

When I use a different set of input rasters GDAL VRT works fine. So I think it has something to do with the input raster types, rather than GDAL or QGIS.

This is the GDAL command:

gdalbuildvrt -resolution highest -r nearest -input_file_list /private/var/vw/buildvrtInputFiles.txt /Users/forVRT/vrt.vrt


texts: Iterating through sub-element `sections` of GeoJSON on a Leaflet map?

texts:
I have the following code which I use to display a marker on my map with a popup containing the feature.properties.name element from my GeoJSON:

var featuresUTMGroundHazards = L.geoJson(dataGuardianUTM, {
    filter: FuncGroundHazards,
    style: setStyle,
    pointToLayer: function (feature, latlng) { return L.marker(latlng, { icon: L.icon({ iconUrl: feature.properties.iconUrl, iconSize: [25, 25], }), }) },
    onEachFeature: function (feature, layer) {
        layer.bindPopup(
            "" + feature.properties.name + ""
            + "
"
            );
    },
});

If they are present in the GeoJSON, I wish to expand this and itterare through any and all available properties.display.sections elements and include the title and text elements from each, within the marker popup.

I have tried the following code but it displays undefined where I would expect some text to be.

var featuresUTMGroundHazards = L.geoJson(dataGuardianUTM, {
    filter: FuncGroundHazards,
    style: setStyle,
    pointToLayer: function (feature, latlng) { return L.marker(latlng, { icon: L.icon({ iconUrl: feature.properties.iconUrl, iconSize: [25, 25], }), }) },
    onEachFeature: function (feature, layer) {
        layer.bindPopup(
            "" + feature.properties.name + ""
            + "
"
            + "
Title: " + feature.properties.display.sections.title
            + "
Text: " + feature.properties.display.sections.text
            + "
"
            );
    },
});

I would have expected that to simply display the first one it found, but it does not.

What do I need to change here in order to itterate through all/any of the sections?

The GeoJSON is structured as:

{
  "geometry": {
    "coordinates": [
      [
        [
          -0.58585870265960693,
          51.223640441894531
        ],
        // etc...

      ]
    ],
    "type": "Polygon"
  },
  "id": "8A7183778EF92F2FB1A84F5C02A7C233ED71B048",
  "properties": {
    "hazardFactor": "40",
    "hazardFactorName": "Warning",
    "fillColor": "#ffbb00",
    "strokeColor": "#b88702",
    "fillOpacity": "0.35",
    "strokeWidth": "1",
    "strokeOpacity": "0.8",
    "detailedCategory": "amenity:university",
    "iconUrl": "https://aa-ne-prod-public-api.altitudeangel.com//v1/map/icon?icon=amenity_university.png",
    "name": "University of Law, Guildford",
    "category": "groundHazard",
    "filters": [
      {
        "name": "Ground Hazards",
        "property": "show",
        "active": true
      }
    ],
    "display": {
      "category": "Ground Hazard",
      "detailedCategory": "University",
      "title": "University of Law, Guildford",
      "sections": [
        {
          "iconUrl": "https://aa-ne-prod-public-api.altitudeangel.com//v1/map/icon?icon=amenity_university.png",
          "title": "Pedestrian Hazard",
          "text": "It is likely that members of the public will congregate in and/or around this area. Please be mindful of their safety."
        },
        {
          "iconUrl": "https://aa-ne-prod-public-api.altitudeangel.com//v1/map/icon?icon=amenity_university.png",
          "title": "Privacy Concerns",
          "text": "This classification of area or point is likely to be additionally sensitive to privacy concerns. You may have greater obligations to disclose your intended flight over, within or near these areas to any controlling authority or land owner, particularly if you plan to fly with surveillance equipment."
        },
        {
          "iconUrl": "https://aa-ne-prod-public-api.altitudeangel.com//v1/map/icon?icon=warning.png",
          "title": "Summary",
          "text": "Yellow zones indicate regions where operation of your drone may raise security, privacy or safety concerns."
        }
      ],
      "actions": []
    }
  },
  "type": "Feature"
},


texts:
I need to create a new sde service with the command line for the first time with mssql on ArcGIS Server 10.0 so that 3rd party software can access it. Direct connect will not work unfortunately.

As far as I can understand I need to do the following:

sdeservice -o create -d SQLSERVER -p spatial.data -i esri_sde -H z:\arcgis\arcsde\sqlexe
manually edit %windir%\system32\drivers\etc\services and %SDEHOME%\etc\services.sde

How could I find out which port numbers should I use?


texts: Buttons that toggle different column visualizations?

texts:
I have a census data table in Carto and I want the end user to be able to toggle through several different demographic visualizations via buttons on the map.

For example one button shows median age, another shows median income, etc. All of these are in the same table just different columns.

Any idea of how I would do this?

I am currently more of a GISer than a coder.


texts:
Do any of you use the Clear Cache function much OR upped your cache from the default 100mb?

We export PDFs from MXDs and it has been reported to us by a few that upping this number or clearing the cache helps generate PDF exports successfully - I dont find much documentation on this in HELP and was wondering if any one had any thoughts on this:

You can find it : ArcMAP (931) > Tools > Options > Data Interoperability > Cache "CLEAR CACHE"

Found a thread on the ESRI forums: http://forums.esri.com/Thread.asp?c=93&f=982&t=235621


texts: Can we create views over a raster table in Postgis2.0?

texts:
Is it possible to create PostgreSQL views over a table that contains raster data?

I have succeeded in creating the view but the entry in raster_columns is not correct (all the raster info is missing). Also, I couldn't open the SQL view using QGIS (most probably because of the raster_columns not being correct).

It's the same situation as this question on gemetry_columns & views but for the raster column I couldn't find the correct syntax to force the typemod during creation of the view.

All that I could find is this:

  raster(ARRAY['8BUI','8BUI','8BUI'], 4326,0.5,-0.5, ...);

but what follows after "..." is a mystery.

Has anybody managed to create a view with correct entries in raster_columns?

A short example will be highly appreciated.


texts: (QGIS) How can I join overlapping polygons, retaining each overlapping pairing with proportions?

texts:
I have two polygon shapefiles--blocks groups (BGs) and neighborhoods. Generally, BGs are smaller than neighborhoods. Sometimes they lie entirely within a neighborhood, but often they overlap two or more.

I can merge these so that the first joining attribute (neighborhood) is retained, but this isn't what I need. How can I join these two files so that the attribute table contains a row for each BG-Neighborhood pairing, with a variable giving the proportion of the BG contained in that neighborhood?


texts: Deleting existing data and thereafter appending updated data to hosted feature layer (layer within group layer) in ArcGIS Online

texts:
I am trying to figure out how to append new data to an existing feature layer in ArcGIS Online. I have tried to solely append the data (in ArcGIS Online), but the problem arises when I receive duplicated data in the attribute table from this action. What I want to do is to remove all original data in the feature layer and thereafter add the new data with its expanded geometrical properties.
Does anyone have any solutions to this problem?


texts: Is it possible to change size of an annotation layer?

texts: So a project got passed on to me where the previous guy annotated his labels into a geodatabase, but made the labels too small. Is there any way I can change the size of the labels? I've tried selecting the labels individually and changing the size but it will not change. I've also tried to go into the properties but it's locked.


texts: How do I adapt existing responsive web site to display Leaflet maps?

texts:
My website works well on different devices but it seems that only works with images, although I would like to show other media, like Leaflet maps and videos.

Just about all pages are laid out on a grid using classes.

Here is the current html code for one image:

Cityscapes have a uniquely vibrant feel to them which is unlike any other form of landscape. In creating ours we employ a variety of techniques, including hand drawing, watercolor and digital.

Downtown Manhattan

Each image has three size versions and depending on the width of the viewer's screen the right one is loaded.

In the past I tried to display videos in the same framework (span) but I was not able, so I settled for screen shots with links to an external page to see the video.

Now I want to be able to display Leaflet maps on the same grid but I tried various html code and it didn't work.

For example, I tried the code below, but regardless of the width of the screen the div with 298px is loading:

Any advice?


texts: Reproject all shapefiles in same folder in same projection with raster?

texts:
I have 4 shapefiles in the same folder names like that (1996.shp, 2000.shp, 2006.shp, 2014.shp) and one raster.

I need to reproject all shapfiles in the same projection (4326) of raster, but I do it for one shapefile that it script :

from osgeo import ogr, osr, gdal

    #tif with projections I want
    tif = gdal.Open("D:/projetp/NDVI.tif")

    #shapefile with the from projection
    driver = ogr.GetDriverByName("ESRI Shapefile")
    dataSource =   driver.Open("C:/projetp/shape/1996.shp", 1)
    layer = dataSource.GetLayer()

    #set spatial reference and transformation
    sourceprj = layer.GetSpatialRef()
    targetprj = osr.SpatialReference(wkt = tif.GetProjection())
    transform = osr.CoordinateTransformation(sourceprj, targetprj)

    to_fill = ogr.GetDriverByName("Esri Shapefile")
    ds = to_fill.CreateDataSource("C:/projetp/shape/1996poly.shp")
    outlayer = ds.CreateLayer('', targetprj, ogr.wkbPolygon)
    outlayer.CreateField(ogr.FieldDefn('id', ogr.OFTInteger))

    #apply transformation

    i = 0
    for feature in layer:
        transformed = feature.GetGeometryRef()
        transformed.Transform(transform)

        geom = ogr.CreateGeometryFromWkb(transformed.ExportToWkb())
        defn = outlayer.GetLayerDefn()
        feat = ogr.Feature(defn)
        feat.SetField('id', i)
        feat.SetGeometry(geom)
        outlayer.CreateFeature(feat)
        i += 1
        feat = None

    ds = None

How can I execute the script for all shapefiles in the folder?

I'm new in Python 3.


texts: QGIS - Moving vector layer relative to world map

texts: I have a shapefiles with coordinates that are supposed to be in EPSG:4326 but seem to be located at a totally wrong place on the planet. I confirm that by adding an osm layer to. Can I somehow move my original layer relatively to the osm layer to the correct place on the planet and then, save it?


texts:
I am working in an application to edit attribute and features, for which i am trying to use WFST. I have given write access to all users and for all workspaces.Also,i have changed the service level to 'Transactional' in WFS settings. Restarted the tomcat and tried inserting a feature WFS_transactioninsert.xml in the Demo request. But feature is not being added to the layer. is there anything i need to configure?

<!--
  YOU PROBABLY DO NOT WANT TO RUN THIS QUERY SINCE
   IT WILL MODIFY YOUR SOURCE DATA FILES

   It will add a simple line to the tasmania_roads dataset.

    -->
   <wfs:Transaction service="WFS" version="1.0.0"
   xmlns:wfs="http://www.opengis.net/wfs"
   xmlns:topp="http://www.openplans.org/topp"
   xmlns:gml="http://www.opengis.net/gml"
   xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
   xsi:schemaLocation="http://www.opengis.net/wfs
   http://schemas.opengis.net/wfs/1.0.0/WFS-transaction.xsd
  http://www.openplans.org/topp
  http://localhost:8080/geoserver/wfs/DescribeFeatureType?typename=topp:tasmania_roads">
  <wfs:Insert>
   <topp:tasmania_roads>
   <topp:the_geom>
    <gml:MultiLineString srsName="http://www.opengis.net/gml/srs/epsg.xml#4326">
      <gml:lineStringMember>
        <gml:LineString>
          <gml:coordinates decimal="." cs="," ts=" ">
           494475.71056415,5433016.8189323 494982.70115662,5435041.95096618
          </gml:coordinates>
        </gml:LineString>
      </gml:lineStringMember>
    </gml:MultiLineString>
     </topp:the_geom>
     <topp:TYPE>alley</topp:TYPE>
    </topp:tasmania_roads>
   </wfs:Insert>
  </wfs:Transaction>


texts:
Projection systems are always troubling and I cam across and issue yesterday. I created a buffer and did a select by location using that buffer, however it was producing the wrong results. It would select points that were outside the buffer and not select ones that were.

It turned out there was an issue with the projections on one of the initial datasets, with a scale factor in the projection parameters being rounded up.

It got me wondering what is the process of re-projection in QGIS and how accurate is it?

I imagine it is based on GDALwarp?

But particularly projecting to the British National Grid the proj seems very inaccurate:
+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs

What are the parameters being used and can I impact them?


texts: Geoserver Routing with OSM data and PostGIS

texts: Does Geoserver have a routing functionality for the OSM data in PostGIS?


texts: Shapefile stopped displaying polygons and gives error: 'No extent could be determined' in QGIS 3.10

texts:
I have been using a shapefile as a main data layer since Jan 2018 in a project, it forms the basis of several duplicate layers, it has been working fine until today.

Now 95% of the polygons are no longer showing, even though they still exist in the attribute table. If I select one and Zoom to Selection I receive an error message Cannot zoom to selected feature(s): No extent could be determined

This happened after I deleted a few polygons and copied some others across from another shapefile. This other shapefile was a direct copy of same layer from a few weeks ago that I was using as a temporary working layer.

As I have this direct copy I can repair things in my project by working back and resetting up the duplicate layers with their different styles. But things have moved on since this copy and I am trying to find a solution that reverses what has happened to this shapefile.

I have used the Vector > Geometry Tools > Check validity... (both QGIS and GEOS) I found one digitising error which was fixed.

I have searched here on GIS Stack Exchange and found a few questions mentioning 'No extent could be determined' but none of them are relevant to this scenario.

Now I do tend to use Geopackages and PostGIS databases for my work. It just happens that this particular project is legacy and old and it hasn't been updated yet.

Any quick and easy way to fix this shapefile?

I am in CRS EPSG: 27700 British National Grid

EDIT: Subsequent to the comments I can confirm that I do also receive the error Cannot use preferred transform between EPSG: 4326 and EPSG: 27700


texts: Total Area of different polygons in a table

texts: I have a single table in which there is a column "Name" and want to find the total area of all those polygons where value of Name is Cardiff. I have more than one polygons with this name.


texts:
This is a very generic map production question in QGIS.
I am trying to add a "You are here" symbol in QGIS layouts for a fire evacuation map like this one below?

There are more than one "You Are Here" points depending on the location. So I was thinking of adding a point table with all "You Are Here" locations and then in Composer have automatically created these symbols.

Is it possible?


texts:
I am now working on a system to perform routing on a network of shipping lanes.
See this dataset that I have as the base: http://geocommons.com/maps/109850
I have postgres and pgrouting up and running and I was able to run basic queries on the above dataset after importing it to postgres/postgis.
However I have encountered some problems for which I seek advice from you.

A lot of lines have a duplicate, so say line A goes from source id 1 to source id 2 and line B goes from source id 2 to source id 1 and they overlap completely.
This can be experienced a lot of times, and frankly, I would be happy to get rid of them so that finally I could have a network with only one polyline per direction.
I discovered that the distance value attributes of the referenced dataset are false.
I noticed this when I checked the dataset with Global Mapper by using the Info tool.
In the window that pops up, you can see the geometry textbox with bounds and length data.
For a line I inspected that displayed a length value of 5300 some kilometres while the max value in the length column of the dataset was 110. I seek for a way to acquire those polyline length values, and assign them to a new attribute column. I have already found posts about how to do this in qgis and arcgis, but neither of them provided the results I wanted to get.
The route dataset is in WGS84, latlon

I hope the lines above make sense, and I appreciate any help!
Thanks in advance!

Tamas


texts: Custom Geographic Transformation on QGIS

texts:
Is there a way to create a custom geographic transformation with specifics parameters in QGIS like in ArcGIS?

Okay, I want to tranform SAD69 system into WGS84. Well, I'm in doubt because in ArcGIS you have to enter both systems and then specify the method of the transformation - wich is the Molodensky method - and its three parameters - defined by the brazilian government.

a1 = 6.378.160 m
f1 = 1/298,25
a2 = 6.378.137 m
f2 = 1/298,257222101 ΔX = − 67,35 m
ΔY = + 3,88 m
ΔZ = − 38,22 m

Thanks in advance,


texts: Extract Raster Value into Polygon Attribute

texts:
I have Raster Data for a State which has Elevation Data in its pixels. I also have Vector layer which has the 300 divisions of State. what I need is, I want to get the average elevation value of particular division from the raster data using ArcGIS.
I tried extract by mask in ArcGIS to extract individual polygon to clip the raster into 300 raster then from individual raster i have collected average elevation but it has taken days together to complete the work. Is there any other easy way to get it?


texts:
I got an error message when I try to execute the command

SELECT  pgr_analyzeGraph('Highway','the_geom',0.0001);

Error message:

function pgr_analyzegraph(unknown, unknown, numeric) does not exist
LINE 1: SELECT pgr_analyzeGraph('Highway','the_geom',0.0001);
               ^
HINT:  No function matches the given name and argument types. You might need to add explicit type casts.

I confirmed that the pgrouting, postgis, postgis_topology is included in the database.

What is the possible cause of this message?


texts: How to Extract Data by Attributes (Clip/Zip/Ship) using Attribute Where Clause - Esri Server Tools

texts:
I want to use the out-of-the-box Esri geoprocessing task to extract data via ArcGIS Server (Server Tools --> Extract Data).  Its pretty straight-forward, but it doesn't seem to support use of a "where clause" as an input parameter for the service.  I could add this in using python, but would love if it just supported it without custom development.

Has anybody seen a way to use ArcGIS Server to extract data using an attribute where clause WITHOUT custom development.

Thanks ahead of time.


texts: Subtracting unique values from different areas within the same raster

texts:
I have a raster of an area with various elevated features, situated on an overall slope. I want to normalize the elevations of each feature by subtracting the lowest elevation at each feature. This is unique to each feature in the raster. The end goal is to plot the elevation distribution of each feature.

My current approach (which is very rough): I digitized each feature into a shapefile for the whole area. I used this shapefile to extract by mask, and imported the resulting grid into MATLAB. I took the extents from each polygon and found the minimum point on the perimeter and subtracted this value from all the values within the extent "box". However, some of these features are close to one another so I end up subtracting values twice from the same feature.

I also ran Zonal Statistics to calculate the minimum elevation at each polygon. There has to be some simple tool I'm missing to subtract a unique value from the areas of the raster within each polygon, but I can't find a solution.


texts: How can a geographic coordinate system such as WGS84 draw on the screen?

texts: Geographic coordinate system represents coordinates on a sphere or ellipsoid. How does GIS software (ArcGIS, QGIS, MapInfo) draw these coordinates on the screen? Does the software apply a default projection method, or does it just consider the lat/long coordinates as planar?


texts: How to add custom parameters in GeoServer cascading WMS

texts: I want to add additional parameters (like time, layer_name) to my cascaded WMS which is using a Flask API (for converting vector data to raster) as source. The parameters can take dynamic value, so to make changes in capabilities doc ( as mentioned here ) won't solve the problem.


texts: Grading layers according to CMYK color codes in Attribute Table

texts:
I have been trying for some time to color grade a geological map according to CMYK values given in the attribute table. I have not yet succeeded. I am certain that is just a simple line of code that is needed and hope some of you can help me.
I am currently using QGIS 3.14. The column that I want to grade to colors after are the "cmykFargekode" column.


texts:
I am working on an assignment which requires me to connect to a database, then import a dbf file into the database using the command db.in.ogr (as can be seen below). However, GRASS simply won't import the dbf file. Is there something wrong with the file perhaps? I know for a fact that the file is there, so the problem must be another.

  (Sun May 27 09:40:58 2012)
  db.in.ogr --verbose dsn=C:\datasource\clc2000\CLC2000LEGEND.DBF
  output=landcover_legend GRASS_INFO_ERROR(4804,1): Unable to open data
  source  ERROR: Input table
   not found or not readable (Sun
  May 27 09:41:00 2012) Command finished (1 sec)


texts: What is a good way to share a GeoDatabase over the internet?

texts:
I've got a large, greater than 10 GB, File GeoDatabase on a remote server. I created this GeoDatabase some time ago and am looking to have local access to it. Now I can't simply burn it to a disk or USB and bring it here for logistic reasons. I need to get it from a remote computer to my local computer via the internet. I can only access this computer via Remote Desktop.

Some things I've already tried:

Zipping the File GeoDatabase and then downloading it. This failed for some reason. When I unzip the GeoDatabase afterwards it has some sort of corruption and doesn't work.
Exported features to Shapefile's. This does work, but it's slow, I lose some metadata and in the end, this is not how I want to accomplish this task.
I tried using some file sharing sites, but either they don't have enough space (Dropbox) or they don't allow large uploads (SkyDrive).

I've got full control over ArcServer on the remote PC and I was really hoping that there was a way to host the GeoDatabase on it and download it that way. I was also thinking of exporting it to an XML, but I've had issues with that in the past and because the GeoDatabase is so big and these processes take so long I want a solution that will work.


texts: SSL Connect Error using HTTP Caller in FME

texts: I am using a HTTP Caller in FME Desktop to connect to a SAP OData service. This is a https service. It works in Google Chrome and it works in FME when I use Fiddler, but it doesn't work without Fiddler. I receive a SSL Connect Error. What should I do next?


texts: gdal_translate called via subprocess.check_output works when called directly, but not from within arcgis

texts:
EDIT: Solved it myself. Beginners error... :/

I am writing a python Toolbox that should use gdal_translate to convert an ASCII Gridded XYZ-File into a lower resolution (and doing afterwards a couple of other things).

For that I have written the following snippet:

filepath=r'X:\some_path\to\file\with_only_ascii_and_no_whitespaces\file.xyz'
try:
    inf=filepath
    outf=filepath + '__dgm20.xyz'
    result = subprocess.check_output(
        [
            'gdal_translate',
            '-of',
            'XYZ',
            '-tr',
            '20',
            '20',
            inf.encode(locale.getpreferredencoding()),
            outf.encode(locale.getpreferredencoding())
        ],
        stderr=subprocess.STDOUT
    )
except subprocess.CalledProcessError as err:
    arcpy.AddWarning(unicode(err))
    arcpy.AddWarning(err.cmd)
    arcpy.AddError(err.output)
else:
    arcpy.AddMessage(filepath + '__dgm20.xyz')
    arcpy.AddMessage(result)

When I run this from the python console in cmd.exe or the ArcGIS python console, everything is fine, my file will be created. But when I run it from within a python toolbox (after retrieving filepath from a parameter and normalizing it with os.path) I get this error:

  ERROR 1: At line 3999971, found 1 tokens. Expected 3 at least

The processed file (filepath) ist totally fine, I think.

So, why is this happening? What is it I can't figure out?


texts: Export ImageCollection to an Asset - Google Earth Engine

texts:
I have the same difficulty in exporting an ImageCollection to an asset in Google Earth Engine.
This align with the previous question on this forum.
Exporting image collection to asset


texts: ArcGIS 10.1 - OD Cost Matrix SOAP Request (Error received)

texts:
I've been using Network Analyst in ArcGIS desktop (10.1), and using succesfully some SOAP API operations (findAddressCandidates, GetNALAyerNames, GetSolverParameters, GetSolverParameters2)

But now want to get OD Cost Matrix (Network Analyst) using SOAP API.

I've been following the steps:
http://resources.arcgis.com/en/help/soap/10.1/#/Overview/01vp0000001p000000/

Basically, I call getSolverParameters2 to get parameters and content of NAClasses for the OD Cost Matrix Layer.

Debugging the code, I can see data returned by getSolverParameters2.

When I call solve operation (defined in WSDL)

I get and HTTP 400 Error (Bad request) and the following message:

"The request sent by the client was syntactically incorrect ()"

NAServerODCostMatrixParams params = (NAServerODCostMatrixParams) port.getSolverParameters2("OD Cost Matrix", "", true);
NAServerSolverResults results = port.solve(params);

(Also, I've been trying fill NAServerODCostMatrix 'by hand', unsuccesfuly )

I trace the SOAP request message and I can see the content. The same results I get when I follow this steps using SOAP UI.

Do you have any idea about this problem and/or to solve it? or what can I try?


texts: Geoserver's WMS GetCapabilities request takes too long to respond

texts:
I have about 300 vector and 30 raster layers in a GeoServer instance (version 2.7.4) and when the client makes a WMS GetCapabilites request, it takes about 3 minutes to receive the response.
The vector layers are stored in a PostGIS database while raster layers are stored as GeoTIFF files.

In the field: Limited SRS list of WMS admin page, I have imported the SRSs that I'm using like:

4326, 3857, 2100

and the result of list in WMS GetCapabilities file is:

<!--Limited list of EPSG projections:-->
  <SRS>EPSG:4326</SRS>
  <SRS>EPSG:3857</SRS>
  <SRS>EPSG:2100</SRS>

Nevertheless, the request takes too long to response. I found a similar question:
Limited CRS options of GeoServer WMS

In the answer, @lubojr says the following: "I just missed the small drop-down list with workspace selection on top of the WMS settings page"

I also selected the corresponding workspace in that dropdown menu but the result was the same.

Any suggestions ?


texts: Remove overlaps by moving overlapping features with a specific, scale-based delta

texts:
I need to use QGIS to move multiple polylines for tracking data.

I want to move any overlapping feature to the left or right with a specific delta based on a scale.

I can do that with ArcMap using parallel line copy, and I have searched how to do that with QGIS or Py-Qgis.

Do you know which lib or tool I must use?

Actually I don't have choice and I make processing on data. But a display solution can be a good idea because there is no data manipulation.


texts: I am looking for a way to create map tiles from OSM data but I would like the resulting map tile images rendered at a 45 degree angle. Is there any way this can be done? I have never see such an example.


texts: Why is the "Best Fit" option is grayed out for Data Driven Pages?

texts:
I'm using Data Driven Pages to display delivery routes associated with hub points. On the extent tab, the "Best Fit" option is grayed out. What could be the cause of this? I think "Best Fit" is exactly what i'm looking for.

ArcMap 10.2 Desktop Advanced.


texts: Applying Symbology To MXD Layers From Symbology Files In ArcGIS Stand Alone ArcPy?

texts:
I have a saved .mxd map in directory called Project.mxd and some layer symbology files .lyrs inside a folder called lyr

import arcpy
mxd = arcpy.mapping.MapDocument(r"C:\Data\Maps\Project.mxd")
lyrPath = "C:\Data\lyrs\"

How can I map the symbology of Layers inside the map to the layer symbology files .lyr on an stand alone ArcPy Script?

For Example if I have a layer called Population in the Project.mxd how can I apply its symbology from a file called pop-2016.lyr and so on?

I already tried to use this snippet

import arcpy
mxd = arcpy.mapping.MapDocument(r"C:\Data\Maps\Project.mxd")
lyrPath = "C:\Data\lyrs\"
df = arcpy.mapping.ListDataFrames(mxd, '')[0]
for lyr in arcpy.mapping.ListLayers(mxd, '', df): # Loop through layers
 if lyr.name == "Population":
    arcpy.ApplySymbologyFromLayer_management("Population", "C:\Data\lyrs\pop-2016.lyr")
mxd.save()
del mxd

but I am getting this error

  ExecuteError: Failed to execute. Parameters are not valid. ERROR
  000732: Input Layer: Dataset Population does not exist or is not supported
  Failed to execute (ApplySymbologyFromLayer).


texts:
On a mac (Mountain Lion), after a major upgrade to GDAL (1.10) and related frameworks from KyngChaos page, QGIS 2.0 doesn't show any MySQL connector:

  Add Vector > Database > ...

The driver doesn't appear neither in the terminal running the "ogr2ogr --formats" command.
I read (link), at the very end of the post, that this feature has been removed from QGIS 2.0. Anyway, it appears to me that it depends on GDAL/OGR: it seems that the framework has been build without MySQL support.


texts: How to show the colorbar from the layer properties of a raster? I would like to compose a PNG image with both raster and a colorbar.


texts: Get the outer ring inner ring for a given polygon using ogr python binding

texts:
I am trying to read the geometry of some shapefiles, now I meet some problem when I try to identify the polygon with holes, since a polygon may(or not) have one hole(or multiple holes).

This is the code I am using:

def _get_features(coord,datalayer):
    bbox=_get_bbox(coord)
    datalayer.SetSpatialFilter(bbox)

    layer_def=datalayer.GetLayerDefn()
    geom_type=layer_def.GetGeomType()

    for feature in datalayer:
        geometry = feature.geometry().Clone()
        if not geometry.Intersect(bbox):
            continue
        geometry = geometry.Intersection(bbox)
        fixed_geom=_get_coordinate(geometry)

        #maybe a polygon with hole
        if isinstance(fixed_geom[0], list):
            if geom_type is wktPolygon:
                # now, how to get the rings? outer, inner, multiple inner?

def _get_coordinate(geometry,result=None):
    if result is None:
        result=[]

    if geometry.GetPointCount()>0:
        for i in range(0,geometry.GetPointCount()):
            point=geometry.GetPoint(i)[:2]
            result.extend(map(_fix_coordinate,point))
    elif geometry.GetGeometryCount()>0:
        for i in range(0,geometry.GetGeometryCount()):
            geom=geometry.GetGeometryRef(i)
            _get_coordinate(geom,result)

    return result

Any suggestion?

BTW, is it possible a hole inside a hole of a polygon?


texts: Calculating the distance of raster cells from a vector within a polygon using QGIS

texts:
I am quite new to QGIS and  am trying to use QGIS 3.6 to create a raster layer of calculated values of the minimum distance of each cells within a polygon (orange) from a vector (the river channel in blue).

However, I need to limit this distance calculation to just the area within the fully displayed orange polygon and don't want the cells within the polygon calculating their minimum distance to vectors (channels) outside the polygon. I have nearly 200 of these polygons each of which represents the catchment area of a specific point and will need to replicate this raster calculation for all of these polygons. Could anyone help shine a light?


texts: Verifying time-invariance for spatial gradients of air pollution?

texts:
I have data from a network of NO2 passive samplers. Measurement are repeated in the same positions in different years.
I need to verify if spatial gradients in air pollution remain constant in time, independently from the absolute NO2 measured concentrations (e.g. in a certain area pollution is always higher than in other areas).
Which statistical/geostatistical analysis should I use?

PS. I am wandering if a possible solution could be to normalize each sampled NO2 concentration by the average concentration of all samplers in each year, then do some kriging and compare in some way different maps...


texts: Find point in Polygon/MultiLineStrings with multiple interior areas, using Shapely and JSON

texts:
I've been working on this for a few days, adding in features, and fixing the problems as they crop up.

The goal is to read in a CSV containing addresses for locates, parse/normalize them (that functionality will come later), geocode them, and compare them to our service territory, which I currently have as both a KML and GeoJSON (I used Ogre to conver the KML). I've uploaded the JSON to PasteBin. I've made a very simplified program that only reads in a point and compares it to the JSON for now.

Our service area is bounded on the outside, and within it, there are several other areas that we do not service, with the exception of two - those two I've added placemarkers, which I assume I could check for. My first thought was to use Shapely's contains function, but after trying, and reading the documentation, I don't know if it would deal with the MultiLineStrings making up my boundaries.

The only way I get a positive hit out of this program is if I put in one of the two Points I made. I've tried inside the main boundary, inside a smaller bounding box (non-serviced area) inside that, and inside the boxes that we do service.

import json, shapely, time

from shapely.geometry import Point, shape#, MultiLineString

jsonFile = 'P:/Territory.json'

point = Point(-86.76642, 40.67597)

with open(jsonFile, 'r') as jsonFile:
    js = json.load(jsonFile)

    for feature in js['features']:
        #polygon = MultiLineString(feature['geometry'])
        polygon = shape(feature['geometry'])
        if polygon.contains(point):
            print("In our territory." + str(feature))
            time.sleep(2)
        #else:
            #print("Not in our territory.")

The MultiLineString was throwing a KeyError: 0 (I can include the entire traceback if needed), so I commented it out and tried polygon. The else is commented out because as written, it checks every feature for a match, and generates a huge list of Nopes. I can fix that later; right now I just want to get the point in polygon working correctly.

It's also worth mentioning that I'm happy to use any Python library and file format to accomplish this. I also have an Arc license if arcpy has anything of use here.


texts: How to Intersect between a vector (administrative areas) and a raster layer (slope)?

texts:
What I've got:

1: A vector-layer with administrative areas of my country

2: A raster-layer with the slope of the terrain of the same country (it is derived from SRTM elevation data with QGIS)

What I want:

A new vector-layer with the 'flat' land of all administrative areas.

How do I achieve this in QGIS?

A) How to extract the area below a certain slope (in the raster-layer)? - The option for rule-based display is not available (maybe this is just there for vector-layers).

B) How to 'intersect' the flat-land-raster-layer(-mask) with my admin-vector-layer?

Update:

A) Actually the "Clip to MinMax" with an appropriate set Min- and Max-value (e.g. Min: 10, Max: 100 - or inverse Min: 0, Max: 10.1) does the job. Now I want to use this as a 'clipping-mask' for my vector-layer, but how?


texts: Differentiate style used to render vertices of lines and polygons from style used to render plain point features?

texts:
Context
I'm using OpenLayers 6 to let a user make some (multiple geometries) drawings (i.e. vector features) on map.
I build a function to render the different geometries.
I also wanted to display the vertices of lines and polygons.
For this latter point, I found this really useful, so I based my function on that, as follows:

setDrawingStyle = function() {
  var imageForPoints = new ol.style.Circle({
    radius: 5,
    fill: new ol.style.Fill({
      color: 'green',
    }),
    stroke: new ol.style.Stroke({
      color: 'green',
      width: 2
    })
  });
  // Insert here the new variable 'imageForVertices' defined hereunder
  return [
    new ol.style.Style({
        image: imageForPoints, // Insert here the new function defined hereunder
        geometry: function(feature) {
          var geomType = feature.getGeometry().getType();
          var geomCoords = feature.getGeometry().getCoordinates();
          // The if statement I'm not totally satisfied with:
          if (geomType == "MultiPoint") {
            var coordinates = geomCoords;
          } else if (geomType == "MultiLineString") {
            var coordinates = geomCoords[0];
          } else if (geomType == "MultiPolygon") {
            var coordinates = geomCoords[0][0];
          }
          return new ol.geom.MultiPoint(coordinates);
        }
    }),
    new ol.style.Style({
      stroke: new ol.style.Stroke({
        color: 'red',
        width: 3
      }),
      fill: new ol.style.Fill({
        color: 'red'
      })
    })
  ];
}

I'm not really happy with the "if" statement that checks for the different types of geometry here, but it's working (for simple geometries).
Notice, I had to write the very first part of this if statement which checks for Point geometry type, otherwise a TypeError: e is undefined is raised.
So now every point, whether it is a vertex or true point, is rendered identically.
But I would like to differentiate the image value for true points from the one used for vertices composing lines and polygons.
What I've tried so far
In order to achieve this, I build a new variable imageForVertices under the existing imageForPoints which is defined as follows:
var imageForVertices = new ol.style.Circle({
  radius: 2,
  fill: new ol.style.Fill({
    color: 'red',
  }),
  stroke: new ol.style.Stroke({
    color: 'red',
    width: 2
  })
});

And I replaced the line
image: imageForPoints,

by this:
image: function(feature) {
  var geomType = feature.getGeometry().getType();
  if (geomType == "MultiPoint") {
    return imageForPoints
  } else {
    return imageForVertices
  }
},

But it throws me an error and I'm going round in circles right now:
TypeError: s.getImageState is not a function

Question
And ideas on how to achieve this?


texts: Falsecolor raster from multiple rasters in QGIS

texts:
How can I combine multiple single band rasters into a falsecolor RGB background in QGIS?
I can, of course, use the multiband color renderer to select specific bands, call them A1,A2,A3 from a single raster file source A to R,G,B channels. But I have the desired bands in different files, so it's R=A1, G=B2, B=C1 that I need.

For the moment, I am creating and maintaining a monster .vrt with A1...An,B1...Bn,C1...Cn munged into one, and then using the multiband color renderer. This requires annoying editing of an XML file to add <Description> tags and to overcome to inability to add different -b switches with -separate in the gdalbuiltvrt call the QGIS Build virtual raster... algorithm uses. I'm hoping for a way where I could assign the right rendering settings to each layer in the QGIS layer tree so that they would individually render as in the 3 rightmost images above, and when all 3 are checkmarked in the layer tree they would add as at left.
I've tried using the Singleband grey renderer for each input layer and colorizing it in Layer rendering settings. But that keeps whites white and only shifts the midgreys. So I created my own color maps from #00000 to #ff0000 (etc.). But now the topmost layer wins, unless I adjust transparency in which case the color gets washed out. So I've been playing with the blending mode (e.g. Darken and Screen as opposed to Normal) but can't seem to get the logic right.
This ought to be simple, so I'm probably barking up the wrong tree!
BTW, the reason for this is exploration of different spectral bands and ratios across a time sequence of Sentinel 2 images to discover what combination of data will best pinpoint forest insect damage.


texts: Joining spatial and not spatial table in oracle

texts:
I have spatial table with rows of data representing polygon(A). I have another table that contains values for polygon field(B). Now I would like to join these two table and create a spatial table. I selected spatial column GEOLOC and data from another table. I want to make view out of this which would be mappable in mapinfo. But

Making table mappable:

While making table mappable:

I tried this with MS-SQL geometry colum and it works but fails with oracle spatial.


texts: FME - Conditionally make parameters visible

texts:
I have an FME workspace in which the user sets several parameters before running.  The first parameter is a multiple choice where the user selects which layers they would like to process.
Depending on which options they choose, different WorkspaceRunner transformers are run.

There are additional parameters that feed into the other Workspace Runners.  However, I don't want them to be visible unless the relevant options are chosen in the first parameter.

Essentially, if someone selects "03 ROW Width" in the first parameter, then I want the associated parameter below ([03] Easement Width (L : R)) to appear.  Likewise, if the user selects "08 Minimum Cover" in the first parameter, then I want the associated parameter below ([08 Depth of Cover (mm)).  If they deselect those options, the associated parameters turn off.

Is this possible?  Does FME support validation in this way?

I'm open to a python solution if this can only be done in a startup script.


texts: Creating multi-value choice list of field values and using selected for select in ArcPy?

texts:
I have created a multi-value choice list following this blogpost:

https://blogs.esri.com/esri/arcgis/2012/07/16/generating-a-multivalue-choice-list/

I changed the validation code to my needs but now I am stuck.
Basically I want the script to use the selected values of the choice list and 1) change them all to a another field value (optional)
2) apply a select analysis on it.

My validation code so far looks like this:

import arcpy
class ToolValidator(object):
"""Class for validating a tool's parameter values and controlling
the behavior of the tool's dialog."""

def __init__(self):
"""Setup arcpy and the list of tool parameters."""
    self.params = arcpy.GetParameterInfo()

def initializeParameters(self):
"""Refine the properties of a tool's parameters.  This method is
called when the tool is opened."""
    return

def updateParameters(self):
"""Modify the values and properties of parameters before internal
validation is performed.  This method is called whenever a parameter
has been changed."""
    if self.params[0].value and self.params[1].value:
        self.params[2].filter.list = sorted({row[0] for row in arcpy.da.SearchCursor(self.params[0].value, self.params[1].value.value) if row[0]})

def updateMessages(self):
"""Modify the messages created by internal validation for each tool
parameter.  This method is called after internal validation."""
    return

My user interface looks like this:

How can I do this?

Update:
The part of the script, where I don't know what to do:

inputFC= arcpy.GetParameterAsText(0) #input Feature Class, type: feature layer
inputField =arcpy.GetParameterAsText(1) #input Field, type: Field
Values = arcpy.GetParameterAsText(2) # field values, type: String
OutputFC = arcpy.GetParameterAsText(3) #output, type: feature layer

where = "\"inputField"\" = Values
arcpy.Select_analysis(inputFC, OutputFC, where)

The problem lies in the parameters of the select_analysis part. I don't know how to access the values the user selected to use them in an SQL statement.


texts: NearestNeighbour with NetTopologySuite with wgs84

texts:
I've got a list of polygons in wgs84 WKB, I read them into NetTopologySuite Igeometries put them into an STRtree and then use the NearestNeighbour function to get the closest polygon to a point (also wgs84).

This doesn't work properly as long lat isn't flat.
I know I need to re-project. But should I re-project my polygons when reading them from the WKB (and how do I do that?) or should I reproject every point in the IItemDistance?

A somewhat simplifed version of my code.

private STRtree<IGeometry> LoadSTRtree(List<byte[]> wkbs)
{
    var sourceLayer = new STRtree<ShapeDetail>();

    foreach (var wkb in wkbs)
    {
        var geometry = WKBReader(NtsGeometryServices.Instance).Read(wkb);
        sourceLayer.Insert(geometry.EnvelopeInternal, geometry);
    }
    sourceLayer.Build();
    return sourceLayer;
}

private IGeometry GetClosest(IGeometry location, STRtree<IGeometry> sources)
{
    return source.NearestNeighbour(location.EnvelopeInternal, location, new MyItemDistance());
}

private class ShapeDetailItemDistance : IItemDistance<Envelope, IGeometry>
{
    public double Distance(IBoundable<Envelope, IGeometry> item1, IBoundable<Envelope, IGeometry> item2)
    {
        return item1.Item.Distance(item2.Item);
    }
}


texts: Using each feature of vector layer as clip mask over separate point layer, storing results in either a Python list or in-memory scratch layer

texts:
I have a point layer (representing cities), and a set of vector layers whose features are hundreds of concentric discs, each of which I'd like to use as a clipping layer to select only the cities that lie within that disc, perform some operations on those cities, and then move to the next disc.

Here are some screenshots of the project. The concentric blue-grey discs are the clipping polygons - 25/50/100/150 mile buffers around military base locations(represented as green diamonds). The red dots are the city features (displayed in graduated sizes by violent crime rate) that I'm trying to use as the input layer :

I was trying to use the QGIS processing module to do this. However all of the documentation I've been able to find on the qgis:clip function shows its output going to a shapefile. For instance the answers to the question Is there a way to call the clip function in pyQGIS from the python console? suggest:

processing.runalg("qgis:clip",inputlayer,overlaylayer,"output_file.shp")
processing.runandload("qgis:clip",inputlayer,overlaylayer,"output_file.shp")

What I want to do instead is iterate through all of the discs in the clipping layer(s), and save the results from using each into a temporary in-memory vector layer instead of a shapefile.

Here's the (non-working) code that I tried to write to do this:

for disc in troopBuffers:
    # Clipping requires the cities layer and troop buffer disc to both be selected
    buffLayer.select(disc.id())
    cityLayer.selectAll()

    clipLayer = QgsVectorLayer("Point", "temporary_points", "memory:tmpClipLayer")
    QgsMapLayerRegistry.instance().addMapLayer(clipLayer)
    processing.runandload("qgis:clip", cityLayer, disc, clipLayer)

    clippedCities = [f for f in clipLayer.getFeatures()]

    for c in clippedCities:
        # do something with each city that lies within this disc ...

    buffLayer.removeSelection()
    QgsMapLayerRegistry.instance().removeMapLayer(clipLayer.id())

... however, when I run this clippedCities appears to always be empty, so something is obviously wrong.

How would I clip like this without actually creating shapefiles on disk for each iteration?

Can I just have processing output the results of each clip to an in-memory "scratch layer", or (even better) just get a list of QgsFeature objects from the cities layer for each disc?

Or am I going to have to use something else besides processing here?


texts: Does QGIS preserve Feature Class Relationships (1to1 or 1to many or Many to many relationships)?

texts: Does QGIS preserve Feature Class Relationships (1to1 or 1to many or Many to many relationships) and can the related attributes be clicked and viewed through the identify feature tool?


texts:
I have three geodatabases (176 files each) with raster data. let's call the first one as G, the second one as D and the third one as S. Now I want to perform a raster calculation so that the G1-D1-S1.... then G2-D2-S2 and so on

I think I would need two submodels (one submodel that will iterate through D files, and the other one which will iterate through S files) and one main model which will iterate through G files and then perform the calculations and then save the results in a new geodatabase .

The problem is I don't know how to do that.


texts:
I'm pretty new to using the Overpass API but we are looking to calculate areas & ways (of certain types) which either contain a particular coordinate or are close to a specified coordinate.

So far we have this, which 'seems' to work:

[out:json];
(
way(around:25, 28.416830, -81.581877)(if: (is_tag("building")) || (is_tag("leisure") && t["leisure"] == "park") || (is_tag("tourism")) || (is_tag("place")));

  is_in(28.416830, -81.581877)->.a;
area.a(if: (is_tag("building")) || (is_tag("leisure") && t["leisure"] == "park") || (is_tag("tourism")) ||
 (is_tag("place")))->.building;
way(pivot.building);
);

out geom;

But we intend to be making our requests from mobile devices so want to try and remove as much from the resultant json as possible e.g. don't include non english name tags.  Similarly we aren't interested in the 'nodes' but the geometry is required.

Is this possible or would we need to do some of our own server based work to remove those elements we don't require?


texts: Deleting required field from ArcGIS Diagrammer?

texts:
while creating a feature class in ArcGIS Diagrammer I set the field id_feature_class as required and as type= text now I want to change text to integer but I can't it's telling me id_feature_class is required

how can I delete a required field in ArcGIS or change the field type from text to integer?


texts: Fitting .svg on layer with georeferencer using QGIS

texts:
I am using QGIS 3.4.

I'm trying to fit an image (as an exercise) using the georeferencer, but I don't understand the result.

The srs I choose seems to be correct, and my coordinates (on several cities) too, but I got huge errors (like Warsaw going in south of Moscow ?) and really don't get it.

Here is the image :

Can you help me to master this tool?


texts:
I am applying RPC's via a vrt file and observing a weird output shown here:

.

The gdalwarp command I am using is:

$ gdalwarp -overwrite -dstalpha -rpc test.vrt test_rpc.tif

Here are the links to VRT file
and the Test file.

The top of the image should be shifted to the left and bottom part of image should be shifted to the right.  I am running gdal version 2.1.3.


texts: QGIS, exporting geometry columns doesn't give decimal degree coordinates?

texts: I am working on a project in qgis with my distance and coordinates units set to degrees and decimal degrees. However the shapefiles I work with have meters listed in their 'units' section (under information option). Therefore when i use 'export geometry columns' I always get the xy coordinates in meters. However I need these in decimal degrees. Is there anyway to get them in the units I need? I can't change the units in each shapefile (I might now have found a way to do that).


texts:
I am having trouble getting geotiff rasters to get correctly placed even with a world file available.  I have tracked down some discussion from 2010 that suggests that if the geotiff already has the registration in the header, QGIS will ignore the world file.  It looks like the problem is that the geotiff header states the location is Unknown.  Since it is stated, QGIS then ignores the world file with the correct information in it.  I thought the purpose of a world file was to that priority over the geotiff header to correct mistakes like this.  If this is the case, shouldn't the code be changed?

So in the meantime, does anyone have a work-around?


texts: Determining distribution followed by Create Random Points tool of ArcGIS Desktop?

texts:
I have created a random sample in ArcGIS using the Create Random Points tool. I know that the tool generates a pseudo-random number generator and in my case is ACM599 with seed number 0(default).

How can I know which distribution it follows?


texts: Specifying a rectangle as two latitude/longitude pairs, one for center and one for size

texts:
I am trying to make sense of the Google S2 library for spherical math (it works like this). I need to get a rectangle for searching for nearby entities in a database, based on a point and a distance (say, 10km, the distances do not have to be exact, I filter the entities later).

In the S2 library there is a function with the signature:

RectFromCenterSize(center, size LatLng)

Both parameters are latitude/longitude pairs.

I am having difficulty understanding what this means (the documentation is unfortunately very poor). The center point I get, but how does the a lat/lng pair specify the size?

Any help would be greatly appreciated.


texts: How to filter dropdown options in QGIS attribute table?

texts:
I have the following layer table
Prim  | Secondary      |Tertiary

Apple | Apple Pen      | Apple Pen White
Apple | Apple Pen      | Apple Pen Red
Apple | Apple Phone    | Apple Phone Black
Apple | Apple Phone    | Apple Phone Silver

Banana| Banana Bread   | Banana Bread Chocolate
Banana| Banana Smoothie| Banana Smoothie Sugar
Banana| Banana Yoghurt | Banana Yoghurt Smooth
Banana| Banana Yoghurt | Banana Yoghurt Fruit

Car   | Car Toyota     | Car Toyota Yaris
Car   | Car Toyota     | Car Toyota Corolla
Car   | Car Ford       | Car Ford Figo

I was able to use value relations and a value map to create a drop-down for 2 levels.
I used a value map for the Primary

I used a value relation widget for the secondary.
I used the expression : "Primary Code"= current_value('Primary')

However, I am struggling to get the third level.
It should be able to limit the choices you can choose when editing the attribute table.
eg.
If I choose Apple in the Primary column then in the second column, I should only be able to select Pen or Phone, if I chose Pen then in the Tertiary column, I should only be able to select White and red  and Black or silver if I chose Phone.


texts: How to put marker on vector layer with lat/lon on QGIS

texts: I have a vector layer, and bunch of lat/lon values. I want to put markers on the vector layer for each lat/lon pair. I have never used QGIS, would you help with it.


texts: Creating a map with selectable layers in Cartodb

texts: I tried to make a map with selectable layers using the tutorial at http://developers.cartodb.com/tutorials/toggle_map_view.html, but was not successful. I really like this map by the LA Times http://graphics.latimes.com/responsivemap-pollution-burdens/ with selectable layers. I tried working off from this http://bl.ocks.org/andrewxhill/10506396, but when I load my map it appears with just the buttons and a black basemap (nothing displayed). My attempt is here https://dl.dropboxusercontent.com/u/253185126/index.html. Any tips or solutions would be helpful.


texts: r stars: calculate ndvi for every time step on stars proxy

texts:
I would like to calculate the ndvi on every time step of a 4d data cube (x, y, bands, time) in on an stars proxy object.
Do you have any ideas of how to achieve this on the proxy object?
My current approach throws the error:

  Error in ds[-MARGIN] : only 0's may be mixed with negative subscripts

Here is the code that produces this error:

tif = system.file("tif/L7_ETMs.tif", package = "stars")
t1 = as.Date("2018-07-31")
x_prox = read_stars(c(tif, tif, tif, tif),
               along = "time",
               RasterIO = list(nXOff = c(1),
                               nYOff = c(1),
                               nXSize = 50,
                               nYSize = 50,
                               bands = c(1:6)),
               proxy = TRUE)
x_prox = stars::st_set_dimensions(.x = x_prox, which = "time",
                                  values = c(t1, t1+1, t1+2, t1+3))
st_dimensions(x_prox)
st_get_dimension_values(x_prox, "time")

calc_ndvi = function(x) (x[2] - x[1])/(x[2] + x[1])
x_prox = st_apply(X = x_prox, MARGIN = c("x", "y", "time"), FUN = calc_ndvi)
x_prox_ndvi = st_as_stars(x_prox)

Changing proxy = FALSE returns correct values.


texts: Running a Python script to extract raster data in QGIS

texts:
I am trying to use QGIS to extract raster statistics for shapefiles (I have shapefiles of species distribution and I want to extract environmental data from within each species' range). This is possible using the Zonal Statistics plugin however I have about 300 shapefiles I need to get data from. I've tried merging the shapefiles into one but this results in massive files that cause QGIS to crash (using shapefiles >15mB causes problems).

Is there a way to automate the procedure using a Python script? I found a post that details how to run ZonalStats using Python (How to calculate raster statistics for polygons? second answer) but can't get it to work. I know absolutely nothing about Python and so I'm not sure if there's something that needs to be done before attempting to run the script or if the problem is something else. A similar script that will run on its own without user input would be ideal.

I'm running QGIS 1.8 on a Linux machine.

Any help would be much appreciated.


texts: Creating model that would take statewide parcel layer and split it into all individual parcels by counties using ArcGIS Desktop?

texts:
Using ArcGIS Desktop 10.5, I want to split the statewide shape into counties to speed up the process.

Is it possible to make a model that creates individual shapefiles for each county with ModelBuilder?  " Sorry I was not specific. I am working with data provided by a client that often has parcels that overun the typical public shapefiles " Edit


texts: Leaflet .js reference layers as variables

texts:
I am advised to reference my layers as variables rather than as functions, so I can order the layers using .bringtofront.

I have tried several combinations of var skull = L.geoJSON(fatalities)etc.)  but everything I have tried simply causes all the layers to vanish from the map.

Here's the map: https://pdxcyclesafetymap.neocities.org/

And here's the .js for one of the layers, as it works now:

var skull = L.icon({
    iconUrl: 'images/blackSkull.svg',
    iconSize: [38, 95],
    popupAnchor: [-10, -30],
    });

function fatalitiesLayer (feature,
layer) {
    layer.bindPopup("
Cyclist Fatality

 "+feature.properties.name+
    "
"+feature.properties.description+
    "
Date: "+feature.properties.date+
    "
Location: "+feature.properties.location+
    "
Cause: "+feature.properties.cause+
    "
");
    layer.setIcon(skull);
    }

L.geoJson(fatalities,{
    onEachFeature: fatalitiesLayer
    }).addTo(map);


texts: GeoTools CRS.getAxisOrder(CRS.getAxisOrder("EPSG:4326")) returns EAST_NORTH (LON_LAT)

texts:
All is in the title, with GeoTools if I do :

AxisOrder order = CRS.getAxisOrder(CRS.getAxisOrder("EPSG:4326"));

It returns EAST_NORTH (= LON_LAT), but EPSG:4326 is defined as ordering latitude first and then longitude, as stated here : http://spatialreference.org/ref/epsg/wgs-84/gml/

So, what is happening here ? Am I doing something wrong ?


texts: Why are Landsat 8 files on AWS compressed cf those from USGS?

texts:
I want to shift from downloading Landsat8 imagery directly from USGS to the AWS service (https://pages.awscloud.com/public-data-sets-landsat.html).

I have some automated processes that are not accepting the TIF files from AWS. Upon closer inspection the file sizes differ substantially.

For example for scene LC80930842016347LGN00 if I download from EarthExplorer as LC80930842016347LGN00.tar.gz and extract all bands except B8 are 116,191KB.

When I download the bands individually from AWS (http://landsat-pds.s3.amazonaws.com/L8/093/084/LC80930842016347LGN00/index.html) all bands except B8 are around 60MB and the BQA band is only about 1MB.

Is there some compression going on in the AWS storage and if so, is there a process to transform these back to original so the Landsat_LDOPE tools will work?


texts: Editing harvested metadata records without data loss during next run in GeoNetwork?

texts: I want to edit harvested metadata records without losing data during the next run in GeoNetwork 2.10.4 because whenever I am updating keywords it gets lost next day.


texts: Select the holes in a vector shapefile in QGIS

texts:
Is there a way to select the holes (gaps) in a vector shapefile in QGIS, like the "fill holes" function which recognizes them indirectly?
I need to get the holes from the original shapefile and work with them.
Is there a way in other software like Arc?


texts: Updating null geometries from attribute using ArcGIS Field Calculator?

texts:
I have a large point feature class and have realised that there are several hundred features with null geometry. Fortunately the coordinates are stored as an attribute. What I need to do is update the geometry of the points with the x and y stores in the attribute table.

After a little search I found an Esri Support Services Blog post on How to Update the Location of a Point Feature and Its XY Fields which includes a python function that will do this for me in the field calculator.

However, it seems that this will only work on features with a geometry. When I run it on the null geometry features it returns a syntax error. presumably this is because there is no shape or there is no shape.getpart(0).

Is there something I could add to this script to allow it to update null geometries, or does anyone have another method?


texts: Methodology for Hansen classification of Global Forest Watch?

texts:
I have just discovered amazing Hansen´s classification about forest loss and re-growth, available at http://www.globalforestwatch.org/, published in Science, 2013 as:
Hansen, M. C., Potapov, P. V, Moore, R., Hancher, M., Turubanova, S. A., & Tyukavina, A. (2013). High-Resolution Global Maps of 21st-Century Forest Cover Change. Science, 342(6160)(15 November), 850–854. doi:DOI:10.1126/science.1244693.

However, I can´t find in this article/on the website the exact methodology how to reproduce such a map so which classification has Hansen used?

The only thing I can find is that the supervised learning algorith was used to identify tree cover, but it is quite a broad term.

If it is possible, I would like to use the same methodology (but apply it on 90th years), so before Hansen´s classification in my selected area.


texts: How can I extract variable information from an OpenLayers 3 object?

texts:
I have tried to make a simple map with OL 3's Geolocation module. I succeeded with a resource wasting method. There must be a simpler way to extract variable information from an OpenLayers 3 object. My code is the following:

var view = new ol.View({
    center: [0, 0],
    zoom: 16
    })

var map = new ol.Map({
    target: 'map',
    layers: [
        new ol.layer.Tile({
            source: new ol.source.OSM()
        })
    ],
    view: view
    });

var geolocation = new ol.Geolocation({
    projection: view.getProjection()
});

geolocation.setTracking(true);

var geolocation2 = new ol.Geolocation({
    projection: view.getProjection()
});

geolocation2.setTracking(true);

var featureSource = new ol.source.Vector({
});

var feature = new ol.Feature();

feature.bindTo('geometry', geolocation, 'position')
    .transform(function() {}, function(coordinates) {
      return coordinates ? new ol.geom.Point(coordinates) : null;
    });

featureSource.addFeature(feature);

var featureLayer = new ol.layer.Vector({
    source: featureSource,
    name: 'Your position'
});

view.bindTo('center', geolocation2, 'position');

map.addLayer(featureLayer);

The code works, but not the way it should operate. It uses two Geolocation instances and because the bindTo method the view periodically resets. I couldn't manage to extract the position variable of the Geolocation object. I've tried Geolocation.getPosition(), Geolocation.get('position') and Geolocation.p.position. Every one of them acted like a void method and returned an undefined value. I use Google Chrome version 35.0.1916.153 on a Debian Wheezy.

To complete the story, if I use only one Geolocation instance and bind it to the feature and the view objects, the point follows the center of the map with every change in the view.


texts: Distance matrix returning incorrect values in QGIS 2.18

texts:
I'm having issues with the distance matrix tool in QGIS 2.18 giving output distances that do not match the distance using the measuring tool. I've checked previous threads and so far can't find an answer to this problem.
I'm working with several thousand points, distributed over many layers, and want to calculate the distance between the points on each layer as a set of matrices. I've been using the linear distance matrix setting (though I have also tried the standard distance matrix). The data processes with no errors reported but when I come to check the results they are too large by a factor of approx 2.275. (I haven't yet checked to see whether this factor is common across the whole domain)
My data is all in CRS WGS84/Pseudo-Mercator [EPSG:3857] because I've been digitizing my data from Google Earth. I've done a crude check that the 'measure line' tool in QGIS is working correctly by calibrating it against an object of known size and trust the results.

Example: distance between two features
Measure line tool: 18.572m - trusted
Distance matrix result: 42.2596m - NOT trusted, definitely too big

From reading previous threads, I know that the distance matrix tool should output the results in the same measurement units as the layer CRS, which for Pseudo-Mercator should be in metres.
If it were a handful of points, then I'd be happy to use the measure line tool, but it becomes unworkable with the size of my data set. I've tried saving and reloading the data sets from a .csv file, which doesn't make any difference. I've also fiddled with the project and layer CRS to see whether that has an effect, but I keep coming back to the same results from the tool.
Any ideas?


texts: Opening ArcGIS Online files without ArcMap gives error 403?

texts:
Although the topic about pkinfo was addressed in the question What is *.pkinfo file? on this occasion it does not matter about "what is a pkinfo file" but rather if there is any solution to the problem of the "error 403" message after following the instructions on the support page from esri https://support.esri.com/en/technical-article/000011109

I tried to open a item.pkinfo file without ArcMap, following the steps detailed in Esri support page, Meaning :

Open  the item.pkinfo file with notepad and copy the URL between  tags

When I go to the specified URL, I get this message:

"You do not have permissions to access this resource or perform this operation. error 403"

Does anybody know how to solve this issue?


texts: GeoServer - page not found http error 404

texts:
I have an html file with the arcgis api for js on it stored in the following location:

C:\Program Files (x86)\GeoServer 2.10.5\webapps\agis_api\js.html

I'm running GeoServer on port 8080 localhost for development and am using the following url to view the page:

http://localhost:8080/geoserver/agis_api/js.html

I'm getting error 404 not found. GeoServer is definitely running. Not sure why the page isn't coming up.


texts: Finding longest line segment within polygon?

texts:
I have been looking around to find a way to calculate the longest line within a polygon.

Determining longest line segment within polygon that passes through its centroid? and Calculating maximum distance within polygon in x-direction (east-west direction) in PostGIS? is what I found on here however these are quite advanced solutions. Also they are pretty old (5 and 6 years).

Is there by now a pre made library or github repository that I can use to do this or is there no other way then to write it myself?

Secondly, if I have to do it myself, could anybody post some kind of pseudo code or plan that I can follow.

I am trying to understand the second link at the moment but i'm having a hard time.


texts: 'Use unfiltered layer' option in QGIS query builder

texts:
When you create a filter, a 'use unfiltered layer' option is displayed in the query builder window. I expected this option to allow the user to toggle between the filtered and unfiltered layer. It doesn't seem to work like that as nothing seems to happen when I select it. Have I missed something?
qGIS 3.14.16-Pi Windows 10 64bits


texts: ArcGIS Python how to get service name from the .msd file?

texts:
How to get the service name from the .msd file on ArcGIS server using Python, the sample below from the Stack Overflow display how to get the layer name from the .msd file.

import zipfile
from xml.etree.cElementTree import iterparse

zz = zipfile.ZipFile('C:\\Temp\\Untitled.msd')

serviceNames = []
for name in zz.namelist():
    if name == "DocumentInfo.xml" or name == "layers/layers.xml":
        pass
    else:
        for _event, elem in iterparse(zz.open(name)):
            if elem.tag == "Name":
                serviceNames.append(elem.text)
                break
zz.close()

Is there any way of retrieving the service name?


texts: Related Tables can be viewed in my webmap but not in Collector

texts:
I have a related table set up in an ArcGIS webmap. I can view the related data (1:M) just fine in the webmap. However, when I try to view the related data in Collector, there is nothing there. Is there a setting I did not set up correctly? Have others faced this problem?

Here are my current settings:

The layer with the related table is editable.
The popup for the related table is configured and enabled.
'Filter Related Types' under Collections Options on the Collector app is off.
The related data was set up using the Create Relationship Class on Desktop and then published.
There are no filters in place nor is visibility limited for the layer.
The Collector app is 'connected' However, it the app should be functional in 'disconnected' mode.

I have an almost identical webmap & Collector app for a separate geographic area working just fine.


texts: Where to find out about QGIS new versions?

texts:
What's the right way to find out about new current QGIS version releases?

I'm on Win64 and use OSGeo4W.
The roadmap gives dates for new releases, confirmed by the main web site header which lists time until next freeze.

However, while the time until next freeze shifted over on schedule on 22 March, the text still says current version is 3.6.0 (3.4.5 LTR). I tried over the weekend, and at the time the OSGeo4W installer still said latest was 3.6.0-1. I reran the installer this morning, and it had silently shifted over to 3.6.1, which I've now installed. I've since found that the release of 3.6.1 is also listed in the GitHub commits list. However, the main website still lists 3.6.0.

All of this makes me think I'm painstakingly sleuthing in back alleyways while notification of release is probably blazing in neon lights if I merely looked in the right place. But what is that right place?

I'm perhaps overeager, but I wasted a day a month ago doing a careful workaround for something that had just been fixed in 3.6.0, I just didn't know it was out since the main website took a while past release date to update and I wasn't really paying attention. And I have a pesky intermittent crash I'm perhaps overoptimistically thinking will miraculously disappear in a new release.


texts: How to Page POST to web page to Print Map with ArcGIS for Server 10.1?

texts:
I am using ArcGIS for Server 10.1. I need pass in a lat/long to a webPage as query parameter and print a map with a TiledMapServiceLayer, DynamicMapServiceLayer, and a 1 mile buffer graphic around the lat/lng.  Any ideas how I can do this?

I tried the ESRI javascript API but problem is that I have to let the map render before I print, otherwise this option would work.  I thought about using the REST API by it self and then mashing the images up together, but I am not sure how I would print the buffer graphic.

I need to be able to call a webpage with a query param generate an image because the current SSRS takes the response for the request (a path to the image) and embeds it in the report.  I know it's not the best way to do this, however, I am just trying to modify something that already exists without reworking everything.


texts: Developing find and replace script tool?

texts:
I'm being tasked for developing a script tool that can find ultimately needs to be able to find any string and replace it with another string in an attribute table. Below is what I have so far, and it seems to be on its way to working except but in a very strict and limited way.

Some ways I want/need to spruce it up include:

Being able to skip over OBJECTID and SHAPE fields so they cannot be changed...I figure this could be achieved by having a parameter that allows several entries, but I'm not sure how that would work. I found some sort of ".split(;)" thing to attach to the targetFields parameter earlier today but it didn't seem to cooperate.
Being able to find and replace partial strings
I would also like the script to be able to handle integers and floats, AND I would like radio buttons in the final tool that would let me select between strings, ints, and floats.

import arcpy

FC = arcpy.GetParameterAsText(0)
oldText = arcpy.GetParameterAsText(1)
replaceText = arcpy.GetParameterAsText(2)
targetFields = arcpy.GetParameterAsText(3)

queryString = '"' + targetFields + '" = ' + "'" + oldText + "'"

with arcpy.da.UpdateCursor(FC, (targetFields,), queryString) as cursor:
    for row in cursor:
        row[0] = replaceText
        cursor.updateRow(row)


texts: Python: how to create a point shapefile from a text file

texts:
I'm writing a Python code to read the points in a polygon shapefile and save them to a point shapefile.

I'm not sure about my approach so if you have a better idea (even if it is totally different from mine) please let me know.

So first I made a text file and stored the points' (x,y) in that .txt file. then I tried to make a point shapefile from the text file but it gave an error.

Here is the code (just the last part):

creat point shape-file from text file
import fileinput
import string
import os
env.overwriteOutput=True
outpath="C:/roadpl"
newfc="newpoint.shp"
arcpy.CreateFeatureclass_management(outpath, newfc, "Point")
infile="C:/roadpl/roadL5.txt"
cursor=arcpy.da.InsertCursor(newfc, ["SHAPE@"])
array=arcpy.Array()
for line in fileinput.input(infile):
    point.X, point.Y = line.split()
    line_array.add(point)
cursor.insertRow([arcpy.Point(array)])
fileinput.close()
del cursor

Here is the error:

Traceback (most recent call last):
  File "C:\roadpl\P_Code_L5", line 49, in <module>
    point.X, point.Y  = line.split()
  File "C:\Program Files (x86)\ArcGIS\Desktop10.1\arcpy\arcpy\arcobjects\_base.py", line 87, in _set
    return setattr(self._arc_object, attr_name, cval(val))
RuntimeError: Point: Input value is not numeric

I made the input file (.txt file) from a polygon shapefile using:

cur = arcpy.SearchCursor("C:/roadpl/L5P.shp", ["SHAPE@"])
f=open("C:/roadpl/roadL5.txt", "w")

for row in cur:
    geom = row.shape
    row = geom.getPart()
    for part in row:
        for point in part:
          print >>f, (point.X, point.Y)
f.close()


texts: Querying OpenStreetMap for tags does not return ways or relations?

texts:
I try to obtain all ways/nodes/relations with a special tag from OpenStreetMap. I use the following code to receive all hotels (in a small town in Austria):

//[out:csv(name, "addr:street", "addr:city")][timeout:25];
area
  [boundary=administrative]
  ["name"="Fügen"];
out body;
// gather results
(
  // query part for: “building=* and name=*”
  node["tourism"="hotel"](area);
  way["building"]["tourism"="hotel"](area);
  relation["building"]["tourism"="hotel"](area);
);
// print results
out body;

You may access the result via Overpass. However, this result gives me only the nodes. Ways like #89729697 also having the tourism = hotel tag are not returned.

How do I also query for the ways and relations having the tag I am looking for?


texts: How do I deselect a single feature in QGIS?

texts:
In QGIS 2.12.0, I have all the features that intersect with another layer selected by using the Select by Location tool. However, I want to deselect some of these features while keeping the rest selected.

I am using a Mac, and when I hold 'command' down it lets me select multiple features. However, I can't find a corresponding deselect tool, except the Deselect Features from All Layers button which doesn't help me.


texts: We are using QField for the front end in data capturing. Currently we are doing a census survey and need to capture fingerprint data, is there a module or plugin available in QGIS and how do we implement this.


texts: Is there a way to populate a spreadsheet of city names with their latitudes and longitudes from Google Maps or some other service?

texts: I have a long list of cities in the United States, and would like to get the coordinates for each city without looking up each one individually. Is there a way to do this?


texts: How to add the Geometry type LineString/ MultiLineString within the function QgsGeometryGeneratorSymbolLayer using PYQGIS

texts:
Within a function of pyqgis I am trying to create a layer according to multiple steps. I want to change the Geometry type option to LineString / MultiLineString in the function so that this will be done automatically. Is there a variable I can add to the function?
shape_sym = QgsGeometryGeneratorSymbolLayer.create({'geometryModifier': 'difference(difference(make_line(start_point($geometry), centroid(offset_curve($geometry, length($geometry)/-10.0)), end_point($geometry)), buffer(start_point($geometry), 10000)), buffer(end_point( $geometry), 10000))', 'outline_color': 'black'})


texts: Converting bulky shapefiles (54 MB) to GeoJSON

texts:
I am trying to convert my shapefile to GeoJSON using the below code snippet:

StringWriter writer = new StringWriter();
String geoJson = null;
FeatureJSON fjson = new FeatureJSON();
fjson.writeFeatureCollection(col, writer);

This works fine for small shapefiles but when I try to convert a file that is of 54 MB the system gets hung.

What should I do differently? Or is there another approach to convert big shapefiles to GeoJSON?


texts: Reading an OGR VRT with geometry from an Oracle OCI data source

texts:
I am having trouble creating/reading OGR VRT files from an Oracle OCI data source.

I am running GDAL 2.2.1 on Windows 7 Professional SP1.

My data source is an Oracle 11g (11.2.0.4.0) 64bit database running ArcSDE (10.2.2 Patch 290200).

I am trying to read spatial and non spatial ArcSDE tables and views. I can read non spatial data without a problem. I can even read the non spatial columns of spatial tables by isolating them using the "SrcSQL" tag. Problems occur only when a spatial column is present.

Below are some simple examples that I have come up with to illustrate/isolate the problem.

This works (non spatial data):

<OGRVRTDataSource>
  <OGRVRTLayer name="TEST_tab">
    <SrcDataSource>OCI:user/pword@(DESCRIPTION=(ADDRESS=(PROTOCOL=tcp)(HOST=host)(PORT=1521))(CONNECT_DATA=(SID=sid)))</SrcDataSource>
    <SrcSQL>SELECT cast('This is text' AS VARCHAR(255)) AS col_text FROM dual</SrcSQL>
  </OGRVRTLayer>
</OGRVRTDataSource>

This doesn't work (spatial data - WKB linestring):

<OGRVRTDataSource>
  <OGRVRTLayer name="TEST_tab">
    <SrcDataSource>OCI:user/pword@(DESCRIPTION=(ADDRESS=(PROTOCOL=tcp)(HOST=host)(PORT=1521))(CONNECT_DATA=(SID=sid)))</SrcDataSource>
    <SrcSQL>SELECT cast('This is text' AS VARCHAR(255)) AS col_text, sde.ST_GeomFromText('LINESTRING(-71.160281 42.258729,-71.160837 42.259113,-71.161144 42.25932)',4269) AS col_geom FROM dual</SrcSQL>
    <GeometryField encoding="WKB" field="COL_GEOM"></GeometryField>
  </OGRVRTLayer>
</OGRVRTDataSource>

Ogrinfo returns the following error - "Unable to identify source field 'COL_GEOM' for geometry."

I have tried specifying "wkbLineString" in a "GeometryType" tag, but the extra tag/specification returns the same error.

After @user30184 @Albert Godfrind for the sdo_geometry cast suggestion, here is a working VRT using the sdo_geometry solution :

<OGRVRTDataSource>
  <OGRVRTLayer name="test">
    <SrcDataSource>OCI:user/pword@(DESCRIPTION=(ADDRESS=(PROTOCOL=tcp)(HOST=host)(PORT=1521))(CONNECT_DATA=(SID=sid)))</SrcDataSource>
    <SrcSQL>SELECT cast('This is text' AS VARCHAR(255)) AS col_text, sdo_geometry(sde.st_asbinary(sde.ST_GeomFromText('LINESTRING(-71.160281 42.258729,-71.160837 42.259113,-71.161144 42.25932)',4269))) AS col_geom FROM dual</SrcSQL>
    <LayerSRS>EPSG:4326</LayerSRS>
  </OGRVRTLayer>
</OGRVRTDataSource>

Note the absence of the "GeometryField" tag. When the tag was included, ogrinfo would return the same "Unable to identify source field 'COL_GEOM' for geometry." error. I am assuming that the error is returned because there is no "SDO_GEOMETRY" encoding option - only "WKT", "WKB" or "PointFromColumns.

The VRT works in QGIS as well!


texts: ArcGIS Desktop 10.5 - adding custom key metadata to raster dataset

texts:
I wonder if there's option to create custom key metadata to raster/image datasets. By default ArcMap gives few key metadata types:
http://desktop.arcgis.com/en/arcmap/latest/manage-data/raster-and-images/raster-dataset-properties.htm#ESRI_SECTION1_D03DC4400F7B4D48A128773A2C75709E

ArcPy gives function SetRasterKeyMetadata() which can only insert default metadata.


texts: When importing shapefile attributes (in a DBF) does the order of columns matter?

texts:
I'm working to import a batch of 250+ shapefiles into PostGIS.  I'm using the shp2pgsql utility to load data.
While doing some basic preparation I've realized that a few are missing 1-2 columns (varies depending on the file).
I'm wondering, during the import process for attributes is it just the column headers that are used to identify the fields or is it the placement in the DBF file?
I have written a few basic C# programs that can add fields to a DBF without a problem.
My question is whether I can add the missing columns at the end of the DBF and if it'll import correctly.
Here is an example:
Template
Column A, Column B, Column C, Column D, Column E, Column F

Missing file (adjusted)
Column A, Column C, Column D, Column F, Column C, Column E

Trying to get a sense of if this will work before embarking on code that rewrites the entire DBF file.  It is a LOT easier to add fields at the end.


texts: I'm trying to use shapeorientation in QGIS, but it doesn't create a new attribute table in my shapefile. When I open the csv file it's empty. Does anyone know if i need to an installation on PYTHON in QGIS with this plugin?


texts: Cannot toggle ''edit'' of shapefile in QuantumGIS 1.8

texts:
My problem is as follows:

I imported a number of shapefiles from a .gpx-file from a gps-device. Each shapefile contains a number of waypoints. For the study I am doing, I need to give all of them unique ID's. I simply selected each shapefile, toggled the editing-switch and changed the name in the attribute table's name-column. Though for some reason, I cannot toggle editing of two of the 20 shapefiles. All of them have the same coordinate system as the project, so I don't understand why I cannot find the solution.

I am still in the rookie-phase of using GIS, so I might have overlooked some obvious solution. I hope someone here can help me.


texts: Classifying line layer based on number of intersection points with polygon layer

texts:
In QGIS 3.16 in Windows 10 I have:

a line layer 'centerline.shp'
a polygon layer 'o31001_AX_Gebaeude_f_12342.shp'

I am trying to classify my line layer based on the number of intersection points with the polygon layer.
I use the Field Calculator of the line layer and wrote this:
CASE
    WHEN aggregate(layer:='o31001_AX_Gebaeude_f_12342', aggregate:='count', expression:="ID", filter:=intersects($geometry, geometry(@parent))) = 0 THEN 5
    WHEN aggregate(layer:='o31001_AX_Gebaeude_f_12342', aggregate:='count', expression:="ID", filter:=intersects($geometry, geometry(@parent))) = 1 THEN 1
    WHEN aggregate(layer:='o31001_AX_Gebaeude_f_12342', aggregate:='count', expression:="ID", filter:=intersects($geometry, geometry(@parent))) > 1 THEN 2
END

I can not get proper results. I think it might be a problem with the:
filter:=intersects($geometry, geometry(@parent))

How can I find if a line has an intersection with the polygon in QGIS? If yes, then how many points do they share.