texts:
Is it possible to interpolate points between multiple points on a drawn line in leaflet? I am creating a site where users can draw running routes and see an elevation chart, but currently only where the users clicks is processed - this makes a sharp looking elevation chart. I was looking to interpolate between points so the elevation chart has more points to work with and therefore looks smoother.
I have been looking at Leaflet.GeometryUtil and interpolateOnLine, but it is a bit above me and after a few failed attempts - I'm not sure how to go about it.
Where the drawing takes place:
map.on(L.Draw.Event.CREATED, function(event) {
    var layer = event.layer;
    console.log("layer", layer.getLatLngs());

    drawnItems.addLayer(layer);
    drawData = layer.toGeoJSON();

    var profile = drawData.geometry.coordinates;

    //reverse coordinates
        function reverseCoordinates(_arr){
            let i;
            let number = []
            let xx
            for(i=0;i<_arr.length;i++){
                xx = _arr[i].reverse()
                number.push(xx)
            }
            return number;
        }
        fetch('https://open.mapquestapi.com/elevation/v1/profile?key=tHXSNAGXRx6LAoBNdgjjLycOhGqJalg7&shapeFormat=raw&latLngCollection='+reverseCoordinates(profile))
          .then(r => r.json())
          .then(data => {
            var latlngs = [];
                for (var i=0; i<data.shapePoints.length; i=i+2) {
                    //console.log(data.shapePoints[i])
                    latlngs.push([data.shapePoints[i+1], data.shapePoints[i], data.elevationProfile[i/2].height]);
                }

            geojson = {
                "type": "Feature",
                "geometry": {
                    "type": "LineString",
                    "coordinates": latlngs
                }
            };


texts: Updating GPS coordinates to overlapped locations using QGIS?

texts:
If I have two many locations overlapped each other, and I want to update their coordinates to make the view clear and be able to spot the locations?

How can I do this easily?

I have more than 512 locations.

I don't want to zoom in to the district to have the view clear, I want to make the points far a little from each other to make the view clear.

Is there a quick way to do that?


texts:
I have successfully created a custom dialog box which uses an array to read in a list of developer names from a table, and populates the selected record based on the row number pulled in from the developer table.

Table is
Developer ID   Developer Name
1              A developer
2              B developer
3              C developer
4              D Developer
5              AAA Developer
6              C2 developer

I have also had success in then populating the selected value from the popup menu back to the table.  At present the developer list contains about 400 names.  As new names are found, they need to be added to the table. In above example- ID 1 to 4 were added initially then 5 and 6 later and they aren't sorted alphabetically.

My question is does anyone know if there any way that I can I get the table from the array in alphabetical order without upsetting the row numbers from the developer table?  If I add say another 50 developers to the bottom of the table then it becomes hard to find them as they are not sorted alphabetically - whilst the first 400 I put in are.

Given that the popup menu reads in the row number and assigns this to the value in the popup menu then writes it based on the value selected in the popup, I am unsure whether this can be done.

iRows1 = TableInfo(developer, TAB_INFO_NROWS)

'redim developername array based on number of rows in table
ReDim sDeveloperName1(iRows1)

For iLoop1 = 1 to iRows1
    Fetch rec iLoop1 from developer
    sDeveloperName1(iLoop1) = developer.developername
Next

'Developer details groupbox
control groupbox
title "Developer Details"
Position 5,206 width 130 height 97
control statictext position 10,216 title "Primary developer"
control popupmenu position 10,228 width 120 ID 20
Title from Variable sDeveloperName1 into iDeveloperID1out

Sub ReadIn_PopUpMenu_Items
    'converts table value into the popupmenu item for primary developer listed
     dim developerID1temp as integer
     developerID1temp=selection.developerID1
     if developerID1temp=0 then
         alter control 20 value 1
    Else
        alter control 20 value developerID1temp
    End if
End

'if 1st value in popup menu is selected then developer is unknown so write the 'value back to the selected record as 0
If iDeveloperID1out=1 then
    Update Selection set col20=0
else
    Update Selection set col20=iDeveloperID1out
End If


texts:
I have a vector line, but it has no z-dimension.

I have some points along the line, which have the elevation included.

Can I add the elevation from the points to the line at that specific point?

QGIS 3.4.4


texts: UK regions is incomplete in Ordnance Survey BoundaryLine data

texts:
This question is about United Kingdom Ordnance Survey(OS) open Mapsets.

I have downloaded bdline_essh_gb from OS Open maps. But I don't see the boundaries for all the regions there - I am only interested in English regions currently. The one missing mainly are Midlands regions.

Does anyone else has the same problem?


texts:
I'm new with PyQGIS.
I need to load a raster layer from gpkg file, and i tried this way:

gpkg_raster_land_layer = "GPKG:" + path_to_gpkg + ":" + layername

r_layer = QgsRasterLayer(gpkg_raster_land_layer, layername, "gdal")
if not r_layer.isValid():
    print("Land layer failed to load!")
else:
    QgsProject.instance().addMapLayer(r_layer)

But it didn't work: Land layer failed to load!
I didn't understand how to proceed to load the raster, can someone help me?


texts: Missing QGIS python plugin database connectors

texts:
I am working with python plugins for QGIS.I developed my plugin and inorder to incooperate the same into QGIS i copied the folder in

C:\Documents and Settings\comp90\.qgis\python\plugins

I opened my QGIS exe and installed the same.It is shown under **Databases-->PG & SL.
My code is as follows:

SUPPORTED_CONNECTORS = {'postgis':'1.5', 'spatialite':'1.6'}
MISSED_CONNECTORS = []

def listDatabases(self):

    actionsDb = {}
    actions = DbConnection.ConnectionManager.getAvailableConnections()

    for a in actions:
        actionsDb[ unicode(a.text()) ] = a

    return actionsDb

def getAvailableConnections(self, conntypes=None):

    if conntypes == None:

        conntypes = ConnectionManager.SUPPORTED_CONNECTORS.keys()
    if not hasattr(conntypes, '__iter__'):
        conntypes = [conntypes]

    connections = []
    for c in conntypes:
        connection = self.getConnection( c )
        connections.extend( connection.getAvailableConnections() )
    return connections

def getAvailableConnections(self):
    connections = []

    settings = QSettings()
    settings.beginGroup( "/%s/connections" % self.getSettingsKey() )
    keys = settings.childGroups()
    for name in keys:
        connections.append( Connection.ConnectionAction(name, self.getTypeName()) )
    settings.endGroup()

    return connections

but the actionsDb is blank.In QGIS, the connection is established. When I try to open the plugin in QGIS ,it does not open because actionsDb is blank. What can be the problem??


texts: wxpython frame run from ArcMap loses focus almost immediately after it's shown

texts:
I'm trying to prepare a Python add-in using "Python Add-In Wizard". I want my add-in to perform some geoprocessing but first it should take data from a user and there's already a problem. I'm displaying a wxpython GUI and it loses focus and hides behind the arcmap window almost immediately after it's shown. I have no idea how to debug it and what could cause the problem. I decided to go this way after getting hands on this example.

My code is here:

addin.py

from dialogs.import_dialog import ImportDialog

class ImportButton(object):
    """Runs the import function."""
    dlg = None

    def __init__(self):
        self.enabled = True
        self.checked = False

    def onClick(self):
        if self.dlg is None:
            self.dlg = ImportDialog()
        else:
            self.dlg.Show()

import_dialog.py

class ImportDialog(wx.Frame):
    def __init__(self):
        wxStyle = wx.CAPTION | wx.CLOSE_BOX | wx.SYSTEM_MENU
        wx.Frame.__init__(self, None, wx.ID_ANY, "some msg", style=wxStyle, size=(350, 200))
        self.panel = wx.Panel(self, wx.ID_ANY)

        self.path_label = wx.StaticText(self.panel, wx.ID_ANY, u"some msg", pos=(8, 8))
        self.input_path = wx.TextCtrl(self.panel, wx.ID_ANY, value="", pos=(8, 36), size=(330, 21))

        self.path_label = wx.StaticText(self.panel, wx.ID_ANY,
                                        u"some msg",
                                        pos=(8, 65))

        self.group_layers = wx.ComboBox(self.panel, wx.ID_ANY, pos=(8, 93), size=(330, 21))

        self.ok_button = wx.Button(self.panel, label="OK", pos=(8, 130))

        self.Bind(wx.EVT_CLOSE, self.OnClose)

        self.Center()
        self.Show(True)

    def OnClose(self, event):
        self.Show(False)


texts: Search API for PSScene4Band/PSScene3Band/PSOrthoTile not return item ID

texts:
I'm trying to search and then download clips of some shots across California, but am running into an issue. It's potentially the same issue as Assets not available via API, but that question doesn't have enough detail for me to know.

Using the API, I cannot get an image IDs for images in california that I can otherwise download with /explorer, for example item 20171218_181808_102f. I can get image IDs when I query the REOrthoTile item type, but the search API returns empty feature arrays for PSOrthoTile, PSScene4Band, or PSScene3Band queries, both analytic and visual. I think my permissions should let me access the other item types since it's in california. I checked my filters to make sure it matches my /explorer query, and I think it's the same.

Can you help me figure out what's going on?

My search request for PSScene4Band returns something like this:

{'_links': {'_first': 'https://api.planet.com/data/v1/searches/b0acb989064a4c5fa5da0b13c86c07a3/results?_page=eyJxdWVyeV9wYXJhbXMiOiB7fSwgInNvcnRfcHJldiI6IGZhbHNlLCAicGFnZV9zaXplIjogMjUwLCAic29ydF9ieSI6ICJwdWJsaXNoZWQiLCAic29ydF9zdGFydCI6IG51bGwsICJzb3J0X2xhc3RfaWQiOiBudWxsLCAic29ydF9kZXNjIjogdHJ1ZX0%3D', '_next': None, '_self': 'https://api.planet.com/data/v1/searches/b0acb989064a4c5fa5da0b13c86c07a3/results?_page=eyJxdWVyeV9wYXJhbXMiOiB7fSwgInNvcnRfcHJldiI6IGZhbHNlLCAicGFnZV9zaXplIjogMjUwLCAic29ydF9ieSI6ICJwdWJsaXNoZWQiLCAic29ydF9zdGFydCI6IG51bGwsICJzb3J0X2xhc3RfaWQiOiBudWxsLCAic29ydF9kZXNjIjogdHJ1ZX0%3D'}, 'features': [], 'type': 'FeatureCollection'}

Searching for REOrthoTile gets me a big features array filled with stuff, including the id key that I'm looking for.

Here's the search code:

https://github.com/earthmancash2/planet-satellite-fetching/blob/master/planet_image_fetching.py#L95-L148

https://github.com/earthmancash2/planet-satellite-fetching/blob/master/PlanetAPIParams.py

Any ideas for what is happening?


texts:
I have DEM file for a small mountain created from TIN interpolation of points. The DEM looks correct on the flat QGIS working area. The moment I try to visualize a 3D of the DEM using threejs, I get infinitely high column without the mountain surface.

I use Balkan Zone 7 as CRS, and my scale and shift option for the threejs plugin are set like this:
-Base size: 100
-Vertical exaggeration: 1.5 (default, I tried '0', but I get flat terrain)
-Vertical shift: 0


texts:
I have a few hundred folders, each with a set of shapefiles that are identical across folders (ie each folder has a 'lot.shp').  They are all supposed to be polygon features, but unfortunately many are line features.

I am trying to convert all line features to polygons using Feature to Polygon.

To do this I have tried the following code which should work through each folder converting the desired shapefiles and saving back into the same folder but with a '2' in the name:

import arcpy, glob, os

filelist = ['buildingpoly','lot', etc etc]

path = "C:/folder_containing_desired_subfolders/*/"

for x in filelist:
    list = glob.glob(path + x + ".shp")
    arcpy.FeatureToPolygon_management(in_features=list, out_feature_class=path + x + "2.shp", cluster_tolerance="", attributes="ATTRIBUTES", label_features="")
    del list

But this gets the following error:

  Runtime error  Traceback (most recent call last):   File "",
  line 10, in    File "c:\program files
  (x86)\arcgis\desktop10.4\arcpy\arcpy\management.py", line 2524, in
  FeatureToPolygon     raise e ExecuteError: ERROR 000210: Cannot create
  output C:/folder_containing_desired_subfolders/*/buildingpoly2.shp
  Failed to execute (FeatureToPolygon).

Any ideas? I'm at a loss on this one.

An alternative that would also work is if I can have all the new polygon features saved in the same directory but this would require each getting a unique name.

EDIT:

Here is code that works.  It works through the folders and converts all polylines to polygons.

import os, sys, arcpy
directory = "C:/working_folder/"

for path, path_names, data_names in arcpy.da.Walk(directory, type="Polyline"):
    for data_name in data_names:
        arcpy.FeatureToPolygon_management(in_features= (os.path.join(path, data_name)), out_feature_class= (os.path.join(path, "two"+data_name)), cluster_tolerance="", attributes="ATTRIBUTES", label_features="")


texts: Finding residential buildings with overpass-turbo?

texts:
I am looking to find all residential buildings with overpass-turbo.eu.
I am trying to use the following query:

way
  ["building"~"residential|house|terrace|detached|apartments"]
  ({{bbox}});
(._;>;);
out;

However, not many residential buildings are found in this way.

Which tags should be used to find the buildings through overpass-turbo, which may be multi-use but also residential (residential and commercial, or residential and retail...)?


texts: ArcGIS JSAPI: How to select graphic from Graphics layer

texts:
In ArcGIS JSAPI there are ways to select Features from FeatureLayer or using querytask like this or this or this. But how to select the Graphic from Graphics Layer?

I am trying to do this in AGS JSAPI which is accomplished by  SelectInteraction in OpenLayers. This lets user select features from map by clicking on them. Base on other interactions of OpenLayers such as (Draw, Modify) AGS JSAPI also provides toolbar which let user Draw or Edit Graphics but lacks in selection.

Is there any similar functionality available to JSAPI as well?


texts: How to install Spatialite-GIS on Windows?

texts:
I was wondering if there is a way to install Spatialite-gis on a Windows & machine?  I have found may posts and links about this, but I have yet to find an installer for spatialite-gis.  This Google Group post links to different sites, but there are not installers available, am I missing something?
https://groups.google.com/forum/#!topic/spatialite-users/5qHN1oUq5oQ

This is the direct link from GAIA to spatialite-gis, but there are no installers.
https://www.gaia-gis.it/fossil/spatialite_gis/index

This is a link to other downloads, including spatialite-gui, but nothing for spatialite-gis
http://www.gaia-gis.it/gaia-sins/windows-bin-amd64/


texts: How to join CSV data to SHP without ArcGIS?

texts:
I'd like to create a chloropleth map displaying quantitative data for an upcoming event. Using TileMill, MapBox.

I am having trouble joining my data to a .shp because the #'s correspond to COUNTY level data. I have the .shp for every county but have had absolutely no luck trying to georeference it or join it to the data without ArcGIS or a comparable program (on a Mac and work won't let me DL any additional programs)!

Or, should I try and convert my .shp to .csv to manually join?


texts: GeometryEngine Offset or GeodesicBuffer - to get minimal corners.net

texts:
We have a a polygon which might have between 8 and 20 points, which we need to draw a boundary around, that is at least 1km around at all points.  The restriction is that we can't increase the the number of vertices dramatically (or better still at all).

If I use GeometryEngine.GeodesicBuffer() to do this, it generates rounded corners at external vertices and hence greatly increases the number of vertices.

So I used GeometryEngine.Offset(). The Offset() has an OffsetType, which allows me to specify Square, and this returns basically the same number of vertices as I started with. However during testing it became obvious that the offset function has issues with the WGS84 spatial reference. Our East-West coordinates are OK, but any coordinates in the North-South plane are extended by 25% (which may be a factor of the latitude our tests are at).

So, the question is,

Can I use GeodesicBuffer in a way that I can stop the corners being rounded, OR do I keep using Offset, but convert/project into a different spatial reference (and if so which one)?

Note: Using the .Net framework of arcgis.


texts: Navigating an Object of class 'SpatialPolygonsDataFrame' to filter on Polygon Area

texts:
I have an Object of class 'SpatialPolygonsDataFrame' with the structure shown below.
What's the secret to getting to the @area slot?  Ultimately, I'm trying to filter this file so that I'm only left with polygons >= 64m^2.  The resultant shapefile will be imported to OCAD (orienteering mapping software), in an attempt to semi-automatically map a class of vegetation.
I thought I had some understanding of S4 classes, but clearly I have a gap.
among many similar failed attempts:
> xfile@polygons$Polygons@Polygons$Polygon
Error: trying to get slot "Polygons" from an object of a basic class ("NULL") with no slots
> xfile@polygons$Polygons@Polygons
Error: trying to get slot "Polygons" from an object of a basic class ("NULL") with no slots
> xfile@polygons@Polygons
Error: trying to get slot "Polygons" from an object of a basic class ("list") with no slots
> xfile@polygons$Polygons
NULL
> xfile@polygons$'Polygons'@Polygons$'Polygon'@area
Error: trying to get slot "Polygons" from an object of a basic class ("NULL") with no slots
> xfile@polygons$'Polygons'@Polygons
Error: trying to get slot "Polygons" from an object of a basic class ("NULL") with no slots
> xfile@polygons$'Polygons'
NULL
> xfile@polygons$area
NULL
> xfile$polygons
NULL
> xfile@polygons
( the whole file )

Formal class 'SpatialPolygonsDataFrame' [package "sp"] with 5 slots
  ..@ data       :'data.frame': 1 obs. of  1 variable:
  .. ..$ layer: num 1
  ..@ polygons   :List of 1
  .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
  .. .. .. ..@ Polygons :List of 191
  .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
  .. .. .. .. .. .. ..@ labpt  : num [1:2] 513021 4319805
  .. .. .. .. .. .. ..@ area   : num 9.64
  .. .. .. .. .. .. ..@ hole   : logi FALSE
  .. .. .. .. .. .. ..@ ringDir: int 1
  .. .. .. .. .. .. ..@ coords : num [1:19, 1:2] 513018 513025 513025 513024 513024 ...
  .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
  .. .. .. .. .. .. ..@ labpt  : num [1:2] 513053 4319805
  .. .. .. .. .. .. ..@ area   : num 1.85
  .. .. .. .. .. .. ..@ hole   : logi FALSE
  .. .. .. .. .. .. ..@ ringDir: int 1
  .. .. .. .. .. .. ..@ coords : num [1:5, 1:2] 513051 513054 513054 513051 513051 ...

......

.. .. .. .. .. [list output truncated]
  .. .. .. ..@ plotOrder: int [1:191] 18 118 101 168 110 148 61 135 166 90 ...
  .. .. .. ..@ labpt    : num [1:2] 513181 4319787
  .. .. .. ..@ ID       : chr "0"
  .. .. .. ..@ area     : num 38206
  ..@ plotOrder  : int 1
  ..@ bbox       : num [1:2, 1:2] 512793 4319369 513380 4319807
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "x" "y"
  .. .. ..$ : chr [1:2] "min" "max"
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
  .. .. ..@ projargs: chr "+proj=utm +zone=13 +datum=WGS84 +units=m +no_defs"


texts: Error QODBC driver not loaded when trying to connect to an MSSQL database from QGIS 3.4.3 (Ubuntu)

texts:
I have installed QGIS Desktop (latest one: 3.4.3 Madeira) along with the unixODBC-dev for SQL Server 2017 (following instructions here). I used sqlcmd from msql-tools and successfully connected to the database server. I opened QGIS Desktop and when I tried to connect to an MSSQL database I received the following message:

Warning: QSqlDatabase: QODBC driver not loaded
Warning: QSqlDatabase: available drivers: QSQLITE QSPATIALITE

I went searching for a driver and I found out that I am supposed to build the driver (!).

Assuming that QGIS Desktop is loading the available drivers from /usr/lib/x86_64-linux-gnu/qt5/plugins/sqldrivers, I wen on and installed Qt5.10.1, built the ODBC driver using these instructions and copied the resulting libqsqlodbc.so to the /usr/lib/x86_64-linux-gnu/qt5/plugins/sqldrivers.
This directory now contains the following:

-rw-r--r--  1 root root  65048 Apr 12  2018 libqsqlite.so
-rwxr-xr-x  1 alas root 123784 Jan  4 20:03 libqsqlodbc.so*
-rw-r--r--  1 root root  71832 Dec 22 06:50 libqsqlspatialite.so

I run again QGIS Desktop, I tried again to connect to the database server and I got again the same error.

My questions are:

a. Is there an easy way to add the QODBC driver to QGIS Desktop?

b. Am I doing something wrong above?

c. IS THERE ANY way to make QGIS Desktop in linux work with a SQL Server???


texts:
I'm trying to learn QGIS. How do I change the Extents (Highlighted in green below).
I tried to make a custom CRS, but I couldn't figure out a way to change these.


texts: QGIS - Change opacity of heat map in selected boundaries

texts:
I've only started using QGIS recently (used MapInfo). I would like to change the opacity of a heat map so that the areas within a political boundary has a different opacity than other areas. Would I use somthing like 'Select Feautres by Area' to select the area and then edit each smaller polygon indivdually?

Is this possible?


texts: Ggplot: Incorrect aspect and map position for northing and easting shapefile

texts:
Stats Canada released new census data in false northing and easting format. The shapefile from 2011 on long/lat format worked fine, but now the map is rotated and distorted.

The code below downloads the shapefiles, extracts them and plots the data.
I tried various "coord_" options without any success. Please help!

library(rgdal)
library(rgeos)
library(scales)
library(ggplot2)
library(dplyr)
library(tidyr)
library(ggmap)
library(maptools)

temp <- tempfile()
download.file("http://www12.statcan.gc.ca/census-
recensement/2011/geo/bound-limit/files-
fichiers/2016/lccs000a16a_e.zip",temp)
unzip(temp)
unlink(temp)

Boundaries2016<-readOGR(getwd(), "lccs000a16a_e")
Boundaries2016<-Boundaries2016[Boundaries2016@data$PRNAME %in%
c("Manitoba","Saskatchewan","Alberta"),]

Boundaries_df<-fortify(Boundaries2016, region="CCSUID")

plotData<-Boundaries_df

p <- ggplot() +
geom_polygon(data = plotData, aes(x = long, y = lat, group = group), color
 = "black", size = 0.25)
 #+coord_map()

p


texts: OSM Country import - missing border-near lakes, rivers

texts:
I'm importing a country (e.g Switzerland) from the cloudmade OSM branch and I'm failing to get the bordercrossing lakes on display. The same issue came up, when I tried to import the geofabrik .dpf branch this morning.

Somehow it makes sense and I don't think it's a bug, but I struggle to get the border crossing stuff to work. How to solve this issue without importing the european-western fileset?

Thanks!
Martin


texts: Publishing ArcGIS geoprocessing service - upload custom python modules?

texts:
I have an ArcGIS geoprocessing service which causes problems when publishing to arcgis server. The python scripts behind the service are structured as a custom python module, in a folder with an 'init.py' file, because the scripts contain many thousands of lines of code and there was a design choice to develop this way. My problem is that the script 'main.py' is published to the server correctly, but all of the other scripts are not. 'main.py' calls all of the other scripts. Therefore when I try to run the service I get ' Import Error - No such module as Scriptsforservice', where main.py is importing.

Is there a way to get arcgis to upload the custom module properly during publishing? I have found a hack which is to copy the module directly into my server's local python libraries, but this is not an acceptable long term solution.


texts: Can ArcGIS ModelBuilder open Excel file for manual editing, then feed those edits back into model?

texts:
I have been tasked with creating a model using ArcGIS ModelBuilder, that initially runs a couple of geoprocessing tasks, before then joining an Excel table to the result of those geoprocesses. However, the model needs to actually open the Excel file first, for the user to manually edit it, before the user then saves their edits and closes the file. It is then expected that the model will join this Excel file to the attribute table of a dataset (still part of the model) before the model continues on with a number of other geoprocessing tasks.

I include a high-level schematic outline below to help illustrate the point (I understand this might not make a lot of sense to start with).

So the model would run, complete a couple of geoprocessing tasks, pause to open an Excel file, wait for the user to edit it, save it and close it, before resuming to join that Excel file to an attribute table and then running a number of other geoprocessing tasks.

I should imagine this would require incorporating a Python script or two into the model. Currently though I am just interested to know whether the concept is actually possible or not.


texts: How to add a float or integer number (i.e 255) to a raster layer?

texts: Is there any code to add a float or integer value with a raster layer? I can easily do mathematical operations between two raster data. But I am unable to do the same between float or integer number and  raster layer. For example if want to execute (256-band3), how can I execute?


texts: Is there documentation anyplace for adding weather to PostGIS/Geoserver?

texts:
We have an internal PostGIS/Geoserver for rendering map tiles for the world based on OpenStreetMaps.  We would like to overlay weather, but would like to do so without referencing external URLs.  The weather data at this URL (http://radar.weather.gov/ridge/RadarImg/NCR/) uses a gif/gfw format, but I am so new to GIS that I am not familiar with how I can take advantage of these for our use.

Does anyone know of any documentation on how to get weather layers in PostGIS/Geoserver, or even better, using the gif/gfw format from above?


texts: Function that runs perfectly fine outside of QGIS keeps crashing QGIS (3.16) Python 3.7

texts:
I'm very new to QGIS, PyQt and even Python so I'm sorry if there is a very obvious mistake I am making, but I can seem to get one of my functions to run from within a plugin I'm writing for QGIS (3.16). When I run it directly from a test main as a script it functions perfectly fine. Ultimately, I would like to run it on a separate thread but for now, it even crashes QGIS when run on the main thread. (it also crashes QGIS when ran from a separate thread, I have tried it with QgsTask.fromFunction and TaskManager as well as PyQt5's ThreadPool and QRunnable subclass).
What I'm wondering about is the following: when I call the exact same function from a test main that I made, from within a virtual environment that is exactly the same as the QGIS Python setup, it runs just fine. I can obviously share some of my code, but I was wondering whether that fact alone pointed to something I'm missing. The line of code that seems to be causing the crash calls geopandas sjoin.
Further relevant information:
Both dataframes involved in the sjoin are created within the function that I want to run, so I doubt it's a matter of scope?
I have a general try/except clause set up that should catch (and currently print and ignore) all exceptions. However, that doesn't seem to prevent the crash.
The structure of the code:
pluginclass.run()
    button.clicked.connect(parent_of_crasher)

parent_of_crasher calls crasher

crasher calls subfunction_of_crasher

In my test main, I just call the function directly from main.
One argument is being passed to the crashing function (a python class object), but it's passed to every function I call from the parent and hasn't created any issues anywhere else, so I'm pretty sure that isn't the problem.
Where it crashes:
#(fixeddrainage and conn_nodes_geo are both geodataframes passed into the subfunction_of_crasher by the 'crasher')
try:
   joined = sjoin(
                fixeddrainage, conn_nodes_geo, how='left', op='intersects', lsuffix='fd', rsuffix='conn')

except Exception as e:
   QgsMessageLog(e)
   raise e from None

#After testing it on a hunch: it also crashes on
    overl = overlay(fixeddrainage, conn_nodes_geo, how="intersection")

I can print and modify both fixeddrainage and conn_nodes_geo without crashing QGIS, which leads me to believe it's not a memory/scope issue (I tried it by renaming columns in place, and doing a groupby). Furthermore, it never prints an exception to the message log, even when I pass instead of raising the exception.
It does seem to only happen with geopandas functions.
I am really at a loss here, so any insight is very welcome! I can provide more detailed information if that would help, but I figured I would describe the outline in case it points to a common mistake!
Edit: geopandas version is 0.8.1 in case that is relevant.
Edit: I have been experimenting further and I have found the following does work:
for index, row in fixeddrainage.iterrows():
    nodes_in_area = conn_nodes_geo.intersection(row.geometry)
    nodes_in_area_index = nodes_in_area[~nodes_in_area.is_empty].index

    fixeddrainage.loc[index, 'area_nodes_m2'] = round(
        conn_nodes_geo.loc[nodes_in_area_index, 'storage_area'].sum(), 0)

While this:
joined = sjoin(
            fixeddrainage, conn_nodes_geo, how='left', op='intersects', lsuffix='fd', rsuffix='conn')

group = joined.groupby([peil_id, 'multipolygon_level'])[storage_area].sum()

fixeddrainage = fixeddrainage.merge(group, how='left', on=[peil_id, 'multipolygon_level'])

However, I would rather not use the iterrows method. Can anyone tell me what the difference is (keeping in mind that the function executes fine outside of qgis) between the two methods?


texts: Using proxy with ArcGIS API for JavaScript?

texts: I am using ArcGIS with JavaScript API and I am trying to get a proxy since some functionality requires specifically to install and set up a proxy. I have looked at the manual but it is a little bit hard for me to follow, I am a new developer and am using github as a server.


texts: 3D line Populating with ArcGIS API for JavaScript

texts:
I'm trying to populate 3d line using ArcGIS API for JavaScript with my ArcGIS online PipeLine url, but somehow I'am not able to populate. I have attached the code.

<script>
      require([
        "esri/Map",
        "esri/views/SceneView",
        "esri/layers/FeatureLayer",
        "esri/widgets/Legend"
      ], function(Map, SceneView, FeatureLayer, Legend) {

  var map = new Map({
         basemap: "topo-vector",

          ground: "world-elevation"
        });

        var view = new SceneView({
          container: "viewDiv",
          map: map,
          zoom: 15,
          center: [103.805699, 1.276566],
           heading: 300,
           tilt: 75
       });
         const options = {
          profile: "quad",
          cap: "round",
          join: "miter",
          width: 5,
          height: 30,
          color: [200, 200, 200],
          profileRotation: "all"
        };

        var renderer = {

          defaultSymbol: {
            type: "polygon-3d",
        symbolLayers: [
                    {
                      type: "path",
                      profile: options.profile,
                      width: options.width,
                      height: options.height,
                      join: options.join,
                      cap: options.cap,
                      anchor: "bottom",
                      profileRotation: options.profileRotation
                    }
                  ]
          }
        };

        const PipelineLayer = new FeatureLayer({
          url:
            "https://services6.arcgis.com/itEfD6Bf5W0iWJuO/arcgis/rest/services/PipeLine1/FeatureServer/0",
          elevationInfo: {
            mode: "relative-to-ground",
            offset: 3
          },
         renderer:renderer
        });
        map.add(PipelineLayer);

      });
    </script>


texts: How can I draw squares in PostGIS layers at different angles?

texts:
I am trying to draw two squares on a 45-degree angle and then merge them using PostGIS to edit. QGIS can draw squares but only on a perpendicular plane (eg, 0, 90, 180), and I need to be able to draw at different angles (eg, 45-degrees).

Thank you!


texts: Converting REST Service to file geodatabase feature class?

texts:
I'm trying to convert a rest service to a feature class but Python is getting hung up in the middle of the process.

Input:

import arcpy
import urllib
import json

# Setup
arcpy.env.overwriteOutput = True
baseURL = "https://maps.somewhere.com/arcgis/rest/services/Bathymetry/MapServer/0"
fields = "*"
outdata = "E:/Tools/Rest/rest.gdb/bathymetry"

# Get record extract limit
urlstring = baseURL + "?f=json"
j = urllib.request.urlopen(urlstring)
js = json.load(j)
maxrc = int(js["maxRecordCount"])
print(("Record extract limit: %s" % maxrc))

# Get object ids of features
where = "1=1"
urlstring = baseURL + "/query?where={}&returnIdsOnly=true&f=json".format(where)
j = urllib.request.urlopen(urlstring)
js = json.load(j)
idfield = js["objectIdFieldName"]
idlist = js["objectIds"]
idlist.sort()
numrec = len(idlist)
print(("Number of target records: %s" % numrec))

# Gather features
print("Gathering records...")
fs = dict()
for i in range(0, numrec, maxrc):
    torec = i + (maxrc - 1)
    if torec > numrec:
        torec = numrec - 1
    fromid = idlist[i]
    toid = idlist[torec]
    where = "{} >= {} and {} <= {}".format(idfield, fromid, idfield, toid)
    print(("  {}".format(where)))
    urlstring = baseURL + "/query?where={}&returnGeometry=true&outFields={}&f=json".format(where,fields)
    fs[i] = arcpy.FeatureSet()
    fs[i].load(urlstring)

# Save features
print("Saving features...")
fslist = []
for key,value in list(fs.items()):
    fslist.append(value)
arcpy.Merge_management(fslist, outdata)
print("Done!")

Output

Record extract limit: 1000
Number of target records: 43001
Gathering records...
OBJECTID >= 1 and OBJECTID <= 1000
OBJECTID >= 1001 and OBJECTID <= 2000
OBJECTID >= 2001 and OBJECTID <= 3000
OBJECTID >= 3001 and OBJECTID <= 4000
OBJECTID >= 4001 and OBJECTID <= 5000
OBJECTID >= 5001 and OBJECTID <= 6000
OBJECTID >= 6001 and OBJECTID <= 7000
OBJECTID >= 7001 and OBJECTID <= 8000
OBJECTID >= 8001 and OBJECTID <= 9000
OBJECTID >= 9001 and OBJECTID <= 10000
OBJECTID >= 10001 and OBJECTID <= 11000
OBJECTID >= 11001 and OBJECTID <= 12000
OBJECTID >= 12001 and OBJECTID <= 13000
OBJECTID >= 13001 and OBJECTID <= 14000
OBJECTID >= 14001 and OBJECTID <= 15000
OBJECTID >= 15001 and OBJECTID <= 16000
OBJECTID >= 16001 and OBJECTID <= 17000
OBJECTID >= 17001 and OBJECTID <= 18000
OBJECTID >= 18001 and OBJECTID <= 19000
OBJECTID >= 19001 and OBJECTID <= 20000
OBJECTID >= 20001 and OBJECTID <= 21000
OBJECTID >= 21001 and OBJECTID <= 22000
Traceback (most recent call last):
  File "<string>", line 42, in <module>
  File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\arcobjects\arcobjects.py", line 421, in load
    return convertArcObjectToPythonObject(self._arc_object.Load(*gp_fixargs(args)))
RuntimeError: RecordSetObject: Cannot open table for Load

Can anyone explain to me what might be going on here?


texts: How to layer statistics from Excel onto a map?

texts:
Would another free app or program be better and easier to do this than mapbox.com ?

I simply want to layer census 2010 data by shading areas onto a map like on theatlanticcities.com/posts/map/ & developmentseed.org.


texts: Similar function of ArcGIS Network Analyst in MapInfo

texts:
I am currently using MapInfo and trying to create a road network base walking catchment map, rather than just a cycle buffer. The finished plan should look similar as follows.

This map was created in ArcGIS with the Network Analyst feature. I wonder whether MapInfo has the similar feature to generate this kind of map? maybe in the latest version v12.


texts: Reading GDB file with GeoPandas throws MemoryError

texts:
I have a gdb file that's 32GB and I'm tasked with extracting the attribute values from each of the layers and converting them to CSV for analysis. Here's my code:
import geopandas as gpd
import fiona as fiona

dataset_path = '...' # for example purposes left out the exact path

layers = fiona.listlayers(dataset_path)
for layername in layers:
    # according to docs settin ignore_geometry=True loads gpd file to pandas dataframe
    geodata = gpd.read_file(dataset_path, driver='fileGDB', layer=layername, ignore_geometry=True)
    geodata.to_csv(project_path + f'{layername}.csv')

This code works but there are some layers that are very large (~5000 columns x 40'000 rows) and for those, it throws MemoryError and everything on my computer crashes because the PC has ran out of RAM.
Is there anything I can do to optimize this besides trying to acquire a more powerful computer? (I have 16GB of RAM).

I already tried to load the file into QGIS and export it from there, but it throws a similar error and crashes anyway.


texts:
We are building up an application to collect map services from different owner. One of our big questions at the moment is how we can merge various map services from different institutions or municipalities into one map service.
The sources are using different software but we can always depend on international standards.

To give you an example:

Every municipality collects it's own building data.
Every municipality publishes it's building data set as a map service
We would like to combine all of these map services into one service
The end user just has to add one map service into his system to have a look at all the buildings of all municipalities
The user should be able to request (identify function) background data from each building

The central application is build on ESRI products, but we can add extra software to merge the map services.

I tried to find a solution with GeoServer, since I do not have a much experience in ESRI. But until now I did not find a proper way to do it.


texts: Displaying data from two UTM zones on one map?

texts:
With at least some USGS quad maps, maps that overlap into adjacent UTM zones continue to use the same false easting (rather than start at zero in the middle of the map).

Is this a required standard convention in the UTM specification or merely done for convenience?


texts: Hiding specific OpenStreetMap basemap tags?

texts:
I have been trying to find a way to hide or suppress specific tags on OSM within my .js code. There is just too much information packed into the OSM basemap and it takes away from what is important for my specific map use. Ideally it would be as easy as a simple call to remove a specific tag (like a picnic table icon or a tree icon). However, I have not found a way to do this yet.

Can this be done, and if so, how can it be done?

If there is not a way to do this how should I go about developing my own custom basemap to meet my needs?


texts: Distorted H3 hexagons - even though projection is matching

texts:
Using H3-py's h3_to_geo_boundary a WKT linestring (via shapley polygons) can be generated.

When later visualizing the hexagons, why are they distorted? How can this be fixed?

To my best knowledge, H3 works in EPSG:4326.
As you can see on the plot, the base-map coordinate system is also set to 4326. And QGIS also interprets the CSV of WKT strings as 4326.

Why do I see the distortion? How can it be fixed?


texts: Cost Path: multiple pairs of sources and destinations

texts:
I have seen many questions related to this problem but very few answers. I am very new to scripting and using ModelBuilder, but I am trying to calculate the cost path for 200 polygons with a single source and a single destination for each polygon.

So far, I have successfully calculated the Cost Distance and Backlink for all these polygons and their respective source, but now I would like to select the right destination, cost distance and backlink to calculate the cost path for each polygon.

All three inputs needed for Cost Path should therefore have the same number incorporated in their name related to the polygon they belong to (eg: Desti_001, CostDis_001, Backl_001) and use these three inputs to calculate cost path.

My questions are:

how do I select these three inputs based on the number is their names?
how do I create a loop to calculate the cost path of all polygons?

I use ArcGIS 10.1


texts:
A have an openlayers map with with two layers: a Tiled image layer (served through zoomify protocol) and a Vector layer. I draw boxes/rectangles on the image then get their geojson coordinates and save them in my DB.

After upgrading from OpenLayers 2 to 3 on my site, all my rectangular features are displayed in the wrong place. It turns out that with my OL2 implementation the Y axis is centered at the bottom of the map and points upward whereas in my OL 3 implementation the Y axis is centered at the top of the map and also points upward.

For backward compatibility reasons I don't want to update the coordinates in the DB. So I want the existing coordinates to show in the right place in my OL 3 image viewer.

How should I define the layers/views/map in OL3 to centre the Y axis of my features at the bottom of the image?

I've tried to reproduce the issue here:

https://jsfiddle.net/zjee9sdw/

You'll see that the box appears above the image. I want this box to appear in the bottom left corner of the image without changing it's input coordinates (or transforming them).


texts: Use Lambert 93 imageOverlay with Leaflet on top of OSM tiles layer

texts:
I'm trying to use a Lambert 93 image overlay on top of an Open Street Map title layer.

var map = L.map('map');
L.tileLayer('http://{s}.tile.osm.org/{z}/{x}/{y}.png', { }).addTo(map);
var coords = {
    x : [-0.734794,-0.418381],
    y : [47.376276,47.572680]
};
var imageBounds = [[coords.y[0], coords.x[0]], [coords.y[1], coords.x[1]]];
L.imageOverlay('image.png', imageBounds, {opacity:0.8}).addTo(map);

The resulting image overlay does not match the underlying OSM map.

It seems that this is because the image.png is in lambert 93 projection.

Can I use a lambert 93 image overlay on top of an OSM tiles layer ?
If not possible, can I convert my lambert 93 image.png to be usable with OSM tiles layers using tools such as gdalwrap ?

What I tried :

1) Proj4Leaflet

I looked at the Proj4Leaflet project which enables projections and provides an imageOverlay with projection example but the projection systems seems to concern the whole map and I want lambert 93 only on the imageOverlay.

This question seems to exclude the possibility of multiple projections
but as it only focuses on tiles layers I'm not sure this is relevant in my case.

2) gdalwrap

gdalwarp -s_srs EPSG:2154 -t_srs EPSG:3857 -of PNG image-lambert.png image-webmercator.png

=> Output driver `PNG' not recognised or does not support

gdalwarp -s_srs EPSG:2154 -t_srs EPSG:3857 image-lambert.png image-webmercator.tiff

Copying color table from image-lambert.png to new file.
ERROR 1: Unable to compute a transformation between pixel/line
and georeferenced coordinates for image-lambert.png.
There is no affine transformation and no GCPs.


texts: How to use data from field papers to QGIS

texts: I work on Field Papers and I need to get the information from the QR Scan code on the field papers on to QGIS. Can you please share how it can be done.


texts: Clip several overlapping polygons at a time - QGIS

texts:
I have a shapefile with lots of overlapping polygons made out of buffering street lines. I need to clip them so that the overlapping part goes to only one feature and gets excluded from the other, which I found I can accomplish with the Clipper plugin, as the image shows:

However, with that extension I would have to select one feature at a time and click on the 'Clipper' plugin button.
Is there any way I can create a script or some sort of action so that QGIS runs this automatically (select one feature, 'clip', select next feature, 'clip' and so on)? If so, could it select the features in a specified order, such as, according to the values in my 'length' column (from highest to lowest)?

If there's no way to do it, any ideas on how to accomplish clipping all these overlaps (thousands of polygons)?


texts: Reclassification of raster surface using ArcGIS for Desktop?

texts:
I have a raster surface with values ranging from 0 -100. I want to create a categorical classification (percentage categories of 5%, e.g. 0-5%, 6-10%...) and reclassify my surface into these categories. I am not sure exactly how to do this. I was wondering whether anyone could guide me through the steps I need to take. I am not sure whether I need to start in the Layer properties>Symbology before using the reclassify tool and if so what exactly do I need to choose in terms of the classification (e.g. natural, equal Interval etc.) and that sort of thing. I read the online help pages but I am still not getting anywhere!

I am using ArcGIS 10.2 for Desktop.


texts: I'm trying to create a geometric network of pipes, sewer structures (points) and manholes (points) inside a of an SDE.  I can create the network in a file geodatabase but when I try creating it inside the SDE i get an error about DBMS tables not being found and an error inside the topology engine detected.  Anyone have any insight? I've attached the error screenshot (I think, this is my first post)


texts: Generating map with polygon aggregation in PostGIS

texts:
Generate the map of federal units from the municipalities map in PostgreSQL with PostGIS.

Example below:

The table has only the geometry information and the "id".

I tried:

SELECT st_union(m.geom)
FROM uf, municipios as m
WHERE ST_Intersects(uf.geom, m.geom)
AND st_contains(uf.geom, m.geom)
GROUP BY uf.gid

But he returns geometries that do not belong to the state.

Note: The operation has the purpose of studies


texts: gdal/python: extracting projection info from hdf file

texts:
I'm trying to do something fairly simple. Take two subdatasets from an hdf file, read them as an array, do some math and write the new array back to a geotiff, using the same projection as the input files.

Because I'll be repeating this process for a list of .hdf files where the projection varies, I cannot hardcode the projection string for my outputs. While the projection varies from hdf file to file, it is the same for all subdatasets in a particular hdf file. So I only need to get the proj from one of the subdatasets. I thought it would be fairly easy to store the projection info into a variable, but the normal routines I use aren't working on the hdf file.

Here are few different things I've tried:

ds = open(hdfFile)
sds_b3 = ds.GetSubDatasets()[4][0]
sds_b4 = ds.GetSubDatasets()[5][0]

# this gives me an empty string:
proj = ds.GetProjection()

# both of these give me "AttributeError: 'str' object has no attribute 'GetProjection":
proj = ds.GetSubDatasets()[4][0].GetProjection()
proj = ds.GetSubDatasets()[4][1].GetProjection()

# this last one works in that it prints the correct projection info, in addition
# to a bunch of other metadata, but I don't know how to extract just the projection
# and store it into a variable:
os.system('gdalinfo %s' % sds_b3)

How can I store the projection information in a variable, like 'proj = dataset.GetProjection()' normally would, so that I can use it when writing the new array to a geotiff? ie:

output = gdal.GetDriverByName("GTiff").Create(outputname, nlines, nrows)
output.SetProjection(proj)

I either need to know the correct way to use GetProjection() when using hdf's or how to isolate the projection part of gdalinfo as a variable.


texts: Does this DXF file contain geolocation data?

texts:
I'm new to autocad and have been sent a file with polylines which I need to plot on a map.
However I cannot see that the file contains any geolocation data to plot it accurately.
The GeoLocation tab does not appear and reading the AutoCad knowledge base, I suspect that the file does not contain any location data.
Can someone download the file and confirm this or if I'm missing something?


texts: Connecting to Azure PostgreSQL using domain user account through QGIS

texts:
We are using QGIS over the top of an Azure PostgreSQL database.  We have it working using SQL user / password combination but security policy dictates this should be a domain user authenticated connection.
Can you connect to an Azure PostgreSQL database using a domain user account through QGIS?


texts:
I have built qgis source on windows 7. Libs have been downloaded from osgeow4 site. I even installed postgress as a standalone.

I am getting this error when the app launches..

  File
  "C:/development/qgis2.9.0/./python/plugins\processing\algs\qgis\postgis_utils.py",
  line 28, in    import psycopg2
  File "C:/development/qgis2.9.0/./python\qgis\utils.py", line 483, in
  _import mod = _builtin_import(name, globals, locals, fromlist, level)
  File "C:\development\OSGeo4W\apps\Python27\Lib\site-packages\psycopg2__init__.py",
  line 71, in
      from psycopg2._psycopg import BINARY, NUMBER, STRING, DATETIME, ROWID
  File "C:/development/qgis2.9.0/./python\qgis\utils.py", line
  483, in _import mod = _builtin_import(name, globals, locals, fromlist, level)
  ImportError: DLL load failed: The specified module could not be found.

I have check all my python paths etc, but am not sure what about the psycopg2 plugin is missing.

Any direction would be great.


texts: Recently we have started migrating application from Javascript API 3.1 to 3.5 and we are finding this jsonp error coming sometimes during adding the layer.When we add the layer we also zoom to the layer .Below is the snapshot of the error I get.


texts: Survey123 Summary Report syntax with dates

texts:
I'm working on a summary report for pump station inspections.  I'm using a repeat tag to fill out the table which will insert a row for each pump station with it's summed inspection data.  I was able to get the syntax correct to sum the data for each question for each individual pump station.  But I'm struggling with how to show the most recent inspection date for each pump station.
You'll see in my syntax I tried to use stats: max,REPDATE,InsDate but that gives me the most recent inspection date of ALL the inspections.  I only want the most recent inspection date per pump station. I've tried shifting things around, but I keep getting errors and when I don't get an error I get the most recent date overall.
I inserted what I have so far for the report template table below.
$<$Summary>

Station No.
Last Inspection Date
Pump 1 Run Hours

${#PumpStationLogs | stats:"count,FACILITYID; max,REPDATE,InsDate; sum,P1_RNHrs,SumP1; sum,P2_RNHrs,SumP2; sum,BIOXIDE,SumBiox; sum,GENERATOR,SumGen; sum, ELCTMTR_R18,Electric; sum, H2OMTR_R18,SumWater":"PSNAME": "" | where:"1=1" | orderByFields:"PSNAME ASC" }${PSNAME}
${InsDate | format:"DD"}
${SumP|round:2}

$</>


texts: Arcpy Copy Feature Classes with Overwrite to specific feature datasets

texts:
Objective - I have a target GDB "schema" with set feature dataset
names and feature classes (no features) within them.  I want to copy
submittal geodatabase feature classes into that target.  I first
reproject the target to match the submittal coordinate system.

What I want to do with python is copy all feature classes from the
submittal over to the target, where if a submittal feature class name
matches one in the target GDB, I want to replace it in the target GDB
in the feature dataset location where the empty feature class was in
the target.  If submittal feature class has no name match in target
GDB, I want to copy it over to the target GDB outside any feature
dataset.

I am a python novice. I used http://support.esri.com/en/knowledgebase/techarticles/detail/40831 as a start and now got code (below) to work to copy all feature classes from the source GDB to the target GDB, deleting out feature classes within feature datasets with matching names between the source and target GDBs, but I can't get it to put the source feature classes in the feature datasets where the matched name feature classes are in the target, it just puts them in no feature dataset in the GDB.

UPDATE: Code does work when I make a copy of my submittal database. This revealed that code didn't work as intended because the reproject I do on my target database feature datasets initially isn't truly bringing in all the spatial reference information.  In order to get this to work the way I want, I believe I have to add code to create new feature datasets in the target GDB with same names as the target feature datasets with spatial reference information from the submittal feature datasets, then copy over feature classes where names submittal names match target names, then delete the original target feature datasets (or overwrite when creating).

Code snippet below updated per comment instruction:

import arcpy, os, string

arcpy.env.overwriteOutput = True

def FeatureClassToFeatureClass_conversion(start_db,end_db):

#Set workspaces
arcpy.env.workspace = start_db
wk2 = end_db
datasetList = arcpy.ListDatasets()

#for feature classes within datasets
for dataset in datasetList:
    fclist = arcpy.ListFeatureClasses("*","",dataset)
    for fc in fclist:
        print "Reading: {0}".format(fc)
        name = arcpy.Describe(fc)
        new_data=name.name
        if arcpy.Exists(wk2 + os.sep + new_data)==True:
            arcpy.Delete_management(wk2 + os.sep + new_data)
            arcpy.Copy_management(fc, wk2 + os.sep + new_data)
            print "Completed copy on {0}".format(new_data)
        else:
            arcpy.Copy_management(fc, wk2 + os.sep + new_data)
            print "Completed copy on {0}".format(new_data)
#Clear memory
del fc

if __name__== "__main__":
    start_db = arcpy.GetParameterAsText(0)
    end_db = arcpy.GetParameterAsText(1)
    FeatureClassToFeatureClass_conversion(start_db,end_db)


texts: Snap/Trace to layer in MapInfo without it being selectable?

texts: Is it possible to snap/trace a layer in MapInfo without it being selectable? I am preparing workspaces for a survey and I don't want the surveyors being able to select the background mapping (they're not too tech-savy) but still be able to trace features from it.


texts:
I want to calculate zonal statistics as table for multiple rasters, but can't figure out how.  What I am trying is:

import arcpy, os
arcpy.env.workspace=r'F:\Sheyenne\Atmospherically Corrected Landsat\Indices\Main\NDWI\main_master'
#zones to calculate stats on
inpoly=r'F:\Sheyenne\SNG_sections.shp'
#where to save output files
destination_path=r'F:\Sheyenne\Atmospherically Corrected Landsat\Indices\Main\NDWI\zonal_stats'
 #list of rasters to calculate statistics on
rasters = arcpy.ListRasters("*water_.tif*")
#start loop
for raster in rasters:
    destination_path = os.path.join(destination_path, raster + ".dbf")
    arcpy.gp.ZonalStatisticsAsTable_sa(inpoly, "OBJECTID",raster, destination_path,"DATA","ALL")

The error returned is :

Traceback (most recent call last):

  File "<ipython-input-10-9522776a92dd>", line 1, in <module>
    runfile('F:/python codes/zonal_stats.py', wdir='F:/python codes')

  File "C:\Users\spotter\AppData\Local\Continuum\Anaconda\lib\site-packages\spyderlib\widgets\externalshell\sitecustomize.py", line 685, in runfile
    execfile(filename, namespace)

  File "C:\Users\spotter\AppData\Local\Continuum\Anaconda\lib\site-packages\spyderlib\widgets\externalshell\sitecustomize.py", line 71, in execfile
    exec(compile(scripttext, filename, 'exec'), glob, loc)

  File "F:/python codes/zonal_stats.py", line 21, in <module>
    arcpy.gp.ZonalStatisticsAsTable_sa(inpoly, "OBJECTID", raster, destination_raster,"DATA","ALL")

  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\geoprocessing\_base.py", line 498, in <lambda>
    return lambda *args: val(*gp_fixargs(args, True))

ExecuteError: ERROR 999999: Error executing function.
Create output table failed
Failed to execute (ZonalStatisticsAsTable).


texts: Converting shapefiles to text (ASCII) files?

texts: How can I convert a shapefile (GIS) to text, or, how can I extract the information in a shapefile?


texts: Batch Geocode Postcodes/Addresses in UTM

texts:
My work has a database of previous jobs (within housing plots) which are meant to have been georefernced with UTM coordinates, but this has not been completed for a number of years. New jobs can use data from old ones, and there is a facility in the database to check for properties within a certain radius of the new job using UTM coordinates.

I have batch georeferenced addresses using Lat Long before using sites like GPS Visualizer, but cannot find a means to georeference in UTM, and converting lat long to UTM seems difficult.

Does anyone know of a site which will help me do this, or know of a site which will convert directly to UTM coordinates?

The UTM Zone is 17N (Winnipeg).


texts: Dealing with raster: creating it in projected CRS and exporting it in geographic CRS

texts:
When creating a raster from a polygon (the CRS of the polygon is a geographic coordinate system), the raster is created in the projected coordinate reference system of the data frame, which is a projected one. I need to do so because I need to specify the Cell size of the raster in meters.

I also need to deliver this raster (previously converted to polygon) in WGS84 (GCS); I am considering to do it so: Changing the Data frame coordinate reference system to WGS84, and, on the result raster, do Data > Export Data and, for the option “Use the same coordinate system as:”, choose “the data frame”.

Is this an accurate way of dealing with that? I mean, aren’t there modifications in areas or coordinates of the delivered shape layer?


texts: Comparing two rasters in scatter plot with different resolutions using ArcMap?

texts:
I have 2 rasters taken from 2 cameras with different resolution over the same field. One shows NDVI values, the other is temperature. I need to make a scatter plot between these rasters but they have different cell sizes due to resolution of the camera and various number of rows/columns.

I tried Resample in Arcmap 10.5. to resize the smaller raster but it will not get the exact cell size and number of rows/col.

How can I compare 2 rasters with different resolutions if resampling function does not work?


texts: Displaying graph plots with mouse click on map in ArcGIS for Desktop?

texts:
I need to build a tool/add-in for ArcGIS 10 which will display a time distribution plot at the location where mouse is clicked. The TOC consists of results of a simulation at different time instants. For the clicked point, the tool should read-in the output values at different time instants from TOC and plot it as a graph and show it as a pop-up.

So far, I have managed to draw the time distribution graphs for a set of points stored in a shapefile and save it to a working folder with a Python script but couldn't manage to get X,Y coordinates for the clicked point and enter it as input and then display the plot at the end. I have really no idea if Python supports mouse-click events like VB does.

If it is not possible with Python, which other language should I turn to (e.g., VB.net)?


texts: Could not install Applications Manager (ArcGIS Server 10)

texts:
I try to install ArcGIS Server version 10.
Everything worked fine with installing the GIS Services.
The Server is running. I can create services and see them in ArcCatalog.
But when I try to install the Web Applications I couldn't select the "Applications Manager".

I try to install it on Windows Server 2008 (64bit) with IIS7.
I double checked the settings and required extensions for IIS and it seems to be correct.
I even changed the IIS seting from 64bit to 32bit.
But the "Application Manager" in the Select Feature Dialog is still not shown.

I appreciate any help with this.

Kai


texts: Windows security warning (ActiveX) when clicking into fields in ArcToolbox tool?

texts:
I have been running ArcGIS 10.5 on my PC (Windows 10 Home) since December with no issues.  This past week I have started to get ActiveX Security Warnings when I click into a parameter field in any ArcToolbox tool.

  Windows Security Warning

  An ActiveX control on this page might be unsafe to interact with other parts of the page.  Do you want to allow this interaction?

It doesn't seem to matter whether I click Yes or No, I can still modify the field and the tool still functions fine.

Answering Yes:

Answering No:

I am working through the answers to Why do ArcToolbox tools give ActiveX error?, although there doesn't appear to be an actual answer on that Q&A. Additionally it's back in the Windows XP/7 days with IE 7/8, and I'm running Win 10 with IE 11 (although that's never used).  If it matters, the browsers we use are Firefox 53.0 (64-bit) and Chrome 58.0.3029.81 (64-bit).

The issue I am having sounds different as it's happening only when I click on the parameter field - not when I open the tool or, as I read elsewhere, ArcMap.  The warning messages on that Q&A also have different text so it may be a completely different problem.

I have opened Internet Explorer (for the first time since getting this PC!) and looked at the Local Zone settings and everything for ActiveX is marked Enabled, including one about running without prompting:

I have also Disabled another one "Automatic Prompting for ActiveX Controls" but that didn't stop the message appearing.

I went through the suggested ActiveX settings from Problem: Issues with ActiveX controls when launching ArcMap, updated all my settings to match (the only one I had to change was the "Automatic Prompting for ActiveX Controls" one I changed in my testing above).

No new software has been installed on my PC in the past week.

What could have caused this to start appearing, and how do I fix it?


texts:
I am working with QGIS 3.8.1 Zanzibar.

I have a layer with Offshore Licenses holders. The attribute table has separate fields for "Participant 1", "Participant 2", ..., "Participant 13".

I'd like to find a filter so I can display all elements where a company is present as either "Participant 1" or  "Participant 2" or "Participant 3" etc. I have used

"particip1" IN ('Company Name') OR "Particip2" IN ('Company Name') OR ...

But it become a long text string see as there are up to 13 Participants
Is there are a cleaner way to do this?


texts: Identifying suitable satellite imagery of the South China Sea

texts:
I've been using the Google Earth Engine API in an attempt to access Landsat 8 imagery for land use change detection in the South China Sea. Initial investigations suggest that Google may have removed imagery of this area from their collection, as when I run the code in the link below all images are clipped to the surrounding mainland coastline of each country.

https://code.earthengine.google.com/b2d2927d5f8abb7c77444cc03dc8ece8

Having checked the following, https://earthexplorer.usgs.gov/ I can confirm that suitable Landsat 8 imagery does exist for this ROI, so my question is:

Is there any documentation where Google specifies why this has been excluded from their Landsat 8 collection? Or have a made a mistake in the code above that somehow excludes the imagery I am looking for?


texts: Counting points within polygon with ArcPy?

texts:
I am trying to script a process where I need to make a selection on a point feature class, count the number of points that fall within each polygon on another feature class, and then update a COUNT field with the new value.  I can preform the selection and creating a new layer from the selection via ArcPy but I cannot figure out how to count the points within the polygon via ArcPy.  I tried

import arcpy
try:
   #select all trip hazards not marked 'Closed'
   arcpy.SelectLayerByAttribute_management ("Sidewalk_Trip_Hazards",   "NEW_SELECTION", "\"Status\" <> \'Closed\'")

   #make layer of selected hazards
   arcpy.MakeFeatureLayer_management("Sidewalk_Trip_Hazards","notClosedHazards")

    #spatially join points to hexbins
    arcpy.SpatialJoin_analysis("notClosedHazards", "LaytonHexBins", "in_memory/SpatialJoin")

    #summarize points based off of HexID
    arcpy.Statistics_analysis("in_memory/SpatialJoin", "in_memory/SummarizedTable", [["HexID", "Count"]])

except:
    print(arcpy.GetMessages())

Is there a way via ArcPy to do a points in polygon count?


texts: Writing Python batch script for Analyze and Compress on ArcSDE Geodatabases?

texts:
I created and modified python scripts based off of ERSI's example (http://resources.arcgis.com/en/help/main/10.2/index.html#//003n000000v7000000), however, what would be great is to be able to batch process these files using python.

I attempted trying ModelBuilder by adding two database connections to the compress tool. (That didn't work.) Obviously, right-clicking on the tool in ArcToolbox brings up the option for batch, but I have no idea how to do this for python.

I can not seem to find anywhere the ability to do this. It's a pain to have ten scripts running the same thing but with different databases.

Anyone know how to do this with python?


texts: How to access the map coordinates show in QGIS status bar?

texts:
Is there a way to access the map coordinates show in qgis status bar? I found only this,

iface.mainWindow().statusBar().currentMessage()

But I know it's not the way.

Or is there any way to capture the coordinates when cursor moves on map?


texts: How can I align perfectly 4 topo sheets?

texts:
I have 4 topo sheets and i want then to align perfectly

Would I have to georeference each image again? Or is there a more simple method?


texts:
I create two raster using QGIS Interpolation tools starting from two vector (Point): one with an accessibility index( between 1 and 5) for each point  and one with population (number of people) for each point.
I need to combine the two raster to have the % of population for each index value.
I've never worked with raster in QGIS.


texts: Converting .cdr to shapefile using ArcMap

texts:
How do I convert a .cdr (Corel) to a shapefile?
I am using ArcMap (ArcEditor) 10.0.  I am brand new to using 10.0 after using 9.3 and brand new to to the ArcEditor level license after using ArcView.
I have never had to work with .cdr before.


texts: dissolving polygon records that only fall completely in another polygon from the same table

texts:
I have a polygon table that consists of 700 records with the subtype of  'OutdoorRecreationalFacility'. this postgres table is supposed to be a table of 1000 foot buffers around certain parcels. However in some instances the parcel geometry exists in a separate record shown in the picture below

The problem is none of these records have addresses so I cannot st_union them by address which would solve my problem. There is no attribute link between the parcels and the parcel buffer. they are only spatially linked. I am trying to come up with a query that would either

select the larger geometry if one geometry falls completely within another geometry or select the larger geometry
for all the geometry records that fall within another geometry record st_union them together

I formulated this query but it just keeps running and running (has an index)

 select case
when st_within(st_makevalid(t1.geom),st_makevalid(t2.geom)) then st_union(st_makevalid(t1.geom),st_makevalid(t2.geom))
else t2.geom end geom
from buffer_unions t1 cross join buffer_unions t2
where t1.subtype = 'OutdoorRecreationalFacility' and t2.subtype  ='OutdoorRecreationalFacility'

My thought process behind this query is if a geom falls within another geom, union them together, if not just return the t2.geom ....


texts: Georeferencing polygons with no spatial reference?

texts:
I have a Google earth image which has drawn polygons on it.  I would like to create a polygon feature class in ArcGIS using this data.  However, the actual image has been manipulated in powerpoint.

Is there a simple method to do this?

The image is a JPEG with 23 polygons drawn on it indicating different areas.

The process is not making much sense to me as the original image has no spatial reference aside from the Google Earth projection.


texts: Diagram Ring Elements Locator (Tough Problem)

texts:
I asked https://stackoverflow.com/questions/35579843/get-only-ring-nodes-from-a-diagram on stack overflow but there was no help at all. So I hope I get help from someone on here :)

I have a diagram that basically forms into a ring, each node is connected to another node. some nodes are part of the ring some are not. As the picture below shows, there is a ring of nodes (0,1,2,7,8,9). I want to get those nodes from the diagram and store them into a list.

Each node object contains details on the nodes it is connected to, so node 0, contains a list of connected nodes which basically has 1, and 9 in it.

I am trying to get the nodes within the ring, but the method I wrote does not work for all ring diagrams.
Here is what I wrote that works on some diagrams

    private bool SetMainRingList(StructureFeature strct, StructureFeature root, List<StructureFeature> mainRing) {
        if ((strct.Equals(root) && mainRing.Contains(strct))) {
            return false;
        }
        var children = strct.GetConnectedStructures();
        if ((children.Contains(root) && mainRing.Contains(strct))) {
            return false;
        }
        mainRing.Add(strct);
        foreach (var structureFeature in children) {
            if (mainRing.Contains(structureFeature)) {
                var strcture = mainRing.Find(x => x.Oid == structureFeature.Oid);
                if (strcture.ParentFeature == null)
                    continue;
                if (strcture.ParentFeature.Equals(root)) {
                    bool skip = false;
                    var crntChildren = strcture.GetConnectedStructures();
                    foreach (var childContainerse in crntChildren) {
                        if (!mainRing.Contains(childContainerse)) {
                            skip = true;
                            break;
                        }
                    }
                    if (!skip)
                        return false;
                }
                continue;
            }
            structureFeature.ParentFeature = strct;
            var leaf = SetMainRingList(structureFeature, root, mainRing);
            var exchangeSite = structureFeature as ExchangeStructure;
            if (leaf && ReferenceEquals(exchangeSite, null)) {
                mainRing.Remove(structureFeature);
            } else {
                return false;
            }
        }
        return true;
    }


texts: How do I overlay vector data over a 3D terrain basemap?

texts:
I would like to create a map that is similar to the one found here: https://www.reddit.com/r/MapPorn/comments/lpo557/geological_map_of_north_america/ .
I'm thinking that you could replicate this with a transparent terrain basemap and add the appropriate layers on top of this. However, I am also open to other suggestions if anyone has better ideas on repeating this style of map.

Edit after everyone's feedback:


texts: Unsplit lines based on length and attributes

texts:
I have a contour dilemma.  Using Microstation, 1-ft contours have been made on a til-by-tile basis in our project.  These DGNs were imported into a GDB, so the polylines have already been merged together.  There are four types of contours in this dataset (index, intermediate, index depression, and intermediate depression).  I have unsplit the lines based on elevation and contour type.  The problem I am having is that Microstation was not consistent in specifying contour type.  So, from tile to tile a contour may have two different types.  I would like to unsplit the lines so that the resulting line retains the attributes of the longer feature.  Does anyone know of a way I can achieve this?  Or, is there a way to prevent this in Microstation?  Much appreciation in advance.

EDIT: I have increased my model buffer to include more model key points.  This removes many of the smaller errors (contours that are within x distance of the processing tile) but larger errors (contours that extend outside of x distance) remain.  Increasing my model buffer to include more tiles also dramatically increases processing time.  Small projects that allow me to increase the model buffer to include all tiles does not produce errors.  This, however, is not possible on large projects.


texts:
I have two different parameters raster data, Which are temperature (X) and precipitation (Y). For temperature (X) and precipitation (Y), each parameters have 12 data files (Tiff format).

I want to perform Pearson's correlation coefficient analysis between two parameters and want output as a raster format. I have tried to do it in Arc GIS but I think in ArcGIS not providing such type computing facility.

So I am looking for a solution in an R based software package. I am new in R.

I keep my rasters in two different file, Below I have attached my folder situation which is shown below. How can i assign here path for the two different raster parameters.


texts: Using Python to list all entries in numpy array that are lower than all adjacent values

texts:
I'm using ArcGIS 10.4.1, arcPy and NumPy. I'm trying to create a list of tuples that identify all entries in a NumPy array (DEM) who's value is lower than all adjacent neighbors (I.e. It's a depression).  Ideally I'll find a way to create a list of tuples with col, row and cell value that I can ultimately sort from highest to lowest cell value.

Any thoughts?

My end goal is to also determine the outlet point (pour point) for each of these sinks, but I'm learning one step at a time.

So far I've got:

import arcpy
from arcpy.sa import *
import numpy as np

lowerLeft = arcpy.Point(inDEM.extent.XMin,inDEM.extent.YMin)
arcpy.env.extent = inDEM
outSink = Sink(FlowDirection(inDEM))

SinkArray = arcpy.RasterToNumPyArray(outSink, lowerLeft)

a,b = np.where(SinkArray >= 1)
x = zip(a,b)

arcpy.Delete_management(outSink)
del SinkArray

Which gives me a list (x) of row/cols where sinks exist, but i'm sure there is a way to just import my initial DEM and identify local minima.


texts: Layer creation failed: Read-only file system

texts:
Ran into this error while trying to create a shapefile layer for training sites in an object segmentation project. This is an issue with QGIS itself, as I tried to create new vector layers on previous projects and got the same error.

I did edit some path variables while trying to connect QGIS to an external python interpreter, and I suspect my issues stem from that based on navigating to Settings > "Options > Environment > Current environment variables (read-only - bold indicates modified at startup)". There was a custom PATH variable that I had set earlier in my ~/.bash_profile and commented out, but clearly it still exists according to Q.

For anyone experienced in customizing environment variables, is there a simple solution to this problem? I'm only just starting to become familiar with customization but I've surely bungled something up along the way.


texts: How do you convert a raster of digital values of ndvi to an index of NDVI?

texts:
I have a raster of the Mean May Normalized Difference Vegetation Index (NDVI) that is a composited measure of VGT daily images, acquired between 2000 and 2009 in May, using the Mean Compositing algorithm (Vancutsem et al., 2007) in order to generate an average composite. The processed VGT-S10 products is in digital numbers (from the vegetation sensor). However, in my model I transformed these values from the digital number to an NDVI index ranging from 0 to 1, the original range of NDVI, using the following equation: NDVI = (0.004∗Digital Number)−0.1.  I would now like to convert the raster so I can use the predict function to model species ditribution models which include this predictor. How would I best do this?

I am using ArcGis 10.2.


texts: There is Change in classification total area after Exporting and reclassifying it from **Google Earth Engine** to ArcMap

texts:
I have done classification of 6 classes in google earth engine on an available forest mask that has an area of about 2493 sq. km. The forest mask basically uploaded to Google Earth Engine and satellite imagery clipped on it for classification purposes. After classification analysis, I have exported the classified image to ArcMap and performed basic analysis e.g(mosaicking, projection to the required coordinate system, and reclassification). After reclassification, I have computed the area of each class and total it which is 2196 sq. km it variate 296 sq km from the original area. I have done the same procedure many times but got the same issue. I don't know where I have done wrong while analyzing. If anyone goes through the same issue or knows how we can check it so that area should be the same as the mask used.
For clarification, I am attaching screenshots. eq3 is the original area and the table is computed one.


texts: Referencing Enterprise Geodatabase feature classes from ArcPy with ArcGIS Pro?

texts:
When using ArcPy with the ArcGIS Desktop (10.4.1) architecture I am able to set a variable to a feature class and make a feature layer from it using the code below:

rock_unit_table = "Database Connections\\TEST PROD [personal].sde\\Test.TEST.DETAILED_SURFACE_COMP"
arcpy.MakeFeatureLayer_management(rock_unit_table, "rockUnit_lyr")

but when I try to do something analogous using ArcPy with ArcGIS Pro 1.3.1:

rock_unit_table = r"Databases\TEST PROD [personal]\Test.TEST.DETAILED_SURFACE_COMP"
arcpy.MakeFeatureLayer_management(rock_unit_table, "rockUnit_lyr")

I get an error (with some anonymizing applied):

  Traceback (most recent call last):   File "N:\path\Test.py", line 324,
  in
      arcpy.MakeFeatureLayer_management(rock_unit_table, "rockUnit_lyr")   File "C:\Program
  Files\ArcGIS\Pro\Resources\ArcPy\arcpy\management.py", line 5607, in
  MakeFeatureLayer
      raise e   File "C:\Program Files\ArcGIS\Pro\Resources\ArcPy\arcpy\management.py", line 5604, in
  MakeFeatureLayer
      retval = convertArcObjectToPythonObject(gp.MakeFeatureLayer_management(*gp_fixargs((in_features,
  out_layer, where_clause, workspace, field_info), True)))   File
  "C:\Program
  Files\ArcGIS\Pro\Resources\ArcPy\arcpy\geoprocessing_base.py", line
  506, in
      return lambda *args: val(*gp_fixargs(args, True)) arcgisscripting.ExecuteError: Failed to execute. Parameters are not
  valid. ERROR 000732: Input Features: Dataset Databases\TEST PROD
  [personal]\Test.TEST.DETAILED_SURFACE_COMP does not
  exist or is not supported Failed to execute (MakeFeatureLayer).

In ArcGIS Pro I can see the feature class under Databases in the Project pane so I think it is my syntax in ArcPy that is astray.

Can anybody suggest what may be astray with the syntax that I am trying to use?


texts: QGIS - look up Table to feature Layer where Table answers questions to Multiple Columns in Feature Layer

texts:
In the image example below:

Feature Line Layer = Transport - District

Lookup Table = TRAFFIC_IMPACTOR_CODE

What is the most efficient way to "connect" (join/relation/query/field calculator) both tables so the DESCRIPTION column of the Lookup Table 'TRAFFIC_IMPACTOR_CODE' is returned for FROM_TRAFFIC_IMPACTOR_CODE and TO_TRAFFIC_IMPACTOR_CODE of  'Transport - District' feature layer?

A simple QGIS JOIN can only return an answer for ONE column in the 'Transport - District' feature layer.

A Relation is specific only for a project, and is useful only when using the Identify Features Tool as far as I can tell.

A Query (Virtual Layer) has the advantage of displaying a subset of the feature layer, how does one Query based on Both tables?

A Field Calculation (Virtual Column) would be nice to have as a persistent set of answers to multiple column relations between the Feature Layer and the Lookup Table.  How does one ask for the DESCRIPTION answer from the Lookup Table for both FROM_TRAFFIC_IMPACTOR_CODE and TO_TRAFFIC_IMPACTOR_CODE columns?


texts: I would like to know if their are any open source map desktop based map viewers something similar to ArcReader. Preferably if this can use the files created from QGIS.


texts: UpdateCursor gives Too few parameters Expected 1?

texts:
Can anyone tell me why my script isn't working?

I can't figure out where I went wrong. I am working with an Access personal geodatabase. Adding a field, then updating that field with yes if the field Surface_Cover equals CONCRETE, or else it returns NO.

I am getting an error message "Too few parameters. Expected 1".

#import tools and set workspace
import os
import arcpy

arcpy.env.workspace = "C:\Users\sutton\Documents\Python\PythonExcercise/PlanoData.mdb"

try:

    #create new field to hold new values
    arcpy.AddField_management("SmokeMain","Indicator","TEXT","10")
    print("Indicator field added to SmokeMain")
    with arcpy.da.UpdateCursor("SmokeMain",("Surface_Cover","Indicator")) as cursor:
        cntr = 1
        for row in cursor:
                #update the indicator field
                if row[0] == "CONCRETE":
                        row[1] = "YES"
                else:
                        row[1] = "NO"
                cursor.updateRow(row)
                print ("Record number " +str(cntr) + " updated")
                cntr = cntr + 1
#error checking
except Exception as e:
    print(e.message)

#script status
ScriptStatus = "Script executed"
print ScriptStatus


texts: Georeferenced vector file to raster file problem

texts:
When I rasterize the vector file it keeps getting an error 6. How do I manually put a -3d,-a, or -burn value?
I am only a beginner.
Feel free to ask for clarifications about my problem.


texts: How can I get back my Entity from a Primitive and change its color?

texts:
I put an Entity (Rectangle) on the map:

var greenRectangle = viewer.entities.add({
    name : 'The Area',
    rectangle : {
        coordinates : rectangleSelector,
        material : Cesium.Color.BLUE.withAlpha(0.2),
        heightReference: Cesium.HeightReference.CLAMP_TO_GROUND,
    }
});

Now I pick it from mouse move:

var handlerLatLon = new Cesium.ScreenSpaceEventHandler(scene.canvas);
handlerLatLon.setInputAction( function(movement) {
    var ray = viewer.camera.getPickRay(movement.endPosition);
    var position = viewer.scene.globe.pick(ray, viewer.scene);
    if (Cesium.defined(position)) {
        var pickedFeature = viewer.scene.pick( movement.endPosition );

        // the code to change Rectangle color

    }

}, Cesium.ScreenSpaceEventType.MOUSE_MOVE );

But pickedFeature is Primitive ( Object { primitive: Object, id: Object } ) and I have no idea how to change its color or take my Rectangle back.

I need to access the Name property too.


texts: Accessing data stored in related tables from a network dataset

texts:
According to ESRI's transit data model, this implementation, and several thesis papers, any non-spatial information (such as public transport schedules, employee info, trips etc) which is not stored in feature classes should be stored in geodatabase tables (I agree).

The problem I have is with accessing the data in those tables. Once the network dataset is created, how can the scheduling data be "loaded" or used for the route analysis? Or any of the other data, for that matter?

The articles I've linked to are from several years ago, and they all used different (custom) applications to retrieve the data from the geodatabase tables for use with the network dataset. They all used relationship classes to relate the tables to the various feature classes in the network dataset. Is there a way to do this using the new Python module in 10.1 (or any other way?)


texts: Newbie PostGIS Geometry and Multilinestring Clarification

texts:
I'm having some trouble figuring out if I've got the right data loaded for PgRouting.

Using the bash script included in PostGIS 2.0, I loaded Tiger2010 data for U.S. State California.

The edges table contains the_geom column, data type geometry.  Using Underdark's example, it seems I need the roads in a Multilinestring format in order to begin generating routes.  Can the data in the edges table be converted into type multilinestring?  The table query is below.
I altered my geocoded location data to data type geography.  PgRouting's shortest_path function needs the data in integer format.  How do I convert the geography type point into an integer that shortest_path can use?

Thanks for your patience.

CREATE TABLE tiger.edges
(
  gid integer NOT NULL DEFAULT nextval('edges_gid_seq'::regclass),
  statefp character varying(2),
  countyfp character varying(3),
  tlid numeric(10,0),
  tfidl numeric(10,0),
  tfidr numeric(10,0),
  mtfcc character varying(5),
  fullname character varying(100),
  smid character varying(22),
  lfromadd character varying(12),
  ltoadd character varying(12),
  rfromadd character varying(12),
  rtoadd character varying(12),
  zipl character varying(5),
  zipr character varying(5),
  featcat character varying(1),
  hydroflg character varying(1),
  railflg character varying(1),
  roadflg character varying(1),
  olfflg character varying(1),
  passflg character varying(1),
  divroad character varying(1),
  exttyp character varying(1),
  ttyp character varying(1),
  deckedroad character varying(1),
  artpath character varying(1),
  persist character varying(1),
  gcseflg character varying(1),
  offsetl character varying(1),
  offsetr character varying(1),
  tnidf numeric(10,0),
  tnidt numeric(10,0),
  the_geom geometry,
  CONSTRAINT edges_pkey PRIMARY KEY (gid),
  CONSTRAINT enforce_dims_the_geom CHECK (st_ndims(the_geom) = 2),
  CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'MULTILINESTRING'::text OR the_geom IS NULL),
  CONSTRAINT enforce_srid_the_geom CHECK (st_srid(the_geom) = 4269)
)
WITH (
  OIDS=FALSE
);
ALTER TABLE tiger.edges OWNER TO postgres;
GRANT ALL ON TABLE tiger.edges TO postgres;
GRANT ALL ON TABLE tiger.edges TO gis_group;

-- Index: tiger.idx_edges_tlid

-- DROP INDEX tiger.idx_edges_tlid;

CREATE INDEX idx_edges_tlid
  ON tiger.edges
  USING btree
  (tlid);

-- Index: tiger.idx_tiger_edges_countyfp

-- DROP INDEX tiger.idx_tiger_edges_countyfp;

CREATE INDEX idx_tiger_edges_countyfp
  ON tiger.edges
  USING btree
  (countyfp);

-- Index: tiger.idx_tiger_edges_tfidl

-- DROP INDEX tiger.idx_tiger_edges_tfidl;

CREATE INDEX idx_tiger_edges_tfidl
  ON tiger.edges
  USING btree
  (tfidl);

-- Index: tiger.idx_tiger_edges_tfidr

-- DROP INDEX tiger.idx_tiger_edges_tfidr;

CREATE INDEX idx_tiger_edges_tfidr
  ON tiger.edges
  USING btree
  (tfidr);

-- Index: tiger.tiger_edges_the_geom_gist

-- DROP INDEX tiger.tiger_edges_the_geom_gist;

CREATE INDEX tiger_edges_the_geom_gist
  ON tiger.edges
  USING gist
  (the_geom);


texts: Dissolve Overlapping Polygons and Combine Attributes

texts:
I have a selection of polygons (specifically circles) in one feature class. Some of these overlap, some do not. The feature class is set up as such.

 OBJECTID   MY_KEY  SHAPE_LENGTH    SHAPE_AREA
 1          1       3141.592654     785398.1634
 2          4       3141.592654     785398.1634
 3          5       3141.592654     785398.1634
 4          8       3141.592654     785398.1634
 5          10      3141.592654     785398.1634

Polygons with MY_KEY 4 and MY_KEY 10 overlap. I want to create a new feature class that combines overlapping polygons and also creates a new field that lets me know what polygons have been combined, so I would end up with something like

SHAPE_LENGTH    SHAPE_AREA  POLYGON_KEY
3141.592654     785398.1634 POLYGON_1
4398.229716     1539380.4   POLYGON_4_10
3141.592654     785398.1634 POLYGON_5
3141.592654     785398.1634 POLYGON_8

I know that dissolve is the tool to use for combining overlapping polygons, but I do not how to achieve separate features for all non-overlapping polygons and obtaining a way to know what polygons have been combined where they do overlap.

How can I achieve this?


texts:
For elevation lines, I would like to have a different style of

full: 1 m
half: 0.5 m
quarter: 0.25 and 0.75 m

lines. Manually, I managed to extract the different lines with the following via field calculator in the attibute table:

Extract 1 m values:                 if(  "ELEV"  = floor( "ELEV" ),"ELEV", '')
Extract 25, 50 and 75 cm values (name: “no1m”): if(  "ELEV"  = floor( "ELEV" ),'', "ELEV" )
Exctract 50 cm values:                  if(  ("no1m" - 0.5) = floor( "no1m" ), "no1m",'')
Extract 0.25 and 0.75 values:               if(  ("no1m" - 0.5) = floor( "no1m" ),'',"no1m")

Now that I have that problem more frequently so I want to setup a QGIS layer style file to fasten the process of styling and labelling the elevation lines.

For the labels, it works to label the full lines by setting up the rule [if(  "ELEV"  = floor( "ELEV" ),"ELEV", '')] and for the half and quarter lines, I use a rule without filter but with label [if(  "ELEV"  = floor( "ELEV" ), '',"ELEV")]. Since I just want the 1 m labels to be larger, it is ok like that.

Now for the style, I can extract all 1 m lines with [if(  "ELEV"  = floor( "ELEV" ),"ELEV", '')] but as soon as I set a rule for the half m lines  [if(  "ELEV"  = floor( "ELEV" ), '',"ELEV")], I dont see any lines. If I only click on one of the two labels, I can see it but together they seem to overrule each other and nothing is shown. How can this be solved? Also, how do I style the quarter lines in the rule based styling?


texts: Using wildcard to represent specific number of characters in ArcPy?

texts:
I have a bunch of feature classes that are labeled in this pattern

CL25_E6_5_4_3

CL266_E6_5_4

CL266_E6_5_4_2

CL4_E6_5

I want to make a list in ArcPy that contains all feature classes with the same number of characters after the E (ex number of characters after *E). So in this example, only feature classes called CL25_E6_5_4_3 and CL266_E6_5_4_2 would be included in the list (because they have the same number of characters after the E) which will eventually be used to union the listed feature classes together.

Is there a wildcard symbol that represents only one character (Ex: CL266_E??????? would represent 7 characters after CL266_E) or is there a way to retrieve a count of characters after a point in the name (Ex CL266_Elen()) or specify number of characters to be in the wildcard (CL266_E[0:6])?

I'm not experienced with python, but am just trying to come up with ideas.


texts: Installing external python packages for QGIS 3?

texts:
While there are tools and methods for QGIS 2 I cannot figure out for 3
I want to install external Python packages such as sqlalchemy etc.
I tried the following things :

Tried to download osgeo4w for QGIS 3.02 as of now there is no advance installation option like it use to be before  http://download.osgeo.org/osgeo4w/osgeo4w-setup-x86_64.exe.

Tried to install setuptools but its not working in the 2nd step. it
returns 1.
https://www.lutraconsulting.co.uk/blog/2016/03/02/installing-third-party-python-modules-in-qgis-windows/

Went to Python of the C:\Program Files\QGIS 3.0\bin and tried pip
but it does not work.


texts: How to subset X number of features from shapefile in QGIS

texts:
I'm looking for a means to take a rather large shapefile and split it out into separate, smaller shapefiles based on number of features.  So for example, if I have a file with 1 million features, subset it out into smaller files of 100k features each.

I know how to grid a shapefile and subset into smaller files, but I don't want to grid it as I want to preserve the original feature geometry (e.g. gridding polygons will ruin the geometry).

Does anyone know how this can be done with QGIS?


texts: Error in write Array in Raster file: ""ValueError: expected array of dim 2""

texts:
I'm using Python in my research and how I don't have large experience in this language, I need Help to understand this error that occur when I`ll write a array as raster file.

data_res_b1_corrected = [np.zeros((col1, row1), np.float32)]

data_res_b1_corrected = (data_res_b1 - (coef_red_nir * (subt_nir_nir_med)))
print(data_res_b1_corrected.shape)
output1 = res_b1.GetDriver().Create('D:/CESAR_PHD/10.Lyzenga/output/Reservoir_RED_corrected.tif', col1, row1, 1, gdal.GDT_Float32)
output1.SetGeoTransform(geoT)
output1.SetProjection(proj.ExportToWkt())
output1.GetRasterBand(1).WriteArray(data_res_b1_corrected)
output1.FlushCache()

The error occur when I`ll run this line:  "output1.GetRasterBand(1).WriteArray(data_res_b1_corrected)"

The console returns: ValueError: expected array of dim 2


texts: Query areas from DynamicMapService Layer?

texts:
I have two dynamicMapServiceLayers. I have to query the area where three dropdowns are there. From a districts dropdown I have to select one district according to that mandals should be loaded to the Mandal dropdown, upon selecting mandal villages will come to villages dropdown. Then after selecting village, the geometry area of the village should zoom.

How to do retrieve that data from one of the layer?


texts: .exe Packaging for Python files with ArcPy?

texts:
I've written a little program that's meant to function as a stand-alone, not as an add-in. I'd like to be able to pack it up as an .exe, rather than ask my users to run the script. It does import ArcPy.

I tried pyinstaller, py2exe, cx_Freeze, several other things ...
I cannot wrap this thing up into an executable. Can't figure it out.

I found cx_Freeze the least ornery thus far, yet I still cannot pack up my Python (2.7) into an executable successfully.
What am I missing?

Here is the error I continue to receive:

import arcgis scripting
ImportError: DLL load failed: A dynamic link library (DLL) initialization routine failed.

(The image is an example of what the error looks like. I'm using ArcGis 10.4)

I see in the ESRI community that someone else had this problem as well, 5 years ago, scroll to bottom of link:
https://community.esri.com/thread/57485

Yet, no answer from the ESRI community.

Here is my code, feel free to chuckle. All the additional mods and included files are from chasing around error messages until I was stumped.

This file, setup.py, is in the same PyCharm Project with my program:

import sys
import os.path
from cx_Freeze import setup, Executable
import arcpy

# source:
# https://www.youtube.com/watch?v=DoHWJV8iVTQ

PYTHON_INSTALL_DIR = os.path.dirname(os.path.dirname(os.__file__))
os.environ['TCL_LIBRARY'] = os.path.join(PYTHON_INSTALL_DIR, 'tcl', 'tcl8.5')
os.environ['TK_LIBRARY'] = os.path.join(PYTHON_INSTALL_DIR, 'tcl', 'tk8.5')

additional_mods = ["numpy.core._methods", "numpy.lib.format"]
# Dependencies
build_exe_options = {"packages": ["os"], "include_files": ["sqlite3.dll", "tcl85.dll", "tclpip85.dll", "tk85.dll"],
                     "includes": additional_mods}

base = None
if sys.platform == "win32":
    base = "Win32GUI"

exe = Executable(script="MyProgramINeed.py", base=base, icon="myiconineed.ico")
setup(   name = "MyProgramINeed",
         version = "0.1",
         description = "description of my program",
         options = {"build_exe": build_exe_options},
         executables = [exe])


texts: Is there an event handler for a Selection being made using ArcPy?

texts: I want to be able to capture an arcpy.SelectByAttribute() or arcpy.SelectByLocation() event being triggered. There doesn't seem to be an extension that is able to handle this.  Is there a way to do this?


texts: How to sieve a raster with 1 and NaN values?

texts:
As a next step to my previous question "How to calculate where raster a is not equal to raster b?", I now want to sieve (using QGIS -> Raster -> Analysis -> Sieve, or gdal_sieve in the shell) a raster which consists of the values 1 and NaN.
The problem is: no matter what parameters I tell gdal_sieve to use, it still returns the original raster, nothing happens.

Above you can see the raster (black pixels) that I want to sieve. Actually, this screenshot was taken immediately after the sieve algorithm ran through. You will find that nothing was sieved.

I guess this is because gdal_sieve cannot handle NaN, so the nearest possible value always is 1 again, if it skips the NaNs.

Because I want to use Zonal Statistics with this raster and some polygons afterwards, it is necessary that this raster consists just of 1 and NaN, because I am interested in the total count of pixels with value 1 in a polygon. If I would use a raster with, say 0 and 1, it would always accumulate these two in my count column.

Does anyone have an idea how to workaround this?


texts:
I am using QGIS 2.18.

In the overview map, is it possible to remove attribute labels?

I just want the labels in the primary map, but not the overlay.

All that yellow in the overlay image is the labels for the attributes.


texts: Mass downloading Google Earth historical imagery

texts:
I've been working with the historical imagery in Google Earth and I'm wondering if there's an easy way to download the geoTIFFs for my study area. It doesn't seem possible from Google Earth itself.

Is there a website where I drag a box around/upload a polygon of my study area and be able to download geoTIFFs of all the historical satellite images for that area in one zipped file?