texts: CartoDB torque chart stop from playing on page load

texts: Is there a way to use a CartoDB torque chart such that the chart does not automatically play upon page load?


texts: Network Analyst: Location Allocation for Model Builder

texts:
I am trying to create a model to answer a question using a network analysis.

Background:

10 possible fire station locations (shapefile)
Calls for demand points (shapefile)
Impedance is 5 minutes from facility to demand

Product: I need to make a model to select 1-3 stations from the list of 10 for the best locations to meet as many calls as possible with a 5 minute drive from the station.

Below is an image of what I have so far, but when I try and run it I get an Error saying:

Error 000735: Input Network Analysis Layer: Value is required
WARNING 030036: A layer with this name exists

Not sure how to fix these two items.


texts:
When I have vector layer with several objects and label them, can I change the text of individual labels? I have tried the labeling tools and the tool change label but the field with the label text is grayed out.

Is it possible to change individual label texts?


texts: How to snap lines to points automatically in qgis using python

texts:
I have two vector shapefiles: One file is points representing manholes and the other is a polyline file representing pipes leading to and from the manholes. I need to snap the lines within a specified distance to the closest manhole(centroid of point manhole).

How can i do it in qgis using python. I have found the solution in sqlite or postgis but i want to do it in qgis.

SELECT
    f.gid as gid,
    ST_Snap(f.Geometry, g.Geometry, 2) as geom
FROM
    pipe as f,
    (SELECT ST_Collect(Geometry) as Geometry
     FROM mh) as g


texts: Getting db manager error from QGIS 3 on Mac?

texts:
I've got a new QGIS 3 install from Kingchaos on osx, but have few issues.
One of them is with db manager and sqlite database.

This is the error I got:

Une erreur est survenue lors de l'éxécution du code Python:

AttributeError: 'sqlite3.Connection' object has no attribute 'enable_load_extension'
Traceback (most recent call last):
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_model.py", line 441, in rowCount
    self._refreshIndex(parent, True)
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_model.py", line 490, in _refreshIndex
    if item.populate():
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_model.py", line 167, in populate
    if not connection.connect():
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_plugins/spatialite/plugin.py", line 80, in connect
    return self.connectToUri(uri)
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_plugins/plugin.py", line 109, in connectToUri
    self.db = self.databasesFactory(self, uri)
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_plugins/spatialite/plugin.py", line 66, in databasesFactory
    return SLDatabase(connection, uri)
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_plugins/spatialite/plugin.py", line 109, in __init__
    Database.__init__(self, connection, uri)
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_plugins/plugin.py", line 218, in __init__
    self.connector = self.connectorsFactory(uri)
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_plugins/spatialite/plugin.py", line 112, in connectorsFactory
    return SpatiaLiteDBConnector(uri)
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins/db_manager/db_plugins/spatialite/connector.py", line 51, in __init__
    self.connection = spatialite_connect(self._connectionInfo())
  File "/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/qgis/utils.py", line 596, in spatialite_connect
    con.enable_load_extension(True)
AttributeError: 'sqlite3.Connection' object has no attribute 'enable_load_extension'

Version de Python : 3.6.4 (v3.6.4:d48ecebad5, Dec 18 2017, 21:07:28) [GCC 4.2.1 (Apple Inc. build 5666) (dot 3)]
Version de QGIS : 3.0.0-Girona Girona, exported

Chemin Python :
/Applications/QGIS 3.app/Contents/MacOS/../Resources/python
/Users/arthur/Library/Application Support/QGIS/QGIS3/profiles/default/python
/Users/arthur/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins
/Applications/QGIS 3.app/Contents/MacOS/../Resources/python/plugins
/Library/Frameworks/GEOS.framework/Versions/3B/Python/3.6/site-packages
/Library/Frameworks/GDAL.framework/Versions/2.2/Python/3.6/site-packages
/Library/Frameworks/Python.framework/Versions/3.6/lib/python36.zip
/Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6
/Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6/lib-dynload
/Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6/site-packages
/Users/arthur/Library/Application Support/QGIS/QGIS3/profiles/default/python
/Users/arthur/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/mmqgis/forms


texts:
Consider this question here. Many servers do not set CORS headers correctly. As I recover here, even esri does not do so, but bing does.

That begs for the question, what are some WMS servers that have their CORS headers properly?


texts: How to remove points on the Openlayers map?

texts:
I am experimenting with openlayers. I am adding points to the map using this example. However, I also want to delete unwanted points from the map. For example, I clicked 6 points on the map and I want to remove the fifth point, but keep 1, 2, 3, 4, and 6 th point on the map. How can I do this in openLayers. So far, I tried this example and this example, but I had no luck so far. Please help me out.

Update: I have the following JsFiddle shows the code to add points to my map. I want to remove one point, but keep the rest of them on the map.

var raster = new ol.layer.Tile({
    source: new ol.source.OSM()
  });

  var map = new ol.Map({
    layers: [raster],
    target: 'map',
    view: new ol.View({
      center: [-11000000, 4600000],
      zoom: 4
    })
  });

  var features = new ol.Collection();
  var featureOverlay = new ol.layer.Vector({
    source: new ol.source.Vector({features: features}),
    style: new ol.style.Style({
      fill: new ol.style.Fill({
        color: 'rgba(255, 255, 255, 0.2)'
      }),
      stroke: new ol.style.Stroke({
        color: '#ffcc33',
        width: 2
      }),
      image: new ol.style.Circle({
        radius: 7,
        fill: new ol.style.Fill({
          color: '#ffcc33'
        })
      })
    })
  });
  featureOverlay.setMap(map);

  var modify = new ol.interaction.Modify({
    features: features,
    // the SHIFT key must be pressed to delete vertices, so
    // that new vertices can be drawn at the same position
    // of existing vertices
    deleteCondition: function(event) {
      return ol.events.condition.shiftKeyOnly(event) &&
          ol.events.condition.singleClick(event);
    }
  });
  map.addInteraction(modify);

  var draw; // global so we can remove it later
  var typeSelect = document.getElementById('type');

  function addInteraction() {
    draw = new ol.interaction.Draw({
      features: features,
      type: /** @type {ol.geom.GeometryType} */ (typeSelect.value)
    });
    map.addInteraction(draw);
  }

  /**
   * Handle change event.
   */
  typeSelect.onchange = function() {
    map.removeInteraction(draw);
    addInteraction();
  };

  addInteraction();

Update 2: I can remove the last element using removeLastPoint() inside ol.interaction.Draw.


texts: Extract multiple admin levels from planet OSM

texts:
I am trying to create a multi-level choropleth map based on OSM data.

I need the outline/boundaries. I would like to do this from the command line if possible.

Here is what I think the steps should be.

Download the planet OSM in XML or PBF format
Extract the administrative boundaries for features with an admin level of 2, 4, 6 or 8. This is my current attempt:

osmconvert planet-latest.osm.pbf -o=planet.o5m

osmfilter planet.o5m --keep-nodes= --keep-ways-relations="boundary=administrative and (admin_level=2 or admin_level=4 or admin_level=6 or admin_level=8)" > planet.osm
Simplify the geometries to help with the file size without losing track of which features are "inside" other features.
Convert the final product into GeoJSON.

I am struggling with #2 the most and it seems like the most important step.

If any other process I should try, please let me know.


texts: Why does this simple Python OGR code create an empty polygon?

texts:
Please excuse my ignorance...I am new to Python OGR.

I am trying to create a simple square polygon and fail...this produces an empty polygon.

Where did I go wrong?

from osgeo import ogr
driver = ogr.GetDriverByName('ESRI Shapefile')
datasource = driver.CreateDataSource('c:/temp/testPolygon.shp')
layer = datasource.CreateLayer('layerName',geom_type=ogr.wkbPolygon)

#create ring and polygon object:
myRing = ogr.Geometry(type=ogr.wkbLinearRing)
myRing.AddPoint(0.0, 0.0)
myRing.AddPoint(0.0, 10.0)
myRing.AddPoint(10.0, 10.0)
myRing.AddPoint(0.0, 10.0)
myRing.AddPoint(0.0, 0.0)#close ring
myPoly = ogr.Geometry(type=ogr.wkbPolygon)
myPoly.AddGeometry
#this produces myPoly.GetArea() of zero!

#create feature polygon:
feature = ogr.Feature( layer.GetLayerDefn() )
feature.SetGeometry(myPoly)
layer.CreateFeature(feature)

Even the following returns an area of zero:
    from osgeo import ogr

#create polygon object:
myRing = ogr.Geometry(type=ogr.wkbLinearRing)
myRing.AddPoint(0.0, 0.0)#LowerLeft
myRing.AddPoint(0.0, 10.0)#UpperLeft
myRing.AddPoint(10.0, 10.0)#UpperRight
myRing.AddPoint(10.0, 0.0)#Lower Right
myRing.AddPoint(0.0, 0.0)#close ring
myPoly = ogr.Geometry(type=ogr.wkbPolygon)
myPoly.AddGeometry
print (myPoly.GetArea() ) #returns zero, not 100...


texts: Adding images to placemarks multiple times

texts: I have a route that I travel that is saved by coordinates to a .csv file for use in GE.  I have images (pictures) that I would like to add to placemarks I have along the route.  Is there a script or macro I can run that would allow me to add the url of the images to the .csv file so that I could have the images appear inside the placemarks balloon without having to add the images individually to the balloon?


texts: Writing no_data_value within jpeg2000[kdu] with gdal_translate

texts:
I would like to specify a NO_DATA_VALUE within a JPEG2000 file (I use the JP2KAK driver). When I call gdal_translate -of JP2KAK -co QUALITY=100 input.vrt ouptut.jp2 , I get an aux.xml file alongside, that holds the value of the nodata value (initially present at band level in the vrt) in some PAM fields.

I am aware that compression messes with nodata values, but I'm using a lossless quality, so it shouldn't be a problem.

I would like to have the jp2 file independent of any other files (in order to be able to send it from service to service without creating a bundle whatsoever).

Can gdal_translate help me in this matter?

[EDIT : gdal version 2.2]


texts: I am using ArcGIS Pro and I have grouped a number of similar layers together creating a 'group layer'. Is it possible to export the group layer as a shapefile or can I only do this with each of the individual layers within the group? I have looked but cant find an answer.


texts: Height in GNSS in acquiring and postprocessing data

texts:
First I am kind of a newbie using GNSS receivers and processing the data (I am an archaeologist). We acquired a couple of emlid reach rs+ receivers. Our first experiment was to set up a base for PPK postprocessing. We acquired 4 hours of data, corrected with data from our local airport data (via CORS). The postprocessing was made using RTKPOST. We obtained a lat/lon/height coordinate. We know that coordinates are in lat/lon WGS84 and that the heights are ellipsoidal. We know that we can use UNAVCO Geoid height calculator to get orthometric heights. Now is when the discussion begins. We set up a series of points and logged data for an hour with the rover. Base was in the previously calculated point and was set up to acquire data with manually defined coordinates, the original ones with ellipsoidal heights. The questions derived from the discussion with my colleagues are:

When setting up the base in the known point, Do I use the ellipsoidal or the orthometric height?

When postprocessing, Do I use the ellipsoidal or orthometric heights?

I know there are two questions here, but I think they are closely related.
Thanks before hand
Gerardo


texts:
I don't know so much GIS specific terminology, so I am sorry for using the wrong words. I'll try to make my point clear by first describing what I want to achieve and then what I tried.

[goal] I am working with tweets that I have stored in a csv file. Now I want to visualise them on a map. I want to sort them according to their date (or some other column in my file). So I want to have a way to display only tweets from a certain day or date or some other property. Also, because there are many tweets for each property, I want them to be clustered. So I can see from the point's colour how many are at the same spot instead of having to zoom in. I need this to be nicely displayed on a static picture.

[what I did] I loaded my .csv file and categorized it according to date. This worked, so now I have a layer with all tweets but I can switch on/off each date. They are not stored in a different layer, but as something else. I don't know how to call this. Then I tried the Heatmap Plugin, which I thought is what I need. But this only seems to work on normal layers. So I cannot visualise my categorized data in different parts, just the full layer. I also tried to categorise the heatmap somehow, but couldn't find out how.

I hope you can help me to achieve my goal. Of course I could create different .csv files for all the properties that I want to sort my data for, but that seems stupid. I think GIS should be able to do that easily.

By the way, I use QGIS 2.10.


texts: Converting KML to Layer in ArcGIS for Desktop?

texts:
I have a KML file that I know for a fact has 400 some points in it. But every single time I try to convert it using the toolbox conversion tool, it keeps ending with this result:

  ERROR 000401:No features were found for processing

I have tried this on two different computers so far, and I have opened it with Google Earth both times to verify that the file does have the points.

Not sure what to do next, does anyone have any other things I could try besides this conversion tool?

I am using ArcCatalog 10.4.


texts: How to create a text box with inverted rounded corners

texts:
I just saw this post on the tex SE site and wondered how I would do the same in ArcGIS.  I created a text box in a layout and drilled down in all of the dialog boxes but couldn't see anything similar.  Sure, it's easy to create a rounded box, but how do you invert the corners?

One solution I had was to manually draw a rectangle with circles on the corners.

Then, I used the subtract tool on the Graphics Operations context menu for each of the corners.

The result was satisfactory:

I could easily put some text in there and it would look fine.  However, if the amount of text changes and I need to stretch the box, the corners stretch as well:

Does anyone know of a way to create such a frame around text in the layout in ArcGIS?

I'd be happy to know if it's possible to do this using Python.


texts: Find and replace of ArcGIS i got message "An invalid SQL statement was used"

texts:
I was replacing text in a field of attribute table, by using "Find and Replace".
After saving and making new map by [Layer properties--Categories--Unique values--VAlue field--add all values--ok]

Now when i tried to replace some text values as before i am getting this msg

"An invalid SQL statement was used" {i am attaching screenshot}

Please help, i want to replace some text like "Calcareous duricrust" with "{Calcareous Duricrust}"

What can i do?


texts:
I am new to spatial analysis and even newer to hydrology. I have a 30arc DEM elevation raster and I was wondering if I can use it as an input to get the approximate presence of rivers in a geographically determined zone.  I know I can use rivers shapefile to do this, but since I am working with cells as my unit observation, it is cleaner get the data from the DEM raster.

raster package in R seems appropriate to this. Nonetheless, it only calculates hill, slope and roughness.


texts: How to add gridlines as !staticlayer to a MapFish3 print template?

texts:
I am trying to add gridlines to a template by using a !staticlayer attribute with the !addBackgroundLayers processor. I am able to output a PDF map, however, no gridlines are visible.

Here is my config.yaml:

templates:

 #===========================================================================
 # Demo Template
 A4 landscape: !template
 #===========================================================================
reportTemplate: template.jrxml
attributes:
  title: !string {}
  subTitle: !string {}
  mapMain: !map
    maxDpi: 300
    width: 780
    height: 390
  overlayLayers: !staticLayers
    default:
      layers:
        - type: "grid"
          gridType: "lines"
          numberOfLines: [4, 4]
          labelColor: rgba(0,0,0,0)
          haloColor: rgba(0,0,0,0)
  scalebar: !scalebar
    width: 240
    height: 40
    default:
      type: line
      unit: m
      orientation: "horizontalLabelsAbove"
      renderAsSvg: true
      backgroundColor: "rgba(0, 0, 0, 0)"
      align: "right"

processors:
- !reportBuilder # compile all reports in current directory
  directory: '.'
- !createMap
  inputMapper: {mapMain: map}
  outputMapper: {mapSubReport: map}
- !addBackgroundLayers
  inputMapper: {
                overlayLayers: staticLayers,
                mapMain: map
               }
- !createScalebar {}

All the other items on my map appear as expected (scalebar, title, etc). I have only been successful in applying gridlines to the map by passing the gridlayer in the requestData json.


texts:
I have a table which contains address-points.

With following sql-statement i have created a 100-meter buffer around each address-point:

  CREATE TABLE BUFFERPOLYGONS AS SELECT a.OBJECTID, SDO_GEOM.SDO_BUFFER(a.shape, m.diminfo, 100,  'unit=meter arc_tolerance=0.05') AS SHAPE
   FROM ADDRESSPOINTS a, user_sdo_geom_metadata m
   WHERE a.ADDRESSPOINTS='MAINADRESS' AND m.TABLE_NAME='ADDRESSES';

After this, the created table BUFFERPOLYGONS contains a lot of circle-polygons and now I want to union them (like dissolving in ArcGIS), so I did this sql-statement:

  CREATE TABLE BUFFERPOLYGONS_DISSOLVE AS SELECT 1 AS OBJECTID, SDO_AGGR_UNION(
    MDSYS.SDOAGGRTYPE(shape, 0.005)) AS SHAPE
    FROM BUFFERPOLYGONS;

But after this union-statement, i get this ORA-13346 error message:

Cause: Invalid definition of an arc. An
arc is defined using three non-collinear coordinates.
Action: Alter the definition of the arc,
or set the SDO_ETYPE or SDO_GTYPE to the line type.

I don´t really know what I have to do now with my data (SDO_GEOMETRY)

any suggestions?


texts:
I am new to the endless world of GIS. I am working on the project where I need to create a "Blocks" which are a combination of some countries.
e.g: I have Block "Balkans" where Bulgaria, Serbia and Romania are added. Each user with a particular role will be able to create blocks. The blocks can be from 2 to 20 countries. The problem that I have is with the performance when I try to combine Bulgaria and Romania cost me 10 seconds.

INSERT INTO regions (boundaries)
VALUES (
    (Select ST_MULTI(ST_UNION(boundaries)) FROM regions WHERE (name='Bulgaria' OR name='Romania') AND level='country')
);

This is the query for inserting both countries.
Here is a Pastebin with the Multipolygon for both countries https://pastebin.com/Hq6EhZd1 if you need to run it and test it.

I tried to add it as BINARY but the query took 10 seconds as well.

INSERT INTO regions (boundaries)
VALUES (
    (SELECT
            ST_AsEWKB(
                ST_MULTI(
                    (SELECT ST_UNION(boundaries) from regions WHERE (name='Bulgaria' OR name='Romania') AND level='country')
                )
            )
    )
);

I am using this version of PostgreSQL and PostGIS

SELECT PostGIS_Full_Version();

Result

POSTGIS="2.4.4 r16526" PGSQL="100" GEOS="3.6.2-CAPI-1.10.2 4d2925d6" PROJ="Rel. 4.9.3, 15 August 2016" GDAL="GDAL 2.1.4, released 2017/06/23" LIBXML="2.9.1" LIBJSON="0.12" LIBPROTOBUF="1.3.0" RASTER

Postgre

SELECT version();

PostgreSQL 10.6 on x86_64-pc-linux-gnu, compiled by gcc (GCC) 4.8.3 20140911 (Red Hat 4.8.3-9), 64-bit

Analyze result is

Insert on regions  (cost=209.86..209.88 rows=1 width=684) (actual time=4.634..4.634 rows=0 loops=1)
  InitPlan 1 (returns $0)
    ->  Aggregate  (cost=209.84..209.86 rows=1 width=32) (actual time=3.883..3.884 rows=1 loops=1)
          ->  Seq Scan on regions regions_1  (cost=0.00..209.84 rows=1 width=22613) (actual time=0.078..0.251 rows=2 loops=1)
                Filter: ((level = 'country'::text) AND ((name = 'Bulgaria'::text) OR (name = 'Romania'::text)))
                Rows Removed by Filter: 591
  ->  Result  (cost=0.00..0.02 rows=1 width=684) (actual time=3.943..3.944 rows=1 loops=1)
Planning time: 0.400 ms
Trigger refresh_site_regions_on_region_change: time=412200.845 calls=1
Execution time: 412210.485 ms


texts: Adding layers to an OpenLayers map from a WPS process automatically after processing

texts: I want to process two layers visible from a an OpenLayers client map interface such as Geonode and then add the result of the output directly to the map automatically after the processing is finished. is there anyone who would have a general description of how to go about it. I have medium level knowledge in OpenLayers.


texts: Creating circular buffers around a point in QGIS with Python

texts:
I have a script that creates a square buffer of 1 degree around a point in QGIS, through the python console but I am unable to find a way to create circular buffers instead. Also, would it be possible to add to the script to add multiple buffers around the same point of different distances, as I currently have 1 degree but would like 5 degrees, 10 degrees... etc. The script I have currently written:
layer = iface.activeLayer()

feats = [ feat for feat in layer.getFeatures() ]

epsg = layer.crs().postgisSrid()

uri = "Polygon?crs=epsg:" + str(epsg) + "&field=id:integer&field=x:real&field=y:real&field=point_id:integer""&index=yes"

mem_layer = QgsVectorLayer(uri,
                           'circular_buffer',
                           'memory')

prov = mem_layer.dataProvider()

for i, feat in enumerate(feats):
    point = feat.geometry().asPoint()
    new_feat = QgsFeature()
    new_feat.setAttributes([i, point[0], point[1], feat.id()])
    tmp_feat = feat.geometry().buffer(1, -1).boundingBox().asWktPolygon()
    new_feat.setGeometry(QgsGeometry.fromWkt(tmp_feat))
    prov.addFeatures([new_feat])

QgsProject.instance().addMapLayer(mem_layer)


texts: What is *.gdb.fm file that appears to be associated with geometric network of ArcFM Water?

texts:
I have been doing some work with a file geodatabase (test.gdb) which has a geometric network that I copied (using ArcCatalog or the Catalog window of ArcMap) from an ArcFM managed water network stored in SQL Server via ArcSDE at ArcGIS 10.2.x.

I have noticed that the test.gdb folder has an accompanying test.gdb.fm file, which appears to have been created at the same time as the geometric network was copied.

What is this *.gdb.fm file used for?

When the Catalog window of ArcMap is used to rename a file geodatabase, I notice that this file does not get renamed, which suggests that ArcCatalog is unaware of it.


texts: Dropping duplicate strings when using group_concat in Virtual Layer

texts:
When using Virtual Layer to group_concat() a string field, the result would contain duplicates string. For example, here's a polygon layer "city_lot" with "land_use" and "owner" fields.

I want to create a virtual layer that merges lots by owner. So here's the query:
select st_union(geometry),"owner",
    group_concat("land_use") as 'land_use',
    group_concat("fid") as 'fid'
from "city_lot"
group by "owner"

And here's the result:

Notice that the 3rd record has two "residence" in the "land_use" field. My question is, instead of being "residence,commerce,residence", would that be possible to make the result be "residence,commerce", merging identical string?


texts: How to run spatial query get fields value of 2 layers?

texts:
Last month I asked vb.net arcobjects spatial query get fields value of desired layer only.

Now there is some change in the requirement, I need to get value of 2 layers after use "Select Features" tool in arcmap,  1 polyline and 1 polygon.

Can I use the same method?

One thing to highlight is the editable variables in attribute table of the 2 layers are the same.

My code is here:

    pMxDoc = My.ArcMap.Application.Document
    pMap = pMxDoc.FocusMap
    pActiveView = pMxDoc.FocusMap

    Dim layerNum = GetIndexNumberFromLayerName(pActiveView, "test_db.DBO.TestTable")
    Dim layerNum2 = GetIndexNumberFromLayerName(pActiveView, "test_db.DBO.TestTable2")

    Dim pFeatureLayer As IFeatureLayer
    Dim pFeatureLayer2 As IFeatureLayer

    pFeatureLayer = pMap.Layer(layerNum)
    pFeatureLayer2 = pMap.Layer(layerNum2)

    Dim pFeatureSelection As IFeatureSelection
    pFeatureSelection = pFeatureLayer
    logger.Info("pFeatureLayer = " + pFeatureLayer.Name)

    Dim pSelectionSet As ISelectionSet
    pSelectionSet = pFeatureSelection.SelectionSet
    logger.Info("pSelection number = " & pSelectionSet.Count)

    Dim pFeatureCursor As IFeatureCursor = Nothing
    pSelectionSet.Search(Nothing, True, pFeatureCursor)

    Dim tempArr As New List(Of String)
    Dim pFeature As IFeature
    pFeature = pFeatureCursor.NextFeature

    Do While (Not pFeature Is Nothing)
        tempArr.Add(pFeature.Value(pFeature.Fields.FindField("Item_ID")).ToString)
        pFeature = pFeatureCursor.NextFeature
    Loop

Sorry of the unclear requirement of the whole process. After performing the spatial query, I will retrieve data related to features of these 2 layers (from database) and allow user to select and zoom to any results from the spatial query.


texts: QGIS Time Manager with Multiple Vector Layers

texts:
I have to animate a map of Europe (31 countries) which contains three separate vector layers:

The first layer is a simple map (with just the political boundaries) of Europe with 31 countries (each country being a distinct feature). This is a colour-coded map with 4 shades of gray, depending on how much energy a country can store during an hour (dark gray indicating the highest amount of stored energy, white the lowest).

The second layer is a collection of 31 circles (one over every country) with the size of each circle being proportional to the amount of energy produced during an hour. Each circle is a separate feature.

The third layer consists of transmission lines drawn between the circles. The width of each line indicates the amount of energy it transmits during an hour. There are 50+ lines in this layer and each line is a distinct feature.

Hence on a given map, these 3 layers (the circle layer superimposed on the storage layer, and the line layer superimposed on the circle layer) represent one picture of the energy stored, produced and transmitted in a particular hourly time-step. Using MATLAB, I can generate these layers for successive time-steps and I would like to animate those pictures for a period of 2 days (hence 48 pictures in total). Can you please guide me as to how this can be done most efficiently? Should I merge the three types of vector layers into one layer, e.g., for the first hour (and append the first time-stamp to them), then do the same for the next three layers for the second hour (and append the second time-stamp), and so on, and then finally merge all the 48 pictures together? Could there be an easier and less time-consuming way to animate those pictures?

I'm sorry if this is a basic question but this is the first time I am using QGIS Time Manager.


texts: How to flatten a PolygonZ onto the XY-plane using pygis


texts: Using polygon layers to cut holes in other layers

texts:
I have a number of polygons in a layer (old rehabilitation areas) and want to be able to selectively bring one to top and use this to cut a hole in the others and fill this with the cutting polygon information (rework rehabilitation). Is this easy to do in QGIS? I know I can use the add ring to put a hole in a single polygon but this uses a drawn polygon (selection of a polygon to use for cutting not available?) and produces an actual hole; similiarly add feature requires the digitising of a polygon. So I require the following:

1) Selectively bring a polygon to the top (is there a bring to front function?)
2) Select the now top polygon to use as a cutter
3) Cut top polygon shape out of underlying polys, this hole to be filled with cutting poly info. Thus leave a layer containing non-overlapping polygons.

UPDATE.  I've found I can get this to work using using the Difference tool.

1) Separate the overlying polygons into their own layer
2) Make base polygon layer editable
3) Use the Difference tool to cut out the overlying polygon shapes in the base polygons
4) Copy and paste the overlying polygons into the base cutouts

Is there a quicker way? Also can this be done easily in ArcGIS? Thanks.


texts: Georeferencing multiple "slightly" overlapping drone orthomosaics

texts:
Situation:
I collected around 2000 drone images over 200 ha. My computer does not have the resources to render one large orthomosaic at once, so I decided to split my dataset into 8 and produce 8 different orthomosaics. I am using WebODM for my image processing.
Because I did not use GCP's, all my 8 orthomosaics do not align properly with one another, and I want to merge them into one large orthomosaic.
Question:
What software/technique can I use to align my 8 orthomosaics so I can merge them into one?
I have tried:
I found this Python library arosics (https://pypi.org/project/arosics/) to co-register overlapping imagery. However, I cant get it to work with my drone orthomosaics. I believe it has to do with the pixel type, but not sure. I can share the code and output here if someone is interested.
I know cloud processing is a solution, but I do not have the funds at the moment to make this a viable solution.


texts: Exporting images attached to a geodatabase using QGIS

texts:
I have an ESRI file geodatabase that was exported from their survey123 app, and the geodatabase has images attached to it.

I do not use ArcGIS, and was wondering how to export them from QGIS with their associated metadata?

Here's a Google Drive link to my gdb:

https://drive.google.com/drive/folders/1KIJMea9bhF0cZrMc-y1uCLjLpzScWukd?usp=sharing

Once open in QGIS, I can see the point features, and below that in my layers panel is a table with the photos names, etc.


texts: Calculate 500m south of coordinate in PostGIS

texts:
I have a point layer punkte_test. From that I want to make a line: first point of the line is the point A from the point layer itself. Second point of the line shall be y(B) = y(A) - 500 m . How can I realize that?
That's my idea:
SELECT ST_MakeLine(geom, concat(ST_X(geom), (ST_Y(geom)-500)))
FROM liechtenstein.punkte_test

I use coordinate system UTM, so I want the second point to be 500 m south of A.


texts: encoding problem with QGIS 1.8 and CSV import

texts:
I have following problem with QGIS 1.8. Lisboa (standalone installer, Win 7 Pro 64bit) - I import CSV - text in CP1250 separated with ";" and using "," as decimal separator.

I correctly set the separated text importer, get correct point layer which has all text entries displayed correctly. I save this as shapefile, set the CP1250 encoding but when I load this layer in QGIS (choose CP1250 again during SHP loading) but get point layer with damaged text entries.

I use this tool for a long time and did not have such problems. Did you have encounter such problems?

thanks


texts: Finding only records with numbers in a character string field?

texts:
I need to isolate those records containing numbers from a field with character string (it's a descriptive field). For example the string might be as follows:

  "A residential development of 12 no. two storey dwellings comprising
  10 no. 4 bedroom semi-detached houses and 2 no. detached 3 bedroom
  houses all with associated on-curtilage parking."

Table name: PLANAPP_iPlan_Merged

Field Name: Description

The numbers can occur anywhere in the string so I can't set Left%/ Right$ parameters. Anyone any ideas?


texts:
I am new to web mapping. I am trying to display roads  on my map, hence I build a server side application with flask, where I am selecting the roads table and saving the geometry column via AsGeoJson in Geojson format.

@app.route('/addstreet', methods=['GET', 'POST'])
def addstreet():
    conn = g.db.cursor()
    conn.execute('''select id,name,sub_type, AsGeoJSON(ST_GeomFromText(geom))  FROM roads''')
    result = conn.fetchall()
    result = json.dumps(result)
    return result

Ok I saved my roads in an array like this:

 $("#roads").click(function(event) {
   //add streets
   $.get(my_url, function( data ) {

   features =[];
   data = JSON.parse(data);
    for (i=0; i< data.length; i++) {
        roads = data[i][3];
        features.push(roads);
    }

    console.log(features);

    var myStyle = {
        "color": "#ff7800",
        "weight": 5,
        "opacity": 0.999
    };

    L.geoJson(features, {
        style: myStyle
    }).addTo(map);

});
});

The error message is gone. But on my map I can't see any roads....


texts: Polygon shapefile has irregularities when viewed in ArcMap?

texts:
i am using ArcGIS 10.4 for Desktop and when i try loading my shapefiles into the arcmaps it shows extra incorrect lines in my polygon shapefile is there anyway to fix this inconsistencies automatically?

I tried using the topology tools and snap tool but it didn't fix the issue.

this is how the data was before being converted to a line


texts:
I have a DEM file that I created a slope map of using gdaldem slope. This works fine and if I view and colorize it, I can now see how steep a hill is but now I'd love to create a shapefile out of that. So I'd like to have polygons with increments of let's say 5% for each level of steepness.
I though I would go ahead and create gdal_contour of that slope-DEM file and using the -p option to create polygons rather than lines, but the problem is, that this takes FOREVER. Even a small sample file of a few square-kilometers (and about 20MB in file-size) is taking a couple of hours.
My command was: gdal_contour -i 5 -amin 30 -p input.tif output.shp So increments of 5% starting at 30% and polygons.
How can I approach my problem?
Additional Info:
For @gabriels anwser. Here's what that looks. So you can basically see each "pixel".


texts:
I display OpenStreetMap with some SLD styles.

This map is not the same OSM.

I want that this map is equal OSM map.

I use GeoServer 2.3.5

Does anyone display OpenStreetMap through MapServer, ArcGIS for Server or Geoserver? - here are some ideas.

Update: The map may be not relevant OSM map. But it must be beautiful.


texts: Summing attribute values from many points geocoded on same location using ArcGIS Desktop?

texts: I have geocoded a batch of addresses and many of these are geocoded on the same location. Now I want to keep just one point on each location but sum the field values of all points from the same location. How I can do that using ArcGIS?


texts: Is editing feature service through the REST API an okay approach?

texts:
I have a database that receives real-time data from alarms throughout the day, or on request.

I'd like to create a web map that displays the "bad" alarms.

My current idea is to create a Python script that queries the database every 15 minutes and gets a unique ID and a lat/long of the bad alarms. I will have a feature service setup and use the REST API (through Python) to add or delete alarms based on the previous query. I'll then have a basic web map set up to display this feature service with other related layers.

Is this a good approach? This data will only be viewed internally. I've never used the REST API before, but it seemed like a straightforward solution.

Thanks for any guidance.


texts: Make diagonal grid in QGIS, divide square into squares

texts: I have square polygon which is aligned under some degree and I want to divide this square to, for example, 25 smaller squares, which must be under the same degree in this square polygon. How this can be achieved?


texts: Install pgRouting on Postgres 9.1 and Win 32

texts:
I referred some documents on GIS-stackexchange and over internet regarding installation of pgRouting on Postgres v9.1 and windows 32.
I am not getting the correct version of pgRouting to install.
Can anyone guide me to get the right pgRouting version for Postgres 9.1 and windows 32.

I tried winnie builds and sanak builds but I am getting errors after running SQL query CREATE EXTENSION pgrouting;. The error are:

Using winnie build for Postgres 9.1 and win 64
ERROR: incompatible library "C:/Program Files/PostgreSQL/9.1/lib/librouting.dll": version mismatch
Using pgRouting for Postgres 9.2 and win 32
ERROR: incompatible library "C:/Program Files/PostgreSQL/9.1/lib/librouting.dll": version mismatch
Using Sanak build:
ERROR: could not load library "C:/Program Files/PostgreSQL/9.1/lib/librouting_bd.dll": The specified module could not be found.


texts:
I have some text labels in Leaflet which I'm building from JSON data using L.divIcon. They're looking like this:

lbl_coord = feature.properties.coord;

lbl_text = L.divIcon({
        className: /*...*/,
        iconAnchor: [0,0],
        popupAnchor: /*...*/,
        html: "<b>" + feature.properties.name + "</b>"
      });

  var label = L.marker(lbl_coord, {icon: lbl_text});
  label.addTo(labelLayer);

Well, this works but the performance is rather bad because each marker is rendered as individual DOM element. It becomes very slow on panning operations. I'd like to load 10 000 labels at some point and am struggling to make it use canvas. I couldn't get the Canvas Markers Plugin to work as I'm experiencing the same issue that is discussed here (The HTMLImageElement provided is in the 'broken' state). However, the workaround proposed there won't work for me (error: Layer is not defined). Besides, there is no license given, so I don't know if I might use this at all for my project. Here's another Plugin but this wouldn't work with divIcon functionality.

At this point, I am wondering why this functionality is not part of Leaflet core. I know >10 000 markers is a lot but performance is already quite bad with 100 markers. I hope somebody's got a solution.


texts:
I have a bunch of MapInfo Tables in raster format (.tab / .grd) that I would like to consolidate into one table.

Is this possible within MapInfo?


texts: How to set an ol.style.text using a source in OpenLayers?

texts:
I have the following code:

var cbcWFS = new ol.source.GeoJSON({
    projection: 'EPSG:3857',
    url: 'http://localhost:7070/geoserver/CuerpoBomberosCurico/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=CuerpoBomberosCurico:CuerpoBomberosCurico&maxFeatures=50&outputFormat=application/json'
  });

var cuerpoBomberosVector = new ol.layer.Vector({
      title: 'Cuerpo Bomberos',
      source: cbcWFS,
      style: new ol.style.Style({
        image: new ol.style.Icon({
          scale: 0.04,
          src: 'https://dl.dropboxusercontent.com/u/27798645/fireman/truck65.svg'
        }),
        text: new ol.style.Text({
          text: 'Hello',
          scale: 1.3,
          fill: new ol.style.Fill({
            color: '#000000'
          }),
          stroke: new ol.style.Stroke({
            color: '#FFFF99',
            width: 3.5
          })
        })
      })
  });

But now, in the part of "Text", I'd like to use the value that it comes from the source.

Any idea to develop it?


texts:
I am a history teacher and I like to make my own maps (I can never find exactly what I want).  I am using QGIS 1.8 and Natural Earth data.

I want to reproject their cross-blended hypso rasters into projections that are more familiar looking to my students, like Lambert Conic Conformal.

I have tried changing the crs for the whole project, and the warp and assign functions from the raster menu, and changing the crs for the raster layer itself.  None of these seem to work properly - the resulting maps are unrecognizable and extremely large (file size).

I feel like I am making some kind of basic mistake but I don't know where to begin with this.  Any help would be appreciated.


texts: Extracting lat lon of points from list of postal codes

texts:
I have scraped a list of postal codes for which I want to extract (lat, lon) data.

extracted postal codes

In Postal Codes Germany, I was able to extract lat lon for the postal code by is there a way I can for a list of postal codes all at a time?

postal code


texts: Replacing characters after certain word in Python parser of ArcGIS Field Calculator?

texts:
I have a table with field which contains path to database. The same thing in this field will be the extension of file.
Here is how it looks like

D:\Проекты\Москва для тестов\Test.gdb\graf_fd

I need to have this

D:\Проекты\Москва для тестов\Test.gdb\

Probably I should replace everything after 'gdb\'. I know that I can type something like delete last 7 charachters, but the name of dataset can be another in other databases as well as the length of string data before 'gdb\'. So are there codes for removing or replacing the characters after certain word?


texts:
I am very new at PostgreSQL / PostGIS. I have got a table with geometry column set in (Polygon,32733), so I would like to update the Area column so as to see the area of each polygon.

I tried this:

ALTER TABLE public."EDIFICIO_2" ADD COLUMN area double precision;
update "EDIFICIO_2" set area = ST_AREA(Polygon,32733);

But the code gave no results.


texts: Add features to a layer where another vector layer does not intersect

texts: I have a vector layer showing the premises of hospitals in a city (purple layer) and a layer of all buildings (not just hospitals). From this I used Clip to find the hospital buildings (pink layer). Some hospital premises are present but the buildings within them aren't (within the red boxes). How can I select the premises where there are no buildings and add them to the pink layer? Essentially finding the polygons where the all buildings layer does not overlap with the purple layer.


texts: Saving QGIS layer to PostGIS using PyQGIS?

texts:
How to save the memory layer from QGIS in PostGIS database in specific schema (public) using PyQGIS or python?

Here some example from my code :

outputs_QGISDIFFERENCE_1=processing.runalg('qgis:difference', layer1,layer2,False,None)
load=QgsVectorLayer(outputs_QGISDIFFERENCE_1['OUTPUT'], 'myname', 'ogr')
QgsMapLayerRegistry.instance().addMapLayer(load)

One way is using line :

load=QgsVectorLayer(outputs_QGISDIFFERENCE_1['OUTPUT'], 'myname', 'ogr')

to get the layer from algorithm.

Or another way is :

load='c:some/path/shapefile.shp'

to get the layer from the path.

But how to save this layer in PostGIS database?


texts: Why are points from a MySQL WMS layer, displayed on a OpenLayers 2 map, reading out a lat/long location that is off by 500 meters,

texts:
I have a MySQL table containing MySQL point geometries (MySQL 5.7). I have configured a GeoServer WMS layer using this data. I cannot set the native srs (not enabled in config dialog), declared srs is set to EPSG:4326 and srs handling is set to force declared. The layer displays well but when zooming close and hovering the mouse over point, in the OpenLayers lat/long readout is typically around 300-500 meters off. Seems worse in the x (longitude) axis. I verified that the stored point have 5 decimal places. Verified by selecting point from table and displaying at text  with AsWKT function.
Any thoughts or suggestions?


texts: When I create flood areas using DEM in ArcGIS, which parameters do I consider to determine DEM's pixel size and resolution?


texts:
How can I deactivate a selected basemap in QGIS 2.18 and QGIS2WEB 2.34.0 and completely dispense with the use of basemaps?

If I see it right, once a basemap has been selected, I can only choose between different basemaps. Can't deactivate or am I wrong?


texts: How to calculate Field using part of filename in ModelBuilder?

texts:
I am new to ModelBuilder and have not been able to figure how to use part of the FC name: For example the first FC in my GDB containing around 100 files is "A1005_GPS_Residence", and I need to calculate the new field SID = "A1005".

I have seen some posts with python and manually using the Field Calculator.

I tried using one expression, shown in the screenshot, but I am not sure if I can combine field calculations with inline substitution or the like !%Name%![1:4].

Screenshot of MB


texts: Start postgresql+postgis without installation on Windows

texts:
I have used Windows7x64 and want start PostgreSQL+POSTGIS without installation process. I copy PostgresSQl\9.1 folder from other laptop to c:\opt. I need create .bat file with enviroment variables and add pg_ctl.exe with any parameters. is it possible start Postgresql server without creating Windows services and using any regedit manipulations?
Sorry if my question is not correct compiled


texts: I am looking for a buildings layer for Ontario, particularly the Chatham-Kent area to properly assess setback for wind turbine development in the region.  Currently, I am utilizing the buildings data provided on Geogratis; however, I would like to know it's accuracy and if other buildings data is available free or for purchase.


texts: How I can reproject a WMS layer in OpenLayers 3?

texts:
I was trying to reproject a WMS layer from EPSG:25830 to EPSG:3857, but it didn't work. I read that it is not possible to reproject a WMS, only a WFS, It is true?

Here is the code:

<!doctype html>

I now created other code that has the CRS, even though I continue to have the same problem:

<!doctype html>

<title>Single Image WMS</title>
<link rel="stylesheet" href="http://openlayers.org/en/v3.15.1/css/ol.css" type="text/css">


texts: OpenLayers saving GeoJSON features in different projection

texts:
I've used ol.interactin.Draw to draw and save features in GeoJSON format. But the drawn features are stored in EPGS:3857 and my preference is WGS 84. For that, I've used this snippet code to apply dataProjection and featureProjection properties but nothing has been changed.

draw = new ol.interaction.Draw({
        source: vectorSource,
        type: selectedType.value
    });

    draw.on('drawend', function (evt) {

        var extraFeature = evt.feature;
        var features = vectorSource.getFeatures();
        var allFeatures = features.concat(extraFeature);

        featuresGeoJSON = parser.writeFeatures(allFeatures, {
            dataProjection: 'EPSG:3857',
            featureProjection: 'EPSG:4326'
        });
});


texts:
Is there a way to have a start and end date in a Carto timeline (instead of the usual start OR end date?

I don't have any programming experience.


texts: GetFeatureInfo is not working on multiple layer

texts:
This is my map URL: http://34.84.56.42:8080/app/
When I click on the map, I'm able not to get the information (GetFeatureInfo) of the overlay layers. I can get only GetFeatureInfo of the top layer.
I just changed the WMS source and layers in the example code.

In the example,

When I click on the map, the information of the two overlapping layers appears.
How can I fix it?

requirejs.config({
    'baseUrl': './lib',
    'paths': {
        'leaflet.wms': '../src/leaflet.wms' //.js'
    }
});

define(['leaflet', 'leaflet.wms'],
function(L, wms) {

var tiledMap = createMap('tiled-map', true);

function createMap(div, tiled) {
    // Map configuration
    var map = L.map(div);
    map.setView([37.568,126.993035], 17);

    // 베이스맵 설정
    var basemaps = {
        'OpenStreetMap': basemap().addTo(map),
        'Blank': blank()
    };

    // Add WMS source/layers
    var source = wms.source(
        "http://34.84.56.42:8080/geoserver/yh_test/wms",
        {
            "format": "image/png",
            "transparent": "true",
            "attribution": "BetaCityCenter",
            "info_format": "text/html",
            "tiled": tiled,
        }
    );

    // 레이어 설정
    var layers = {
        'yh_test:200814_JungGu_JongroGu_Buld(epsg5186)': source.getLayer("200814_JungGu_JongroGu_Buld(epsg5186)").addTo(map),
        'yh_test:200702_Japaneese house_GUGA': source.getLayer("200702_Japaneese house_GUGA").addTo(map),
        'yh_test:200702_Korean house_GUGA': source.getLayer("200702_Korean house_GUGA").addTo(map),
        'yh_test:guga_yh': source.getLayer("guga_yh").addTo(map),
        'yh_test:Sewoon_area(epsg5186)': source.getLayer("Sewoon_area(epsg5186)").addTo(map),
        'yh_test:3D_scan_UOS_process': source.getLayer("3D_scan_UOS_process").addTo(map),
        'yh_test:Process_area': source.getLayer("Process_area").addTo(map),
    };

    // 레이어 선택창 생성
    L.control.layers(basemaps, layers).addTo(map);
}

// 오픈스트리트맵으로 베이스맵을 설정
function basemap() {
    var attr = 'OpenStreetMap';
    return L.tileLayer("https://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png", {
        attribution: attr
    });
}

// 베이스맵 없는 상태
function blank() {
    var layer = new L.Layer();
    layer.onAdd = layer.onRemove = function() {};
    return layer;
}
});


texts: R GIS layers not lining up (spplot/plot)

texts:
I'm trying to make a choropleth in R using a polygon file, with a points layer of local care homes on top. When I plot each layer individually they're fine, however when I put them both together (using spplot for the polygons and plot for the points) the coordinates not longer line up, in that various points which should be in the area are now off the map (in the sea!). Is it possible that using an spplot/plot combo would cause the problem?

I'm confident there are no errors in the original files.

library(rgdal)
library(maps)
library(maptools)
library(RColorBrewer)

setwd(....)
shape <- readOGR(".","LSOA Boundaries")
setwd(".....")
values <- read.csv("lsoamappingvalues.csv")
plotme <- merge(shape,values, by.x="LSOA11CD",by.y="LSOA")

setwd(".....")
carehomes <- read.csv("care home points.csv")

x <- carehomes$LONG
y <- carehomes$LAT
homecords <- data.frame(x=as.numeric(x),y=as.numeric(y))
coordinates(homecords) <- c('x','y')

spplot(plotme, "POP",
col.regions = colorRampPalette(brewer.pal(9, "Greens"))(100), col=NA)
plot(homecords, add = TRUE, col='red', pch=19)


texts: I have a excel and a shapefile.I want to check whether the attributes present in the excel matches the one in the shapefile using FME


texts: Data being excluded from table joins (QGIS)

texts:
I've got a map with spatial data that I want to merge with some non-spatial tables. The tables will become spatial based on unique numerical identifiers contained in a column shared with the spatial data.

I import the tables as .csv files and join it to the relevant layer without problem. However, I noticed that when I then query the combined attributes some have not made it into the joined table. This leaves me with mystery gaps in my data that are not easy to discover in most cases.

Has anyone had this experience before? Using QGIS 2.18.11


texts: Accessing MODIS subdatasets from the command line?

texts:
FYI I'm new to command-line geospatial operations. I'm working with the MODIS BRDF-Albedo Quality product (MCD43A2) and am trying to:

1) do a simple gdal_calc.py operation to mask certain QA values in 100s of .hdf4 files

2) use these masks in the 'mask' option  in gdal_translate to create .tiff files.

I can access the subdatasets of one hdf file in the command line with:

gdalinfo -sd 1 MCD43A2.A2013225.h12v09.005.2013234130154.hdf

^ This prints the information, but I need to be able to reference this subdataset in the subdataset_band1.hdf portion of the gdal_calc expression below:

gdal_calc.py -A subdataset_band1.hdf --outfile=subdataset_band1_mask.hdf --calc=“1*(A<3)”

How can I reference or save a single subdataset to access in the input file of the calc command?

Note: I'm trying to avoid first creating .tiffs of the QA, then creating the masks, then applying them in a gdal_translate with other .hdf files -- I plan to apply .hdf masks to other .hdf files in step 2 above... suggestions for more efficient ways to do this welcome!


texts: Create Geographic Reference for X,Y Coordinates

texts: I have a set of x,y data that do not have any sort of geographic reference.  I have lat/longs (WGS 84) for a small set of those points.  How can I calculate lat/long for the remaining points so I can get them geographically placed?  Thanks!


texts: Data frames not showing up in Extent Indicator Property Box

texts:
I have a data frame showing an overview of a city. I also have multiple  data frames showing zoomed in areas. I want to show the extent indicators for my zoomed in areas on my overview data frame. When I got to the Extent indicator tab for the overview I don't have any data frames to add. However when I go to the extent indicator tab for the zoomed in view I can select the overview data frame. This seems a little backwards to me. Has anyone else run into this issue or have any suggestions?
I'm using Arcmap 10.1

Here is a screen shot of what I have set up and where the blow ups are.

All data frames have the same projection and and set to a fixed extent. I just added a new data frame with a larger extent than my overview data frame. The over view extent shows up in the new data frame but none of the other zoomed in views do. The overview is set at a 1:3000 scale and the zoomed in data frames range from 1:100-1:500. For some reason the zoomed in data frames think they can show the extent of the overview. Unfortunately I didn't set up the mxd in the beginning, is there a setting or something that I have missed that might cause this?


texts: Exporting Network Analyst layer objects (orders, routes, etc) created in Python after solve?

texts:
I'm trying to export Orders dataset after solving a Vehicle Routing Problem (VRP) into a shapefile using Python only. The problem is that apparently arcpy without mxd open doesn't recognise my VRP layer. Below is the code I'm working on:

outNALayer = arcpy.na.MakeVehicleRoutingProblemLayer (inNetworkDataset,
                                                      outNALayerName,
                                                      "Time",
                                                      "Length",
                                                      "Minutes",
                                                      "Meters",
                                                      "",
                                                      "",
                                                      "High",
                                                      "Low",
                                                      "ALLOW_DEAD_ENDS_ONLY",
                                                      ["Oneway","Speed"],
                                                      "NO_HIERARCHY",
                                                      "",
                                                      "NO_LINES")
outNALayer = outNALayer.getOutput(0)

#Get the names of all the sublayers within the VRP layer.
subLayerNames = arcpy.na.GetNAClassNames(outNALayer)
#Stores the layer names that we will use later
ordersLayerName = subLayerNames["Orders"]
depotsLayerName = subLayerNames["Depots"]
routesLayerName = subLayerNames["Routes"]
pairsLayerName = subLayerNames["OrderPairs"]

arcpy.na.AddLocations(outNALayer, ordersLayerName, inOrders, orderFieldMappings,"")

#List name of attributes in Depot layer
depotFieldMappings = arcpy.na.NAClassFieldMappings(outNALayer, depotsLayerName)
arcpy.na.AddLocations(outNALayer, depotsLayerName, inDepots, depotFieldMappings, "")

#List name of attributes in Route layer
routeFieldMappings = arcpy.na.NAClassFieldMappings(outNALayer, routesLayerName)
arcpy.na.AddLocations(outNALayer, routesLayerName, inRoutes,  routeFieldMappings, "")

#List name of attributes in Route layer
pairFieldMappings = arcpy.na.NAClassFieldMappings(outNALayer, pairsLayerName)
arcpy.na.AddLocations(outNALayer, pairsLayerName, inPairs,  pairFieldMappings, "")

arcpy.na.Solve(outNALayer, "SKIP")
#arcpy.management.SaveToLayerFile(outNALayer,outLayerFile,"ABSOLUTE")
arcpy.CopyFeatures_management(outNALayerName + r"\Orders", inOrders + '_' + str(iterator))

The error message is that the Orders does not exist or is not supported. I have tried saving the layer file on the disk and accessing Orders from that layer but get the same error. Is there a solution to access the actual sublayers within solved layer just in python. I think Select function is available in model builder only.
Note: the above code runs just fine when is run in the mxd python window.

The error is at the last line. The error message is:

  Traceback (most recent call last):   File
  "C:\Users\zarem\OneDrive\Total_Transport_Home\Code\Vehicle Routing
  Problem_Only - Copy.py", line 143, in
      arcpy.CopyFeatures_management(outNALayerName + r"\Orders", inOrders + '_' + str(iterator))   File "C:\Program Files
  (x86)\ArcGIS\Desktop10.3\ArcPy\arcpy\management.py", line 2435, in
  CopyFeatures
      raise e ExecuteError: Failed to execute. Parameters are not valid. ERROR 000732: Input Features: Dataset VRP_Lyr_Grp_74\Orders does not
  exist or is not supported


texts: How do I project a point stored as Lat / Lon to a state system?

texts:
I have some points stored as WGS84 (EPSG:4326) data. I want to write a SELECT statement to return the projected values, in EPSG:2811, and unprojected values.

Here is my table SQL:

CREATE TABLE points (name varchar, point geometry);

Here is the SQL used to inserted the points:

INSERT INTO points(name, point) VALUES(2, ST_GeomFromText('POINT(<lon> <lat>)', 4326));

This gives me a readable view of the imported data for a quick QA check:

SELECT name, ST_AsText(point) FROM points;

How do I add the transformed coordinates to my SELECT statement? I think I need to use ST_TRANSFORM but beyond that I'm lost.

Install details:

POSTGIS="2.1.5 r13152" GEOS="3.4.2-CAPI-1.8.2 r3924" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.11.1, released 2014/09/24" LIBXML="2.7.8" LIBJSON="UNKNOWN" TOPOLOGY RASTER


texts: Is it possible to detect TileWMS layer by click on map?

texts:
Here is a map with TileWMS.

Is it possible to detect if a user clicked on a red square?

How can I detect by mouse clicking on the map to get the TileWMS layer?


texts:
I want to create a query to select all ways and their nodes that exist within a bounding box using postGIS. The bounding box shall includes all details as osmosis "--bounding-box" command will retrieves.

Is there any way to do that?


texts: Is there a way to test if a feature class is registered with geodatabase?

texts:
I've got feature classes inside an Oracle 11g ArcSDE geodatabase.

I'd like to perform the arcpy method RegisterWithGeodatabase_management on all those feature classes.
But before that, I'd like to test if the feature class is already registered, as the registration method will return an error in that case.

How can I test if a feature class is registered with geodatabase or not?


texts: Can the Identify popup be made relocatable by the end-user?

texts: I've noticed that the Identify widget's popup always tends to go right over the feature that is being identified and cannot be moved.  Is there a way to enable functionality to move the Identify (InfoWindow) popup around to a more favorable position?


texts: Using 3D Models in Cityengine, Cesium and 3D scenes

texts:
I recently started working with 3D models, and the following question may be a bit too abstract? If I need to build a 3D city with buildings and other features, what software should I use as a database, editing and displaying the result as a web map. I may be comparing apples to oranges, but maybe I can get a general overview of the available softwares and their possibilities and limitations.
I have previously only heard of Esri solutions, such as CityEngine and 3D Scenes. I have an ArcGIS Pro with data interoperability extension, so I can use FME if needed. What is the difference between CityEngine and 3D scenes, how are the files exported to a web service? Do I have to use the ArcGIS Online for representation, or can I use ArcGIS API for JavaScript so I can display them in other apps? What can I specifically do with ArcGIS API for JavaScript as it pertains to 3D modelling and web maps?
There is also an open source solution called Cesium. I have never used it, but it seems to have similar functions. How is different from Esri?


texts: Checking whether shapefile has overlapping polygons?

texts:
Using QGIS, what is the quickest way to check whether a shapefile has polygons that overlap each other?

I guess I could import it into PostGIS, then see whether:

SUM(ST_Area(geom))

is the same as:

SUM(ST_Area(ST_Union(geom)))

But that doesn't feel like it's likely to be the most efficient way to do it, especially since my shapefile has about 40k features.

Any ideas?


texts: Why does QGIS display an empty dialog box at the top and force me to quit?

texts: When doing basic joins I get this empty dialog box that locks me out of doing anything, forcing me to quit. it's not a crash, but I cant click out of it or anything


texts: Publishing raster data in OpenLayers applications

texts:
I need tips for publishing raster data in my OpenLayers application.I published shapefiles through my postgis db (PgAdmin) to geoserver, not sure how to do that for raster? I readed that is not cleaver idea [in this answer][1]https://gis.stackexchange.com/a/181575/160696
I have map of municipality tiff format (1,6 GB) and i want to import that map as Layer to include in LayerSwitsher window.So my question is how? I had already published my raster data in db (not on Geoserver), but now I'm thinking what to do? Maybe to try directly?


texts: Set definition expression from click event

texts:
I'm trying to set a definition expression according to a user-selected value from a drop down menu, but I can't get the map to update the feature layers according to the expression. The user is able to click the drop down and select the values, but nothing happens when they click one.

HTML:

<select id="level" name="level">
  <option value="" selected="selected">Select Building Level</option>
  <option value="1">1</option>
  <option value="2">2</option>
  <option value="3">3</option>
</select>

JS:

//creates and adds feature layers to map
var generalIssuePoints = new FeatureLayer(".../FeatureServer/0", {
  mode: FeatureLayer.MODE_ONDEMAND,
  outFields: ["*"]
});
var firePullPoints = new FeatureLayer(".../FeatureServer/1",{
  mode: FeatureLayer.MODE_ONDEMAND,
  outFields: ["*"]
});

map.addLayers([generalIssuePoints, firePullPoints]);

//handles drop-down change event
on(dom.byId("level"), "change", function(e) {
  var level = e.target.value;
  var definitionExpression = "LEVEL_ = " + level;

  //sets definition expression
  generalIssuePoints.setDefinitionExpression(definitionExpression);
  firePullPoints.setDefinitionExpression(definitionExpression);
});


texts: How to load OSM file from local machine to the server using osm2pgsql?

texts:
I wanted to upload the OSM file from my local machine to server using osm2pgsql, but the data is not imported.

I have just created the database s1(no table), do I need to create the table also, if yes how?

My osm2pgsql query

osm2pgsql -c -d osm -U s1 --password xx -k -H yy  /home/g/Downloads/niedersachsen-latest.osm.pbf pbf -r -S /home/g/default.style -s

I am getting an error

Osm2pgsql failed due to ERROR: Unknown file format '-S'.


texts: Defining Color band to recognize from sing layer raster in R

texts:
I have a 2 part question, the first part was easy & is answered but I have one issue :

1) Load map in R and convert to Raster : Done

2) Tell R to only show a defined color band : NOT done

As you can see in my example map, it is a map of Namibia, and it is a single layer black and white map, with one blue dot.  The RGB code of this blue dot is:  Red == 14 , Green == 6, and Blue == 209

What I would like to do is tell R to Isolate this dot, and to leave out the rest, so in essence, when I plot the map, I have only the blue dot.

Here is my current code, and here is an example map.

rm(list=ls())

library(raster)

library(maps)

library(magick)

library(rgdal)

# working directory

setwd("Working here")
# load first image
Species_id<-1
Species_i <- image_read(paste("Species (", Species_id, ").png", sep=""))

Species_i_rotated<-image_rotate(Species_i,0)
image_write(Species_i_rotated, path = "Temporary.png", format = "png")

# create raster
Converted <- raster("Temporary.png", band=1)

#### This is where my problem needs solving
#This defines the threshold to classify where the dot on the image is present as 1
# While removing all the unwanted features

Converted<-Converted == 14 & Converted==6 & Converted==209
Converted[Converted==0]<-NA

plot(Converted) # I would like this to show only the blue dot.


texts: Key authentication module Geoserver - authkey

texts:
I have a problem with this Geoserver plugin. I have done the following:

downloaded the plugin from
http://gridlock.opengeo.org/geoserver/2.2.x/community-latest/
extracted archive to /var/lib/tomcat6/webapps/geoserver/WEB-INF/lib/
restarted tomcat6 uncommented required line from /geoserver_data/security/authkeys.properties restarted tomcat6 again
and still I can access all WMS services without authkey in the URL

What am I doing wrong? I am probably forgetting something and can’t figure what. Please help? THX


texts:
I am currently working on a 3D mapping project. And I have created some KML files from ArcGIS Online, and I have edited them in Google Earth to adjust the 3D attributes. I am also designing a webpage with JavaScript, CSS, and HTML to display the results with 3D features. I have brought in the Google Earth API, and it properly displays what I want it to, but it is also very tricky to get anything to work in terms of creating functions to interact with those layers, rather I am just left to embedding them in the API. I know that this type of stuff is doable with ArcGIS Online and the scenes provided there, but I am also running into a lot of difficulty with uploading the layers/scenes from ArcGIS Pro to ArcGIS Online, and I would like to be able to have mroe control over what capabilities I am putting in there (which using JS API helps me to accomplish).

I really would like to stick to ArcGIS rather than Google, as I have more experience with using the ArcGIS JS API, and I was wondering if there is a 3D equivalent to ArcGIS JS API currently available, I have not been able to find anything yet, but I may be overlooking something out there, or searching for the wrong things.

Please help!


texts:
If I have a vector dataset with XYZ data and convert it with QGIS from EPSG:31468 to EPSG:4979, for example, and the heights (z) are based on EPSG:31468.

Will these values also be converted to the EPSG:4979 height system (i.e. based on the WGS84 ellipsoid) during the conversion?

My source data comes from the first EPSG while I want to display the converted data in Cesiumjs, which only supports the second EPSG, and the heights are relevant.


texts:
A lot of career websites allow you to select jobs from within a radius of your chosen location.

Since cities and towns can be considered irregular in shape, sprawling in the definition (do you discount outer suburbs etc.).

How is the centre of a city defined? And is it the same as the location used as described above?


texts: Groups the pixels into the nearest large class

texts:
I have clustered an image with 3 classes using the K_menas, but I would like the isolated pixel groups to be incorporated into the nearest large class. Is there a quick way to do this in the K_means parameter tuning, or do you need to change the classification algorithm? I have read that some use ISODATA, maybe it can be useful in my case, but I don't know.
with rasterio.open(raster)as ds:
    arr = arr[0,:, :]

K = KMeans(n_clusters=n_clusters, init='k-means++', n_init=10, max_iter=400, tol=0.0001,
           verbose=False, random_state=None, copy_x=True, n_jobs=1)
K.fit(arr)
X_clustered = K.labels_
point = X_clustered.shape[0]
X_clustered = X_clustered.reshape(point,1)


texts: Creating training areas for parametric classifcation

texts:
I am creating a signature file for a supervised classification for a subset of a LANDSAT image. I know what classes I would like to create, and have identified ideal areas for training. However, I am running into a problem in some of these areas in which the covariance matrix cannot be inverted, thus I am not able to run a parametric classification, such as Maximum Likelihood. I understand that the covariance matrix must be linear for this to happen and is a function of the pixels themselves and I cannot edit this. I am using layer-stacked LANDSAT 8 imagery with 14 bands (7 from spring imagery, 7 from winter imagery to detect leaf-on and leaf-off deciduous areas).

Is there a way to determine areas that would lend themselves to fit into a linear covariance matrix? I ran across literature that suggests that an increase in bands makes it more difficult to achieve this linear satisfaction for covariance matrix inversion. I am using ERDAS 2014. It has the option to run a non-parametric (such as parallelepiped) along with a parametric (max likelihood), but I would like to stick to maximum likelihood. Would larger training areas fix this? Or more smaller areas merged together? It seems that when adding new training areas it is random as to whether the cov matrix will be invertible or not. Thanks in advance!


texts: How to print an OpenLayers map with its popups open using Mapfish printing module?

texts: I'm using GeoExt and Mapfish printing module for geoserver. There is a functionality, in the application I'm developing, that allows users to create popups on the map. I'd like that printing functionality allows popups to be printed with the map. Is that possible with Mapfish ?


texts: Getting characteristic high points from DEM

texts:
I use QGIS 3.14 and I have a DEM file which I'd like to use to get a point layer with most characteristic high points from that DEM. So every top of the mountain, hill, hillock, low etc should be extracted from that DEM to a point layer with an attribute "height" of that point. Can you recommend me any solution? Even source code would be great.
By "characteristic" i mean:

the tops of mountains, hills and mounds
lowest point at every valley, basin
highest point at every watershed (might be hard to realise but i have also vector layer with rivers..)
at roads intersection (same as above, might be hard but i have that vector layer.)

The second rule is that per every 1 km² it should be around 3-6 points, no more.


texts: Changing font size of annotation class in ArcGIS Desktop?

texts:
I recently started using annotation feature classes to see what capabilities they offer, so I'm a newbie to them.

Previously I would add my labels in layout view when all is good and done with a map, but when a view needs to be changed, it would mess up labeling.  So that is what brought me to try the annotation feature classes.

I've created an annotation feature class with a few classes to it (ex runway, taxiway, highway, street, etc.).  I set the text symbol size initially and when I find I don't like that size, I go into ArcCatalog to change the size for that specific annotation class.  My thought is it would carry over to all text within that class.  However, it does not do that for me.    Any annotation that has been created goes under an "old" label.  For example, when I switched my Street annotation from size 12 to size 8 in ArcCatalog, any current Street annotation in the map went onto another annotation class that was automatically called "Street_Old".  So to change this text that I had already done, I have to highlight and change it back to the Street class to change the font size.

How do I change the font size so it changes automatically for all labels associated with that class and does not create this "_Old" class?

I am using ArcGIS Desktop 10.2.2 with a Basic License.


texts: How to make highlighted polygons when search on Leaflet map from qgis2web plugin?

texts: I´m using the version 3.12 and I´m trying to make that  when I search on the exported map  on "qgis2web" using Leaflet, it highlights the polygons that correspond with the searching. I´m novice using this plugin, so please be quite specific.


texts: Splitting legend in QGIS in to 2 columns

texts:
I am new to QGIS and I am making a map of India showing some quantity. I added the legend to my map layout, I have 10 symbols in my legend and it is very long , so I wish to split it into 2 columns. I tried it using the column option in legend properties. But when I increased the column number to 2 and apply the changes into my legend, only half part of the legend is seen. The second part is get cropped. Anybody please help to split the legend into 2 and show both halfs. I will add the screenshot below.

The changes when I tried to split into 2 ,I highlighted how I tried to change into 2 columns


texts: All UK census geographies in one file irrespective of scale

texts:
Is it possible to get all UK census geographies in one file?

For example, I am currently looking at a table of data titled

Life satisfaction by country, region, unitary authority and county, April 2011 to March 2012

Which turns out to mean unitary authorities or counties, or local authorities in Scotland.  Nevertheless each entry has its own code (e.g. E06000006).  Is there a single shapefile that will match all of these without having to download all the geographies separately?


texts: Getting ExecutionEngineException when using GetCount in loop

texts:
I've created a small ArcObjects app that searches for all matching files in a directory tree and returns the matched files with the number of rows in each file. When performing this loop on a small number of files it seems to process fine. However, when running it over a larger number of files (81), it throws a System.ExecutionEngineException error that I'm having a hard time tracking down. Here is the (allegedly) relevant code.

FYI: The Execute method is thanks to one of Kirk's excellent answers in another question about debugging, I believe.

Following the stack trace of

Multinet Merging Application.exe!Multinet_Merging_Application.shp2gdb_form.Execute(ESRI.ArcGIS.Geoprocessor.IGPProcess proc = {ESRI.ArcGIS.DataManagementTools.GetCount}) Line 460 + 0x10 bytes  C#
Multinet Merging Application.exe!Multinet_Merging_Application.shp2gdb_form.RunCount(string[] inFile = {string[81]}) Line 448 + 0x8 bytes  C#
Multinet Merging Application.exe!Multinet_Merging_Application.shp2gdb_form.SearchAndReport(int option = 1) Line 106 + 0x3a bytes  C#
Multinet Merging Application.exe!Multinet_Merging_Application.shp2gdb_form.searchBtn_Click(object sender = {Text = "Search and Report"}, System.EventArgs e = {X = 125 Y = 13 Button = Left}) Line 326 + 0xd bytes  C#

Here are the functions--

 private void searchBtn_Click(object sender, EventArgs e)
        {
            if (infolderselected == true && filetypeselected == true)
            {
                string dir = folderTextBox.Text;
                string[] oDir = Directory.GetDirectories(dir);
                if (oDir.Length > 500)
                {
                    if (MessageBox.Show("There are more than 500 directories under this root folder. Are you sure you wish to perform this operation?", "Warning", MessageBoxButtons.YesNo) == DialogResult.No)
                    {
                        searchrep.Text = searchrep.Text+ "\nSearch canceled!";
                        return;
                    }
                }
        }

private void SearchAndReport(int option)
    {
        string searchstr = type + shpdbf;
        string dir = folderTextBox.Text;
        try
        {
                searchrep.ForeColor = Color.Black;
                setOutputText(searchrep.Text + "\n" + sb.ToString());
                setOutputText(searchrep.Text + "\n--Performing Count. Please wait...");
                setOutputText(searchrep.Text + "\n--Total rows in listed files: " + RunCount(oFiles).ToString());
                setOutputText(searchrep.Text + "\n---------- Search and Report Complete at: " + DateTime.Now.ToString() + " ----------\n");
                searchBtn.Enabled = true;

        }
        catch (UnauthorizedAccessException error)
        {
            MessageBox.Show(error.Message + "\nThis process cannot continue.", "Error");
            setOutputText(searchrep.Text + "\nSearch failed due to error.");
        }
    }
private static int RunCount(string[] inFile)
    {
        finalcount = 0;
        GetCount count = new GetCount();
        foreach (string file in inFile)
        {
            count.in_rows = file;
            Execute(count);
            finalcount += count.row_count;
        }
        return finalcount;
    }
public static void Execute(IGPProcess proc)
    {
        Geoprocessor gp = new Geoprocessor();
        gp.AddOutputsToMap = true;
        gp.OverwriteOutput = true;
        gp.RegisterGeoProcessorEvents((IGeoProcessorEvents)new GPEvents());
        IGeoProcessorResult2 result = gp.Execute(proc, null) as IGeoProcessorResult2;
        IGPMessages msgs = result.GetResultMessages();
        for (int i = 0; i < msgs.Count; i++)
            Debug.Print("{0} {1}", msgs.GetMessage(i).Description, msgs.GetMessage(i).Type);
    }

Apologies, that seems excessive. The gist of the process is that I am attempting to use the Get Count GP tool to return the number of rows from all the files returned in the search. It runs peachy until the number of files returned gets to a certain point, and then it goes kaput.
Any ideas?


texts: Making subset of Sentinel-2 with rasterio?

texts:
I am looking for a solution to my problem with rasterio. I want to download just a specific part of the image from AWS based on lat lon coordinates. Rasterio offers a great tool (window) although it proivdes just pixel values and not geographycal coordiantes (latlon).

Is there any solution for my problem?

Here is my code:

red = 's3://sentinel-s2-l2a/tiles/34/T/DS/2017/10/12/0/R10m/B04.jp2'
nir = 's3://sentinel-s2-l2a/tiles/34/T/DS/2017/10/12/0/R10m/B08.jp2'

bands=[red, nir]
def ndvi_calc(bands):
    window = rasterio.windows.Window(100, 100, 100, 100) # How can I set here lat/lon values?
    red_band=bands[0]
    nir_band=bands[1]
    with Session:
        with rasterio.open(red) as red_band:
            red=red_band.read(1, window=window)
        with rasterio.open(nir) as nir_band:
            nir=nir_band.read(1, window=window)

    ndvi = (nir.astype(float) - red.astype(float)) / (nir+red)


texts: Delete duplicate records from feature class in Python

texts:
I saw some examples with for loop to find duplicates.

I had 1.5 million records in feature class, what would be fastest way to find duplicates and delete that particular rows.

Based on three fields, I'm trying to find duplicates[ID,S_NAME,TIME] using ArcPy Python 2.7.

EDIT-1:
Delete identical getting error with sde connection-

arcpy.DeleteIdentical_management (r"C:\FolderConn\DEV\master.sde\DatasetName\FeatureClassName", ["ID","S_NAME","TIME"])`


texts:
I am (very) new to GIS in general. I would like to convert a satellite image of a province into raster, assigning each pixel its own data (i.e., a number and a name) based on its html color code. For example, shades of red and/or blue would be assigned the data "id=1; name=Residential" and shades of light green would be "id=2; name=grassland" and green to dark green would be "id=3; name=forest".

I'm thinking of taking screen shots (zoomed in as much as possible) and then combining them to capture the whole province. Next, I'll georeference the image.

How could I go about this next?


texts: Is it possible to change or reconfigure a domain owner on SDE?

texts:
We have had to reconfigure a few things within a geodatabase managed in SDE. This inadvertently has caused some issues updating the domains within the gdb, to the point where not even the DBO can update a domain, only the originator. The DBO and owner are actually the same user/person now. Does anyone know of a method to change the domain owner to the DBO, rather than the listed user that originated the domain? We are trying to get around the message...

"Failed to update one of the domains. Must be the owner to perform this operation."


texts: matching coordinate systems by transforming and projecting

texts:
matching coordinate systems of all my data (transforming, projecting)

My data uses the GCS_Bessel_1841 coordinate system. Would there be any inaccuracies if i would transform it to DHDN 3-Degree Gauss Zone 3? Or should i rather match all my other data to the DEMs coordinate system?

How would i do that?


texts: Options to improve conversion raster to polygon?

texts:
Is there a way to convert a raster.tif into a single polygon?

# Loading raster.tif
DEM <- raster ("C: /Users//raster.tif")
# We generate boundary of the raster
poly_DEM <- DEM> -Inf
polygon <- rasterToPolygons (poly_DEM, dissolve = TRUE)***

I have tried rasterToPolygons if it works, but for longer rasters it takes too long, I think it is because it has to dissolve many pixels.

In summary, I just want to extract the edge of a raster in a single polygon.


texts: Wrong with adding time series in the same chart in Google Earth Engine

texts:
Recently I've been trying to print time series for different regions and show them in the same plot. Everything goes fine if I exclude the 4th region into the time series (1-3 time series are fine). But when I add the 4th region, time series points at other locations are gone. The output plot only shows the time series of the last region (only remains 2 points of other regions).

So the problem is: how to add the time series of the 4th region into the plot?

Excluding the 4th region in the feature collection would work properly, but I don't want to omit it.

Here's the link to the code.

var timeField = 'system:time_start';
// Define test regions.
var regions = ee.FeatureCollection([    // Test Regions are in Hunan, China.
  ee.Feature(    // Region 1
    ee.Geometry.Rectangle(113.201806473291, 27.6867096858736, 113.251537728778, 27.7276180415286), {label: 'Region1'}),
  ee.Feature(  // Region 2
    ee.Geometry.Rectangle(113.152404668942, 27.6389241071109, 113.177553060666, 27.6780758597365), {label: 'HengJiangCun'}),
  ee.Feature(  // Region3
    ee.Geometry.Rectangle(113.196381164761, 27.7769042331009, 113.216382570871, 27.7999150021881), {label: 'ChenJiaBa'}),
  ee.Feature(  // Region4 Adding this region would omit most of the points
    ee.Geometry.Rectangle(113.038857564920, 27.9581880215077, 113.042355165475, 27.9719283228899), {label: 'ZhaoShanXiang'}),
]);

// Add base maps.
Map.addLayer(regions, {}, 'Test Regions'); // Layer: Test Regions.
print(regions)

// Load Sentinel-2 TOA reflectance data.
var masked = S2
    .filterDate('2017-01-01', '2017-12-31')
    .filterBounds(regions);
print(masked);

// Function to calculate NDVI in S2 images.
function S2NDVI(image){
  var bands = image.select('B4','B8');
  var NDVI = bands.expression('(b(1)-b(0))/(b(0)+b(1))').select(['B8'],['NDVI']);
  //return NDVI;
  return image.addBands(NDVI);
}

// Calculate S2 NDVI.
var S2VI = masked.map(S2NDVI);
print(S2VI);

// Generating 25 percentile time series chart.
var reducer = ee.Reducer.percentile([25]);

// Time series chart settings.
var TSChart = ui.Chart.image.seriesByRegion(
    S2VI.select('NDVI'), regions, reducer, 'NDVI', 20, 'system:time_start','label')
        .setChartType('ScatterChart')
        .setOptions({
          title: 'CI-Red Edge time series at 4 regions',
          vAxis: {title: 'CI Red Edge'},
          lineWidth: 1,
          pointSize: 4,
          series: {
            0: {color: '008744'}, // Region 1
            1: {color: '0057e7'}, // Region 2
            2: {color: 'd62d20'},  // Region 3
            3: {color: 'ffa700'}, //Region 4
}});

print(TSChart);
Map.addLayer(S2VI.select('NDVI'),{},'Sentinel VI Collection');


texts: creating a view of multiple tables using Postgresql

texts:
This might be the stupidest question for postgresql user. Please do understand that i just started to use postgresql.

what i am trying to do is joining 3 tables using same key field (which is polyong_id)

keyfield of table 1 = polygon_id, and it has additional fields that i want to keep
keyfield of table 2 = polygon_id, and it has venue_id field that i want to append
keyfield of table 3 = polyong_id, and it has feature_type field that i want to append

so what i did is

Create View BULD as
Select t1.shape, t1.polygon_id, t1.name, t1.height, t1.ground_clearance,
       t1.iso_country_code, t2.venue_id, t3.feature_type, t3.main_feature_type
From
       BULD_base as t1, Venue as t2, Feature as t3
Where
       t1.polygon_id = t2.polygon_id and t1.polygon_id = t3.polygon_id

BULD's total number of records is about 60K and the result returned around 150 records.
I changed the AND to OR and it gave me same result.

I may create a view between t1 and t2, then
and create another view between the table from the result above and t3 for the final result.

but I'd like to make a neat and smart query for this.

can anyone help?


texts: Correcting addresses with Python or Model Builder

texts:
I am trying to make a script in Python (or Model Builder in ArcGIS) for correcting addresses.

I have a point shapefile with the addresses (it has X and Y, and fields for "number", "street name", and a unique code for locality), and
a line shapefile with the name of the streets and also the unique locality code (we assume that there are no errors in the streets shapefile).

Example of errors I have to find:

in the same locality I can't have two points with the same number on
the same street
no major leaps between numbers (for example the number 32 between 12
and 16)
only even numbers on a side of the street and only odd numbers on the
other side
isolated numbers on the same street (ex: 2-3 km away from the next
point on the same street)
etc.

I'll have to divide it into many sub-problems and try to fix them separately.

Do you have any suggestions about how I should proceed?


texts:
I have been provided with a number of blank shapefile to use as a schema for Work As Executed data. I have a survey in DWG format, from which I export data as a CSV file.

Is it possible to import this directly into an existing shapefile layer?  I can only see that you can create a new layer, but would like to avoid this step if possible.  Ideally I would be able to map some of the attributes also, but even the geometry would be useful.


texts: Understanding ordnance grid reference such as TQ38-7152-9751?

texts: I have a list of OS grid references, such as TQ38-7152-9751, which I assume is at a resolution level below 1km, but cannot find a grid for this or a way of converting to eastings and northings?


texts: How to select 200 random points within a 100m radius in Postgis

texts:
I'm currently learning to use PostgreSQL and PostGIS. I would like to know how to write the SQL that would allow me to select 200 random points within a radius of 500m from a user defined centre point. Currently I have the database containing all my records, with an individual row for ID, latitude and longitude. The co-ords are in WGS84. The 'points' refer to shops.

Once I have the 200 randomly selected points, I would then like to edit a further two coloumns ("Status" and "Time") to show the selected points are closed and the time which they closed.

What would you suggest is the simplest way to perform this?


texts: What's the most common way to join photos with geographic data?

texts:
Increasingly spatial data will contain a column that has URL's pointing to images hosted on Flickr, Dropbox or some other cloud based storage system.

What's the most common/efficient way to join photos with geographic data? Are there any fast ways to do it in bulk, or do most people just enter them manually?


texts:
I try to add a function to the qt event loop via QTimer.singleShot(...) from within a QgsTask. However the function is never executed. I wrote a short script that isolates this issue:

import time
from functools import partial

class TestTask(QgsTask):

    def __init__(self):
        super().__init__('test task', QgsTask.CanCancel)

    def run(self):
        time.sleep(1)
        QTimer.singleShot(10, partial(QgsMessageLog.logMessage, 'this message will NOT show up in the message log'))
        time.sleep(1)
        return True

def start_task():
    task = TestTask()
    QgsApplication.taskManager().addTask(task)

    return task

t = start_task()
QTimer.singleShot(10, partial(QgsMessageLog.logMessage, 'this message will show up in the message log'))

My guess is this happens due to the fact that threads do not have their own event loops, but I could be mistaken. I've looked around to find a way to add an event loop to a QgsTask / to add a function to the main event loop but could not find anything.

How can I schedule the execution of a function in a QgsTask without relying on QTimer?


texts: Replacing shapefiles in MXDs using only Python options available prior to ArcPy?

texts:
I am trying to come up with a python script that will go through a designated folder, locate all MXD files within the folder, and then check each MXD to see if it contains a specified shapefile, which should be replaced with another shapefile if it is found. I am working in ArcGIS 9.3 so arcpy.mapping is not available to me.

I have heard that ArcObjects might help but I don't know that much about it and am not sure if it would be feasible with my modest python experience.

Any Suggestions?


texts: Lapply function on raster iterating through SPDF polygons in R

texts:
I have a raster DTM and a SPDF. I'd need to mask a different elevation interval of the DTM within each overlaid polygon.

I tried a for loop to iterate through the polygons, but as one might expect this process requires a huge amount of memory (not to mention the hours) and after 2 or 3 iterations the loop stops for errors

 # open DEM and SPDF
 dem <- raster('dem.tif')
 SPDF <- readOGR(dsn='foo', layer='SPDF')

 # open elevation df
 df <- read_excel('elevation_intervals.xls') # the df is structured as the SPDF@data

 # loop
 for(c in 1:nrow(SPDF)) {
   dem_masked <- mask(dem, SPDF[c,])
   min_elev <- df[c,]$elevation-150
   max_elev <- df[c,]$elevation+150
   reclass <- data.frame(from=as.numeric(c(minValue(dem_masked), min_elev, max_elev)),
                         to=as.numeric(c(min_elev, max_elev, maxValue(dem_masked))),
                         becomes=as.numeric(c(NA, 1, NA)))
   dem_masked <- reclassify(dem_masked, reclass)
   assign(paste('DEM_', c, sep=''), dem_masked)
   rm(dem_masked) # I tried to gain memory
   gc() # idem
 }

This is the error message

 Error in matrix(unlist(ini), ncol = 2, byrow = TRUE) :
    'data' must be of a vector type, was 'NULL'

and probably relates to memory issues

I guess that using an *apply function would be much effective but I can't find a way how write properly the function since the iteration relates to the argument of mask and not to the x, i.e. the DEM

Otherwise, any other suggestions?


texts:
Does anyone know how the instanceData property of a Geocortex request is encoded?

The entire request is JSON- key-value pairs, arrays, quoted strings, etc. - but the instanceData value is encoded. It looks similar to Base64 until you look closely and see slashes (/) and plus signs (+) throughout the string. I tried various online Base64 decoders but those I tried have been unable to decode the string.

Here's a very small piece of what I found with Fiddler. The entire request is 2+ MB.

"instanceData":"Rrfqn+xQO2IEFTMoYXBfnZQU9FlQf8fCqfRUvCMjc78zthft4xWTh/sIaKW0ySzVnpbSDAog293K+XX32jR/20bxw51t


texts: Meter Coordinates to decimal Coordinates

texts:
I have an excel file in which I have XY coordinates, example on the attached immage
that is the form
XCoord           |      YCoord
15822.426100     |      -2747650.026000

I need them in the form similar to this
XCoord           |      YCoord
-24.7236869      |      31.1995142

Please help, Im using ArcGIS 10.2


texts: How to load point data in JSON format from an URL into Qgis using the python console?

texts:
Via an URL: http://www.odaa.dk/api/action/datastore_search?resource_id=c3097987-c394-4092-ad1d-ad86a81dbf37 I have what appears to be a file in JSON format that contains points which I want into QGIS directly from the URL, preferably with the Python console. Within the JSON file columns for longitude and latitude are present.

Do any of you know how this can be done?


texts:
There is one map and one dataset in my account. The dataset is synced to a Google Sheet to which my students add data.

All the layers of the map are derived from the one synced dataset.
One of the layers is an animated timeline. The timeline worked fine before syncing, but now that it is synced, the data type of the "date" column is set as "string", and I'm not able to change it.

I found this answer on Stack Exchange, but the solution it refers to is a bit beyond my SQL ability (setting up a dynamic query I could probably manage, but the linked was reference to dates in Postgres is opaque to me at this point.):
How do you change the data type of a column for a synced table?

Is there a way to maintain a "date" data format in a synced table?


texts: Import pure (regular) JSON file to Leaflet

texts:
I would like to add pure JSON file to Leaflet. Saying "pure" I mean plain JSON structure that is not GeoJSON compliant.

There are lot of examples advicing how to add GeoJSON files:

Leaflet load JSON data
Load polygon JSON to Leaflet
Adding Multiple JSON Layers to Leaflet Map
Feature Layer from JSON file

My JSON data looks like this:

  {
 "Sheet1": {
"-430524279882": {
  "MDU_ID": -430524279882,
  "Status": "Invoiced",
  "THP": 12,
  "Latitude": 52.41635,
  "Longitude": -1.55227
    },
   "-430532279893": {
  "MDU_ID": -430532279893,
  "Status": "Invoiced",
  "THP": 12,
  "Latitude": 52.41635,
  "Longitude": -1.55227
},
"-430622279798": {
  "MDU_ID": -430622279798,
  "Status": "Invoiced",
  "THP": 12,
  "Latitude": 52.41555,
  "Longitude": -1.55155
},
"-430519279894": {
  "MDU_ID": -430519279894,
  "Status": "Invoiced",
  "THP": 12,
  "Latitude": 52.41635,
  "Longitude": -1.55227
},
"-430558279922": {
  "MDU_ID": -430558279922,
  "Status": "DO NOT DESIGN - SDN",
  "THP": null,
  "Latitude": 52.41635,
  "Longitude": -1.55227
}
...

It comes straight from Google Sheets. Unfortunately this format is not recognized by L.geoJson function:

  job2 = L.geoJson(data2, {  //layer2 Virgin Media start
    style: style2,
    pointToLayer: function(feature, latlng) {
            feature.properties.myKey = '<b>'+ feature.properties.Owner +
  '</b>; ' + feature.properties.Address
            label = String(feature.properties.Owner)
            if (feature.properties.Post_Survey_Home_Count >=100)
        return L.circleMarker(latlng, {

How can I use this data in a similar way as GeoJSON?


texts: Unicode Decode Error for shapefile in Python

texts:
I have Python script for converting shapefiles to Geojson. My code works well if I use some usa data, but when I want to use some other encoding, it gives me this error:

  Traceback (most recent call last):   File
  "C:\Python33\lib\tkinter__init__.py", line 1475, in call
      return self.func(*args)
  File "C:/Users/name/PycharmProjects/desktopApp/core.py", line 105, in
  shp2geojson
      for sr in reader.shapeRecords():
  File "C:\Python33\lib\shapefile.py", line 553, in shapeRecords
      for rec in zip(self.shapes(), self.records())]
  File "C:\Python33\lib\shapefile.py", line 525, in records
      r = self.__record()   File "C:\Python33\lib\shapefile.py", line 501, in __record
      value = u(value)   File "C:\Python33\lib\shapefile.py", line 60, in u
      return v.decode('utf-8')

  UnicodeDecodeError: 'utf-8' codec can't decode byte 0xe1 in position
  10: invalid continuation byte

How to fix this problem?

  import shapefile
  import json

  def shp2geojson():

        filename = askopenfilename()
        reader = shapefile.Reader(filename)
        fields = reader.fields[1:]
        field_names = [field[0] for field in fields]
        buffer = []

        for sr in reader.shapeRecords():
            atr = dict(zip(field_names, sr.record))
            geom = sr.shape.__geo_interface__
            buffer.append(dict(type="Feature", \
                               geometry=geom, properties=atr))

        geojson = open('shp.json', "w")
        geojson.write(dumps({"type": "FeatureCollection", \
                             "features": buffer}, indent=2) + "\n")
        geojson.close()


texts: Keep marker constant on Google Maps Javascript

texts:
I want to be able to keep a marker constantly in the center of Google Maps on my mobile application.

I have tried using bindTo functionality of Google Maps, which binds the marker to the center of the map, but it seems a bit clumsy.

I have tried doing it using HTML and CSS also, but Google Maps does not allow it.

The same as in Uber app.


texts: How to mark a stored vector layer as exported

texts:
Through PyQGIS using QgsVectorFileWriter.writeAsVectorFormat() we can store a vector layer as .gpkg, is great!
However, once we close the project and reopen it again, the layer name is showed but their associated gpkg data are lost.
QGIS provide a way to permanently store a memory layer on disk, but it seems this is a manual process through the QGIS tool environment (context menu on every layer).
I'd like to export the layer using writeAsVectorFormat() but as next stemp mark it as exported, so it would be available in the future even if I close the project.
Therefore, the 'memory icon' would disappear too.
Does anybody know it this is feasible using PyQGIS?