texts: Performing Full Text Search of Esri geodatabase?

texts:
How would you go about searching the contents of features, the values in the attribute columns, in an Esri geodatabase? (preference for file-gdb, but SDE would be useful too.)

I can envisage workflow where one dumps the whole gdb using Export XML Workspace and searches the result for strings (yuk!), or loops through all the tables in a python version of Select by Attributes. While both of these would work, they are not appealing.

The goal is to be able to efficiently answer a question like: which of these 50 feature classes has something to do with "Keno"?


texts: Counting longest distance between two lines using QGIS?

texts: How do I count longest distance between two lines as shown at the picture in every bump using QGIS?


texts: Where are the radiation models for Japan catastrophe?

texts:
I'm interested in seeing authoritative model results showing how radioactive contaminants from Japan might circulate through wind and water.

The AGU blog points out this bad model, but in skimming the sites they reference I didn't see any other models.  I see no mention of radiation at NOAA's Global Forecast System. Where are the reputable models?


texts:
Can anyone help explain me how to use the ContourGenerate method of GDAL.

The official documentation states two sets of arguments for the method, with one being

ContourGenerate ( **args, ***kwargs )

And the other being:

ContourGenerate(Band srcBand, double contourInterval, double contourBase, int fixedLevelCount, int useNoData, double noDataValue, Layer dstLayer, int idField, int elevField, GDALProgressFunc callback=0, void * callback_data=None)

I do not understand how this is supposed to work.

I tried working my way by trial and error method, by intentionally leaving out the arguments and filling them as per the demands of error message. The method kinda worked, by asking me the arguments of second definition (in the same order as given above). This worked up until I entered elevField = -1, and then it just gives an error message:

  TypeError: not a sequence

I tried googling the error but to no avail. I also tried including the callback=0, callback_data=None part in the argument list, but the same error occurs.

Can anybody help me out?


texts: GUID to open Geoprocessing Enviroment Settings Window through C#?

texts:
I am trying to let the user click the enviroment button at the bottom of my custom form and have the Enviroment Window Pop Up. I cannot seem to locate the Enviroment Settings Window GUID anywhere?

    void environmentsButton_Click(object sender, EventArgs e)
    {
        try
        {
        UID uid = new UIDClass();

        uid.Value = "esriGeoprocessing.GPEnviromentSettingsWindow";
        Debug.WriteLine("Setting");
        ESRI.ArcGIS.Framework.IDockableWindow win = Program.ActiveDockableWindowManager.GetDockableWindow(uid);
        Debug.WriteLine("Showing");
        win.Show(true);
        }
        catch (Exception ex)
        {
            Debug.WriteLine(ex.Message);
        }
    }


texts: Not able to get the PostGIS Layer Viewer running in pgAdmin or command line

texts:
I've installed PostGIS Layer Viewer as a plugin for pgAdmin. I've downloaded postgis_viewer.py from the official repository on GitHub! and followed the install instructions described on Web site!. I've also checked the installation prerequisites and made sure that all ini-files are set accordance with my own directories layout. After running the plugin in pgAdmin the application stalls. And by stalling I mean that no layer becomes visible.

The pgAdmin version I am using is version 1.20.0 on Windows7 64bit OS. The QGIS version I am using is 2.8.1.

The question I have is how to get the plugin running?

NB! I also tried running the application from the command line but the results were the same.


texts:
I have a multilayer image with hundreds of bands and I would like to extract them all automatically to single raster files, one file for every layer.

I'm looking for an automated process. I've got aware of the Orfeo Toolbox split image tool, but I don't manage to make it work.

Any suggestion?


texts: Resampling while drawing rasterfiles in QGIS 3?

texts:
I can not find a solution to resample my rasterfiles in QGIS 3. With Arcmap usually I use "resample while drawing" with the method cubic or bilinear interpolation for continuous data to smooth my drawing (like the second image). The cell size is quite big so it makes no sense for my task (heatmap) to visualize the rastercells themselves.

Original rasterdata:

Resampled while drawing:


texts: How to fix ImportErrors when loading Python plugins in QGIS?

texts:
I have recently installed QGIS 1.8 in Windows 8. Although QGIS seems to be working fine, I cannot load plugins after their installation. Instead I keep getting this message:

Traceback (most recent call last):
File "C:/OSGeo4W/apps/qgis/./python\qgis\utils.py", line 143, in loadPlugin
__import__(packageName)
File "C:/OSGeo4W/apps/qgis/./python\qgis\utils.py", line 309, in _import
mod = _builtin_import(name, globals, locals, fromlist, level)
ImportError: No module named linear_ref

Any ideas?


texts:
I do not understand the function behaviour. I am giving the function a point B that is North of point A so I would expect something around 90 degrees but I get 263 degrees instead. In another example, I keep the same point A and I am giving the function a point B that is on the South of point A and this time I do get a correct value of 279 degrees.

The function definition states that it is using the horizontal of the vector between A and B. But I don't understand the part which says:

  Angle is computed clockwise from down-to-up: on the clock: 12=0;
  3=PI/2; 6=PI; 9=3PI/2.

Here are the examples, maybe I am doing something wrong with the SQL
point b NE of point A:

SELECT degrees(
                ST_Azimuth(
                    ST_Transform(ST_GeomFromText('POINT(127.702 38.169)',4326),32652),
                    ST_Transform(ST_GeomFromText('POINT(126.728 38.178)',4326),32652)
                           )
                )

This returns 271 degrees whilst it should return ~45 degrees

Someone care to explain? Thanks in advance!


texts: How to CORS enable ArcGIS Server 10.2.1 to Access REST Services without Using proxy.ashx?

texts:
Form the ENABLE-CORS.ORG I see that editing web.config is enough to CORS enable IIS7

<?xml version="1.0" encoding="utf-8"?>
 <configuration>
  <system.webServer>
   <httpProtocol>
     <customHeaders>
       <add name="Access-Control-Allow-Origin" value="*" />
     </customHeaders>
   </httpProtocol>
  </system.webServer>
 </configuration>

Anything else I need to do, so my ArcGIS REST end points are accessible from ArcGIS JavaScript API web application without using proxy page:

esri.config.defaults.io.alwaysUseProxy = false;
esri.config.defaults.io.corsEnabledServers.push("server.organization.com");


texts: How do I use gdal_proximity.py from the QGIS python console?

texts: How do I run gdal_proximity.py from the QGIS python console?


texts: Join attributes of polygons to another polygon feature based on location QGIS

texts: I have two large polygon data-sets. One highlights the building footprints for a county and the other highlights the parcel boundaries for the same county. Some parcel polygon features have multiple building footprint polygons that are within the boundary of the parcel polygon. I want to take each individual building footprint polygon feature and join the parcel ID number to each building footprint feature that lies within the larger perspective parcel polygon feature it's located in. The actual exported feature I want in the end are the building footprint features with the parcel ID # in the attributes.


texts: Export multiple Images in one task - GEE Python API

texts:
Context: I process a set of image using the Google Earth Engine Python API. I want to apply the process on a set of image and save the output (in my case export to GEE Assets).
What I have done: Looping through the image set and export each image one by one. That works well.
Problem: since I would like to apply my process over 60k images, I face the  number of task restriction (3000).

Traceback (most recent call last):
ee.ee_exception.EEException: Too many tasks already in the queue (3000). Please wait for some of them to complete.

What I found: the geetools modules. However, by looking at the code, it also seems to iterate over the image collection and exports each image one by one.
Question : Is it possible to export multiple images in one task ? What are the solutions ?
Here the code:
col = ee.ImageCollection("COPERNICUS/S2") \
        .filterDate("2018-01-01", "2018-12-31") \
        .filterBounds(geometry)

# Get names of all image as a python list
image_names = get_name_collection(col).getInfo()
total = len(image_names)

# Create new image collection if don't existe
folder = "users/ab43536/masks_4_methods"
if folder not in [elt["id"] for elt in ee.data.getList({"id": "users/ab43536"})]:
    ee.data.createAsset({'type': "ImageCollection"}, folder)

# Apply process to all images
for i, name in enumerate(image_names):
    print("{:4d}/{} = {:05.2f}%   Image {}".format(i, total, i / total * 100, name))
    # The long process
    mask = computeCloudMasking(name)

    # Export image one by one
    export_image(mask, folder, getGeometryImage(ee.Image(name)),
                 name.split('/')[-1], num=i, total=total)

Where export_image is:
def export_image(image, asset_id="users/ab43536/", roi=None, name=None, num=None, total=None):
    """ Export one image to GEE Asset
    Arguments
        :param image: image to export
        :param roi=None:  specify the roi, default compute from image dimension
        :param name=None: name of the image
    """
    if roi == None: roi = getGeometryImage(image)
    if name == None: name = image.id().getInfo()
    description = "Default export"
    if num != None and total != None:
        description = "Image {} on {} equal {:05.2f} pourcent".format(num, total, num / total * 100)
    # print(description)
    assetId = asset_id + name

    # Create a task : export the result as image asset
    task = ee.batch.Export.image.toAsset(image=image.clip(roi),
                                         description=description,
                                         assetId=assetId,
                                         scale=30,
                                         region=roi.coordinates().getInfo(),
                                         )
    task.start()


texts:
I am interested in making reproducible PDF figures in QGIS3, and therefore want to work from the console/PyQGIS.

Using some previous helpful advice from PDF Output with Python the problem is successfully reduced to setting up the QgsPrintLayout. I am not sure how to do this, below is my best attempt so far:

#setup the background classes for managing the project
canvas =iface.mapCanvas()
manager = project.layoutManager()
#new layout
layout = QgsPrintLayout(project)
layout.initializeDefaults()
layout.setName('console')
itemMap = QgsLayoutItemMap(layout)
itemMap.setExtent(canvas.extent())
layout.addLayoutItem(itemMap)
#add the layout to the manager
manager.addLayout(layout)

The "console" layout is successfully added to QGIS 3 and I can open it to see if the content is satisfactory. The map does not render and its extent is nan,nan,nan,nan.

How do I setup the extent and size of the map object in the layout?

SOLVED. The following code produces a map without issues.

#get a reference to the layout manager
manager = project.layoutManager()
#make a new print layout object
layout = QgsPrintLayout(project)
#needs to call this according to API documentaiton
layout.initializeDefaults()
#cosmetic
layout.setName('console')
#add layout to manager
manager.addLayout(layout)
#create a map item to add
itemMap = QgsLayoutItemMap.create(layout)
#using ndawson's answer below, do this before setting extent
itemMap.attemptResize(QgsLayoutSize(6,4, QgsUnitTypes.LayoutInches))
#set an extent
itemMap.setExtent(canvas.extent())
#add the map to the layout
layout.addLayoutItem(itemMap)


texts:
I'm trying to calculate the average house price within different sized buffers using MapInfo.

I have thousands of properties on my map as points (in one layer) and around 60 different sized buffers. Individually I can click each buffer using the boundary select, and then open the statistics tab and manually write down the number of points within the buffer and the average price.

But it's extremely slow, is there an export tool that will allow me to export the statistics tab figures into a table, for all 60 buffers?


texts: I installed OSGeo4W using a Windows installer... it contains GDAL 1.7.3.  I can't get GDAL to recognize PDF as a file type using commands in MSYS shell.  Does this version of GDAL not support PDFs?  If not, what do I need to do to add support to PDFs for GDAL?


texts: Area Attribute in Feature Classes Gets Renamed to full FeatureClassName.AREA

texts:
I have a number of feature classes in an Oracle-backed Enterprise Geodatabase which have an attribute called "AREA" which is type double. All of these feature classes are part of the database user/schema called "AMIS". The problem, is that the attribute is automatically getting renamed to the AMIS.FeatureClassName.AREA.  For example, in the feature class called CAMP_AREA, the AREA field is getting called "AMIS.CAMP_AREA.AREA". I have tried deleting the field and re-adding "AREA" again, and it just defaults to "AMIS.CAMP_AREA.AREA".

I've checked the list of Oracle Reserved Words, and AREA is not one of them.
I have actually also seen this behaviour on a Microsoft SQL Server-backed Geodatabase with the same results.
What is going on here?


texts: Generate dynamic XY coordinate labels for dataframe corners in ArcGIS layout

texts: I need to label the four corners of a dataframe as shown in the image. Can this be done with dynamic text scripts?


texts: ST_MakeLine returns Unsupported geometries not rendered in geometry viewer

texts:
I am trying to run a ST_MakeLine query that will return the shortest distance between given coordinates
The code is as follows:
SELECT ST_MakeLine(route.geom)  as test
FROM (SELECT *
      FROM pgr_fromAtoB('roads2_noded',73.1789,33.5183,73.1768,33.5191)
      ORDER BY seq
     ) As route

The result looks like this:

it's clear that the result does not show null. But when I view it in the Geometry Viewer, it gives the following message:

Can someone tell me what this message means and how can I find a workaround to this problem?
The geometry type of my roads shapefile is MULTILINESTRING.
St_AsText (route.geom) shows the following result:
SELECT ST_AsText( ST_MakeLine(route.geom)) FROM (SELECT *FROM pgr_fromAtoB('roads2_noded',73.1789,33.5183,73.1768,33.5191)
      ORDER BY seq
     ) As route

Using simple Select from pgr_fromAtoB gives:
SELECT * FROM pgr_fromAtoB('roads2_noded',73.1789,33.5183,73.1768,33.5191) ORDER BY seq

The roads2_noded is the roads shapefile I am working on shown below:


texts: How do you pull up the entire attribute table in Esri Collector app

texts: I am new to the Esri Collector app. How do I pull up the attribute table using the Esri Collector app on a tablet? I have searched everywhere on line and Esri but can't seem to find a simple answer to the question.


texts:
I have a set of layers from PostGIS that I have converted to offline. I am attempting to edit them via GPS.

If the layer is point data, then I can add in that data without issue. However when I go to edit a polygon or line, it will not allow be to add a feature based on GPS location.

The geometry i'm using is MultiLine and MultiPolygon - I am about to try and edit this, but is there anything else causing issue here?

Secondly - QGIS docs say I can convert a WFS to offline - however whenever I try to convert my WFS (served through GeoServer) it just crashes.


texts: Line network integration analysis in QGIS

texts:
I want to analyse how integrated a line (e.g. road, stream) is within a network.
I want to perform the "integration" type of analysis identified by Javadi (not the connectivity or choice analysis) i.e. the integration analysis performed by Depthmap (I've been unable to get the analysis to work in Depthmap due to presumably data issues).

I think I just need a plugin, but I don't know which one does this analysis.

Figure illustrating the difference between connectivity, integration and choice for a sample network. Taken from: Javadi, A., et al. (2017) “Hippocampal and prefrontal processing of network topology to simulate the future”. Nature Communications.


texts: Go from multiple rasters to R compatible table

texts:
I have a set of geocoded, pixel-aligned rasters which I would like to use in a model in R. I would like to produce a table of data that can be imported to R, in which there are as many columns as rasters, and as many rows as points at which data is extracted from the rasters. Say I have 3 rasters: r1, r2, r3. I would end up with a table that has r1, r2, r3 as the column headers. The rows would be point1, point2, point3 - and there would be as many rows as points in the layer used to define which points to extract data from. Forgive the bad table:

------r1---r2---r3

point1-1----4---2

point2-4----3---1

How can I do this? I've tried some Arc tools but they have failed or hung. I can use QGIS, Arc or Idrisi.

Is the Raster package for R a good approach?


texts: Filtering for all buildings with given house number on given streets using OverPass Turbo?

texts:
I already created the following script to filter for all living_streets in a known postal_code.

area[postal_code="12345"]->.a;
(way(area.a)["highway"="living_street"];
 rel(area.a)["highway"="living_street"])->.w;
(.w;.w >;)->.w;
.w out;
(._;>;);
out skel;

Now I need all buildings with a known housenumber.
My problem is that I found no way to get the buildings on the streets/ ways.


texts: QGIS: Change attribute values from NULL to 0 for all attributes using Python

texts: I have an attribute table where I have occurrences of NULL spread throughout in different attributes and features. How can I change NULL to 0 using Python for all attributes?


texts:
I work with GDAL 2.2.1, released 2017/06/23

I had some troubles when trying to resample GeoTiff images to make different data sources match with extent and resolution. I found an answer here: GDAL: Raster Intersection

But I identified a bug that wasn't notified here : when I use float values for the extent -te, for example with longitude / latitude values, the result image is only filled bith NoData values !

For instance, this gives ND image, input and output are EPSG:4326:

gdalwarp -te 7.5 43.5 8.5 44.5 -tr 0.005 0.005 -r bilinear input.tif output.tif

While this gives a correct image with values:

gdalwarp -te 7 43 9 45 -tr 0.005 0.005 -r bilinear input.tif output.tif

Is it a known bug of GDAL that anyone else has experienced?


texts: wrapped C/C++ object of type QgsGeometry has been deleted ubuntu and mac and segmentation fault

texts:
I have a problem with a specific method in qgis 2.8, 2.12 and 2.12.1 on linux and mac (I have not idea if in windows this problem persists)

The example:

from qgis.core import *
from qgis.utils import iface
from ceamos.utils import *
from ceamos.config import *

def run_script(iface):

    #remove all layers
    clearAll(iface)

    lyr = QgsVectorLayer(project_path + '/qgis_data/area.shp',
                          'laflorida',
                          vectorLayerType.ogr)

    assert lyr.isValid()

    QgsMapLayerRegistry.instance().addMapLayers([lyr])
    geom = lyr.getFeatures().next().geometry()
    print geom # <qgis._core.QgsGeometry object at 0x7f3a08715510>

    print geom.contains(QgsPoint(0,0)) # RuntimeError: wrapped C/C++ object of type QgsGeometry has been deleted

Update 25-DIC-2015

I have changed the code to

from qgis.core import *
from qgis.utils import iface
from ceamos.utils import *
from ceamos.config import *

def run_script(iface):

    #remove all layers
    clearAll(iface)

    lyr = QgsVectorLayer(project_path + '/qgis_data/area.shp',
                          'laflorida',
                          vectorLayerType.ogr)

    assert lyr.isValid()

    QgsMapLayerRegistry.instance().addMapLayers([lyr])
    #geom = lyr.getFeatures().next().geometry()
    geom = QgsGeometry(lyr.getFeatures().next().geometry())
    print geom # <qgis._core.QgsGeometry object at 0x7f3a08715510>

    # RuntimeError: wrapped C/C++ object of type QgsGeometry has been deleted
    #print geom.contains(QgsPoint(0,0)) #Segmentation fault

    feat = lyr.getFeatures().next() #Segmentation fault
    #print feat.geometry().contains(QgsPoint(0,0))

But in the line feat = lyr.getFeatures().next() I got a segmentation fault error (in mac and ubuntu 14.04)


texts: Calculating distance between point and MultiPolygon Geoseries object in Python

texts:
I want to calculate the distance of a point from the Multiploygon geoseries object created using cascaded_union function of GeoPandas. Following is my code:

import numpy as np
import pandas as pd
import geopandas as gpd
from shapely.ops import cascaded_union
from matplotlib import pyplot as plt
from shapely.geometry import LineString, Point, MultiPoint

fp = "\\TM_WORLD_BORDERS-0.3\\TM_WORLD_BORDERS-0.3.shp"
data = gpd.read_file(fp)

mycountries = ['Canada','Russia']

req_data = data[data.NAME.isin(mycountries)]
req_data = req_data.reset_index()
req_data = req_data.drop('index',axis=1)

polygons = [req_data['geometry'][0],req_data['geometry'][1]]
boundary = gpd.GeoSeries(cascaded_union(polygons))
boundary.plot(color = 'black')
plt.show()

I require that when I calculate the distance of two points (as shown on the map: red dots) from the MultiPolygon, it should consider Canada for the one in the left side and Russia for the right one. Also, the distance in miles/km.

I tried using distance function but the result which I am getting is not correct. Moreover, if I consider point somewhere in the Indian Ocean, it gives a result as 0 which is again not what I require.

Shapefile can be downloaded from the TM_WORLD_BORDERS-0.3.zip file available at http://thematicmapping.org/downloads/world_borders.php


texts:
I need a few customizations over grass and would use qgis gui for accessing grass. I am building qgis with grass enabled (WITH_GRASS=TRUE). While I am doing that I get the following.

Linking C executable ../../../output/lib/qgis/grass/modules/qgis.g.info
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRGetLinearUnits@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `VSIFree@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRExportToWkt@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRGetAttrValue@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `SetCSVFilenameHook@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRExportToProj4@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRIsProjected@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRNewSpatialReference@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRExportToPrettyWkt@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRMorphFromESRI@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRIsGeographic@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRMorphToESRI@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRImportFromProj4@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRGetUTMZone@GDAL_1.8'
/usr/lib/grass64/lib/libgrass_gproj.so: undefined reference to `OSRDestroySpatialReference@GDAL_1.8'
collect2: ld returned 1 exit status

***
GDAL 1.9.2, released 2012/10/08
GRASS version: GRASS GIS 6.4.2
Python: 2.7.3
Platform: Ubuntu 12.04.1 LTS / x86_64
***
The C compiler identification is GNU
The CXX compiler identification is GNU
Check for working C compiler: /usr/bin/gcc
Check for working C compiler: /usr/bin/gcc -- works
Detecting C compiler ABI info
Detecting C compiler ABI info - done
Check for working CXX compiler: /usr/bin/c++
Check for working CXX compiler: /usr/bin/c++ -- works
Detecting CXX compiler ABI info
Detecting CXX compiler ABI info - done
QGIS version: 1.8.0 QGIS (10000)
Found GRASS: /usr/lib/grass64 (6.4.2)
Touch support disabled
Looking for openpty
Looking for openpty - not found
Looking for openpty
Looking for openpty - found
Found Proj: /usr/local/lib/libproj.so
Found GEOS: /usr/lib/libgeos_c.so
Found GDAL: /usr/local/lib/libgdal.so (1.9.2)
Found Expat: /usr/lib/x86_64-linux-gnu/libexpat.so
Found Spatialindex: /usr/lib/libspatialindex.so
Found Qwt: /usr/lib/libqwt-qt4.so (5.2.2)
Found Sqlite3: /usr/lib/x86_64-linux-gnu/libsqlite3.so
Found PostgreSQL: /usr/lib/libpq.so
Looking for Q_WS_X11
Looking for Q_WS_X11 - found
Looking for Q_WS_WIN
Looking for Q_WS_WIN - not found.
Looking for Q_WS_QWS
Looking for Q_WS_QWS - not found.
Looking for Q_WS_MAC
Looking for Q_WS_MAC - not found.
Found Qt4: /opt/QtSDK/Desktop/Qt/4.8.1/gcc/bin/qmake (found suitable version "4.8.1", required is "4.4.0")
Pedantic compiler settings enabled
Found PythonInterp: /usr/bin/python (found version "2.7.3")
Found Python executable: /usr/bin/python
Found Python version: 2.7.3
Found Python library: /usr/lib/libpython2.7.so
Found SIP version: 4.13.2
Found PyQt4 version: 4.9.1


texts: How to scatter all points at identical locations in CartoDB?

texts:
I am working on a project in which data points are often at the same location (e.g. Rome). This poses a problem because the data points then overlap and cannot be clicked on.

I found some SQL at http://bl.ocks.org/andrewxhill/8695515 that scattered the points very well when I applied it to one location. However, my final dataset could have dozens of identical locations, and I don't want to enter the same SQL dozens of times for each identical location.

Is there a way to modify this SQL (or something similar) so that it scatters all points at the same location (i.e. it identifies the locations at which there are multiple data points, and then scatters the points at each of those locations)?


texts: Lat Long and timestamp to hyperlink transect videos (YouTube)

texts:
I have a down looking GOPRO video of a shallow reef (from the surface) that has been synced with GPS point data (XY, Lat Longs) from the GPS track. Have uploaded a video to YouTube which has a "share" hyperlink that can be set to any second of the video. I have XY data for every second of the video. I would like to have the point GPS data hyperlinked to open the YouTube video at the correct reference second. I think that Google Earth might be the best way to share and interact with a georeferenced video feed (so would like to build a KML file, unless there is a better way?).
How do I approach this?


texts: Add save/commit button to custom feature from

texts:
I would like to add a save/commit button to a custom QGIS feature from. This is very helpful when there is a linked table dispayed on the same form - the link will work only once the main table has been saved.

The form should look something like this:

Which function should I call for the save/refresh button?

After some research, I found the following solution (which is quite clumsy!)

from PyQt4.QtCore import *
from PyQt4.QtGui import *

myDialog = None

def formOpen(dialog, QGIS_vector_layer, featureid):
    global myDialog
    myDialog = dialog

    global ly
    ly= QGIS_vector_layer #this is a QGIS layer object, not a layer ID

    # Connect a signal
    Refresh_button = dialog.findChild(QPushButton,"Cmd_refresh")
    Refresh_button.clicked.connect(refresh)

def refresh():

    myDialog.accept() #step 1 dump data into editing buffer
    ly.commitChanges() #step 2 save (but this closes the editing mode!)
    ly.startEditing() #step 3 re-open for editing

This works, but it is still not the best method. The form is not emptied from the data entred, which means that each time we hit either refresh or OK button, a new record will be added. But that would be another problem.


texts: How to set my ESPG to geography POINT in PostGIS?

texts:
I'm trying to create a table into my PostGIS 2.3

CREATE TABLE destinations (
    id          serial                      PRIMARY KEY,
    sector      integer                     REFERENCES sectors(id),
    name        varchar(64)                 NOT NULL,
    location    geography(POINT, 25831)     NOT NULL,
    authcode    varchar(16)                 NOT NULL
);

I specify the SRS ESPG:25831 that I'm using to the location field, but I get this error:

ERROR:  Only lon/lat coordinate systems are supported in geography.
LINE 6: location    geography(POINT,25831)      NOT NULL,

I've checked the spatial_ref_sys table and 25831 is there. I also have enabled the PostGIS extension. In fact, my query works if I remove the SRS:

location    geography(POINT)     NOT NULL,

Any idea of what am I doing wrong?


texts: How to save several asc files to geotiff

texts:
In QGIS (2.18), I have several raster files (.asc-format) that have no projection information.
The images display well at EPSG:5677. My aim is to store all the files in a new folder, using the GeoTIFF format and assigning EPSG:5677.

I have no python experience, however some python coding could help me here, I guess. I read something similar but shapefile-related here:
Exporting several files at same time in QGIS?

Of course, this doesn't work with my raster files. So here goes my question: how can I export a bunch of raster files, using the original file names and assigning the desired coordinate reference system without having to export every single file manually?


texts: pyQGIS - How to make a reproducible example with dummy rasters?

texts:
I'm having troubles with a pyQGIS script. I'm trying to create a large number of maps showing NDVI in a certain area over the time. The script opens a raster file, applies a style, changes title and exports map.

I made a loop to do that and it doesn't work, the canvas doesn't refresh and many maps are from the same raster file. I Use solutions like How to refresh MapCanvas while PyQGIS-script is running? or Map on canvas renders after image saving [duplicate] with no success.

But, I'm looking a way to create dummy rasters/vectors prior to recreates pyQGIS problems. R have excellent guides, and this also can be applied to spatial data.

So, is there a way to create dummy rasters with random numbers to recreates problems like this?

By the way, the script is something like the code below. I'm looking a way to replaces files item into dummy rasters.

from qgis.core import QgsProject
from qgis.gui import *
from PyQt4.QtCore import QSettings, QTranslator, qVersion, QCoreApplication,QFileInfo,QSizeF,QTimer
from PyQt4.QtGui import QAction, QIcon, QPrinter,QPainter
import os
import glob

files = glob.glob(r'Path/to/monthly/NDVI/*.tif')

aaa = [some list]

n = 287 # 12 year

for i in range(0,n):
    project = QgsProject.instance()
    project.read(QFileInfo('some/Project.qgs'))
    layer1 = QgsRasterLayer(files[i],'NDVI')
    uri = '/path/to/estilo.qml'
    layer1.loadNamedStyle(uri)
    root = QgsProject.instance().layerTreeRoot()
    group = root.findGroup('Raster')
    QgsMapLayerRegistry.instance().addMapLayer(layer1, False)
    group.insertChildNode(-1, QgsLayerTreeLayer(layer1))
    iface.mapCanvas().refresh()
    iface.mapCanvas().refreshAllLayers()
    composers = iface.activeComposers()
    composers[0].items()[2].setText('Test' + str(files[i]))
    c = composers[0].composition()
    c.items()[8].refreshDataDefinedProperty()
    c.items()[8].updateBoundingRect()
    c.refreshItems()

    def ExP():
        imgOut = c.printPageAsRaster (0)
        imgOut.save(os.path.join(r'/path/to/export/',aaa[i]+'.png'),"PNG")
        return

    QTimer.singleShot(2000, ExP)

    def ReM():
        QgsMapLayerRegistry.instance().removeMapLayer( layer1.id() )
        return

    QTimer.singleShot(2000, ReM)

    QgsProject.instance().clear()


texts: Iterating over a raster layer in ArcGIS using Python

texts:
I have a raster layer for drive times.  I have another raster layer with my bounding geography.  I want to tell ArcGIS:

For every cell in geography:
    if the corresponding cell in drive time is not null:
          geography cell = drive time
    else:
          geography cell = 120

How do you iterate over cells in a raster layer in ArcGIS?


texts: Ways to halo map text over a varying background?

texts:
I am ambivalent with using a halo effect. I am currently working on a small scale map of Africa showing the different types of vegetation and there are some crowded areas.

In terms of consistency, if I am going to introduce a halo for labels, I should add this effect to all labels of the same type (e.g. countries). However, it does just not look so nice in areas with single backgrounds.

What is the best practice?


texts: How to burn geometry to raster without shapefile, just from list of coordinates?

texts:
I have a list of coordinates points(x,y) of some multipolygon like this:

[[x1,y1][x2,y2][x3,y3] ... [xN,yN]]

Also can convert between lat/lon and UTM coord systems.

When we use gdal.Raesterize like on this answer:

gdal.RasterizeLayer(raster_dataset, [1], shape_layer, None, None, [1], ['ALL_TOUCHED=TRUE'])

we get third argument in this method shape_layer from .shp file:

shape_datasource = ogr.Open(file_name.shp)
shape_layer = shape_datasource.GetLayer()

So, how can I use gdal.RasterizeLayer() with my list of coordinates to make mask for image instead shape_layer?


texts:
I am trying to use the ST_Intersection function to intersect points by county and find the counts by county, but I am having no luck.

SELECT counties.name,
FROM counties INNER JOIN routes
ON ST_Intersects(counties.geom, routes.geog1)
group by counties.name;

The error is telling me:  relation "counties" does not exist.

Can anyone help me fix my query. Should I be using the st_contains function?

*2nd day working with posgres/postgis so bear with me please!


texts: One CSV with field that delineates CRS for QGIS?

texts:
I have a CSV with points with X,y coordinates. One of the fields is "Horizontal Datum," which says what the projection is. Here is a sample of the sheet:

QGIS lets me import this csv file, but when I try to clip it into another layer, it says it has trouble with the CRS. Should I make this into multiple layers with different projections?


texts: Adding exit command to ArcGIS Python script tool?

texts:
I have a GIS script tool that creates a copy of an mxd, sets a definition query and exports a set of data driven pages based on the definition query. After it's finished exporting it deletes the mxd copy, some variables and a bunch of database connections.

Is there some way to code an escape code which would cancel the data driven pages exports wherever they are but instead of cancelling the whole script, it would just jump to the end and still delete the mxd copy etc. The block of code looks something like this:

    mxd = arcpy.mapping.MapDocument(mxd_master)
    mxd_folder = os.path.dirname(mxd_master)
    mxd_copy = mxd_folder + "\\WORKING.mxd"
    mxd.saveACopy(mxd_copy)
    mxd = arcpy.mapping.MapDocument(mxd_copy)

    ddp = mxd.dataDrivenPages
    for page_num in range (1, ddp.pageCount + 1):
        ddp.currentPageID = page_num
        row = ddp.pageRow
        mapsheet = row.getValue(ddp.pageNameField.name)
        arcpy.mapping.ExportToPDF(pdf export parameters....)

    del mxd, mxd_master
    arcpy.Delete_management(mxd_copy)

Basically I want to be able to hit the "Q" key (or whatever key) and the pdf export loop would break and the del/Delete_management block would still run. I don't think it's possible since it's running the script through a tool but I figured I'd see if it was as it would stop me from having to manually delete things if I decide to stop the export early.

I am using ArcGIS Desktop 10.4.


texts:
In our project we need to use rasters offline on QField. Our tablets are a very limited software and using uncompressed rasters is too slow to be practical.

First, we try to convert our tif files to ecw, but didn't recognice the file (shows an error like the file is missing). Then we try to use compressed tiffs, with overviews:

gdal_translate -b 1 -b 2 -b 3 -co COMPRESS=JPEG -co PHOTOMETRIC=YCBCR -co TILED=YES ortho.tif ortho_ycbcr_ov.tif
gdaladdo --config COMPRESS_OVERVIEW JPEG --config PHOTOMETRIC_OVERVIEW YCBCR --config INTERLEAVE_OVERVIEW PIXEL -r average ortho_ycbcr_ov.tif 2 4 8 16

The same result. Then we try without overviews and then withour YCBCR with the same results.

The only option who works is to use overviews with no compression. The files are huge, but the performances is pretty decent.

¿Any success with compressed formats?

EDIT:
Now, I've realised that the overviews are generated with the command you see above. Who are compressed as JPEG and use YCBCR photometric command, but QField recognised it without trouble. That's very confusing.


texts: Calculating costs of a line geometry based on a different intersecting geometry

texts:
I have 2 tables, a route table with line geometry and a duct table also with line geometry. Basically I need to calculate the costs of the different routes based on their intersection with the duct table
 route table
 id integer
 geom geometry

 duct table
 id integer
 type integer
 geom geometry

To explain what I needed as a query, it would be something like this:
select route.id ,
case when duct.type = 1 then st_length(route.geom) * 0.3 --intended here was only the intersecting part not the full length
when duct.type = 2 then then st_length(route.geom) * 0.5
else null end as cost
from route
left join duct  on st_intersects(route.geom,duct.geom)

There are a few problems with the above query:
1 - since the join is based on every line that touches, it multiplies the number of routes for every duct that it touches. A route can intersect multiple ducts. I need the cost to reflect the different duct types it passes through, in a single row
2 - the cost is purely based on the intersection between both layers. So if the route table only intersects for 2 meters with the duct table but is 100 meters long, I'm getting a cost for the full 100 meters instead of the 2 meters.
3 - sometimes the duct layer has multiple geometries overlapping each other, I need to only use the cost for the one with the higher value for type and ignore the others
How to get the proper cost based exclusively on when the lines are completely intersecting and not just touching?


texts: How do I determine the subsidence of land in Manila, Jakarta and Bangkok?

texts: Can I see how these cities have lost elevation using digital elevation models? Is there a digital elevation model which data are updated? Or are there other methods?


texts: Adding Basemap and changing projection - still not aligned

texts:
I am using the basemap in ArcMap10 to create a shapefile of hand digitised polygons based on building outlines.

However, when I change the basemap projection to Brit Nat Grid, it does not line up with the OS building polygons that are also Brit Nat Grid. No matter what course of transform/project/tool I use, I cannot get the two to concur


texts:
In an ESRI shapefile, a polygon can have multiple parts.  Each part is a closed ring.  Rings that wrap clockwise are... um... rings.  Rings that wrap counterclockwise are holes in (some upstream) rings.
Suppose I want a polygon to represent the land portions of: a) a state, b) with a lake in the state, and c) with an island on the lake.
Would it be common to represent this as a single polygon with three parts: a) a clockwise ring for the state, b) a counterclockwise ring for the lake, and c) another clockwise ring for the island?
Similarly, a state with two lakes each with one island would have 5 parts.
If so, do I have to figure out, own my own, which parts are clockwise and which are counterclockwise, and which parts are contained in other parts; or is there some more direct way of figuring this out?
If not, what is the proper representation?
I am processing the shapefile on my own (in C#) to learn how it's done.  I'm not really looking for links to libraries that do everything for me.


texts:
I already have the information for US states and Counties. I have the information in the JSON files. I want the information for ALL US city and Zip Code information.

I can't really navigate the US consensus sites so if anyone has the information, I'd be really thankful.

I found this site where I had to pay for Zip codes for different regions. Other places I have found the info for major cities i.e 10 odd cities only. I am looking for the entire US.


texts: Leaflet map - custom markers from CSV and referring to CSV / GeoJSON marker properties?

texts:
I run a non-profit online resource called the Megalithic Portal and have created a Leaflet map that uses Leaflet.geoCSV to load in points from a CSV file
https://github.com/joker-x/Leaflet.geoCSV

To see ultimately what I'm hoping to create, I've made a map with several thousand megalithic sites, based on data from the Megalithic Portal

Here's the full map
http://www.megalithic.co.uk/leaflet_megalith_map-old.html

So far so good but I'm struggling with the JavaScript to get any further

I've created a simpler test map here for debugging purposes:

http://www.megalithic.co.uk/leaflet_megalith_map.html

I can load in the CSV points, and loop around the marker properties as in the code below (var clave) which outputs the CSV data in the marker popups.

What I'm struggling with is how do I refer to the properties individually , e.g.

popup += 'Colour=' + feature.properties.Colour;

(which is wrong as it shows as 'undefined' in the popup)

I think the CSV loads in as an extension of the GeoJSON format
http://leafletjs.com/examples/geojson.html

The second thing I'd like to do is to display custom markers set by a column in the CSV.

eg in the properties I have Icon: tp6
 Icon: tg60 etc

I would like to create a different marker for each point based on this property -
Do I need to do this as part of the onEachFeature and loop around all the markers to set the property?

I've tried the following (full code below) but it doesn't work so I must have the syntax wrong:

iconUrl: 'images/mapic/' + feature.properties.Icon + '.gif'

var bankias = L.geoCsv(null, {
    onEachFeature: function (feature, layer) {
        var popup = 'SID='+ feature.properties['SID'] +'
';
        popup += 'Colour='+ feature.properties.Colour +'
';
        for (var clave in feature.properties) {
            var title = bankias.getPropertyTitle(clave);
            var csvtext = feature.properties[clave];
            //if (title != 'SID' && title != 'Colour' && title != 'Icon' && title != 'CatID' && title != 'TypeID') {
                popup += ''+title+': '+csvtext+'
';
            //}
        }
        layer.bindPopup(popup);
    },
    pointToLayer: function (feature, latlng) {
        return L.marker(latlng, {
            icon:L.icon({
                //iconUrl: 'images/mapic/' + feature.properties.Icon + '.gif', // WRONG
                iconUrl: 'leaflet/images/marcador-bankia.png',
                shadowUrl: 'leaflet/images/marker-shadow.png',
                iconSize: [25,41],
                shadowSize:   [41, 41],
                shadowAnchor: [13, 20]
            })
        });
    },
    firstLineTitles: true
});


texts: How do I relink or open a project in Mapbox Studio after changing the path?

texts: I changed some folder names and lost the ability to open my projects in Mapbox Studio. Changing the folder names back did not help. How do you open or link a project to Mapbox Studio?


texts: Get Multiple Shortest Path using PgRouting in One Query

texts:
I am trying to get multiple shortest path using A* Algorithm available in PgRouting.

I referred Solution given by UnderDark ( Link of the Stakoverflow Post ) . But it is not giving me proper Output.

MY Query (Simple Query - Without Multiple Path Detection)

SELECT * FROM shortest_path_astar('
            SELECT gid as id,
                     source::integer,
                     target::integer,
                     length(the_geom)::double precision as cost,
                     x1, y1, x2, y2
                    FROM connectivity',
            110, 111, false, false);

Output of this Query:

vetex_id  |  edge_id |  cost
110       | 1        | 0.01235
111       | -1       | 0

Now Query as suggested by Underdark on Given link:

 SELECT
source,
target,
(SELECT SUM(cost) FROM
    (  --- Complete Above Query Here
             )
    ) AS foo
) AS cost
FROM connectivity;

Output of this Query:

source  |  target |  cost
9       | 8       | 0.01235
...     |...      | ...
(More 211 records )

Note: In my test table (connectivity) I have 212 Records

So this query is giving me all the records as output without any path relativity.

I do not know where I am doing wrong. any help / suggestion will be helpful for me.


texts: Converting vector layer to raster with cell size same as other raster using QGIS?

texts:
I am a student and was taught by my professor in ArcMap but he wants us to repeat the process in QGIS. And i dont know how to set cell size in rasterize (vector to raster)in QGIS just like in ArcMap.

How can it be done?


texts: unable to upgrade gdal-bin from 1.9 to 1.10 on ubuntu 12.04

texts:
I am trying to install import files to cartodb installed on local machine. But import fails with these errors like

Error 2001 unable to load data
Error 3009 SHP normalization error

Also the ruby console says the ogr2ogr error. So tried to update gdal-bin 1.10.x. Current version of gdal-bin is 1.9.2. But during the sudo-apt get install gdal-bin, these error appears

gdal-bin : Depends: libgdal1h (>= 1.10.0-1~precise1) but it is not going to be installed

Tried to add gdal 1.10 ppa to system. Also removed older version. But not successful.


texts: Using a selected layer with processing in QGIS 3

texts:
How do I create a Processing script which asks the user for a layer and then conducts an operation on it?

import processing
from qgis.PyQt.QtCore import QCoreApplication
from qgis.core import (QgsProcessingAlgorithm,
                   QgsProcessing,
                   QgsProcessingParameterFeatureSource,
                   QgsProcessingParameterFeatureSink,
                   QgsVectorLayer)

class testAlg(QgsProcessingAlgorithm):
    HEXES = 'HEXES'
    OUTPUT = 'OUTPUT'

    def tr(self, text):
        return QCoreApplication.translate('testalg', text)

    def createInstance(self):
        return type(self)()

    def group(self):
        return self.tr('Test')

    def groupId(self):
        return 'test'

     def __init__(self):
        super().__init__()

    def initAlgorithm(self, config=None):
        self.addParameter(QgsProcessingParameterFeatureSource(self.HEXES,
                                                          self.tr('Input layer'),
                                                          [QgsProcessing.TypeVectorAnyGeometry]))

        self.addParameter(QgsProcessingParameterFeatureSink(self.OUTPUT,
                                                        self.tr('Output'),
                                                        QgsProcessing.TypeVectorPolygon))

    def name(self):
        return 'testalg'

    def displayName(self):
        return self.tr('Test Algorithm')

    def processAlgorithm(self, parameters, context, feedback):
        hex_source1 = self.parameterAsSource(parameters, self.HEXES, context)
        hex_source2 = QgsVectorLayer('/Users/username/temp/hexes.shp', 'hex_source2', 'ogr')

        feedback.pushConsoleInfo(str(type(hex_source1))) # Returns QgsProcessingFeatureSource
        feedback.pushConsoleInfo(str(type(hex_source2))) # Returns QgsVectorLayer

        processing.run('qgis:fieldcalculator',
                   {'INPUT':hex_source1, # Dies here as it can't accept a QgsProcessingFeatureSource
                    'FIELD_NAME':'test',
                    'FIELD_TYPE':0,
                    'FIELD_LENGTH':10,
                    'FIELD_PRECISION':6,
                    'NEW_FIELD':True,
                    'FORMULA':'( "c" / sum(  "c" )) * 100',
                    'OUTPUT':'/Users/username/temp/hexes2.shp'},
                   feedback=feedback,
                   context=context)

        return {}

This code fails at the processing.run command as it can't accept a QgsProcessingFeatureSource (hex_source1) as an input (actual error is Incorrect parameter value for INPUT). Using hex_source2 works fine as that returns a QgsVectorLayer which is apparently fine as an input.

I have also tried QgsProcessingParameterVectorLayer which also seems to give a QgsProcessingFeatureSource result.

Should I be converting a QgsProcessingFeatureSource to a QgsVectorLayer? If so, how? Or is there a simpler way to load a layer for processing that I'm missing?


texts: How to use a polygon as a reverse cookie cutter on features in a Postgresql database

texts:
I want to be able to use a polygon in a database to remove all features and parts of features that falls within the polygon. Not as part of a query transferring the database to another database. I want to remove the features and parts of features inside of the database I am working in. Could someone point me to a example of how to do this with a postgis query?
I want to be able to erase the parts of features that overlap the polygon.


texts: NoData values recognised in raster in ArcGIS Desktop but not following symbology rules?

texts:
I have a GeoTIFF raster in ArcGIS Desktop, which can be downloaded from the following two links:

http://www.rtwilson.com/downloads/MODIS_MOD04_16-06-2006_Subset.tfw
http://www.rtwilson.com/downloads/MODIS_MOD04_16-06-2006_Subset.tif

If you load these into ArcGIS Desktop you will see large uniform areas. Clicking on these areas with the Info tool shows the pixel value to be No Data. However, when you choose to display No Data values as transparent (or as any colour) in the Symbology preferences for the layer, the No Data values do not follow that rule.

Interestingly, the No Data values outside of the main area of the image do follow the rule (set the colour to red to see that), but the ones inside the image don't.

I have really struggled to try and fix this problem - and can't work out how to do it at all.

Does anyone have any ideas?


texts: How to install GDAL for python on windows and use for spatial purpose?

texts:
Possible Duplicate:
  How to install GDAL with Python on windows?

I am trying to copy GDAL osgeo to my python library file to run my python program for gis function. following is code.

import osgeo.ogr

shapefile = osgeo.ogr.Open("tl_2009_us_state.shp")
numLayers = shapefile.GetLayerCount()

print ("Shapefile contains %d layers" % numLayers)
print

for layerNum in range(numLayers):
    layer = shapefile.GetLayer(layerNum)
    spatialRef = layer.GetSpatialRef().ExportToProj4()
    numFeatures = layer.GetFeatureCount()
    print ("Layer %d has spatial reference %s" % (layerNum, spatialRef))
    print ("Layer %d has %d features:" % (layerNum, numFeatures))
    print

for featureNum in range(numFeatures):
    feature = layer.GetFeature(featureNum)
    featureName = feature.GetField("NAME")

    print ("Feature %d has name %s" % (featureNum, featureName))

But u able to find missing .py file as some error. Can any one say step is write or any other steps for installing.


texts: Label expression that parses a multi-sentence label so that every 5 words show up on a new line

texts:
I have long notes that I have to put on the map. My method is to add the notes to the attribute table and display them as labels. I'm having a hard time writing a label expression that parses the sentences so that every 5 words show up on a new line.

This solution I found on the internet only parses the label into two lines, I can't figure out how to make it work for multiple lines:

left([Labels_for] ,instr(50, [Labels_for] ," ")) & vbnewline & right([Labels_for] ,len([Labels_for] )-(instr(20, [Labels_for] ," ")))


texts: Overlay raster above uploaded shpfile but behind background?

texts: I have a .shp file of the world oceans (hydropolys.shp) and create a map with it in QGIS, so the background of the image is the land. Does anyone have any suggestions on how I can overlay my raster output on top of the oceans, but not overlain on the land? So I essentially want my image overlain on top of the uploaded .shp file but underneath the background.


texts:
I downloaded QGIS, GeoServer, OpenLayers, Ext, GeoExt, qgis-web-client and osgeo, and I have map shapefiles. I would like create a web site like the Uster site.

I need help, can you provide tutorial links?


texts: Reverse Engineering Movable Labels in ArcMap?

texts:
I have an .mxd where the labels for most features are fully editable and movable. However, they have NOT been converted to notation. They are stored in shapefiles that have no symbology but also have "Label Features" checked off. When a shapefile is checked off, the labels turn off as well. We can't figure out how this is possible but we'd like to be able to recreate this label style since it works really well.

How can I recreate this?

I am not new to GIS, please trust that I understand what a shapefile is.

I am using ArcMap 10.5.1.


texts: QGIS Desktop WFS GetCapabilities fails with "Unexpected end of file"

texts:
When I am trying to add WFS layers with QGIS Desktop it works: GetCapabilities request is sent and I can add any layer I want. However when another person using other machine and QGIS Desktop tries to do the same he gets unexpected end of file. If he tries to GetCapabilities via browser it works fine. From the image he sent me I see that he has selected same WFS version as I am.
Update: installed QGIS network logger and this is what I see on other person's QGIS

Error text says: "server is unsuitable for that kind of operation"
While on my QGIS it says


texts: Attribute table Does Not Show Any Lat Long Field in ArcGIS for Desktop?

texts:
I had import shapefiles data which is a point features from GPS receiver. The points are well on pointed on the ArcGIS display, but when I open each attributes table, there is no field show the latitude and longitude of the points while the points were shown at display in a correct location abse on my georeferenced map. Why? The attributes only show the ObjectID, Shape, and Comments.

I am using ArcGIS 10.3 for Desktop.


texts: How to show an array of 2d points/polygons in Polymaps (Cartesian coordinate system)

texts:
I need to display these 2d points/polygons (Cartesian coordinate system) on Polymaps with GeoJson layer.
I do not use the other layers.

Do I need to convert xy coordinates to the lat lon? If so, how?
How to use Cartesian coordinate system in Polymaps?


texts: Help needed with incorporating SQL query into ArcGIS python script

texts:
I'm new to python and have been trying to learn how to use it within ArcGIS. I have primarily been using the free Penn State courseware available here - http://open.ems.psu.edu/courseware

I have understood most of the content, however, I am struggling to understand how the instructors have incorporated the SQL query into this script:

(The comments and course notes related to this scrip can be viewed here - https://www.e-education.psu.edu/geog485/node/136)

import arcpy

fc = arcpy.GetParameterAsText(0)
affectedField = arcpy.GetParameterAsText(1)
oldValue = arcpy.GetParameterAsText(2)
newValue = arcpy.GetParameterAsText(3)

queryString = '"' + affectedField + '" = '+ "'" + oldValue + "'"

rows = arcpy.UpdateCursor(fc, queryString)
row = rows.next()

while row:
row.setValue(affectedField, newValue)
rows.updateRow(row)
row = rows.next()

del row, rows

The above script is used to update the attribute table in ArcGIS. In line 8 an SQL query is used to identify a particular name that I want to change (I used a work related attribute table here and wanted to change the name of a city from 'Aberdeen' to 'Newcastle'). However, having very little experience with SQL I am not sure what the instructors mean or how to apply it to my own task

Could anyone please help me to understand what the SQL query in line 8 means and how to use it to complete my task?

Thanks,

Chris


texts:
As far as I can tell, the option to "Create point on surface for each part" in the QGIS Centroid tool is not returning the correct geometry when activated in QGIS (version 3.8.1-Zanzibar/GDAL v 2.2.3 on Debian)

If we take a simple "banana-shaped" polygon:

and run the tool (Vector --> Geometry Tools --> Centroids) as follows:

we get the following output, which represents the centroid, but not the point on surface.

We can cross check this output using ogr2ogr as follows:

ogr2ogr -dialect SQLite -sql "SELECT ST_Centroid(geometry) FROM tmp" tmp_centroid.shp tmp.shp

which returns the identical centroid:

whereas the following command forces the point on the surface of the polygon:

ogr2ogr -dialect SQLite -sql "SELECT ST_PointOnSurface(geometry) FROM tmp" tmp_ptonsurf.shp tmp.shp

is this a bug or is there something else that should be run/enabled in order to get the point on surface using the QGIS tool?


texts: Extract very accurate lat lon coordinates from a raster of a given resolution

texts:
I am a beginner in raster object, and I have a particular problem.

Suppose that I have a raster which is focus on Mexico.

Besides, I have created a SpatialPoint Object of coordinates that are accurate of different locations that stand for let's say houses.

I would like to extract the houses within the positive values of the raster. How should I do ?

class       : RasterLayer
dimensions  : 72, 122, 8784  (nrow, ncol, ncell)
resolution  : 0.25, 0.25  (x, y)
extent      : -117.125, -86.625, 14.625, 32.625  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
data source : in memory
names       : layer
values      : 54.00137, 126.3991  (min, max)

Formal class 'SpatialPoints' [package "sp"] with 3 slots
  ..@ coords     : num [1:459382, 1:2] -103.4 -98.2 -114.8 -110.9 -110 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:459382] "1026" "1027" "1028" "1029" ...
  .. .. ..$ : chr [1:2] "Longitude" "Latitude"
  ..@ bbox       : num [1:2, 1:2] -117.1 14.7 -86.7 32.7
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "Longitude" "Latitude"
  .. .. ..$ : chr [1:2] "min" "max"
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
  .. .. ..@ projargs: chr "+init=epsg:4326 +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"


texts:
I have a shapefile with coordinates listed in Easting and Northing in the attribute table. Is it possible to create two new fields called lat and long and populate them with values corresponding to easting and northing using ogr2ogr?

The idea is to try and get this to work without having to use a specialist GIS software and of course automation.


texts: Unable to select interpolation attributes from my shapefile?

texts: My point shapefile has four fields, three of which are numerical: longitude, latitude, and homicides. However, the homicides field does not appear for me to select it as an interpolation attribute (the field name has no spaces).


texts: Extracting *.tif agricultural area from shapefile?

texts:
I am attempting to extract from an image using Esri shapefile that is in the zip file RR_KZKD2012

In ArcGis the spatial reference is Projected Coordinate System: WGS_1984_UTM_Zone_43N.

Then am using arcgis extract by mask to clip my Mg1.tif using  RR_KZKD2012

The output tif file is not extracted properly. I do not understand why see example:

am using different scenarios but not work..

I am sharing the link below...please have a look and see if you can see where it is going wrong.

https://drive.google.com/file/d/0By0CgkNTboOXMDZmal9xMFV2Rjg/view

I have many more to do.


texts: Viewing entire Python script of process ( input settings and environment settings) in ArcMap

texts:
When you "copy as python snippet" from the results of a process, the script contains all of the input settings but does not contain the environment settings. Why is this?
Is it possible to also view the code of the environment settings of a process that you run? meaning the entire script of the process.


texts: Setting "Always run in foreground" within ArcPy code?

texts:
Is it possible to set the "Always run in foreground" option within ArcPy code? I haven't found anything on this yet?

It would be nice to be able to specify within a script, instead of through the Geoprocessing Options window.


texts: Coordinate problem with GPCC netcdf data and QGIS

texts:
I downloaded Global Precipitation Climatology Centre (GPCC) monthly precipitation data in netCDF format from here: https://www.esrl.noaa.gov/psd/data/gridded/data.gpcc.html

According to the site, they should be in WGS-84 projection:

However, when I open the file in QGIS, the latitude values all show up as negative.

How do I fix this using QGIS?


texts:
In QGIS, I want to color a polygon with a hatching in gradient mode. It should look like in the picture.

In this special case, it was realized with a green-to-yellow fill in the background, combined with a green hatching (line pattern fill) in the foreground.

But this is not possible when a gradient-hatching should be applied over different objects with different styles/colors. Any idea how this could be realized?

Here is the situation I want to use with gradient hatching:

The dark-blue hatching should be in gradient mode, but the background is in different colors like red, light-blue, white, another hatching in yellow etc. So the solution from the first example is not possible here. I want to have the dark-blue hatching with strong colors in the middle of the polygon, fading out towards the margins of the polygon-area.


texts: Join points to a line with minimal distance

texts:
I have a point cloud (x,y) in a random order.

Say sorted([(1, 1), (2, 2), (0, 0), (4, 4), (3, 3)]). Is there an easy way to connect them to a simple line like ((0,0),(1,1),(2,2),(3,3),(4,4)) doing this "smart" reordering (shortest path)?


texts: What platform to use for building mapping application using existing ArcGIS Rest Services?

texts:
I have clients that have an existing ArcGIS Server.  Map layers are published as rest services and I am going to build a mapping application that displays their datasets.  Functionality should include users' ability to click on layers and retrieve attribute info.

I've been using Google's JS API v3 for most of my projects.  I am also using leaflet and CartoDB for other projects.  I am trying to stay away from ESRI's JS API for server as it seems a little cumbersome.  I'd prefer to use something like Google Maps or Leaflet for this project.

Any suggestions as to what platform I should use to construct this application?

Will it make my life easier to just use ESRI's JS API for server?
Any major downfalls to loading the REST services in Leaflet or Google ?


texts: Removing space between objects of a shapefile

texts:
I'm trying to visualise the results french presidential elections of 2017 by departements. I found a shapefile of France with departements here

When I try to export this file in TopoJSON format from mapshaper.org it has a lot of space between mainland France and some of the departements.

Is there a way I could make it look like this

If someone has a shapefile already in this format I'd be very thankful if they could share it with me.


texts: Convert N-sided polygon to nearest fitting rectangle?

texts:
I'm working with OSM data and many of my objects are simple 4-polygons (parallelograms), which can be expressed by 4 unique corner points. But the full OSM data has the shape represented by more than 4 points, with multiple points per line segment (see image). This means that the shapes are expressed with redundant points, 16+1 in the case of this polygon:

For an arbitrary polygon, how can I find either

1 the minimum bounding (rotated) rectangle?

2 or the best fitting rectangle?

Example GeoJSON:

{"type":"FeatureCollection","features":[{"type":"Feature","properties":{"osm_id":"1269601","type":"multipolygon","leisure":"pitch","sport":"soccer"},"geometry":{"type":"Polygon","coordinates":[[[6.6131123,46.5124914],[6.6129421,46.5125385],[6.6127998,46.5125783],[6.6126291,46.512626],[6.6125308,46.5124593],[6.6127016,46.5124121],[6.6128452,46.5123724],[6.6130153,46.5123244],[6.6131123,46.5124914]]]}}]}

http://bl.ocks.org/d/465e77bb03acfd7085b2f9a1d2dd08d5

What I've tried:

ogr2ogr's simplify flag and the simplify() function in Shapely, which doesn't force my polygons to be 4 sided, and is also smoothening things I don't want (like corners). It's my understanding that simplify won't work if the points aren't near enough to each other.
Shapely's buffer() function
Bounding Box, which gives me the smallest non-rotated rectangle that fits my polygon.


texts: Transforming Lat/long data to UTM --2 different UTMS

texts:
I am trying to transform Lat/Long data from Italy (UTM 32 and 33) to UTM. I am struggling a bit with this.

The code I am trying to use is one I found in a previous post, I tried to modify it to solve the problem of several UTMs but does not seem to work.

  df <- data.frame(x = Data_group_malesCLEAN$Latitude, y =
  Data_group_malesCLEAN$Longitude, ID = Data_group_malesCLEAN$Group.1)

  long2UTM <- function(long) {   ## Function to get the UTM zone for a
  given longitude   (floor((long + 180)/6) %% 60) + 1 }

  LongLatToUTM <- function(df){   coordinates(df) <- ~x+y
  proj4string(df) <- CRS("+proj=longlat +datum=WGS84")

  zone <- long2UTM(df$x)

  ## Change CRS of the spatial data frame and convert to data frame

  res1 <- spTransform(df, CRS(paste0("+proj=utm +zone=", zone,
  "+datum=WGS84")))   return(as.data.frame(res1)) }

  res1

I also tried to do it with Qgis... but I am in the same problem... It is only possible for one UTM,... any ideas?
Best


texts: Creating density map of roads for each district using ArcGIS Desktop?

texts:
I have to illustrate a map with the road density of each district.

I have the city layer which comes with the bordered districts and another layer with all the roads.

How can I create a density map where road density is displayed for each district independently (within each district's/polygon border)?

I have ArcGIS 10.3.1.


texts: How to update/delete data from a datagridview?

texts:
I wanna be able to update/delete data from a datagridview, and when I do that of course the data in my PostgreSQL database should be edited too.
So I've thought about giving the user the ability to select a layer, choose an attribute, choose a logical operator (=, < or >) and then type the criteria he wants for his search.
So for example he could make a search like this one :

  Layer : Activity

  Attributes : Id

  Logical Operator : =

  Criteria : 15

And when he clicks on Search, the result should display on the datagridview
The thing here is, if I want to update a row or a number of rows, how can I be able to do that? Same thing for deleting?
Is it possible to edit data in the datagridview and then I don't know do like an update button and save the changes accordingly but how am I gonna do that ?

EDIT

So now I have this idea : I have two forms :

  Form 1 : Activity contains :

    the datagridview linked to PostgreSQL database (and the DB is linked to ArcCatalog)
    a button named Update on which I can click and show the Form2 where I can edit my data (I pass the values from the datagrid view to my textboxes fields in Form2)
    a button named Delete on which I can click and delete the IFeature.

  Form 2 : UpdateActivity contains :

    the fields of the selected Ifeature, with the values already passed
    A Validate button on which I can click and update the data after I finish editing, also the form should close and then I see the updated data in the datagridview in form 1

      P.S : How can I call the datagridview in the second form?

Here's some code :
On Activity Form:

// This is the Update button that pop up the Form2, and it's where I pass the values
private void btnUpdate_Click(object sender, EventArgs e)
    {
        updf = new UpdateForm();

        //Loading data from datagridview to UpdateForm
        updf.txtIdActivity.Text = this.TableAttributaire.CurrentRow.Cells[0].Value.ToString();
        updf.txtNameActivity.Text = this.TableAttributaire.CurrentRow.Cells[1].Value.ToString();
        updf.txtAddressActivity.Text = this.TableAttributaire.CurrentRow.Cells[2].Value.ToString();
        //Update Data
        updf.ShowDialog();

     }

On Activity Form:

//Clicking on this button should I create another form like I did with the update ?
//Or just delete straight from the datagridview with a messagebox?
private void btnDelete_Click(object sender, EventArgs e)
    {
        //Haven't done it yet?
    }

On UpdateActivity Form:

//This is the validate click button on which I can click and save and update the data
private void Validate_Click(object sender, EventArgs e)
    {
        IFeatureClass pFeatClass;
        IQueryFilter pQF = new QueryFilter();
        // assume pFeatClass is an IFeatureClass object
        // and pQF is a queryfilter that allows to select the features in question.
        IFeatureCursor pCur = pFeatClass.Update(pQF, true);
        IFeature pFt = pCur.NextFeature();
        pFt.set_Value(indx, pValue); // I don't know what to do here ??
        pCur.UpdateFeature(pFt);
//But where is the datagridview in all this? How can I call it in this second form and update tha data ?
//Also like I said before I should be able to close this form after I validate and refresh the datagridview in Activity Form
    }

Updating/Deleting from the map:

  How can I select an object on my map in Activity Form, and display the Update Form accordingly and be able to edit the data?
  Same thing for deleting?


texts: Is it possible to save frequently used expressions in QGIS?

texts: Is there any way to save frequently used expressions in the selection or field calculator dialogs? Obviously I can cut and paste them into a file where I might save them, but I mean within the context of QGIS like you can do in ArcMap.


texts:
I have been using QGIS 1.8 Lisboa for some time now and would like to be using OSM for my base maps. I am aware this  a very old version now and doesn't seem to support OSM anymore.
I have upgraded to QGIS 3.16 and the OSM plugin works fine.
However in this version the raster Clip by extent on canvas tool window has a different layout which is very awkward to put coordinates into (copied and pasted from Google Maps) compared to QGIS 1.8 (has separate boxes for lat and long)
Can I recreate this interface in 3.16 (or newer)?
The alternative is to be able to load Google Maps into QGIS but am I right in thinking this no longer works?


texts: Using Python, how to get SQL Server queried points in the map?

texts:
Since days I'm struggling with this problem, did a lot of searching - also in this forum - but did not find the answer I need.

I have a layer 'Project' in QGIS that is connected to a spatial table in SQL Server with the same name. The table contains more than 30000 projects (done since 2009) with their position.
If we want to do a new project at a new address with let's say coordinates (lon=8.002, lat=35.453) in WGS84, we want to know which projects nearby we have done already. So, select the projects around this address within a distance of say 0.1°.
I prefer to select these projects using SQL in the database itself, I think that's the quickest way.

In Python, I make a connection with the SQL Server database and run the spatial query that selects the projects.

The code up to now that runs without error messages in Python is:

pr = QgsProject.instance().mapLayersByName('Project')[0]
self.iface.setActiveLayer(pr)
server = 'MyServer'
database = 'MyDatabase'
conn = pyodbc.connect("DRIVER={ODBC Driver 17 for SQL Server};SERVER="+server+";DATABASE="+database+";Trusted_Connection=yes")
cursor = conn.cursor()
cursor.execute("SELECT * FROM Project WHERE Position.STDistance(GEOGRAPHY::Point(1.002, 51.453, 4326)) <= 0.1")
row=cursor.fetchall

I'm looking for the correct Python code to show the points in layer 'Project' that meet the query result.


texts: Using interactivity from multiple layers produced in TileMill in a single map

texts: I have a MapBox.js map with several layers here .  I would like to be able to use TileMill's teaser and full interactivity on a few of the layers, but I have only been able to accomplish this for a single layer.  How can I have interactivty for multiple layers in this map?  The code for the map is here. I saw an example from MapBox.js but I don't know exactly how to implement it into my HTML code.


texts: pgsql2shp error when query is introduced

texts:
I'm trying to export a table from PostGIS using pgsql2shp. The basic query to export the entire table works fine using the following format:

"C:\Program Files\PostgreSQL\9.1\bin\pgsql2shp.exe" -f "D:\path\shapefile.shp" -h
localhost -p 5432 -u postgres -P xxxx mapdata planning.dm_planning-apps

However when I add a query (which works fine in PGAdmin III) the command fails, due to the hyphen in the table name.

Working query in PGAdmin:

SELECT * FROM planning."dm_planning-apps" WHERE "dm_planning-apps"."Modified_Reference"
ILIKE '%np%'

Included in command:

"C:\Program Files\PostgreSQL\9.1\bin\pgsql2shp.exe" -f "D:\path\shapefile.shp"  -h
localhost -p 5432 -u postgres -P xxxx mapdata "SELECT * FROM planning."dm_planning-apps"
WHERE "dm_planning-apps"."Modified_Reference" ILIKE '%np%'"

Produces a "syntax error at or near "-" and highlights the hyphen in "...FROM planning.dm_planning-apps WHERE...". The double quotes obviously aren't doing the same job they were within PDAdmin where the query worked fine.

Can anyone suggest a way to make this query work without renaming tables?


texts:
I have a layer made with ArcGISDynamicMapServiceLayer like this:

return esriLoader.loadModules([
  'esri/layers/ArcGISDynamicMapServiceLayer',
  'esri/layers/ImageParameters',
  'esri/InfoTemplate'
]).then(([ArcGISDynamicMapServiceLayer, ImageParameters, InfoTemplate]) => {

  const infoTemplate = new InfoTemplate();
  infoTemplate.setContent(this.handleContentChange);
  infoTemplate.setTitle(this.props.intl.formatMessage({ id: 'layers.popup.title' }));
  const imageParams = new ImageParameters();
  imageParams.layerIds = [layerCode];
  imageParams.layerOption = ImageParameters.LAYER_OPTION_SHOW;
  imageParams.transparent = true;

  const layer = new ArcGISDynamicMapServiceLayer(layerUrl, {
    id: layerId,
    opacity: 0.5,
    imageParameters: imageParams,
    infoTemplates: { [layerCode]: { infoTemplate: infoTemplate } }
  });

  esriMap.addLayer(layer)

I want to set the extent of the layer that I'm rendering to fit what it's rendered. The problem that I see is the featureExtent and initialExtent comes from the service, and it's the same for each service. But some layers are in different parts of the map, some are big, some are very small, so I need a way to zoom to what it's showed in the map. How can I set the extent on those service to fit it?

For example, with: esriMap.setExtent(layer.initialExtent, true) doesn't work, because in the services

I have this extents:


texts: Integrate AT_GIS_GRID (NTv2 Transformation) in QGIS 2.8 & 2.10

texts:
Compared to previous versions in QGIS 2.8/2.10 no *.gsb-Files seem to be referenced in srs.db (for NTv2 Transformations) - did the system of using NTv2 change after QGIS 2.2 ? Wanted to integrate AT_GIS_GRID in 2.10 like described here: http://isticktoit.net/?p=371
QGIS 2.8/2.10 (Win64) do not allow to edit srs.db and p1 not able to save character-values like AT_GIS_GRID.gsb


texts:
I am trying to transform a GeoJSON object of a certain transformation (EPSG:28992) into my current map coordinate system (EPSG:4326). I have been trying to do a couple of things but none seem to work. I am currently trying to use proj4.js. A  minimal example of my script can be found here:

http://jsfiddle.net/TimLucas/k69ujyev/5/

I am trying to find a method how I can transform my GeoJSON (Layer2) layer into the projection used in my map, as defined per the styles it should show up as green points. It somehow seems to have worked for the other layer (Layer1, red points). I tried to follow the various examples on openlayers.org (which suggest to use proj4.js), but this does not work for some reason. Can you please help me with what I am doing wrong?

Edit: I forgot to mention, the other projection is one that OL supports by default, so the problem here is using external projection systems...

Edit2: Found out that ArcMAP does not export this coord system very well. Changed it to the correct coords by hand but still does not show up, any suggestions? And does anyone have any suggestions on the easiest way to convert ESRI JSON to GeoJSON?


texts: how to display map center coordinate using OpenLayers

texts:
As far as I know OpenLayers can display our mouse coordinate position. However, I do not know how to display map center coordinate?
Can anyone help me?


texts: How to add multiple stop points for routing using ArcGIS in Android

texts:
I have a project where in I have a source and a destination, but in between there can be any number of stop points. I need to display a route from source to destination covering all the stop points.

I started finding route between two points and is working fine. But the moment I add a extra point in between it gives a error saying

  Location "Location 3" in "Stops" is unlocated. Location "Location 2" in "Stops" is on a soft-restricted network element. Cannot preserve last location "Location 3" in "Stops" because it is invalid. "Stops" does not contain valid input for any route.

I tried googling for this error but couldn't find any results.

I thought location 3 is not a valid location but when I try to find a route between source and location 3 it works fine. I am not sure what I am doing wrong. Below I have attached function which gives route result.

private void QueryDirections(final ArrayList<Point> stopPoints) {

dialog = ProgressDialog.show(RoutingSample.this, "Routing Sample",
        "Calculating route...", true);
Thread t = new Thread() {
    @Override
    public void run() {
        try {
            // Start building up routing parameters
            NAFeaturesAsFeature stops = new NAFeaturesAsFeature();
            ArrayList<StopGraphic> stopGraphics = new ArrayList<StopGraphic>();
            for(Point pt : stopPoints){
                stopGraphics.add(new StopGraphic(pt));
            }
            Graphic[] stopGraphicsPrim = new Graphic[stopGraphics.size()];
            stopGraphicsPrim = stopGraphics.toArray(stopGraphicsPrim);
            stops.setSpatialReference(wm);
            stops.setFeatures(stopGraphicsPrim);
            stops.setCompressedRequest(true);
            RouteParameters rp = mRouteTask
                    .retrieveDefaultRouteTaskParameters();
            rp.setStops(stops);
            rp.setReturnDirections(true);
            rp.setFindBestSequence(true);
            rp.setPreserveFirstStop(true);
            rp.setOutSpatialReference(wm);

            // I get error here. When routeTast.solve is called. Please help

            mResults = mRouteTask.solve(rp);
            mHandler.post(mUpdateResults);
        } catch (Exception e) {
            mException = e;
            e.printStackTrace();
            dialog.dismiss();
            mHandler.post(new Runnable() {
                @Override
                public void run() {
                    Toast.makeText(getApplicationContext(),"directions cannot be found ",Toast.LENGTH_SHORT).show();
                }
            });
        }
    }
};
// Start the operation
t.start();

}

I was referring to a similar post Finding shortest route between multiple stops using ArcGIS for Android? but I wasn't able to get any solution from this.


texts: Using Select By Attribute on shapefile in ArcPy?

texts:
I am trying to take an existing shapefile and narrow down the attributes using the 'select by attribute'  The problem I am trying to solve is:  Find parcels that have an acreage >= 100.  The code I have is not narrowing down the search or showing up in arcmap.  Below is the code I have so far.

import arcpy

#Get Parcel layer
countyParcels = "C:\\Users\\Krystle\\Desktop\\WCGIS\\Geog485         \\final\\TarrantCoParcels\Parcels2013.shp"
nameField = "TAXPIN"

try:
    #Make a feature layer of parcels2013
    arcpy.MakeFeatureLayer_management(countyParcels, "AllParcelsLayer")

except:
    print "Could not create feature layer"

try:
    #Create a seach cursor
    rows = arcpy.SearchCursor("AllParcelsLayer")
    #Call SearchCursor.next() to read the first row
    rows = rows.next()

    #Start a loop that will exit when there are no more rows available
    while row:
        arcpy.SelectLayerByAttribute_management("AllParcelsLayer", '"CALCULATED" >=100.00')
        #Print vaules in the current row
        print row.getValue(nameField)
        #Call SearchCursor.next() to move to the next row
        row = rows.next()

    #Clean up cursor and row objects
    del row
    del rows
    #Write the selected features to a new featureclass
    arcpy.CopyFeatures_management("AllParcelsLayer", "parcels_100plus")

except:
    print arcpy.GetMessages()


texts: How can I combine USGS GeoTIFF quadrangles into a single GeoTIFF using `rasterio`?

texts:
There are several similar questions about this but most as about QGIS and are non-specific about the data involved. I have a specific data source/format and I need to use Python, and preferably rasterio.merge to do the heavy lifting.
The files in question are USGS 2020 quadrangles, downloaded from https://ngmdb.usgs.gov/topoview/viewer/#4/40.01/-100.06 in GeoTiff format
I need to combine 4 of them (with a common corner) into a larger, seamless GeoTiff. I have succeeded at this in the past using rasterio.merge similar to the following
import rasterio
from rasterio.merge import merge
from rasterio.plot import show

quadrangles = [rasterio.open(f, 'r') for f in files]
combined, new_transform = merge(quadrangles)
show(combined, cmap='terrain')

But in this case it results in one image (complete with legend, scale, header, footer) stacked on top of another. When I try this with four images, two of them seem to get completely lost.


texts: Displaying single band reflectance from Sentinel in GEE

texts:
I have been able to display NDVI value in my project but I need some single certain band reflectance value for my data. I have tried in showing blue band (B2) but when I check the inspector tab, the number represent DN and not reflectance.
How do I display single band reflectance value in earth engine? Been searching and did not find the answer yet.
Below is my script:
var sentinelns = ee.ImageCollection('COPERNICUS/S2_SR')
                   .filterDate('2019-01-01', '2019-12-30')
                .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 95))
                .map(maskS2clouds)
                .median()
                .clip(geometry);

function maskS2clouds(image) {
  var qa = image.select('QA60');
  var cloudBitMask = 1 << 10;
  var cirrusBitMask = 1 << 11;
  var mask = qa.bitwiseAnd(cloudBitMask).eq(0)
      .and(qa.bitwiseAnd(cirrusBitMask).eq(0));
  return image.updateMask(mask).divide(1);
}

var SingleBandBlue = sentinelns.select(['B2']);
var BlueParam = { min: 0, max: 1000,};

Map.addLayer(SingleBandBlue, BlueParam, 'Sentinel B2');


texts: Postcode mapping (QGIS): how to transform OS csv data to an polygon map (.shp)?

texts:
I am trying to obtain a full England postcode map (polygons) out of free (and reliable) data out there.

The English National Office of Statistics does share a huge list of postcodes and some information I am not acquainted with and I was wondering if there was a way to use QGIS and this info to create a workable map to plot some points?

  This dataset contains information from the Ordnance Survey, consisting
  of a list of all UK postcodes, together with their locations as
  Easting/Northing and as latitude/longitude. Each postcode links to the
  local authority district and LSOA that contains the postcode
  'centroid'.

As per @Joseph advice, I downloaded an OS boundary map, then uploaded it on QGIS. All good.
Then I downloaded the OS Code-Point® Open, which is a .csv file. I uploaded one of the csv "ab" on QGIS (where X = Northing and Y = Easting; CRS = 27700).  I could see some points (see below).

Next I executed the command "Join attributes by location" (see below), but I didn't get the expected results (a map with the postcodes delimitation) - again, see below.


texts:
I have two instances of Geoserver running. One on a local machine (version 2.14.0) on Oracle 1.8.0_161. The other on a server (2.13.2) with Oracle  1.8.0_181.

On the server version when I look at Tile Caching in the layer it is setting the Grid subset bounds dynamically. On the local version when I fist set up the layer it is also dynamic but as soon as I save the settings it switches to fixed, as screengrab below. I want it to be dynamic, how do I configure this?

I'm trying to debug an error which I think is caused by the Grid set not being dynamic


texts: Non-zero error with ImportFromEPSG() in GDAL/OSR/Python

texts:
I'm currently getting an error when testing ImportFromEPSG() on my GDAL.

import osr
testSR = osr.SpatialReference()
testSR.ImportFromEPSG(4326)
#this returns value of 6
print testSR.ExportToPrettyWkt()
#this returns blank/empty string

I have GDAL_DATA in os.evirons set to C:\Program Files\GDAL\gdal-data, which is where gcs.csv is located. The proj.dll file is in C:\Program Files\GDAL. Running 'GDAL_DATA' in os.environs returns TRUE.

Any suggestions?

I consulted this past link but it did not lead to a solution for me.


texts:
I need to have the tiles of a small city as PNGs. The style should be the same as in https://www.mapbox.com/editor/? Can I export the map there for a certain bounding box to tiles (or mbtiles)? I have setup a project but export to mbtiles or tiles does not seem to be available.

Or is it (as a workaround) possible to export just the mapbox default style in order to use this style with TileMill?

I have seen the OSM Bright style, but I am not sure: Is this the same style as the default mapbox style?


texts: Labeling using different expression if attribute is null using ArcGIS for Desktop?

texts:
I am using an expression to label features by concatenating two attributes in ArcMap. It works like this:

<attribute_1>+" "+<attribute_2>

Unfortunately when either attribute is Null the whole label is left blank. How do I change my expression so that if one attribute is Null, the other appears on its own?


texts: QGIS Raster subtraction giving 0.005 value off from correct value

texts:
I have two synthetic DEMs I have created to test a method of detecting changes.  I am using the Raster Calculator to subtract the first DEM from the second DEM.  The formula is simple, just
"first_DEM - second_DEM".
I have set my Z value to 100.0000 for 90% of the DEMs.  The other 10% is 50.0000 for first_DEM and 40.0000 for second_DEM.  The results should be 0 over 90% of the DEM and 10.000 over 10%, which should show up when I calculate the histogram of pixel value vs. frequency in the resulting raster properties.  Instead, I get 90% at 0.005 and 10% at 9.995.  The two DEMs are exactly the same size at the same geographic location with the only differences in values being the 10.000 change over 10%.

Why is the calculator adding 0.005 to the majority of my pixel values, and subtracting 0.005 from the other 10%, when I'm only asking it to subtract the two DEMs?


texts: Installing ArcGIS Pro and ArcMap at same time

texts:
I want to install ArcGIS Pro and ArcMap at the same time.

Are there any problems while I use them in one laptop?


texts: Show multiple points and lines in Garmin GPS Map

texts:
I have two shapefiles that contain (1) 30 point features and (2) an extensive road network (line features) that connects them.

I would like to show all data at once in the mapview of my Monterra Device to be able to follow the road network and reach the closest waypoints in a given region. I followed all export steps as e.g. indicated here:

http://www.gpszone.ca/blog/2014/03/convert-kml-to-garmin-tutorial/

I successfully imported the waypoints and if I go to the waypoint manager I can see the 30 waypoints one by one. Than I'm able to open a given waypoint by clicking the little map symbol in the edit view.

However for my purpose this is unsatisfying because I want to see all features at once so that I am always able to see the closest waypoints next to me (and later also see the road network to get there). Garmin Software seemed to me quit rudimentary but there should be a way to do that right?


texts: Aggregate data from several polygons into one in QGIS

texts:
I am using QGIS 3.10.1, and I am trying to aggregate Census data (polygons) according to another set of polygons (neighbourhood areas that I have drawn). While I am able to clip the Census data to each neighbourhood, the census data polygon boundaries still show within the neighbourhood polygons, and is still separate data.
Does anyone know how I can take all the census polygon data within a neighbourhood boundary and combine it in to one dataset showing averages / sum etc?


texts: Dividing area into polygons from road network using ArcMap

texts:
I'm trying to divide a closed boundary into sub zones. These sub zones will be of shapes based on a road network. I have tried split tool but it is not working for this.
What tool should I use?
I am using ArcMap 10.8.


texts: I have added 'montserrat' into my OSGeo4W64 folder but I'm not sure how/where to add it so I can use it in my project.


texts: Manipulating the shapefile and made a zoom effect of city districts in QGIS

texts:
I need TopoJSON file for use in Power BI with a map where the city districts are zoomed/scaled and placed outside the country in the map.
Here is a map I have made in ArcMap Desktop  https://webguru.no/wp-content/uploads/2021/02/zoomOslobydelar.png by using the  "Scale" tool in ArcGIS.
Is it possible to manipulate the shapefile in QGIS and copy city districts to the right side of the map like the image (url) above?


texts:
Here's a quick summary about what I'm trying to do: I have 3 tables in Postgres, 'a' and 'b', each have a Polygon column, and 'c' has a Point column. What I'm trying to do here is to get the geometries intersections between 'a', 'b' and 'c', and to display such geometries on an OpenLayers vector layer.

I already know how to display any kind of geometry from a String in OpenLayers, but I'm having troubles with the PostGIS' ST_Intersection function, I'm doing this:

SELECT ST_Intersection(a.geom, b.geom) as inter from a, b;

where a.geom and b.geom are both the geometry columns, and I get this error message:

NOTICE:  TopologyException: found non-noded intersection between 515172 2.14408e+06, 497067 2.13373e+06 and 501321 2.13546e+06, 471202 2.14843e+06 500621 2.13576e+06
ERROR:  GEOS Intersection() threw an error!

Also I tried to express the resultant geometry as text using ST_AsText like this:

SELECT ST_AsText(ST_Intersection(a.geom, b.geom)) as inter from a, b;

but it send me this error message:

HINT: No function matches the given name and argument types. You might need to add explicit type casts.

I don't know what I'm doing wrong, I just want to get the Polygons' WKT to display it on OpenLayers, here's how I display a geometry from a WKT:

                    var in_options = {
                        'internalProjection': new OpenLayers.Projection("EPSG:4326"),
                        'externalProjection': new OpenLayers.Projection("EPSG:4326")
                    };

                    var fea= new OpenLayers.Format.WKT(in_options).read(data); //data is the string with the WKT
                    vectorLayer.addFeatures([fea]); //this piece of code works great
                    map.zoomToExtent(bounds);

UPDATE: I tried the next:

SELECT ST_Intersection(a.geom, b.geom) as intersect_ab FROM a INNER JOIN b ON
ST_Intersection(a,b) WHERE ST_Overlaps(a.geom, b.geom)
AND ST_isvalid(a.geom)='t' AND ST_isvalid(b.geom)='t';

but I get the next error message:

ERROR: Function st_intersection(a,b) does not exist.
HINT: No function matches the given name and argument types. You might need to add explicit type casts.

I added the isvalid to verify only valid polygons are being evaluated, but it's telling the error is in the ST_Intersection(a,b), both a, b and c have the same SRID so I'm really confused, sorry if I'm asking too much, but I'm quite new with PostGIS so I hope I'm not bothering you a lot. Thanks.


texts:
I'm using leaflet 1.0.3 with latest plugin leaflet-locatecontrol. When I start to locate myself, it takes arround 5 seconds. My position is being watched constantly and even if I move right after the last position update it takes 5 seconds again to reposition my map. I've tried using parameters such as timeout or maximumAge, but these have no effect. I've searched for any numbers in leaflet.js to replace these with lower values, but that didn't help either.

A position update interval of 500ms would be awesome for my application. From what I've read, Leaflet uses Browser API to process geolocations. But I managed to locate myself in OpenLayers with my desired update interval. Is there any way to hack Leaflet's update Interval?

var lc = L.control.locate({
        position: 'topright',
        strings: {
            enableHighAccuracy: true,
            setView: true,              //no effect
            keepCurrentZoomLevel: true, //no effect
            watch:false,                //no effect
            timeout: 1000,              //no effect
            maximumAge: 1000            //no effect

        }
    }).addTo(map)


texts:
I need to work in WGS 84 so I am trying to project my Shapefile into WGS 84.

Despite changing the default CRS to WGS 84, when I mouse over my shapefile I am still getting a different coordinate system (please see image)

Now, if I write click on my layer and go to properties, and choose WGS 84, the shapefile disappears from the map. So for whatever reason it doesn't seem to want to be this CRS.

In ArcGIS you can use the toolbox to define a transformation - is there equivalent in QGIS?

Changing the CRS makes my image disappear and that's the problem so ideally I need some suggestions beyond the usual steps (Reprojecting vector layer in QGIS?) as to what might be happening


texts: Creating KML that is synced to CartoDB data and styles?

texts:
The one part of CartoDB I am struggling with is connecting my table or visualization to Google earth. I see that you can export an un-styled KMZ from a table.

Has anyone figured out a way to keep a synced or network lined KMZ/KML file with polygons from a CartoDB table or visualization?


texts:
I have some z/x/y.png map tiles. While in the old days tiles were 256x256 pixels, nowadays it is not unusual to have images referencing the same extents in higher resolution (sometimes called "retina" in marketing speak).

This means the tile 1/2/3.png would represent the same geographical area in both cases but rendered in different resolutions.

When I tried to convert my directory tree to a mbtiles file with mb-util, the resulting file displays nothing in QGIS. Does the mbtiles format support higher resolution tiles and can I turn such files into a mbtiles file with mb-util?