texts: Explicitly close a OGR result object from a call to ExecuteSQL

texts:
How can I explicitly close/release the result of an ExecuteSQL statement when querying a Spatialite database with OGR via Python?  I have a simple query to return the SRID of a dataset.  This works as expected.  However, subsequent calls to ds.ExecuteSQL fail with the error 'SQL logic error or missing database' unless I iterate all the way through the result rows.

For example:

ogr.UseExceptions()

# This query returns a single row
sql = 'select distinct srid(geometry) from foo'
result = ds.ExecuteSQL(sql)
row = result.next()
epsg = row.GetField(0)

# This call fails
ds.ExecuteSQL('drop table bar')

Completing the iteration avoids the error:

_ = [r for r in result]

This is fine for simple cases, but not very explicit.  I've tried calling result.Dereference(), row=None; result=None but these do not help.  What is it that I'm missing here?

UPDATE
The exception is only raised when I enable ogr.UseExceptions().  Otherwise the error passes silently and the drop table statement has no affect.


texts: ArcGIS 3D Analyst Line Interpolator not working

texts:
I have a stack (8) of DEM layers. The 3D Analyst line interpolation tool is only extracting values from the first layer in the stack. I specify different layers in the 3D Analyst window butit still only extracts the values from the first layer on the ArcMap screen.
Do I have to add and delete layers to get the values for each layer?
Is there some way I can make the drop down window active?


texts: gdal_translate not using correct types when resampling

texts:
I'm using gdal_translate to downsample some data using averaging. The input data uses an Int data type, and I want to store the float values of the pixel averaging.

However, even when I force gdal_translate to use Float32, the values in the output file have been cast to int. The command I'm using is:

gdal_translate -tr 0.0075 0.0075 -r average -ot Float32 -co COMPRESS=LZW -q in.tif out.tif

When I use the same options for gdalwarp, I get the expected result (float values).

gdalwarp -tr 0.0075 0.0075 -r average -ot Float32 -co COMPRESS=LZW -q in.tif out.tif

Am I missing something, or is this a bug in gdal_translate?


texts:
I have a line layer with a rule-based symbology. For every rule I use a marker line.

Can I hide single labels in my layers window? I don´t want two of the same labels. And I don´t want the markerline symbol in my layers window? Could I hide it, too?

Thanks


texts: Is it possible to highlight polygons coming from a WMS?

texts:
Please have a look at this site: https://kso.etjanster.lantmateriet.se/?lang=en# and click the binoculars in the left sidebar and search for e.g. "ESKILSTUNA KARLKULLA 1:2" you will see a multipolygon, one property with two sites, getting highlighted. I thought this was only possible to do with a WFS, but is it also possible to do with a WMS (it looks like they use WMS...?) - and if so, how is it done?

Or do you have any ideas how they have done it? Combining both WFS and WMS on top of each other?

I'm using GeoServer to cascade this remote WMS boundary layer and then display it on a Leaflet map.


texts: Changing isNullable property for field using ArcPy?

texts:
Once in a while I encounter a field in my featureclasses that have their isNullable property set to No.

Is there a simple way using python to correct this?

I came up with this code, but unfortunately it's not doing anything.

print "Modify the RESERVE_ID field to allow NULLS"
desc = arcpy.Describe("OGMA_Final")
fields = desc.fields
for field in fields:
   print str(field.name)
   if field.name == "RESERVE_ID":
      print "Setting property for field"
      field.isNullable = True
print "Finished"


texts: Combining graduated styling with dots/stripes using QGIS?

texts:
I am using QGIS and I am making a heat map, but want to mark the ones that have missing data with spots.

Is it possible to add dots/stripes to my graduated colour styling?

As far as I can see it is only possible to use simple fill.


texts: How to handle Point Data with date ranges

texts: I am running into an issue with the time sliders in Leaflet. It appears that most all of the sliders developed for Leaflet operate on the assumption that each record in your data layer has a specific time/date. Unfortunately, my data (historical by nature) have date ranges instead of a specific date (i.e. 1787 - 1810 as opposed to 8 August 1787). I was envisioning a time slider that could add the point to the map when the slider crosses into the items date range and then removes the point when the slider leaves the date range. Does my question make sense and if it does, how can I best address this issue?


texts: qgis3: how to edit project save paths using pyQGIS?

texts: how to programmatically change save paths from absolute to relative?


texts:
How can I simply draw a line between two lat long coordinates values on a map using ESRI ArcGIS JS API?! This most basic task is taking a gruesome amount of effort.

I am able to draw a line with the below, but it's not reading my latitude and longitude coordinate values as such, instead it's taking the point to some arbitrary destination.

require([
      "esri/Map",
      "esri/views/MapView",
      "esri/Graphic",
      "esri/layers/GraphicsLayer"
    ], function(Map, MapView, Graphic, GraphicsLayer) {

  var map = new Map({
    basemap: "topo-vector"
  });

  var view = new MapView({
    container: "viewDiv",
    map: map,
    center: [-118.80500,34.02700],
    zoom: 13
  });

  var graphicsLayer = new GraphicsLayer();
  map.add(graphicsLayer);

  // Create a point
   var point = {
   type: "point",
   longitude: -118.80657463861,
   latitude: 34.0005930608889
  };

  var simpleMarkerSymbol = {
   type: "simple-marker",
   color: [226, 119, 40],  // orange
   outline: {
     color: [255, 255, 255], // white
     width: 1
   }
  };

  var pointGraphic = new Graphic({
   geometry: point,
   symbol: simpleMarkerSymbol
  });

  graphicsLayer.add(pointGraphic);

  // Create a line geometry
  var simpleLineSymbol = {
     type: "simple-line",
     color: [226, 119, 40], // orange
     width: 2
   };

   var polyline = {
     type: "polyline",
     paths: [
       [35.9411934679851, -80.96135253906438],
       [-118.808878330345, 34.0016642996246]
     ]
   };

   var polylineGraphic = new Graphic({
     geometry: polyline,
     symbol: simpleLineSymbol
   })

   graphicsLayer.add(polylineGraphic);

  var polygonGraphic = new Graphic({
   geometry: polygon,
   symbol: simpleFillSymbol
  });

  graphicsLayer.add(polygonGraphic);

});

This block specifically, I simply just want to feed in two sets of lat and long coordinates values.... i.e. CoordinatesA and CoordinatesB (as the startpoint and endpoint, why is this so difficult)

   var polyline = {
     type: "polyline",
     paths: [
       [35.9411934679851, -80.96135253906438],
       [-118.808878330345, 34.0016642996246]
     ]
   };

CodePen https://codepen.io/camerondevcam/pen/PowjVqq?editors=0010


texts: Analogous functions for Erdas Autosync in QGIS?

texts:
We have two overlapping rasters with some shifting occurring as below image:

and we need to eliminate the shifting in the overlap zone, so the rasters align like in the below image:

In Erdas Imagine we can use the Autosync function to make this correction, are there analogous functions in QGIS?


texts: Name of children of toolbar in QGIS python and how to disable?

texts:
There are name of toolbar you can see in :
Finding name of QGIS toolbar in Python?

But, I want to know the name of children of toolbar. For example, the children of mapNavToolBar are touch zoom and pan; pan map; pan map to selection; zoom in; zoom out; etc. Then if we know the name, how to disable in python?

For example to disable menuBar you can type the script :

iface.mainWindow().menuBar().setVisible(False)

But, for the child?

I want to show not all mapNavToolBar, Layers Panel, and menuBar.


texts: Create different Server Database Connection file for ArcGIS Server

texts:
I would like to register database with ArcGIS Server. I know it is done through ArcGIS Server Manager under Site / Data Stores / Register: Database.

I am also aware there can be different connection files (.sde or .odc) imported for Publisher database connection and Server database connection. I understand this option is primarily meant for access of different databases.

In my case I do need to access same database by both connections but Publisher database connection must be set through external IP (as desktop can only connect through external IP to DB server), while Server database connection through internal IP (ArcGIS Server can connect to DB Server only internally, can't see external IP). I hope till this point I get it right.

Now how do I create Server database connection (.sde file)? Is the only way to get .sde file (with internal IP) to install ArcGIS Desktop on the server machine? Is not that an awful thing to do to be able to create a tiny .sde file? Am I missing something?


texts:
With Grid References, what can the range be for the third and sixth figures (the tenths of the area squares) be? 1-10, 0-10, 1-9 or 0-9?

So if I had two grid squares like this, where X marks the spot:

|-----|-----|
|     X     | 1
|     |     |
|-----|-----|
   1     2

– would X be at GR019017 or GR020017?


texts:
I have seen this type of question posed a few times on this site. The solutions that I have gleamed from these questions has helped me. But I have not been able to get my script to produce a result that I would like.
I have a layer in my QGIS 3.10 map called "warn", the actual file path and name is "C:\temp\warn1.shp". I am writing a Toolbox to perform a series of processes on this "warn"-file, with the intention that this file will be called something different in the future.
Some Examples that I have used:

Getting layer by name in PyQGIS?
Activating layer by its name in PyQGIS
How get layer by name in PyQGIS 3.0?

I set up the parameters for the Toolbox as described in the QGIS Template:
self.addParameter(
            QgsProcessingParameterFeatureSource(
                'INPUT',
                self.tr('Area'),
                [QgsProcessing.TypeVectorPolygon]
            )
        )

I call this input in the main script:
Polygon = parameters['INPUT']

When I try and confirm the name of the input using a print statement, a value is returned, but the name contains what appears to be an identifier:
print (Polygon)

Returns:

warn1_940166d7_22bc_47bf_b76c_d75462802047

When I use the Polygon variable in the definition of the layer:
layer = QgsProject.instance().mapLayersByName(Polygon)[0]

I get an "index out of range" error. Which is expected, as the name does not exist in the map.
How can I return the name of the layer, as it appears in Layer View (Table of Contents)?
In this case I would like "warn" returned.


texts: LineLocate a sequence of points, preserving order

texts:
I have a linestring, which represents a bus route, and a sequence of points in order, representing the coordinates of each stop of that route. It kinda looks like this (pardon my poor drawing skills):

Basically, I want to get the total distance from the start of the route to each stop. I know this can be done by doing something like this:

SELECT stops.order,ST_LineLocatePoint(routes.line,stops.point)*ST_LineLength(routes.line)
FROM routes,stops
WHERE routes.route_id = stops.route_id

But the problem is that in the LineString, those two segments that are really close in the drawing (those that contain 5 and 8) actually overlap, so with that solution, stop 5 may be wrongly located between 8 and 9, or stop 8 between 5 and 6. But as I have the order of the stops, I suppose there must be a way to get the distances without getting those errors.
Maybe it cannot be done with a PostGIS function or query.

I have also thought about using ST_Line_Substring to "cut" the linestring, as I advance through each stop in order, but I still have the problem when 5 is located after 8.


texts: what is the format of Geometry data type of SQLServer 2008

texts:
I want to know, In which format the Geographic data stored in geometry data type of SQLServer 2008.it looks like a long string,like 0xE6100010466DFEOA4089663524....... for one polygon. My second question is
How can i get set of google latitudes/longitudes from this data format ?


texts: Right raster zonal statistics for NDVI at county scale?

texts:
I am currently doing a modeling project at county scale in which I have disease rate as dependent variable and environment variables (such as temperature, NDVI and evapo-transpiration) as the independent variables.

In order to tell how good my model is, I need to do regression analysis for which I need a particular value of each variable at county level. For this purpose I used zonal statistics tool in ArcGIS and took the mean value of NDVI, temperature and evapo-transpiration.

As a county scale is too coarse, so taking 'Mean' is the right statistics or not?
Also in case of counties having a lot of cloud cover, zonal statistics only considers available pixels (may be too few) to compute. What is the right way to deal with such an issue?

I will really appreciate your solutions and if you have any links/research papers to support this.


texts: Summarizing population by neighborhood using PostGIS

texts:
I'm a novice with PostGIS and am currently working on a lab assignment where I have two geometric tables - one containing demographic data by blockgroup (blockgrp2010) while the other identifies various neighborhoods (nbo_hood). I'm trying to write a query that lists the population for each neighborhood. For now, I'm assuming that if the centroid of the blockgroup is inside the neighborhood, then all the population lies within that neighborhood. I drafted the query below:

select nbo_hood.name, blockgrp2010.pop10
from blockgrp2010
inner join nbo_hood on ST_Within ((ST_Centroid (blockgrp2010.geom)), nbo_hood.geom)
group by nbo_hood.name, nbo_hood.geom, blockgrp2010.pop10;

This returns a table listing the neighborhoods and corresponding population for each blockgroup. However, many of the neighborhoods are listed under multiple rows with a different population assigned (obviously because there are multiple blockgroups within the neighborhood).

How do I modify the above query (or should I start from scratch) to list each neighborhood only once with the population sum of all blockgroups?

I've tried variations of the SUM function with the SELECT statement but nothing as worked.


texts: Determining invalid parameters for geoprocessing service

texts:
I have a loop that feeds the parameters to execute syncronymous geoprocessing services.
When the loop runs it gives me the results when it receives valid params, and I get a generic error for the invalid params. The problem is that I need to know which params were invalid. It is very time-consuming to check manually which params are invalid when the loop runs over 1000 times. Is there way to customize the error message to include the params that caused the error?

this.gp = new Geoprocessor("http://xxx");
this.own(on(this.gp, 'execute-complete', lang.hitch(this, this.onExecutethis)));
        this.own(on(this.gp, 'error', lang.hitch(this, this.onError)));

onError: function (error) {
-----
}

-----
-----
for (i = 1; i < data.length; i++) {
    var params = {
              "var1": data[i][3],
              "var2": data[i][2]
            };
     this.gp.execute(params);
}

onExecutethis:function (results){
----
----
}


texts:
I'm modelling a bus route as a series of waypoints.  Most of these waypoints have just a latitude and a longitude and represent points along a road.  However, when the bus gets to a bus stop, the waypoint has an arrival and departure time associated with it.  This information could potentially be relevant to plotting this route, for example, if I was to animate the progress of the bus.

Since I am going to be exposing these models through an API, I want to keep to standards by using the GeoJSON format.  This will be my first time using GeoJSON.  I know that the first three coordinates of a GeoJSON point are [longitude, latitude, elevation].  According to Tom Macwright, GeoJSON supports multi-dimensional coordinates.  However, I couldn't find any information on any standards beyond those three.

I plan to represent the times in ticks (JavaScript ticks, i.e.  number of milliseconds since midnight of January 1, 1970).  I would like to use the following standard:

{
    "type": "Point",
    "coordinates":  [longitude, latitude, null, arrivalTime, departureTime]
}

Will this format collide with any GeoJSON standard?
Elevation data is meaningless to me.  Do I need the "null" in order to conform to GeoJSON?


texts: Is there any way to symbolize a user-supplied input parameter in ModelBuilder?

texts:
I have a model built in ESRI ModelBuilder that is run from its tool dialog, where the user supplies the input feature layer that then gets processed in the model.  I need the input feature layer to get a certain symbology when the model is run.  I've tried symbolizing it by right-clicking the input layer in ModelBuilder, going to its Properties, and setting the layer symbology using a layer file. I then save & run the model, and the input feature layer is still not symbolized - it just stays the same.

Am I doing it wrong?

I can get output layers in the model to symbolize using layer files just fine with no issues.

Is there no way to symbolize inputs... only outputs?

I'm also open to using Python to symbolize it, but not sure if this is possible?


texts: Adding dataset to ArcMap with ArcObjects and C#?

texts:
I am trying to add feature dataset to ArcMap as a layer

I  tried some code but addlayer needs iLayer but IFeatureDataset can't casting as ilayer

How Add dataset to ArcMap?

       ESRI.ArcGIS.ArcMapUI.IMxDocument mxd = ArcMap.Application.Document as
                                               ESRI.ArcGIS.ArcMapUI.IMxDocument;
        ESRI.ArcGIS.Carto.IMap map = mxd.FocusMap;
        ESRI.ArcGIS.Geodatabase.IWorkspaceFactory wsf = new
                    ESRI.ArcGIS.DataSourcesGDB.FileGDBWorkspaceFactory();
        ESRI.ArcGIS.Geodatabase.IWorkspace ws = wsf.OpenFromFile(@"D:\ARC OBJECTS\TemplateData.gdb", ArcMap.Application.hWnd);
        ESRI.ArcGIS.Geodatabase.IFeatureWorkspace fws = ws as ESRI.ArcGIS.Geodatabase.IFeatureWorkspace;
        ESRI.ArcGIS.Geodatabase.IFeatureDataset fd = fws.OpenFeatureDataset("World");


texts:
I wonder if there is a better way to join attributes based on location.
I have a vector layer with polygons representing lines.
Layer A is the layer where I want to take the correct attributes and Layer B is the layer with the merged/dissolved polygons.
They're split into multi parts and I merged/dissolved them based on 'type' attribute (see image tables).
After merging/dissolving the polygons I want to re-assign the attributes back, based on the location.
I've tried the tool join attributes by location selecting 'contains' and 'one to one' join. The result displayed here is correct but with other features it looks like that it's taking the attributes of the adjacent polygon (which has different type attribute) instead of the
the polygon that falls in the same area and same attribute type.
So I have to check them manually which is a lot of time consuming if you have hundreds of features.

I have also tried using PostGIS, which I'm new, using the following spatial query:
SELECT table_A
FROM
table_A,
table_B
WHERE
st_contains(table_B.geom, table_A.geom);
Using this I get only a few features joined. What am I missing here?
Is there any other way to do this without spending too much time checking the features manually? Using python inside or outside QGIS or Arcmap?


texts: How to always disable the “add saved file to map” checkbox in Save as ...?

texts:
I often use the "Save as.." dialog and just want the exported geojson file there. But he always enables the checkbox "add saved file to map". I really dislike it that way and it does not make any sense to perma-enable it.

Is there a way to permanently disable it?


texts: Why does Openlayers plugin crash QGIS 1.7.4?

texts:
When I open the openlayers plugin my qgis 1.7.4 crashes.

I use Win7.

I'm using the CRS Google Mercator and QGIS crashes.

However, I've install the QGIS 1.9.90 and the plugin works in this version... but... There are always a but... the Google Satellite images don't adjust with my shapefiles. It's a strange behaviour...

Any suggestion to solve this?


texts: Field of view (FOV) of the Camera - ArcGIS Java runtime SDK

texts:
I'm trying to set the 3D Camera of the SceneView in the Java SDK and I cannot find a "field of view" / "angle of view" member.
I was wondering how the camera's field of view is being set, since I need a way to change it according to the horizontal and vertical field of views.


texts: "fill in" .vrt with missing locations in a mosaic

texts:
I have ~2000 NASA SRTM .hgt files.  I wanted to join them together in a single mosaic, and to save disk space I used GDAL to create a .VRT rather than a very large GeoTiff.  This is read by Mapnik.

What I would like to do is "fill in" the areas in the mosiac that aren't covered by one of the original .hgt files.  In other words, I want to pretend as if those uncovered areas are actually covered by a .hgt that has a constant value for all data points.  Thus, the resulting .VRT will be able to report a value for those areas.  (or, if I end up making a .tif instead of .vrt, the .tif would have a value for each of those areas).

Is this possible?

(to elaborate: the reason I want to do this is because I am generating tiles to cover the entire world - this will be an overlay, not base tiles.  My data does not cover everywhere, and for the areas that it doesn't cover I want to be able to style it in a particular color.  Normally, I could use Mapnik's 'background-color' attribute, and it would color any background area in the map.  However, one of my styles for the data areas I want to cover is transparent, and setting the background-color would erroneously fill in the transparent areas.  Hence, my desire to modify the data set)


texts: Assign Values to Breaks in Graduated Colors Symbology - ArcGIS Pro

texts:
In ArcMap you could assign values to each break within a layer's graduated color symbology using classBreakValues and classBreakLabels (A code sample of this can be seen on this page http://desktop.arcgis.com/en/arcmap/10.3/analyze/arcpy-mapping/graduatedcolorssymbology-class.htm)

However, I have now transitioned over to ArcGIS Pro and since the ArcPy between Map and Pro are different I can't figure out how to assign values to each break within a layer's graduated color symbology.

Is it possible to assign values to each break within a layer's graduated color symbology in ArcGIS Pro?

Update to Question:

I now have a piece of code that will assign a value of 2 to the first break within a layer's graduated color symbology and then tell the following breaks to increase by 10. However, I still can't figure out how to personally assign a random value to each break.

Here is my code:

import arcpy, os, sys, locale
p = arcpy.mp.ArcGISProject("Current")
m = p.listMaps("Map")[0]
l = m.listLayers("ProvEPOSigNet_ZIPs")[0]
sym = l.symbology
sym.renderer.classificationField = "SUM_USER_Patient_Count"
sym.renderer.breakCount = 3

breakVal = 2

for brk in sym.renderer.classBreaks:
  brk.upperBound = breakVal
  brk.label = "\u2264" + str(locale.format("%d", breakVal, grouping=True))

  breakVal +=10

l.symbology = sym


texts:
I would like to have label visible for both cryterias in my layer.
One case works fine whereas the other one doesn't.
In the code everything is exactly the same...
    function categories_Area5chamber_14(feature, value, size, resolution, labelText,
                   labelFont, labelFill, bufferColor, bufferWidth,
                   placement) {
            switch(value.toString()) {case 'BT':
                return [ new ol.style.Style({
    stroke: new ol.style.Stroke({color: 'rgba(35,35,35,1.0)', lineDash: null, lineCap: 'butt',
    lineJoin: 'miter', width: 0}),fill: new ol.style.Fill({color: 'rgba(0,0,0,1.0)'}),
    text: createTextStyle(feature, resolution, labelText, labelFont,
                          labelFill, placement, bufferColor,
                          bufferWidth)
    })];
                break;
    case 'Voneus':
                return [ new ol.style.Style({
    stroke: new ol.style.Stroke({color: 'rgba(35,35,35,1.0)', lineDash: null, lineCap: 'butt',
    lineJoin: 'miter', width: 0}),fill: new ol.style.Fill({color: 'rgba(252,0,0,1.0)'}),
    text: createTextStyle(feature, resolution, labelText, labelFont,
                          labelFill, placement, bufferColor,
                          bufferWidth)
    })];
                break;}};

the createTextStyle function looks exactly the same!
Next in the styling:
  var style_Area5chamber_14 = function(feature, resolution){
  var context = {
    feature: feature,
    variables: {}
};
var value = feature.get("Chamber Type");
var labelText = "";
size = 0;
var labelFont = "15px \'MS Shell Dlg 2\', sans-serif";
var labelFill = "#787878";
var bufferColor = "";
var bufferWidth = 0;
var textAlign = "center";
var offsetX = 0;
var offsetY = 0;
 var placement = 'point';
  if (feature.get("Chamber No") !== null) {
    labelText = String(feature.get("Chamber No")); labelFont = "15px \'MS Shell Dlg 2\', sans-serif";
 }
  var style = categories_Area5chamber_14(feature, value, size, resolution, labelText,
                      labelFont, labelFill, bufferColor,
                      bufferWidth, placement);
   return style;

I have got all the variables defined, which unfortunately are displayed with one case only. How about the other case?
For this occasion I'vee made the alterations:
  else if (feature.get("X") !== null) {
    labelText = String(feature.get("X")); labelFont = "15px \'MS Shell Dlg 2\', sans-serif";
  }

but still without the result...
Does anyone know how to solve it?


texts: Create WMS service from existing image tiles

texts: I have image tiles (jpg image files with txt georeference info files; see picture) stored in catalog (about 12GB, size one image tile ~ 275*466 meters). It possible create WMS service of the tiles to load it into QGIS and its use by multiple users in a network. As I understand it can be done with MapServer, GeoServer, PostGIS. Is there a manual how I should load this image tiles in MapServer (or GeoServer, PostGIS)?


texts:
I saw Images not draping onto DEM in ArcScene? but it is different from my issue. I have an image with projection:

NAD_1983_StatePlane_Oregon_North_FIPS_3601_Feet.

I have a DEM and set the base height to itself and I am getting a good DEM representation of the terrain. I added the image and set the base height to the DEM and it doesn't drape on top of DEM, it appears to be below DEM surface.


texts: Using OpenLayers' Refresh strategy for GeoJSON

texts:
I am trying to use the Refresh Strategy reading a GeoJSON file in OpenLayers, however it doesn't work (but it works perfectly using for instance .Fixed() strategy), here is my code:

    var acediot_map = {};
    var map;

    function initMap() {
        map = new OpenLayers.Map("info");
        var mapnik = new OpenLayers.Layer.OSM();
        map.addLayer(mapnik);

        map.addControl(new OpenLayers.Control.ScaleLine());
        map.addControl(new OpenLayers.Control.OverviewMap());
        map.addControl(new OpenLayers.Control.MousePosition());
        map.addControl(new OpenLayers.Control.LayerSwitcher());

        var myStyles = new OpenLayers.StyleMap({
            "default": new OpenLayers.Style({
                pointRadius: 10,
                fillColor: "${getColor}",
                strokeWidth: 1,
                strokeColor: "#800026",
                fillOpacity: 1.0
            },
            {
                context: {
                    getColor : function (feature) {
                        return feature.attributes.severity == 'LOW' ? '#FFFF00' :
                            feature.attributes.severity == 'MID'    ? '#FF8C00' :
                            feature.attributes.severity == 'HIGH'   ? '#FF0000' :
                                                                      '#7CFC00' ;
                    }
                }
            })
        });

        var myInterval = 5000;

        acediot_map.report_layer = new OpenLayers.Layer.Vector("GeoJSON", {
                projection: new OpenLayers.Projection("EPSG:4326"),
                strategies: [
                    new OpenLayers.Strategy.Refresh({interval: myInterval})],
                protocol: new OpenLayers.Protocol.HTTP({
                    url: "myGeoReport.json",
                    format: new OpenLayers.Format.GeoJSON()
                }),
                styleMap: myStyles
        });

        acediot_map.location_layer = new OpenLayers.Layer.Vector("LocationLayer", {
            projection: new OpenLayers.Projection("EPSG:4326"),
            strategies: [new OpenLayers.Strategy.Fixed()],
            protocol: new OpenLayers.Protocol.HTTP({
                url: "myGeoLocation.json",
                format: new OpenLayers.Format.GeoJSON()
            })
        });

        overlays = [acediot_map.report_layer, acediot_map.location_layer];
        map.addLayers(overlays);

        map.setCenter(
            new OpenLayers.LonLat(9.2257418, 45.4788448).transform(
                new OpenLayers.Projection("EPSG:4326"),
                map.getProjectionObject()
            ), 13
        );
}

And here is an example of my GeoJSON file

{"type":"FeatureCollection",
  "features":[
    { "type": "Feature",
      "geometry": {
        "coordinates": [
          9.2257418,
          45.4788448
        ],
        "type": "Point"
      },
      "properties": {
        "accuracy": 20.0,
        "description": "Prima prova ",
        "notes": "",
        "picture": "",
        "severity": "LOW",
        "team": {
          "id": 2,
          "name": "Eugenio Gianniti",
          "org": "Politecnico di Milano"
        },
        "timestamp": 1411920949.921636
      }
    }
  ]
}


texts: Point to raster conversion works as a tool in arcmap but fails when attempted with identical parameters in a script

texts:
I'm working on a script that will automate a stream delineation protocol given a DEM and a culvert layer. I was having difficulty troubleshooting problems with the point to raster conversion tool so I ran the tool in arcmap (successfully) and copied the python snippet from the results. I then closed arcmap to avoid lock problems and tried running the script (shown with my original code commented out) with the snippet. It produced the same error I encountered with my code.

Code:

    # Environment Set-Up and parameter initialization
import arcpy
from arcpy import *
from arcpy.sa import *
from arcpy.cartography import *
arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension ("Spatial")

    #Provide the full filepath for you working geodatabase here. ##DATASET STARTING NAME IS "LIDAR_DEM", MAKE A COPY OF YOUR ORIGINAL AND CALL IT THAT##
WSV = r"C:\Users\rdcoopercaroselli\RowynProjects\LiDAR_Work\AutomationTest\Automation_try2.gdb"
arcpy.env.workspace = WSV

    #Provide the culvert datset here
CulvertFeature= WSV+'/Culverts'

    #Get the cell size, pixel type, desired operating extent, and coordinate system of the LiDAR DEM. Set up DEM as Snap Raster for raster calculations
DEM = Raster (WSV+"/LIDAR_DEM")
CellSize = GetRasterProperties_management (WSV+"/LIDAR_DEM", "CELLSIZEX")

ptype    = DEM.pixelType
code     = {'U1':'1_BIT','U2':'2_BIT','U4':'4_BIT','U8':'8_BIT_UNSIGNED','S8':'8_BIT_SIGNED','U16':'16_BIT_UNSIGNED', 'S16':'16_BIT_SIGNED',
           'U32':'32_BIT_UNSIGNED','S32':'32_BIT_SIGNED','F32':'32_BIT_FLOAT','F64':'64_BIT'}
ptype    = code[ptype]

arcpy.env.extent     = DEM.extent

arcpy.env.snapRaster = WSV+"/LIDAR_DEM"

     # Convert Culvert Shapefile to Raster with same cell size as DEM
#CreateRasterDataset_management (WSV, "CulvertRaster", CellSize, ptype, DEM,1)
#PointToRaster_conversion(CulvertFeature,"ELEV",WSV+"/CulvertRaster","MOST_FREQUENT","NONE",CellSize)
#CulRast = Raster (WSV+"/CulvertRaster")

# Replace a layer/table view name with a path to a dataset (which can be a layer file) or create the layer/table view within the script
# The following inputs are layers or table views: "Culverts"
PointToRaster_conversion("Culverts","ELEV","C:/Users/rdcoopercaroselli/RowynProjects/LiDAR_Work/AutomationTest/Automation_try2.gdb/CulRast","MOST_FREQUENT",
                         "NONE","C:/Users/rdcoopercaroselli/RowynProjects/LiDAR_Work/AutomationTest/Automation_try2.gdb/LIDAR_DEM")

Error message:

Traceback (most recent call last):
  File "C:\Users\rdcoopercaroselli\RowynProjects\LiDAR_Work\AutomationTest\R05_DRAINAGE_DELINEATOR_Backup.py", line 40, in <module>
    "NONE","C:/Users/rdcoopercaroselli/RowynProjects/LiDAR_Work/AutomationTest/Automation_try2.gdb/LIDAR_DEM")
  File "C:\Program Files (x86)\ArcGIS\Desktop10.1\arcpy\arcpy\conversion.py", line 2393, in PointToRaster
    raise e
ExecuteError: ERROR 999999: Error executing function.
The table was not found. [CulRast]
The table was not found. [fras_bnd_CulRast]
No spatial reference exists.
Failed to execute (PointToRaster).

This error message alternately appears with and without a long series of
"A column was specified that does not exist." messages


texts: Getting elevation from polyline feature using QGIS?

texts:
I have imported contour lines from AutoCAD DXF file. I've 3D geometry. Í would like to add a new column to the attribute table of the contour lines. In th field calculator I user the following expression:

z(point_n($geometry, 1))

The values in the new fields are -2,147,483,648, but if I click in a field I get NULL.

If I use the x or y function I get the correct x or y coordinate.


texts: Joining field to a polygon shapefile from another overlapping polygon in QGIS

texts:
I have two polygon shapefiles (parcel and zoning). Parcel shapefile is overlapped by zoning shapefile. Parcel layer does not contain information about zone and I want to add that information to each parcel from zoning shapefile.

My idea is to create centroids from parcel shapefiles and after that to join zoning shapefile by location. On that way, I will add zoning attributes in my centroid point file. When I finish that I will join by location that point file and my parcel file. Unnecessary fields could be deleted after that.

I need to know is there any faster way to finish this job in QGIS?

This is a screenshot of a part of my parcel shapefile:

This is a screenshot of a part of my zoning shapefile:

Just to note that this is a very large data, about 800000 parcels.


texts: How to generate missing points in GPS trajectories?

texts:
I have some GPS trajectories that include gps locations recorded on variable frequencies (every 1 min ,5 min , 20 min,etc.)using a geographic projection WGS84. In my work, i need a fine-grained GPS collection (e.g. 30 sec) to be able to recognize accurately users' activities grom their GPS traces.

I thought then to generate the missing GPS points between two known points (for instance, when the points frequency is every 5 minutes, i need to generate the 8 missing points that should be recorded the frequency was every 30 second). Note that i have to estimate the speed of the missing GPS points in my interpolation .
Could you suggest  some mathematical estimation or an existing tool that can answer my issue?


texts:
Couldn't find a similar question for Linux, just Windows.

I'm running Linux Mint 19.3.

Traceback (most recent call last):
  File "", line 1, in
  File "/usr/lib/python3/dist-packages/qgis/__init__.py", line 70, in
    from qgis.PyQt import QtCore
  File "/usr/lib/python3/dist-packages/qgis/PyQt/QtCore.py", line 24, in
    from PyQt5.QtCore import *
ImportError: /usr/lib/x86_64-linux-gnu/libQt5Core.so.5: version `Qt_5.14' not found (required by /home/jrickards/.local/lib/python3.6/site-packages/PyQt5/QtCore.abi3.so)

Python version:
3.6.9 (default, Apr 18 2020, 01:56:04)
[GCC 8.4.0]

QGIS version:
3.10.5-A Coruña 'A Coruña', exported

Python path:
['/usr/share/qgis/python', '/home/jrickards/.local/share/QGIS/QGIS3/profiles/default/python', '/home/jrickards/.local/share/QGIS/QGIS3/profiles/default/python/plugins', '/usr/share/qgis/python/plugins', '/usr/lib/python36.zip', '/usr/lib/python3.6', '/usr/lib/python3.6/lib-dynload', '/home/jrickards/.local/lib/python3.6/site-packages', '/usr/local/lib/python3.6/dist-packages', '/usr/lib/python3/dist-packages', '/usr/lib/python3.6/dist-packages']

How do I resolve this?


texts: What are valid values for "data defined" labeling settings?

texts:
I'm using QGIS 1.8.
I'd like to customize label features on an existing shapefile using the new label engine and it's "data defined settings" option. I'd like to know what the attribute table field types should be for each of the settings: All 7 font options, all 2 buffer options, and all 6 position options. Should they all be numeric fields, or some text, or both? I want to create these fields in the attribute table correctly so that I can map them to the label field settings. I can't find information on the proper field types anywhere.
Thanks for your thoughts.


texts: Draw Polygon according to an external GeoJSON file (on a server)

texts:
I'm programming a map with OpenLayers (4.6.5) and I have a GeoJSON file stored on a server and I want to use it. This GeoJSON file represents Polygons. Before I was working with Points and all was working fine, so I don't get why it's not the same thing with Polygons.

Here is my code :

<!DOCTYPE html>

      <meta charset="utf-8" />
      <meta name="viewport" content="width=device-width, initial-scale=1.0">
      <title>Carte Quartiers</title>
      <link rel="stylesheet" href="ol/ol.css">
      <link rel="stylesheet" href="ol/ol3-layerswitcher.css">

      <link rel="stylesheet" href="style.css">
      <link rel="stylesheet" href="https://stackpath.bootstrapcdn.com/bootstrap/4.1.0/css/bootstrap.min.css" integrity="sha384-9gVQ4dYFwwWSjIDZnLEWnxCjeSWFphJiwGPXr1jddIhOegiu1FwO5qRGvFXOdJZ4" crossorigin="anonymous">
      <link rel="stylesheet" href="https://use.fontawesome.com/releases/v5.0.12/css/all.css" integrity="sha384-G0fIWCsCzJIMAVNQPfjH08cyYaUtMwjJwqiRKxxE/rx96Uroj1BtIQ6MLJuheaO9" crossorigin="anonymous">

What am I doing wrong ? I have found some examples but the GeoJSON code was inside the script ...


texts: C#/.Net Library that can draw a line across the date line

texts:
I want to be have an embedded map in my application that can handle drawing lines across the date line. I have tried out dotspatial and sharpmap and neither have been able to handle this.

Edit: Here is an early example of the kind of information that would need to draw over the date line http://imgur.com/B3Fg7R5 Please excuse the very low res map product.

The software needs to support loading custom map products via geotiffs or something similar.

This is a c# winforms application.


texts: How to draw SVG images stored in PostGIS with Mapserver?

texts:
I have a table in PostGIS with two columns. In first are stored SVG images as text. In second is geometry (Point in epsg 4326).
Is it possible to draw SVGs from PostGIS on position given in the geometry column with Mapserver?


texts: Export 3D layers automatically in CityEngine to FBX?

texts: How to automatically export all layers in a scene in CityEngine to corresponding FBX? Is there a batch processing or a programming approach for exporting them to FBX?


texts: replaceDataSource issue when iterating through a loop for multiple map documents

texts:
Trying to update data sources across multiple map documents. The first iteration works but fails on the second no mater the order.

The input csv is the output from another script that lists unique broken datasources across the map documents. I then fill in a second column for new links.

If I run on an individual map document it works fine even for the one that failed when second in the loop.

I tested on two map documents that have the same broken data links and still it fails on the second one.

In this instance I am only working with data on SDE.

import arcpy, csv, os
import arcpy.mapping as m

# assess whether shapefile, geodatabase, or sde
# and get parameters for replaceDataSource()
def updateDataSource(datapath, broken_lyr, workspace_path=None, workspace_type=None, dataset_name=None):
    if datapath.endswith(".shp"):
        workspace_path = datapath.rsplit("\\", 1)[0]
        workspace_type = "SHAPEFILE_WORKSPACE"
        dataset_name = datapath.rsplit("\\", 1)[-1][0:-4]
    elif datapath.find(".sde") != -1:
        workspace_path = datapath.rsplit(".sde", 1)[0] + ".sde"
        workspace_type = "SDE_WORKSPACE"
        dataset_name = datapath.rsplit("\\")[-1]
    elif datapath.find(".gdb") != -1:
        workspace_path = datapath.rsplit(".gdb", 1)[0] + ".gdb"
        workspace_type = "FILEGDB_WORKSPACE"
        dataset_name = datapath.rsplit(".gdb\\")[-1]
    else:
        print "Neither SHP, GDB, or SDE"

    if workspace_path:
        broken_lyr.replaceDataSource(workspace_path, workspace_type, dataset_name)

input_csv = r"C:\Users\glen.bambrick\Documents\csv\broken_links_20180523.csv"
input_folder = r"C:\Users\glen.bambrick\Documents\mxd"

# store data from csv in here
# key: broken data link
# value: new data link
broken_dict = {}

with open(input_csv, 'rb') as read_csv:
    reader = csv.reader(read_csv)
    next(reader, None)
    for row in reader:
        broken_dict[row[0]] = row[1]

# list map documents from folder provided
mxd_files = [mxd_file for mxd_file in os.listdir(input_folder) if mxd_file.endswith(".mxd")]

# iterate through map documents
for mapdoc in mxd_files:
    print mapdoc
    mxd = m.MapDocument("{0}\\{1}".format(input_folder, mapdoc))
    # get a list of broken links
    list_broken = m.ListBrokenDataSources(mxd)
    if len(list_broken) == 0:
        continue
    else:
        for broken in list_broken:
            # if the link is in the dictionary then update it
            if broken_dict[broken.dataSource]:
                updateDataSource(broken_dict[broken.dataSource], broken)

    mxd.save()

This is the error I get

map1.mxd
map2.mxd
Traceback (most recent call last):
  File "C:/Users/glen.bambrick/Desktop/test.py", line 45, in <module>
    updateDataSource(broken_dict[broken.dataSource], broken)
  File "C:/Users/glen.bambrick/Desktop/test.py", line 19, in updateDataSource
    broken_lyr.replaceDataSource(workspace_path, workspace_type, dataset_name)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.5\ArcPy\arcpy\utils.py", line 182, in fn_
    return fn(*args, **kw)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.5\ArcPy\arcpy\_mapping.py", line 682, in replaceDataSource
    return convertArcObjectToPythonObject(self._arc_object.replaceDataSource(*gp_fixargs((workspace_path, workspace_type, dataset_name, validate), True)))
ValueError: Layer: Unexpected error

As said already it works perfect on the first iteration and when run one at a time but I cannot get it to run through multiple like it should.
Any ideas?


texts:
I am using ArcGIS 10.2.1. I have a set of clustering points. I am trying to create a zone for these points. Is there a tool or script to group the points in an area and convert them to a polygon?

I don't have the POINT TO POLYGON tool available


texts: How to get linestring length while drawing?

texts:
I've got a control to draw lines:

new OpenLayers.Control.DrawFeature(vecLayer, OpenLayers.Handler.Path, {
            callbacks : {
            //    "done" : doneHandler,
                "point" : function(aPoint){
                    console.log(aPoint);
                }
            }
        })

and now i can get all points from linestring i add to the map. But what i need is to get length of the feature that is not finished (the reason i can't use done callback).

How can i make this? I'm completly stuck, any help!


texts: Update Query on Query layer in ArcGIS Pro

texts:
While updating query on query layer in ArcMap is nightmare process, it seems with ArcGIS Pro this could possibly be better thanks to python CIM Access.
I have tried following code and it works in the way I am able to get / replace SQL query on layer:
aprx = arcpy.mp.ArcGISProject('current')
am = aprx.activeMap
lyr = am.listLayers()[0]  # Picks first layer in active map
lyr_cim = lyr.getDefinition('V2') # Reads CIM definition
query = lyr_cim.featureTable.dataConnection.sqlQuery # Reads SQL Query
query += ' WHERE GEOM.STX > 0' # Example of query update
lyr_cim.featureTable.dataConnection.sqlQuery = query # Set query in definition
lyr.setDefinition(lyr_cim) # Update definition on layer

The problem is at that moment all features from layer disappear and I get red exclamation symbol next to layer name with tool tip 'Click to repair layer's data source.'
There is nothing wrong with the query or data source, ArcGIS Pro just thinks for some reason there is. If you click through user interface of layer query properties, the updated query is there, and is correct. You just must pretend edit (e.g. delete space and add space) to be able to manually click on 'Validate' button, and 'Next', and 'Finish' button in the next window.
Any idea how to tell to ArcGIS Pro programmatically that the layer is alright?
Side note: This is a workaround to workaround to workaround in ArcGIS Pro to be able to programmatically filter features for export to ArcGIS Online. Filtering by extent on layer does not work and sends all data to AOL, filtering by definition query fails too with error - seems publishing can't deal with definition query containing underlying DB methods like STEnvelope, STCentroid etc. Ultimately, what I want is to find some way to be able to clip layers to current map extent on the fly for publishing to ArcGIS Online


texts: Can't plot polygons with R: Invalid number of points in LinearRing found 3 - must be 0 or >= 4?

texts:
I want to plot a map of my region, with colour filled polygons using R. I have downloaded municipal (lines) and province (polygons) boundaries. I've used SAGA to merge municipal and province shapes to produce new boundaries (polygons) to read and plot with R and ggplot2.

I succesfully plotted the original province (polygons) but failed to plot new boundaries (polygons) with this error message:

  Error en createPolygonsComment(p) :    IllegalArgumentException:
  Invalid number of points in LinearRing found 3 - must be 0 or >= 4

Maybe I made some mistake when merging original polygons and lines to create the new boundaries. Shoud I check shp attributes?

You can find original boundaries here and newly created ones here

EDIT: Added link to data file lineas_limites_municipales_etrs89.dbf

And the original R code I used to plot original boundaries, not working for new ones (see error message above):

require("rgdal")
require("maptools")
require("ggplot2")
require("plyr")

# Reading municipal boundaries

 esp = readOGR(dsn=".", layer="lineas_limite_municipales_etrs89")

 muni=subset(esp, esp$PROV1 == "46" | esp$PROV1 == "12" | esp$PROV1 == "3")
 muni@data$id = rownames(muni@data)
 muni.points = fortify(muni,region="id")
 muni.df = join(muni.points, muni@data, by="id")

# Reading province boundaries

prov = readOGR(dsn=".", layer="poligonos_provincia_etrs89")

pr=subset(prov, prov$CODINE == "46" | prov$CODINE == "12" | prov$CODINE == "03" )
pr@data$id = rownames(pr@data)
pr.points = fortify(pr, region="id")
pr.df = join(pr.points, pr@data, by="id")

ggplot(muni.df) + aes(long,lat,group=group) + geom_path(color="blue") +
coord_equal()+ geom_path(data=pr.df,   aes(x=long,y=lat,group=group),color="red", size=0.5)


texts:
Is there a way to add text/annotation to my map in QGIS as a layer?
I have a map of vegetation types an point data the text (field names/flurnamne) should be between these two layers. It there must not be a background or box around it and should not cover the points. It would be ideal if I could arrange the text on a path or line as the shortest name must stand for the largest area.

Is there a way to add annotation just a text, no box, background or point?

At the moment I just added text in print composer and set transparency, but results are not very satisfying

EDIT: The text should be under the points. I have a map with the vegetation types and the point data and want to add the fieldnames/flurnamen of certain area. It would ideal if spacing etc. could be adjusted so that the text covers the whole area, e.g. name_1 in on the line.


texts:
How do I handle map legend with tens of classes using QGIS? I searched for instructional videos and previous post relating to the question but could not find any helpful material. The legend text in the map below are partly invisible, I would like to have all the legend text visible.

N.B: I thought of having the legend part as a second page on the print composer   but on importing the map, the legend part was not printed.


texts: How to assign diameter size to valve feature class based on waterline size?

texts:
I have two shapefiles, valves (points) and waterlines (polyline). The waterlines have a diameter attribute that is populated with the correct size. The valves have an attribute for diameter size as well but some of the values are showing 0 for the size. I need to update the 0 sized valves with the correct size from the waterline it is on.

The valves are in close proximity to the waterline if not directly on it. I don't have access to the "Near" function because of the level of my license. I am using ArcGIS 10.2 for Desktop with a Standard license.

What I have done so far is select valves by 1ft location to specific waterlines with specific diameter size (ie selected valves that are within 1ft of 16 inch waterline). From this point I'm not sure what I should do next or even if this is the correct path to the solution.


texts: Stream network distance between specific paired points

texts: I have a stream network dataset and am attempting to do a distance analysis between paired points located in two different layers. Ex. I have 500 points in layer A with given coordinates, and then 500 points with different coordinates in layer B. I plan to base the distance-pairing on a shared ID. Having looked at other questions related to this, I have often seen the OD Cost Matrix suggested as a solution, or tools that give linear distance. I don't think the OD Cost Matrix will work for me- I want the stream network distance between the points that share an ID, not the distance between Point A and all other points.


texts: Create a rectangular buffer from point shapefile using QGIS/GRASS/SAGA GIS/R

texts: I have a points shapefile and I want to create a rectangular buffer (22 x 16 km) around each point. My goal is to create a 16x22 km grid where each of the points is the centroid of the grid. Seed the attached image for an example of what I want to do.


texts: Embed WMS server login info into MapInfo tab file

texts:
I need to create a .tab file that opens a restricted/protected WMS server without prompting for credentials every time.

Currently when I open it, a prompt for the user name and password comes up. I don't want to give this to all the users on the system but have it embedded in the .tab/.xml file...ideally scrambled.

The .tab is

!table
!version 750
!charset WindowsLatin1

Definition Table
File "qld_sisp_imagery_wms.xml"
Type "WMS"
CoordSys Earth Projection 8, 116, "m", 147, 0, 0.9996, 500000, 10000000
begin_metadata
"\IsReadOnly" = "FALSE"
"\MapInfo" = ""
"\MapInfo\TableID" = "8b6009d6-c314-4bd4-9223-d7d290ed84e4"
end_metadata

the .xml is

<MI_WMS_MapDefinition version="1.0" wms_version="1.3.0"><HTTP GetCapabilities="https://auth-spatial-XXX/ImageServer/WMSServer?" GetMap="https://auth-spatial-XXX/ImageServer/WMSServer"/><Exception format="XML"/><SRS>EPSG:28355</SRS><Image bgcolor="16777215" transparent="true" format="image/png"/><Layer><Name>LatestStateProgram_SISPUsers</Name><Style></Style><LatLonBoundingBox minx="-180" miny="-90" maxx="180" maxy="90"/><BoundingBox SRS="EPSG:28355" minx="-481881.729993" miny="6734668.137079" maxx="1225547.023262" maxy="8981699.047294"/></Layer></MI_WMS_MapDefinition>

QGIS enables this in it's convert to configuration process when you add a WMS feed.

Is there a similar MI process?


texts: Copy multiple feature classes in FME without merging

texts:
I am trying to create an FME workspace that reads a series of feature classes from an Oracle SDE database and copies them to another Oracle SDE database, one by one. The key thing is I do not want to merge them into a single feature class. Do I need to have multiple writers, or can I use some kind of "dynamic" setting to maintain separate tables? I'm not finding very much documentation about this type of workflow. I found this article:

Pass a List of Tables or Layers to an FME Reader using Python Scripted Parameters

but I would prefer to have all of my sources contained within the workspace.


texts: How do I store GeoJSON format in Firebase so it can be used for Leaflet GeoJSON

texts:
I am new with Firebase, Leaflet, and GeoJSON ..

I have a question, how to store data with GeoJSON format in Firebase, so it can be retrieved and used for leaflets GeoJSON to display marker(s). I know there is GeoFire library that allows me to store and query a set of keys based on their geographic location, but I want to use GeoJSON for some reason.

as i know if Firebase store the file as an object and GeoJSON as array, but is there any way? Please let me know.


texts: Expression to extract part of an address

texts:
I have a column of addresses and want to separate into several parts.
Ex. house number, street name, city, zip. For example lets say I have a cell with this string: '111 Example St, City, State 89786 zipcode, country'. This is how I'd like to separate it:
house_numb | street_name | city | state | zip | country |

'111' | 'Example St' | 'City' | 'State' | '89786' | 'US'

I was able to use a left regex expression to create the house number. Since the city, zip, state, and country are the same I was able to populate those columns as well. But I cannot seem to get the street name in the string. What expression can be used to extract the street name?
I tried regexp_substr( "Text", ';(\\d*)\}' ) and replaced the ';' with ' ' and '}' with ',' and did not get the street name between the first space and the first comma.


texts:
Following an answer to my previous question I am planning to work with NASA ASTER GDEM global dataset. I am trying to download it according to instructions given here for a couple of days now.

Meanwhile can somebody answer: does this datum contain objects like trees and buildings?


texts:
This is the second part of my edited two part question (first part here).

So I used the following script:

// define the var
var Catchment = /* color: 98ff00 */geometry;
var landcover = ee.Image('users/roynahas/ESACCI-LC-L4-LCCS-Map-300m-P5Y-2010-v161_RECLASS').select('b1');
// Clip the image to the polygon geometry and add it to the map
var landcover_clip = landcover.clip(Catchment);
var sld_intervals =
'<RasterSymbolizer>' +
 '<ColorMap  type="intervals" extended="false" >' +
 '<ColorMapEntry color="#FFFF00" quantity="1" label="Agriculture"/>' +
 '<ColorMapEntry color="#00FF00" quantity="2" label="Grassland and Shrubland"/>' +
 '<ColorMapEntry color="#008000" quantity="3" label="Forest"/>' +
 '<ColorMapEntry color="#00FFFF" quantity="4" label="Flooded"/>' +
 '<ColorMapEntry color="#FF00FF" quantity="5" label="Urban areas"/>' +
 '<ColorMapEntry color="#808080" quantity="6" label="Bare areas"/>' +
 '<ColorMapEntry color="#0000FF" quantity="7" label="Water"/>' +
 '<ColorMapEntry color="#FFFFFF" quantity="8" label="Permanent snow and ice"/>' +
 '</ColorMap>' +
'</RasterSymbolizer>';
Map.addLayer(landcover_clip.sldStyle(sld_intervals), {}, 'IGBP classification styled');
// Print out the frequency of landcover occurrence for the polygon.
var frequency = landcover.reduceRegion({
  reducer:ee.Reducer.frequencyHistogram(),
  geometry:Catchment,
  scale:300
});
print('landcover frequency', frequency.get('b1'));

To get the following console output:

My question : How can I change the labels in the console output in Google Earth Engine?

Instead of "1", "2", "3", "4", "5", "6", "7", and "8", I want to have "Agriculture", "Grass and Shrubland", "Forest", "Flooded", "Urban areas", "Bare areas", "Water", and "Permanent snow and ice".


texts: I have a set of point cloud data in QGIS. After triangulating them and forming a Raster temporary file, i would like to find the 2D surface area of this file. IS there any tool to do it? Or do i need to use meshing methods for the point cloud file?


texts: Generate polygons for each map extent in ArcGIS DDP layout

texts:
I have a set of maps I am creating with Data Driven Pages. Each map has a scale and rotation set by the feature class that is driving the map creation.  I would like to create a polygon that is the extent of each page that is created - not the polygon from the DDP index feature class. Is this possible?

I am playing around with the Fishnet, as my page size is the same, but not sure I can put the rotation into that.


texts: ArcPy: Using lblclass expression to change font style and color

texts:
I´m trying to label features of a dataset using ArcPy and the lblclass function.
While commands like "VBNewline" and "Round" work pretty well, I´m not able to change colors or font style. An expression that works in ArcMap doesn´t work in the script.
Could anyone help me out please?
I´m using ArcGIS 10.0

Code that works:

 #Adding Labels
    layer = arcpy.mapping.ListLayers(mxd, "")[0]
    if layer.supports("LABELCLASSES"):
        for lblclass in layer.labelClasses:
            lblclass.className = "Flaeche"
            lblclass.expression = "[FID] & VBNewLine & Round([AREA_mm2],2)"
        lblclass.showClassLabels = True
    layer.showLabels = True
    arcpy.RefreshActiveView()

Code that works in ArcMap but not in ArcPy

#Adding Labels
layer = arcpy.mapping.ListLayers(mxd, "")[0]
if layer.supports("LABELCLASSES"):
    for lblclass in layer.labelClasses:
        lblclass.className = "Flaeche"
        lblclass.expression = ""<CLR red='255' green='0' blue='0'>" & [FID] & VBNewLine & Round([AREA_mm2],2) & "</CLR>""
    lblclass.showClassLabels = True
layer.showLabels = True
arcpy.RefreshActiveView()


texts: Not all points overlain on raster map are showing up with R

texts:
Using R, I loaded a raster file (.tif) and overlaid points on the map. I have 46 points but only about half of them appear on the map. I can't figure out why other points are not showing up. Here are the codes. I attached an image of all points and the R plotting results.

'db' (data fame) contains latitude and longitude in decimal degrees. The locations of those points actually showed up on the map are correct. And I can extract cell values of 'all points' from the raster file. It's just that the plot doesn't show all points in the data frame.

## Load points

sPoints=data.frame(
  x=as.numeric(db[1:46,4]),
  y=as.numeric(db[1:46,3])
)
samplingpoints=SpatialPoints(sPoints)

## Create SpatialPointsDataFrame

spFrame=data.frame(ID=1:46)
str(spFrame)

SamplingPointsTable=SpatialPointsDataFrame(samplingpoints, spFrame)

## Project points

proj4string(SamplingPointsTable) <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
SamplingPointsTable <- spTransform(SamplingPointsTable,CRS("+proj=utm +zone=xx +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"))

## Plot sampling points on a map

plot(map.tif)
title(main="Map Title", sub="Sub Title",cex.main=0.8, cex.lab=0.8, cex.axis=0.8, cex.sub=0.8)
plot(SamplingPointsTable, add=T, pch=20, col=c("transparent", "black", size = 1))

Is this only because of scaling (points are too large for the R plot map)?

The projections of the map and the points are matching.


texts: PostGIS "Overlay" Style Union, not "Dissolve" Style

texts:
I'm trying to use PostGIS to do a 'union' of polygons. When I Googled "PostGIS Union", I found the ST_Union function, which performs what I would call a 'dissolve,' not a Union. What I actually want to do is what GRASS refers to as an "Overlay" with the "OR" operator, or QGIS/ArcGIS just simply call a "Union." For example:

From what I can tell, in PostGIS it appears I need to use both ST_Intersection and ST_SymDifference together to get the results I am seeking. I have had some success with the following syntax, but it's terribly inefficient. Is there a better/more efficient way to do what I'm hoping to accomplish?

INSERT INTO "CombinedResults" ("geom")
SELECT ST_SymDifference(
        "test1".GEOM
        ,"test2".GEOM
        )
    FROM "test1","test2";

INSERT INTO "CombinedResults" ("geom")
SELECT ST_Intersection(
        "test1".GEOM
        ,"test2".GEOM
        )
    FROM "test1","test2"


texts: OpenLayers & GeoServer: How to force TILEMATRIX prefix

texts:
I set up a WMTS in GeoServer-GWC. I generated the tiles with the gridset EPSG:4326.

In my OpenLayers client, I have following configuration

    layers: [
      new ol.layer.Tile({
        source: new ol.source.OSM(),
        opacity: 0.7
      }),
      new ol.layer.Tile({
        opacity: 0.7,
        source: new ol.source.WMTS({
          url: 'https://geoserver.***.com.au/geoserver/gwc/service/wmts',
          layer: 'layer',
          matrixSet: 'EPSG:4326',
          format: 'image/png',
          projection: projection,
          tileGrid: new ol.tilegrid.WMTS({
            origin: ol.extent.getTopLeft(projectionExtent),
            resolutions: resolutions,
            matrixIds: matrixIds
          }),
          style: 'default',
          wrapX: true
        })
      })
    ],

No tiles are displayed in the browser and the requests fails as they look like this:

https://geoserver.***.com.au/geoserver/gwc/service/wmts?layer=layer&style=default&tilematrixset=EPSG:4326&Service=WMTS&Request=GetTile&Version=1.0.0&Format=image/png&TileMatrix=6&TileCol=58&TileRow=18

Notice the TileMatrix parameters. It should be EPSG:4326:5 instead of just 5.

Any clue where I have to configure the prefix?


texts: Edits, When Adding to Geodatabase, Revert to Default in ArcMap

texts:
ArcMap 10.3.1

I have several shapefiles that I am attempting to edit in the following way.

1) Right click on shapefile -> properties -> add hyperlink capability to a field in the attribute table -> hide first two fields of data -> rename the alias of one field.

2) Change the color of the "point" icon to be unique for each shapefile.

Then I create a new File Geodatabase -> name it -> right click -> import ->  Feature Class (multiple) -> select all the shapefiles -> select my directory and save.

The problem is that when I open a blank ArcMap and connect to the new file geodatabase all the shapefiles maintain none of the edits I had made to their icon or properties.

I am attempting to send the geodatabase to a client so they can open it and be able to access the hyperlink and see the desired fields without having to go through and hide the superfluous field again or without having to tell the hyperlink field to accept hyperlinks again.

Is there a way to also export your changes to the data when creating a database?


texts: How can the users of my map add their favorite places

texts:
I want to create a map where the users of my website can add their favorite touristic places by clicking on the map and add some attributes.
I don't mind if it's done by Leaflet openlayers or any technology

Edit : I'm still learning but i think i found something about roles and permissions in geoserver that allow only registred users to make changes and edit the map

Here's an example with Mapserver : https://youtu.be/kxshRcjjC_8?t=2m1s


texts: Is there a way to programatically force RGB (3 bands) in a raster with arcpy?

texts:
I'd like to include conversion to a 3 band RGB raster in a script to prepare multiple DEMs for image processing. The input files are rasters extracted from a raster mosaic.

The manual option is to use the Export Data option in the TOC, but I have to repeat this for 15-25 different rasters so I'd like to automate it.

Each input raster is a single band that I'd like to force RGB like the image below (maunal method).


texts: Getting time-zones from latitude and longitude in raster using R?

texts:
I have a raster with longitudes and latitudes for the US with meteorological data (for about 100k grid points). I would like to get the time zone for each point in order to transform the time to local time. Here is how my raster look like:

> raster
class       : SpatialPointsDataFrame
features    : 95025
extent      : -179.9979, 179.9918, 1, 85.33336  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0
variables   : 1
names       : narr.a_221_20141002_1200_000.sub
min values  :                                0
max values  :                              100

This problem is somewhat similar to https://stackoverflow.com/questions/23414340/convert-to-local-time-zone-using-latitude-and-longitude, except that I have lots of points. The suggested solution there does not work as the Google API accepts only limited number of requests per day and I do not think that the brute force solution would be the most efficient.

I am an absolute beginner in geographical data, but familiar with R, so I need a solution which works in R.


texts: osm2po vertices table, what does ref_count stands for?

texts:
i am using osm2po to import osm file, into postgres/gis.

after de commenting in osm2po config file:

#postp.2.class = de.cm.osm2po.plugins.PgVertexWriter

i get a nice vertices table, with turn restrictions.

i also get a ref_count column for each vertex id.

i am trying to understand what it stands for.

initially i thought its how many times this vertex exist in the table, but quickly i discovered its untrue.

i will appreciate any help/insights as to ref_count column and its meaning.

specifically, i am trying to understand, given a target vertex, how many options for turning i have. (how many vertex and which are accessible via this vertex)

example:
if im on vertex id 111, -> can access vertex id 112 and 113

i thought ref_count might help with that. (just the number of possible vertex's will be useful for me, and yes i can count it easily myself, but im curious about ref_count)

thanks for your help!


texts: Map Tile server for backend Java PDF renderer

texts:
We have a Java Application deployed on Tomcat.  There is a PostgreSQL database with PostGIS extension enabled.  We have both geometries and geo-points in the database.

We have a backend thread that generates PDF reports asynchronously (when requested).  We now want to add maps to the PDFs (based on the geometries) with geo-points plotted.

How can we get access to the necessary map tiles for the PDF from Java?  There is always Internet access from Tomcat, but, obviously, no browser or JavaScript interpreter available.

Added
We have found something that suggests installing our own OSM Map Tile Server and using that.  Are there any other options that do not require us to install our own Map Tile Server?


texts: How to save a intersection between layers in a shp file?

texts:
I have calculated the intersection between polygons in different layers, I got it, and now, I need to save the intersection in a new shp file. I have tried to use QgsVectorFileWriter(), however, I haven´t had good results. It´s difficult to prove because the layers are enormous and everytime I try to run my script takes around 15 minutes. So I want to know if someone knows how can I save this file in a shp file?

This is my code,

# -*- coding:utf-8 -*-
import timeit
start = timeit.default_timer()

layer1=None
for lyr in QgsMapLayerRegistry.instance().mapLayers().values():
    if lyr.name() == "ENTIDAD_1":
    #if lyr.source() == "/path/to/shapefile.shp":
        layer1 = lyr
        break

layer2=None
for lyr in QgsMapLayerRegistry.instance().mapLayers().values():
    if lyr.name() == "ENTIDAD":
    #if lyr.source() == "/path/to/shapefile.shp":
        layer2 = lyr
        break

print layer1.name()
print layer2.name()
selections=[] #declares it is a list
intersecciones=[]

interseccion=None
for f in layer1.getFeatures():
    for a in layer2.getFeatures():
        if a.geometry().intersects(f.geometry()):
            intersection = a.geometry().intersection(f.geometry())# I need to save in a shp this intersection
            print intersection.exportToWkt()
            intersecciones.append([f[2],a[2]])

stop = timeit.default_timer()
print stop - start


texts: Adding sequence numbering based on attributes in table?

texts:
Is there a way to auto populate a sequence of numbers based on another attribute in a table?

For example, if field "Name" has a sequence of numbers from 1 to 10, and I want to number all the rows in the table that associated to "Name = 1" from 1 to 10, and then start this process over when "Name = 2" etc (see attached graphic).

Can this be done within the Field Calculator, with a plugin, or does it need to be coded?


texts: I am using geoserver 2.11.2 and tomcat 8.5. I have some french accents in shapefile. The french accents looks fine when viewing the attributes in QGIS. But it appears as special characters when viewing in geoserver. Is there any configuration i need to do ?


texts: OpenLayers plugin for QGIS problem with maps loading

texts: I am new in GIS. I use QGIS 2.0.1 with OpenLayers plugin downloaded from official repository. When I click to Add OpenStreetMap layer the map is not showing and load process is never ending. What is the problem and how can I solve it ?


texts:
I use 2 markers : Icon and Tooltip.
The problem is that given the sames Lon & Lat they are not at the same position.
The Icon position is right, the tooltip position not.
This is my js code for the 2 markers
//Tooltip
var marker = new L.marker([getUrlParameter("Lat0"), getUrlParameter("Lon0")], {
  opacity: 0,
});
marker.bindTooltip(getUrlParameter("D0") + "km", {
  permanent: true,
  offset: [0, 0],
});
marker.addTo(map);
// Icon
var Marker = L.icon({
  iconUrl: "Point.png",

  iconSize: [27, 32], // size of the icon
  iconAnchor: [13, 32], // point of the icon which will correspond to marker's location
  popupAnchor: [0, -30], // point from which the popup should open relative to the iconAnchor
});
var Lon = getUrlParameter("Lon");
var Lat = getUrlParameter("Lat");
var starts = new L.LayerGroup();
L.marker([Lat, Lon], { icon: Marker })
  .bindPopup(Lat + "<br>" + Lon)
  .addTo(starts);

and the result : the tooltip is on the left or on the right depending on the point position on the screen but shifted in both cases (Lat and Lon are the sames)

Edit : the tooltip is binded to the L.marker.
When the marker is on the right of the middle of the screen the tooltip is on the left of the marker (last sample) by 15 pixels (and 25 pixels above), so I've set an offset of [15,25], it's now all right. But When the marker is on the left of the middle of the screen the tooltip is on the right of the marker (before last sample) by 15 pixels (and 25 pixels above), so with an offset of [15,25] it is 30 pixels to the right ! The Xoffset should be +15 or -15 depending on its position on the screen : how to do that ?


texts:
When I use ArcScene 10.4 I can easily export 3D maps (VRML in *.wrl) which I then convert in FME to a *.prc file that can be imported in Adobe Acrobat PDF files (to create interactive 3D maps, readers can navigate, zoom, change colors, etc.)

However when I'm working with ArcGIS Pro, I can't find out how to do this. ArcGIS Pro works perfectly for visualizing 3D and it has many 3D tools and options, but the export of 3D features I haven't found (yet).

How do I export a 3D map from ArcGIS Pro?

(File types *.wrl, *.u3d or *.prc are all fine)


texts: How do I draw multiple copies of a rectangle in QGIS?

texts:
I am trying to create a map of our parish cemetery. I almost have a file with the centre point of each headstone plotted and now want to add a rectangle to each point to lay out the grave positions of both existing and future plots.

I understand that GRASS has a module v.mkgrid but this appears to want to create a uniform grid with uniform spacing where my real world cemetery does not have graves in straight lines and equal distances apart.

A previous answer to a similar question suggested the CAD Tools plugin or Rectangles ovals digitising plugin neither of which I can find using Plugin Manager

I am using QGIS 1.8 Lisboa on a PC.

Two question then:-

Does QGIS have what I need to do this ?
How do I find the 'missing' plugins ?


texts: How can I merge 8-bit RGB tiff files using QGIS

texts:
When I try to merge my tiff files without translating (converting format), my color space becomes grey.
If I try to translate my tiff, I am getting this error:
Error : band 1 has no color table

I want to retain my color space after merging.

This is the gdalinfo of one of tiff image:

Driver: GTiff/GeoTIFF
Files: Muna_20160905_WV3_R2C1.tif
Size is 15957, 10018
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
Origin = (13.160409776871793,11.860700450062284)
Pixel Size = (0.000003210135867,-0.000003210074723)
Metadata:
  AREA_OR_POINT=Area
  TIFFTAG_COPYRIGHT=Image Copyright 2016 DigitalGlobe Inc
  TIFFTAG_RESOLUTIONUNIT=1 (unitless)
  TIFFTAG_XRESOLUTION=1
  TIFFTAG_YRESOLUTION=1
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  13.1604098,  11.8607005) ( 13d 9'37.48"E, 11d51'38.52"N)
Lower Left  (  13.1604098,  11.8285419) ( 13d 9'37.48"E, 11d49'42.75"N)
Upper Right (  13.2116339,  11.8607005) ( 13d12'41.88"E, 11d51'38.52"N)
Lower Right (  13.2116339,  11.8285419) ( 13d12'41.88"E, 11d49'42.75"N)
Center      (  13.1860218,  11.8446212) ( 13d11' 9.68"E, 11d50'40.64"N)
Band 1 Block=512x512 Type=Byte, ColorInterp=Red
Band 2 Block=512x512 Type=Byte, ColorInterp=Green
Band 3 Block=512x512 Type=Byte, ColorInterp=Blue

If I perform identify, I get this:

Muna_20160905_WV3_R2C1.tif TIFF 15957x10018 15957x10018+0+0 8-bit sRGB 503.3MB 0.000u 0:00.010

The provider also tells in description:

These are "8-bit TIFs which ‘tiles’ the image"

Any help?


texts: comparing parcel ids from two tables with arcpy.da.SearchCursor

texts:
My goal is to compare parcel IDs from two tables, one is a standalone shapefile and the other is a feature class inside a geodatabase (SQL Server). Both tables field name is "PARID" and the values are text in the table.

I would like this code to look for differences in list 1 from list 2 and then write the different parcels in a clean formatted text file.

however currently the text file output looks like this:

I am facing few issues here, if I don't set the parcelList.append to (str(value)) I get the error:

TypeError: unhashable type: 'list'

but when I append the value as a string the list is only 1 item that contains all different values I compared, rather than having each parcelID its own item in the list.

My code is below:

import arcpy
#import os
import pprint

arcpy.env.workspace = r"C:\Users\path\tosde.sde"
legacyTable = "Database Connections\\database\_Legacy_Info"
NewParcelTable = r"\newParcels.shp"
legacyParcelList = []
NewParcelList = []
differentPCN = []
sql = "PARID LIKE '38%'"

print "generating list of existing legacy PCN's..."
with arcpy.da.SearchCursor(legacyTable,["Parid"],sql) as cursor:
    uniquePCN = sorted({row[0] for row in cursor})
    legacyParcelList.append(str(uniquePCN))
    #print uniquePCN
print "generating list of new PCN's..."
with arcpy.da.SearchCursor(NewParcelTable,["Parid"],sql) as cursor:
    uniquePCN2 = sorted({row[0] for row in cursor})
    NewParcelList.append(str(uniquePCN2))
    #print uniquePCN2

difference = set(NewParcelList).difference(legacyParcelList)
for item in difference:
    differentPCN.append(item)
print "the values below are in the New Parcel Shapefile list but not inside the Legacy Parceel List: " +'\n'
print "Writing output to text file..."
for item in differentPCN:
    pprint.pprint(item)
txtFile = r"C:\Users\differentPCNs.txt"
x = 0
try:
    for item in differentPCN:
        x +=  1
        f = open(txtFile,"w")
        f.write((str(x)) + ". "+ (item) + "\n")

finally:
    f.close()
print "file output is complete you can find it here: " + txtFile

I believe the issue is with the original search cursor and appending to new lists but not sure how to fix this.


texts: How can one modify map service legend (specifically raster layer)?

texts:
I've got the raster layer in my map service, which in fact is created automatically by asynchronous GP service as its output result ("View results with a map service" option). The raster is symbolized with graduated colors and then the whole map is exported along with its legend to PDF using printing service.

The thing is that the symbol for the legend, according to rest output and flex app layer list widget, looks like three separated html lines with Hi and Low values as labels next to the first and the last square (yes, symbol for the raster actually consist of three squares). It looks ridiculous in the TOC and on the printout.

Layout in PDF is even worse. It doesn't show the lowest square, probably due to transparency options in map (PNG32 is fine though).

The question is: Can I change its look somehow? I'd like to rename or even get rid of labels (I need it in Polish) and I need one rectangle instead of three split boxes.

I'd like to customize the TOC in the flex app to show "normal" raster symbol at least. Recently I even wanted desperately to remove symbol for this layer from layerlist widget in flexapp but I haven't figured out how to do that so far.


texts:
I'm looking to move 119 GB of data, mxd's, and other files from my computers C drive to "documents" (which I just found out is the only location that gets backed up to our server).

How can I do this without needing to re-set all data sources for all layers within all of my mxds? Is there a way?

It would only be the beginning 2 locations of the pathname that would need to be changed because I will keep all folder locations and titles the same within the "documents" location.


texts: Is it possible to append the same type of shapefile captured on two different Trimble Juno 3b's with ArcPad 10?

texts:
I will soon have access to two Trimble Juno 3b's which will both have ArcPad 10. I have done some practice with one Juno 3b and I am very confident with my ability to use them. All I need to know is that if I have the same type of shapefiles, in this case points, and load them into each on of the Trimble Juno's, can I combine them?

I plan to hand one to another person and use the other myself to collect points on different areas of a field. Can I run some sort of tool or even a Python script in ArcMap 10.2.2 to combine both shapefiles so that it creates a single shapefile with the same attributes and spatial extent? I will need them to be combined so that I can use them all in Kriging interpolation.

I believe the union tool is a possibility, but can it be used with points collected in this manner?


texts:
Is it possible to make a style for the line, with arrows in the direction of digitized lines?

No coding after digitization. This is possible when I explode lines to a segments, and I made azimuth field, and then rotate the arrow simbol with azimuth, but I do not want to explode lines.

This screenshot is from another GIS program.


texts: How to convert DTED Level 1 into SRTM DTED

texts:
We have a GIS tool that handles only SRTM DTED level 1 terrain files. SRTM level 1 terrain files are limited to 60 deg lat. So I need to fill up the gap from 60 lat up to 90 deg lat.
I have DTED level 1 data covering the area 60 deg lat and upwards that could be used for that.
Is there a way to convert that DTED level 1 data into SRTM data?
Very much appreciated,

Yvan


texts: Rendering overlapping polygons without holes in multipolygon in Leaflet

texts:
Is there a way to render overlapping polygons without holes in multipolygon in Leaflet?
I found Leaflet has an open issue.
codesandbox.io
I am asking because I have a bound polygon that contains inner polygons.

I would need to drag them all together, but it doesn't work in canvas mode. It could work if I can render this as a multipolygon. This is a follow-up question on my previous question and cross-post from SO.


texts: Open source GIS engine that is similiar to ArcGIS Engine?

texts: If I was to develop a mapping application to be a standalone windows app, is there an open source GIS engine that is somewhat comparable to ArcGIS Engine?


texts: I am new to SQLite and when I try to import a raster image in a SQlite database that I created using Spatialite_gui, I get this message rasterlite_load: command not found after I write the importing line in cmd.exe. What is the issue and how do I fix it?


texts: Making categorized symbols mutually exclusive using PyQGIS

texts: Hello, I am wondering if there is a way in PyQGIS to automate the process of selecting/deselecting the check marks for each categorized symbol in the included picture. I am writing a plugin and would like to be able to have the option of viewing one of the symbols at a time instead of seeing all of them at once. Thanks


texts: Difference between classification methods with Laspy

texts:
I need to add classification to create a LAS file. I used the code below but it was miss-classified. Is there any problem with the code? And what is the difference between:

out.file classification

and

 outfile.Raw_Classification

      hdr = laspy.header.Header()
      outfile = laspy.file.File("output.las", mode="w", header=hdr)
      allx = scanne[:,0]
      ally = scanne[:,1]
      allz = scanne[:,2]

      xmin =np.min(allx)
      ymin =np.min(ally)
      zmin =np.min(allz)

      outfile.header.offset = [xmin,ymin,zmin]
      outfile.header.scale = [0.001,0.001,0.001]
      outfile.x = allx
      outfile.y = ally
      outfile.z = allz
      outfile.Raw_Classification  = scanne[:,3]
      print(hdr)
      outfile.close()


texts: Finding maximum distance point can be from polygon using PostGIS?

texts:
I have a hydrography dataset containing navigable waterways: ponds, lakes, rivers, and streams. I need to find the maximum distance ANY location within the state can be from a navigable waterway. I'm familiar with the st_geometry predicates but I'm not sure how to attack this problem.

I do not have a table of predefined point geometries. The ANY location is the set of X,Y real numbers contained within the state.

I have a collection of polygons and lines, all contained within a much larger region (another polygon) and need to find the X,Y of the most remote point and the distance from it to the nearest waterway.

How do I proceed with figuring out the maximum distance any point can lie from any navigable waterway?

The dataset is large and the geographic area covers an entire state. I don't expect an instant answer and this will not be part of an on-the-fly solution. For development purposes, I have a subset of the dataset covering only a few hundred acres with limited hydrographic features

Postgres 9.6.1 / PostGIS 2.3.1


texts: Getting an sf boundary of a raster or stars object

texts:
I am trying to figure out the best approach to convert an raster or stars object into its boundary returned as an sf object. I think I have a pretty decent handle on how to do it with raster (shown below):
library(sf)
#> Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 6.3.1
library(stars)
#> Loading required package: abind
library(raster)
#> Loading required package: sp

## path
test <- system.file("external/test.grd", package = "raster")

## Using raster
x <- raster(test)
r <- x > -Inf
pp <- rasterToPolygons(r, dissolve = TRUE)
#> Loading required namespace: rgeos
boundary_x <- st_as_sf(pp)
plot(boundary_x)

Here using stars i am forced to rely on st_union for that final aggregation step. For this example, it is pretty quick but for larger examples I first have to create many little polygons (using st_as_sf.stars then use st_union) which takes a long time. Is there a more efficient way of doing this?
## Using stars
y <- read_stars(test)
boundary_y <- st_as_sf(y, merge = TRUE) %>%
  st_union()
plot(boundary_y)


texts: How to create Polygons of a specific size?

texts:
I'm trying to create polygons of specific size in square meters that represent the areas surrounded by roads. The roads are stored in PostGIS as points / lines. If required, I can convert the data to other formats for processing in other softwares like ArcGIS etc.

Because of the randomness of the roads, the polygons could differ in size but their size should not be less than the specified size.

See the below picture. For illustration, I have set different colors to the possible polygons of equal size fulfilling minimum size constraint.

Can somebody please mention how can I achieve this? As I am beginner so sample query or code would greatly help. I have basic know how of PostGIS, QGIS and ArcGIS but I am willing to use any other software if it solves my problem.

If the street network looks like below, polygons like given in red needs to be created (only as an example, in fact all the intersections between streets fulfilling the minimum area constraint should produce polygons like these).

Let suppose we want to create polygons of minimum 2 square kilometer area and let suppose that the area made up of 3 white patches shown in blue is of 2 square kilometer then the code should create a polygon consisting these 3 white patches. You can see that polygon 3 is larger than the polygon 1 but fulfilling the minimum size constraint, similarly polygon 2 has four white patches to fulfill the minimum size constraint.


texts: What ArcGIS geoprocessing tools support exporting tabular data to CSV?

texts:
I have found that most if not all geoprocessing tools that allow tabular output are limited to geodatabase tables and DBF files for supported output formats. Are there any that also support CSV output? Numerous tools supports CSV input, such as Table to Table, Copy Rows, etc.

The only exception I have seen is a script tool in the Spatial Statistics toolbox called Export Feature Attribute to ASCII, but this only supports feature classes, not tables, and it chokes on non-ASCII characters:

<type 'exceptions.UnicodeEncodeError'>: 'ascii' codec can't encode character u'\xa0' in position 213: ordinal not in range(128)`

Edit: Thanks for the alternatives in the answers so far. I actually do use some of those those already myself, but just to be clear my question is specifically about the built-in geoprocessing tools in ArcGIS 10. Alternatives are welcome but I am more specifically wondering why the capability to export to CSV exists interactively in ArcMap but apparently not through the GP framework (without Python scripting).

Incidentally, I trawled around the ArcGIS Ideas site a bit and did come upon a few existing ideas that may be of interest to some:

Add Text Option to Table to Table
Model Builder only Tool - export to csv
Export a table directly to an MS Excel file (not specifically about CSV but has a lot of relevant discussion


texts: QGIS: Converting a { "DisplayValue": "Value being displayed", "Value": "value" } string format from a GeoJSON file

texts:
Currently we are trying to to convert a JSON file, which has an attribute table that is too chaotic and needs to be changed. There is certain string format in a few columns that keeps popping up and we want to covert this into one go.

It always looks like this { "DisplayValue": "Value being displayed", "Value": "value" }.

Note: Value being displayed" and value are placeholders, and yes there are encapsulated into "".
What we suspected is that these things were supposed to be part of a value map...kinda because there a few columns, that are also supposed to be a date.

For example { "DisplayValue": "2020-03-05T16:37:27.7Z", "Value": "2020-03-05T16:37:27.7Z" }

Also we have something like this: { "diepteNauwkeurigheid": "Onbekend", "dieptePeil": "1 m" } that is something that has the same string format, but this does not seem to be a value map.

What actually needs to be extracted is the value, without the special marks.

I have bit of trouble to convert this in one go. Sure I can make some replace and right left commands but can anyone tell me how do this in one go?

Edit: TO clarify the files we have are GeoJSON. And no, it is not just  one single attribute we want, we want multiple columns to have the same kind of string format converted. Here is screenshot, some data is cut because of privacy.

Edit 2: To Illustrate further I have for example the procedure to translate the the theme column into the appriate string:

CASE
WHEN theme= '{ "DisplayValue": "Elektriciteit", "Value": "electricity" }'
THEN
'electricity'

WHEN theme= '{ "DisplayValue": "Telecommunicatie", "Value": "telecommunications" }'
THEN
'telecommunications'

WHEN theme ='{ "DisplayValue": "Water", "Value": "water" }'
THEN
'water'

WHEN theme ='{ "DisplayValue": "Riolering", "Value": "sewer" }'
THEN
'sewer'

Else
NULL

END

Now this is just for translating the text of ONE Column specifically, I am looking for something more specific that I can apply to all columns.

So it should be something like this, this is just to illustrate what must happen, actual value and Displayed Value are placeholders.

CASE WHEN COLUMN = '{ "DisplayValue": "Displayed Value", "Value": " actual value" }'
THEN
actual value

I hope this helps in helping me find a solution.

Edit 3: Leo M gave me a good solution but there is sadly a problem that has arrived. I get this error message, not all values also seem to be translated.

After the okay button has be clicked my changes has been undone.

Edit 5: This is also the error message in the Field Calculator itself, yes it is on a different column but I had to test it out:

Edit 6: Thanks to LeoM the sollution has been found, apparently because I had Null Value's, the expresion did not work so I made a part of a case: For example

Case when ElementTyp IS NOT NULL THEN
string_to_array( "elementTyp" ,'"')[7]
ELSE
ElementTyp = NULL
END

You can replace Null with anything.


texts: How to remove a domain from a field using ArcObjects

texts:
I need to delete a domain in my database. Since the domain is still associated with a field, I first need to delete the association. After reading the documentation for IWorkspaceDomains.DeleteDomain, it say to use the IClassSchemaEdit.AlterDomain.

So I try using this code to remove the association

    Dim workspaceDomains As IWorkspaceDomains = CType(workspace, IWorkspaceDomains)
    Dim feature As IFeatureClass
    Dim schemaLock As ISchemaLock = Nothing
    Dim classSchemaEdit As IClassSchemaEdit

    feature = CType(workspace, IFeatureWorkspace).OpenFeatureClass("FeatureClassName")

    Try
        schemaLock = CType(feature, ISchemaLock)

        schemaLock.ChangeSchemaLock(esriSchemaLock.esriExclusiveSchemaLock)
        classSchemaEdit = CType(feature, IClassSchemaEdit)

        classSchemaEdit.AlterDomain("FieldName", Nothing)
    Finally
        If (Not schemaLock Is Nothing) Then
            schemaLock.ChangeSchemaLock(esriSchemaLock.esriSharedSchemaLock)
        End If
    End Try

When it comes to executing this line classSchemaEdit.AlterDomain("FieldName", Nothing), I get this error :

  A first chance exception of type 'System.Runtime.InteropServices.COMException' occurred in DLLFile.dll

  Additional information: Error HRESULT E_FAIL has been returned from a call to a COM component.

So I when back and check the documentation for IClassSchemaEdit.AlterDomain to see how to remove the association, but I found nothing.

What is it that I'm doing wrong. I know I could just copy the data to a new field and delete the original, but this mean that I would have to modify a toolbox, but I would rather avoid it.


texts: Error when training data for Deep Learning

texts:
I am trying to train a model and classify it using the "Classify Pixels Using Deep Learning" tool. I exported my training sample data for deep learning, but when I add that .emd-file to the tool, I get an error that reads "Cannot access file %s." I saved the .emd-file to the main folder of my project. I have write access to this location and tried saving and closing the software.

I am classifying vegetation, soil, and water classes. I used the rectangle training sample polygon shape. Would this impact the classification if I am depending on classifying pixels and textures, but not shapes?

This is with ArcGIS Pro 2.6.
Has anyone encountered this before? I am not entirely sure what I am doing wrong.


texts: Using wildcard for field names of ArcPy Search Cursor?

texts:
How would you wildcard the field_names for a SearchCursor so that it only searched through rows for columns in a table that begin with the word "Layer"? For example:

FCList = []
FCS = [u'Building', u'Bridges', u'RecreationArea']
LCC = "G:\\CL.gdb\\Layer_CatCode"
with arcpy.da.SearchCursor(LCC, ['Layer*']) as cursor:
    for row in cursor:
        if row[0] in FCS and row[0] not in FCList:
            FCList.append(row[0])

I have tried to construct a sql statement for the field_names with no success, and I have even tried the asterisk option for field_names so that is just searches all columns but the FCList comes out empty. However, this may be due to the fact that when using the asterisk option I'm not sure what to put in the brackets for row["*"].


texts: Imported TIF has dark gray background instead of white

texts:
using QGIS 2.6

I import a TIF file, and expect the import to show an image that has black lines on a white background (with some other colors here and there), as the native windows viewer will see a white background.   The file, 49MB, imports, but

The TIF image has a dark gray background instead of a white background

I can't tell why.   Any hints?

The same thing happens if I import the TIF image into georeference, which is how I discovered this in the first place.

I created/exported the TIF image using PDFCreator (from a print of a PDF to a TIF, using that program).   Within PDFCreator, I chose a TIF export that set the color map to 4096 colors, as the default created a TIF file that was twice as large.

On testing with the default size (16 million colors or some such)....
That file imported and showed the background as white.

And I can't tell why, for that, either.   Any hints?

========
In response to the comment/question by @BradHards:

The TIFFTAG_PHOTOMETRIC is not present in either file.

The gdalinfo for the large file and small file are identical, but for the presence of two extra lines on the file where the white background turns to dark gray.  Each band sets NBITS=4 for Image Structure Metadata: eg:

Band 1 Block=5100x68 Type=Byte, ColorInterp=Red
  Image Structure Metadata:
    NBITS=4
Band 2 Block=5100x68 Type=Byte, ColorInterp=Green
  Image Structure Metadata:
    NBITS=4
Band 3 Block=5100x68 Type=Byte, ColorInterp=Blue
  Image Structure Metadata:
    NBITS=4

Something about the NBITS=4 state turns the background dark gray.  Alternatively, there is some other problem with the PDFCreator export.

These are the full GDALInfo:

Small with dark gray background:

Driver: GTiff/GeoTIFF
Files: small file.tif
Size is 5100, 6600
Coordinate System is `'
Metadata:
  TIFFTAG_DATETIME=2015:01:28 16:41:04
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  TIFFTAG_SOFTWARE=GPL Ghostscript 9.07
  TIFFTAG_XRESOLUTION=600
  TIFFTAG_YRESOLUTION=600
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 6600.0)
Upper Right ( 5100.0,    0.0)
Lower Right ( 5100.0, 6600.0)
Center      ( 2550.0, 3300.0)
Band 1 Block=5100x68 Type=Byte, ColorInterp=Red
  Image Structure Metadata:
    NBITS=4
Band 2 Block=5100x68 Type=Byte, ColorInterp=Green
  Image Structure Metadata:
    NBITS=4
Band 3 Block=5100x68 Type=Byte, ColorInterp=Blue
  Image Structure Metadata:
    NBITS=4

Large with white background

Driver: GTiff/GeoTIFF
Files: larger file.tif
Size is 5100, 6600
Coordinate System is `'
Metadata:
  TIFFTAG_DATETIME=2015:01:28 17:07:16
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  TIFFTAG_SOFTWARE=GPL Ghostscript 9.07
  TIFFTAG_XRESOLUTION=600
  TIFFTAG_YRESOLUTION=600
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 6600.0)
Upper Right ( 5100.0,    0.0)
Lower Right ( 5100.0, 6600.0)
Center      ( 2550.0, 3300.0)
Band 1 Block=5100x68 Type=Byte, ColorInterp=Red
Band 2 Block=5100x68 Type=Byte, ColorInterp=Green
Band 3 Block=5100x68 Type=Byte, ColorInterp=Blue


texts:
I'm blocked with creation of an iterator in my script.
I want coordinates of selected points.
If I do that, it works, it print all points coordinates

##point=vector
layer = processing.getObject(point)
it = layer.getFeatures()
for i in it:
    x = i.geometry().asPoint().x()
    y = i.geometry().asPoint().y()
    print "X Coord :  %f" %x
    print "Y Coord :  %f" %y
    print

So after I add my selection. But it crash.

##point=vector
layer = processing.getObject(point)
nbr = QgsExpression( "ID LIKE '%-1'" )
itnbr = layer.getFeatures( QgsFeatureRequest( nbr ) )
for i in itnbr:
    x = i.geometry().asPoint().x()
    y = i.geometry().asPoint().y()
    print "X Coord :  %f" %x
    print "Y Coord :  %f" %y
    print

I have an OGR error message : invalid data source

What causes that?


texts: Mysterious Lat Long conversion to XY problem using QGIS

texts:
I am trying to import a CSV file with Lat Long values I geocoded in Google Earth into QGIS, and I'm getting some very strange results.

A well-formatted decimal-degree Lat/Long pair (like 30.402716, -88.867475) gets translated into a tiny tiny XY pair like (-0.0008,0.0003). The points show up basically in the right position relative to each other, but these decimal values are way off relative to my other layers, why is this decimal conversion happening and why is it so wrong?

I'm using the "Create Layer from Delimited Text" dialogue, and I set X to Long, and Y to Lat, and selected the CRS to Google Maps Global Mercator EPSG:900913


texts: Switch between heatmap and single symbol renderer at different zoom level in QGIS?

texts:
On a point layer in QGIS, I would like to display the points as normal single symbol points up to 1:5000 and switch to heatmap map between 1:5000 and 1:1000.

Is it possible? I can't seem to find this option in the Property > Style panel of the layer.
The problem I have at the moment is that the heatmap layer is almost invisible on a wide scale. I would like to avoid doubling layers.