texts: Fatal Error when PyQGIS searches for available providers: QgsApplication.initQgis()

texts:
I am setting the prefix path in Python using the code below. When running the code, a pop-up window stating "python.exe has stopped working" shows up. I've gone line by line in the code within  OsGeo4W shell and it seems that 'QgsApplication.initQgis()' is crashing Python.

from qgis.core import *

QgsApplication.setPrefixPath("C:\\OSGeo4W\\apps\\qgis\\", True)
QgsApplication.initQgis()
print QgsApplication.showSettings()

When printing the current settings I get the following paths:


texts: How to let leaflet layer events to propagate to layer below

texts:
I am working on a leaflet map with three layers: a base tile-layer, a GeoJson layer and a label tile-layer. For display, I want the label layer to be on top, but for interaction I would like to let click events on the map to propagate to the GeoJson layer, to show popups.

Here is an adaption of the leaflet geojson example that illustrates my problem:

http://jsfiddle.net/3h6yytg3/

Is there any way to show the GeoJson popups with a setup similar to this?


texts: Can't install the PyQt5 - UI to py conversion issues because of multiple python versions

texts:
I have both Python 2.7 (which comes with ArcGIS) and 3.7 (which comes QGIS 3.4) - in normal command prompt its showing Python 3.7.0 and path from qgis 3.4.
I confirmed version and path running following commands

python

>>> import sys
>>> print(sys.path)

so I tried to install PyQt5 by the following command

pip install PyQt5

but the command prompt throws the

  DEPRECATION : Python 2.7 will reach ...

error.

So as command prompt suggested I tried to upgrade the pip version by using following commands

python -m pip install --upgrade pip

but again I got the different error like

  already up-to-date

I don't know what to do.  Is there any problem with Version (2.7 or 3.7 issues) else am I doing mistake on installing/setting up with environment variables.

Note: I'm in the process of developing a tool which uses PyProj, GDAL, so I used Qt Designer for front-end. trying to convert UI to py. I'm a beginner so messed-up with these things.


texts: Is there a snap to grid function available in sf for R similar to ST_SnapToGrid(geometry geomA, float size)in PostGIS? It is mentioned as a solution to non-noded intersection problem for PostGIS and am looking for a solution in R. Using other snap options with just the polygon lead to geometry errors (even with small tolerance).


texts: Fixing QGIS generated SLD to avoid repeating polygon labelling

texts:
I have generated SLD's in QGIS that I use in GeoServer. Although the labeling shows fine when previewing the layer in GeoServer, when bringing the WMS in Leaflet then labeling repetition occurs. How do I center labeling in the QGIS generated SLD and avoid repeating labeling?

This is how the leaflet map looks like:

Desired outcome:

And here is the QGIS generated SLD:

<?xml version="1.0" encoding="UTF-8"?>
<StyledLayerDescriptor xmlns="http://www.opengis.net/sld" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:schemaLocation="http://www.opengis.net/sld http://schemas.opengis.net/sld/1.1.0/StyledLayerDescriptor.xsd" version="1.1.0" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:ogc="http://www.opengis.net/ogc" xmlns:se="http://www.opengis.net/se">
  <NamedLayer>
    <se:Name>Uis_Erven</se:Name>
    <UserStyle>
      <se:Name>Uis_Erven</se:Name>
      <se:FeatureTypeStyle>
        <se:Rule>
          <se:Name>Single symbol</se:Name>
          <se:PolygonSymbolizer>
            <se:Fill>
              <se:SvgParameter name="fill">#729b6f</se:SvgParameter>
              <se:SvgParameter name="fill-opacity">0</se:SvgParameter>
            </se:Fill>
            <se:Stroke>
              <se:SvgParameter name="stroke">#232323</se:SvgParameter>
              <se:SvgParameter name="stroke-width">1</se:SvgParameter>
              <se:SvgParameter name="stroke-linejoin">bevel</se:SvgParameter>
            </se:Stroke>
          </se:PolygonSymbolizer>
        </se:Rule>
        <se:Rule>
          <se:TextSymbolizer>
            <se:Label>
              <ogc:PropertyName>erf_num</ogc:PropertyName>
            </se:Label>
            <se:Font>
              <se:SvgParameter name="font-family">Arial</se:SvgParameter>
              <se:SvgParameter name="font-size">8</se:SvgParameter>
              <se:SvgParameter name="font-weight">bold</se:SvgParameter>
            </se:Font>
            <se:LabelPlacement>
              <se:PointPlacement>
                <se:AnchorPoint>
                  <se:AnchorPointX>0</se:AnchorPointX>
                  <se:AnchorPointY>0.5</se:AnchorPointY>
                </se:AnchorPoint>
              </se:PointPlacement>
            </se:LabelPlacement>
            <se:Fill>
              <se:SvgParameter name="fill">#000000</se:SvgParameter>
            </se:Fill>
            <se:VendorOption name="maxDisplacement">1</se:VendorOption>
          </se:TextSymbolizer>
        </se:Rule>
      </se:FeatureTypeStyle>
    </UserStyle>
  </NamedLayer>
</StyledLayerDescriptor>


texts: Configuring PyCharm to work with QGIS libraries

texts:
I have been following various posts and tutorials including Anita Graser's posts Anita Graser and https://gis.stackexchange.com/questions/302774/  as I try to set up my PyCharm installation to access and work with the QGIS libraries in my OSGeo4W64 QGIS 3.14 installation
I have created the pygis.cmd and pycharm.cmd files that are required and they appear to do what they should. My pycharm.cmd successfully launches PyCharm
> @echo off
SET OSGEO4W_ROOT=C:\OSGeo4W64
call "%OSGEO4W_ROOT%"\bin\o4w_env.bat
call "%OSGEO4W_ROOT%"\bin\qt5_env.bat
call "%OSGEO4W_ROOT%"\bin\py3_env.bat
call "%OSGEO4W_ROOT%"\apps\grass\grass78\etc\env.bat
@echo off
rem path %PATH%;%OSGEO4W_ROOT%\apps\qgis
path %PATH%;%OSGEO4W_ROOT%\apps\qgis\bin
path %PATH%;%OSGEO4W_ROOT%\apps\grass\grass78\lib
path %PATH%;C:\OSGeo4W64\apps\Qt5\bin
rem path %PATH%;C:\OSGeo4W64\apps\qgis\python\plugins
path %PATH%;C:\OSGeo4W64\apps\Python37\Scripts
rem path %PATH%;C:\OSGeo4W\apps\Python37

set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\qgis\python
set PYTHONHOME=%OSGEO4W_ROOT%\apps\Python37

set PATH=C:\Program Files\Git\bin;%PATH%
set QT_PLUGIN_PATH=%OSGEO4W_ROOT%\apps\qgis\qtplugins;%OSGEO4W_ROOT%\apps\Qt5\plugins

start "PyCharm aware of QGIS" /B "C:\Users\b8060507\AppData\Local\JetBrains\PyCharm Community Edition 2019.3.3\bin\pycharm64.exe"

However, despite following https://gis.stackexchange.com/questions/302774/ I continue to get this error message popping up

I am sure that I have set my Qt Plugin Path.correctly which post 302774 advises to include:
set QT_PLUGIN_PATH=%OSGEO4W_ROOT%\apps\qgis\qtplugins;%OSGEO4W_ROOT%\apps\qt5\plugins

in the pycharm.cmd file
Can anyone offer some advice on where I might be going wrong? I'm hoping 'Ben W' amongst others may see this question as his advice helped others in the 302774 post but it's not working for me in 2020. So perhaps there are other anomalies in the set up now required for later versions of QGIS?


texts: Creating table of point in polygon results using PostGIS

texts:
I have one PostGIS table of polygons like this:

index | lat | lon | geom

where geom is a polygon d-type
ex. POLYGON ((x1,y1,x2,y2.....))

CREATE TABLE buildings2 (
    index integer PRIMARY KEY,
    lat numeric,
    lon numeric,
    geometry geometry(Polygon,4326),
    poly geometry(Polygon,4326)
);

-- Indices -------------------------------------------------------

CREATE UNIQUE INDEX buildings2_pkey ON buildings2(index int4_ops);

And another of points like this:

index | point_id | lat | lon | geom

ex. geom = SRID=4326;POINT(32.123 12.123), lat = 12.123, lon = 32.123

CREATE TABLE pings (
    index SERIAL PRIMARY KEY,
    ad_id text,
    lat numeric,
    lon numeric,
    time text,
    hour integer,
    geom geometry(Point,4326)
);

-- Indices -------------------------------------------------------

CREATE UNIQUE INDEX pings_pkey ON pings(index int4_ops);

I want to create a new table that has all the points that are in polygons of the polygon table. My points file is large (100M rows), buildings(170), so:

What query am I supposed to use to create this table? I have tried contains, intersect, and within, but I am not sure which is the best for my problem, and all the stack exchange answers I have looked at provide a different solution.


texts: Is there a way to syncronize the current view extent in QGIS with Google Earth?

texts:
I already have the openlayers plugin but it's not working has I need it and I wanted to have the same view on both QGIS and Google Earth. Any ideas?

Update:

I'm using the openlayers plugin which enables me to add the google iamgery as a layer in the TOC. But it takes a long time to load and messes up the scale i'm working with. So i would like to have qgis with my layers, and simultaneously have google earth displaying the same extent view.


texts: Merged vector layers in two seperate rows in attribute table QGIS

texts:
I merged two shapefiles together, with the same geometry:

rain runoff data in a watershed
population living in that watershed

But when I open the attribute table, the information merged is presented in two rows, such as in the picture below. Any insight on why I can't have all the info in the same line for this watershed?


texts: Image classification: separating pine forest, lowland and upland hardwood

texts: I want to classify 20 feet-resolution 4-band leaf-off image into lowland hardwood, upland hardwood, pine trees, and other (bare ground, water, roads etc). I am thinking to use randomForest package in R and classify the image in 3 classes (hardwood, pine trees, and other) and then use stream buffers to separate lowland and upland hardwood. Would it improve classification to add layers with NDVI, DEM or/and topographic wetness index? Which other layers should I try?


texts: Can transferring a shapefile between local and network drive corrupt the attribute table?

texts:
Recently I ran into a troubling issue with ArcGIS (10.2.2) shapefiles. I transferred a shapefile with approximately 10,000 records from my local drive up onto a shared network drive at my office, using the basic function in ArcMap where you right-click the layer in the table of contents, select "Data," and then select "Export Data."

I then discovered that the copy of the file on the network drive had its attribute table significantly altered. Specifically, while the table had the same number of records as before and the file could be used normally, one particular record from the pre-corrupted version had somehow propagated itself throughout much of the rest of the table, overwriting nearly 25% of it. Thus, a previously diverse subset of records had all been overwritten with identical values from another record.

This is fairly disconcerting, because I have always been under the impression that "corrupted" shapefiles simply will not open or will raise the classic "Number of shapes does not match number of table records" error.

Has anyone encountered a similar problem before and can confirm that something like this can be caused by upload/download of shapefiles between local and network drives? Are there other potential causes for an attribute table being altered in this unique way? What are some precautions against something like this?

Standard ArcGIS tools like Feature Compare & Table Compare successfully detect the issue, but it would be nice to have enough confidence in file transfer processes to avoid having to incorporate such relatively paranoid checks into my workflow.


texts: Sharing feature classes between different geodatabases?

texts:
We're starting to work on different projects. All of them have their own data and classes but some of them are common (ex. parcels, traffic, road definition, etc.).

How can I share feature classes across projects?

I found only these approaches:

I could take from the source geodatabase but I cannot create relationships or make a copy in each geodatabase but I have the risk that some of them will not be updated (in any case I have to do manually).

My environment is based on ArcGIS and AutoCAD and it could be that we do an upgrade buying ArcGIS server.


texts: file geodatabase api on osx mavericks via brewed gdal

texts:
I've got a new machine with a fresh install of OSX 10.9, and am having trouble getting esri's file geodatabase api to work correctly with my brewed gdal.

This is a known issue, noted here:
http://forums.arcgis.com/threads/95958-OS-X-Mavericks

I've tried John Tull's noted workaround:

# install gcc
brew install --enable-cxx gcc49

# noted fgdb/gdal dependency:
brew install libdap --cc="gcc-4.9"

# gdal
brew install gdal --enable-unsupported --with-postgresql -v --cc="gcc-4.9"

# this points the fgdb libraries at gcc4.9 somehow
install_name_tool -change /usr/lib/libstdc++.6.dylib /usr/local/Cellar/gcc49/4.9-20131103/lib/gcc/x86_64-apple-darwin13.0.0/4.9.0/libstdc++.dylib /usr/local/opt/FileGDB_API/lib/libFileGDBAPI.dylib
install_name_tool -change /usr/lib/libstdc++.6.dylib /usr/local/Cellar/gcc49/4.9-20131103/lib/gcc/x86_64-apple-darwin13.0.0/4.9.0/libstdc++.dylib /usr/local/opt/FileGDB_API/lib/libfgdbunixrtl.dylib

This helps - ogr2info --formats now shows FileGDB as a valid format.
However, running ogr2info on a valid filegdb gives this error:

ogrinfo(45386,0x7fff74779310) malloc: *** error for object 0x10b49dde0: pointer being freed was not allocated
*** set a breakpoint in malloc_error_break to debug

Any idea what I'm doing wrong? I tried removing all of the gdal dependencies and letting homebrew reinstall them when calling brew install gdal with the --c"gcc-49" flag but that didn't change anything.

edit - I now see that this issue already has a thread here Homebrew OGR with FGDB support on OSX 10.9


texts: Can't execute plugin made by Plugin Builder

texts:
I have read this post:
Can't compile plugin builder plugin in QGIS 2.4 (Windows)

But I couldn't resolve my problem.

I have already seen this post, looking for help. I can execute:

pyrcc4 -o resources_rc.py resources.qrc

But when I want to open my plugin made by Plugin Builder, it displays an error, saying that it is "broken".


texts: How to calculate location at given time from GPX track

texts:
I have a number of sound files with timestamps (HH:MM:SS), made when walking a transect. I need to know where I was when each file was saved. To do this I have a GPX track of the route I took - I can extract the time for each point along the GPX track however the times don't match exactly with the sound files - the seconds aren't always the same.

These are my sound file time stamps:

Time       File
19:32:53   1
19:34:26   2
19:40:43   3
19:40:47   4
19:49:21   5
etc

These are the times of the GPX points:

Time       Latitude         Longitude
19:30:29   51.4463364054   -2.4129919242
19:37:44   51.4464039635   -2.4128513597
19:37:59   51.4464663249   -2.4127552193
19:38:01   51.4464837592   -2.4127533752
19:55:08   51.4465564303   -2.4127591588
etc

So in the above example I need to estimate where I was at e.g. 19:34:26.

I could achieve this by removing the seconds however I would lose a lot of accuracy and would leave me with redundant coordinates.

Is it possible to do this in GIS? Or might it be best to plot this as a graph? Either way, how would I do this?

I'm on Linux and ideally would like to use Calc/Excel, QGIS or R.


texts: Sample values from a grid by points in RQGIS

texts: I need to sample values from a grid by points. In RSGA, there is a function called "pick.from.saga.grid" to sample raster data set using points. Is there something similar in RQGIS?


texts:
I've encountered two issues when trying to display some aerial imagery data in QGIS and haven't managed to work out what the problem is. Frustratingly the issues only affect some tiles, and there doesn't appear to be anything obviously different between tiles that display fine and tiles that don't.

The two issues are:

Some tiles do not display at all in QGIS e.g. sx6346.jpg
Some tiles display, but have colourful horizontal lines across them e.g. sx6345.jpg.

In the screengrab above tile sx6346 is not displayed; tile sx6345 displays extraneous colourful lines. Adjacent tiles included to show that they display fine.

Regarding the extraneous colourful line issue there seems to be a threshold, which is different for each image, where they display fine at scales above a certain value. For example sx6345 displays fine at scale 1:19,845 and above, but has colourful lines across it at all scales below 1:19,845.

This threshold seems to be somehow linked to the size of the actual QGIS application window (if I make the QGIS window bigger, the image displays fine at all scales above 1:19,856).

QGIS displays the following errors for the above two issues:

Unsupported Data Source: C:\Data\aerial_imagery\sx\sx6346.jpg is not a supported raster data source libjpeg: Bogus Huffman table definition
Warning 1: libjpeg: Corrupt JPEG data: 285 extraneous bytes before marker 0xda [note: this error doesn't appear when adding the image to QGIS, but it does if you try and create a virtual raster catalog that includes the affected images (Raster > Miscellaneous > Build Virtual Raster Catalog...)]

The same images display fine in ArcGIS, MapInfo and in a range of image applications including GIMP, Paint and Paint Shop Pro. The above screen grabs are from the latest version of QGIS (2.16.1) but I’ve also tested them in multiple previous versions (including the current LTR 2.14.5) and unfortunately the behaviour is the same in all versions of QGIS.

The two error messages displayed in QGIS suggest problems with the source JPEG files, but given they work in every other package I've tried it appears that QGIS is being more rigorous in the code it's using to display JPEG's than everything else.

Using GDAL to do a JPEG to JPEG conversion resolves the problem but given we can't easily identify the affected files we'd have to process them all, which isn't practical given the number of files we have.

The two affected tiles can be downloaded from Dropbox.

Really, I'm asking if anyone has encountered anything like this before?

And do these issues sound like bugs in QGIS given the affected images work in every other package I've tried?


texts: Openlayers - How to put street view icon on google layer?

texts:
I have a basic OL map with google layer.

http://jsfiddle.net/Alophind/bDzg7/

I've read in similiar questions that the following css disable the google street view icon:

img {max-width: 100%;}

So I added (just for case) this css into my code:

img {max-width: none;}

And I still have no google street view icon.

What else can I do ?


texts: Extract specific altitude curve from raster

texts:
Is there any way one can extract a specific curve from a raster DEM qith QGIS? I know I can extract all curves, specifying the desired altitude interval between curves. But I was interested in eextracting a specific curve (i.e., I'd like to state something like "I want to extract a curve or a polygon corresponding to the 200 m values in the raster".

Any suggestions?


texts: How does one control the rendering order of geometries from a postgis geometrycollection?

texts:
I have been using PostGIS for a while now but have only been using simple geometries (points,polylines and polygons).
For a particular project of mine I am considering using the geometrycollection format. But I need some help understanding the concept. From what I have read, a geometrycollection stores points,polylines and polygons in a single table. When the geometries are displayed, how is the rendering controlled? If I want points and polylines to always be above polygons how would I define this?

Thanks for any comments,

Robert


texts:
I am trying to read the RGB values at certain locations on my
map which is a .tif-file I brought into QGIS.
Goal: match species presence points on my map to certain RGB values on the .tif-file.
from PIL import Image
img = PIL.Image.open(iface.layerTreeView().selectedLayers())
col = img.convert("RGB")
val = col.getpixel((0,0)) #tried to read value at center
print (val)

And getting the error:
AttributeError: 'list' object has no attribute 'read'


texts: What is the preferable value for setting tolerance for selecting a feature

texts:
I know it's really a basic question. But I am really confused about this. For selecting a point feature by mouse click, we need to set a tolerance value to create the envelope. And here is the problem lies. in some feature classes points are really so close, for some feature classes, map unit is in feet and for some it is in meter. And the tolerance should vary for these case. For searching feature I personally use IFeatureClass.Search method.

Does there any better option to search a feature from a feature layer? Or does there any way to pick a preferable tolerance?


texts: Distinguishing corridors from fragments in a habitat map

texts:
I was wondering if someone has delt with this kind of issue before:

I want to delineate fragments based on neighboring pixels of the same class (forest) with a distance threshold of 3 pixels (which is biologically significant for my specie of interest).

My concern is, see the attached image for an example, that sometimes these fragments are actually corridors, and often corridors and actual fragments are grouped into a same fragment because of their proximity.

I'm wondering if there is a way to distinguish the corridors from the fragments based on shape, number of surrounding pixels, etc..?

For instance, in the following box, possible corridors are indicated by the red boxes, and fragments by the green ones.

I have access to QGIS and R.


texts: Change feature class from ESRI native geometry to SQL Server geography

texts:
Original question:
I have found a lot of documentation on importing and using query layers from SQL Server into ArcCatalog/ArcMap/etc. However, when you change the projection of a feature class or use the editor tool in ArcMap (say to break up or split a polygon) then try to export it back into SQL Server, you cannot visualize the geometry when querying the table in sql server anymore.

I was wondering if anyone knew why this happened, and if there was a way to fix this either before you export into sql server or within sql server itself.

***************update*****************

I now know why I cannot see the spatial results in sql server; when exporting it is giving my feature class ESRI's default native geometry type, which SQL Server cannot visualize without conversion. From what I can find online, the "Migrate Storage" processing tool in ArcMap 10.3 should be able to fix this issue. However, after I ran the tool (and got a "succeeded" message in the results) it still does not show the shape column as a "geography" data type.


texts: How to check "simplification sensibility" of a geometry?

texts:
It is expected that a "geometry simplification" is a good practice to preserve stable versions of a geometry, avoiding  minor and non-relevant changes (side effect of geometry-editors and other "false updates").

The functions ST_Simplify, ST_SimplifyPreserveTopology or ST_SimplifyVW seems good candidates. It is important

to preserve  original geometry when it not needs simplification;
to control the "need less simplification" by a parameter, like tolerance.

So, the hypothesis is that we can check sensibility by changing a tolerance parameter.

... But it is not working, so the question is "How to check sensibility?" in that context.

I am testing ST_Simplify(geom,0.00001)=geom and is false.. But even ST_AsGeoJSON(geom,3)=ST_AsGeoJSON(ST_Simplify(geom,0.00000001),3)  is false!!

EXAMPLE:

 origianal  = ... [-60.931,2.64],[-60.929,2.638],[-60.925,2.635],[-60.919,2.63],[-60.917,2.628],...

 simplified = ... [-60.931,2.64],[-60.925,2.635],[-60.919,2.63],[-60.917,2.628],...

Usage context

I need to preserve stable versions into a git repository of GeoJSON files.

PS: if a "non-invasive simplification algorithm" not exists (?), them, a workaround,  a simple way to normalize geometry for git stable versions is to apply something as ST_buffer(ST_buffer(geom,0.000001),-0.000001).


texts: Get Lat Long coordinates from LiDAR .LAS 'Selected Point Attributes'

texts:
The third number for 'position' is accurate altitude, but how can I convert the first two to Lat/Long?

Selected Point Attributes:
position = 2237977.620 1354143.960 1460.720
intensity = 4
returnNumber = 1
numberOfReturns = 3
pointSourceId = 208
classification = 4
color = 0 1280 4864
diff with previous = 66.539
vector diff = (45.880 43.390 20.970)

Metadata text file at link below : (USGS_LPC_AL_JeffersonCo_2013_Jefferson_1167_LAS_2015_meta.txt)
https://shrib.com/?v=nc#Luciana7lm4mQv


texts:
I am new to Python and am trying to run a simple QA/QC on a particular set of fields.

I basically want to verify that fields x and y add up to z, I may have to incorporate rounding since a small difference doesn't make a difference and will throw off my results.  This is what I have currently and it is not processing:

def Qcontrol(check):
  if (!Forest_Acres! + !NonForest_Acres!) = !ACRES!:
     return 0
  else:
    return 1

Qcontrol(!check!)

all fields are double. **Also need to compensate for Nulls in Forest and Non_forest...


texts: Listing total number of features into an ArcGIS Online feature pop-up when the service consists of grouped features

texts:
I'm trying to use the count function to list the total number of features in a feature pop-up.  I can do this through configure pop-up/add expression Count($layer) if the feature service consists of separate features, but I cannot if the features are grouped.  If the features are grouped when publishing the service, the "Globals" tab in the add expression window doesn't list the same functions, most importantly the $layer functionality.
Is it still possible to perform the count?
I'd like to avoid ungrouping and republishing all of my services, if possible.  In the image below, you'll see the add expression window for the same feature.  On the left, the feature was published on it's own--it will look this way if published with multiple features, as well.  On the right, it was published as a group.


texts: Javascript syntax error on print map, only when feature count > 8000 (updated)

texts:
I am using arcgis jsapi v3.10, and esri/dijit/Print for my map printing.

I am encountering an very weird bug where the browser return:

SyntaxError: Unexpected token <
at Object.b.fromJson (init.js:231)
at Object.json (init.js:188)
at t (init.js:193)
at c (init.js:79)
at d (init.js:79)
at b.Deferred.resolve.callback (init.js:80)
at init.js:196
at k (init.js:200)
at r (init.js:200)
at f.resolve (init.js:202)

I've done multiple test while eliminating as much variables as possible. I print nothing except a point feature layer with simple maker symbol. As I slowly increases the feature count, I end up seeing the syntax error when the count reaches ~8200. I've also tried to use point from other layers just in case there are some data corruption in certain layer, but the results are similar.

One additional piece of information is that the 'print-start' event has been triggered and fully completed, whereas the 'print-complete' event never triggers due to the js error.

The print map function can do everything correctly including legends, scale, multiple layer, basemap e.t.c. when the feature count is low. Has anyone encounter similar issue? I suspect there might be some sort of feature count limit, but I couldn't find any related ESRI documentation so far.

EDIT:
Ok I've done some further digging. I basically add a log line into init.js to see what is it that triggers the syntax error. turns out that function expects a JSON response and I got a html page instead. The html page is a response from the print API and it says something like "error : WebMapAsJSON cannot be empty". I immediately check the outgoing request and verified that the WebMapAsJSON field is populated correctly. After more research online, it seems like that map server has its own maximum post content size limit and I've exceeded it.

I've followed some instructions and change maxPostSize=0 to disable limitation for post size in C:\Program Files\ArcGIS\Server\framework\runtime\tomcat\conf\server.xml. However, another error pop-up after I retry the test case with the new settings. This time it is in the server error log "Internal Server Error. Error handling service request : Cannot write the request to the server:"

Is there some other settings I am still missing?


texts: Some output cells show up as noData in ArcMap Block Statistics?

texts:
I am using ArcMap 10.2.2, and I need to aggregate categorical data from 30m up to 90m and 900m resolution for species distribution modeling. The algorithms available in resample (nearest neighbor, bilear, etc.) are not appropriate for my purposes, and the Aggregate tool does not have a "majority" option. Thus, I think I have to use the Block Statistics tool. However, I am getting a strange, recurrent error where some cells in the output raster show up as noData. Has anyone else seen this problem? I got the same effects when working with both fgdb rasters and geotiffs. The same apparently random cells show up as noData over multiple runs. The input data does not have any noData cells in the areas where I'm seeing the problem.
As a work around, I can do this in R using the aggregate function in the raster package. Their algorithm is inexplicably slow, though, and it's more convenient for me to work in Arc if possible. Here are the messages from the results panel, just in case they're helpful.

Messages
Executing: BlockStatistics NLCD_2011_WV_WGS_1984.tif C:\Users\mosteele\Desktop\aggregation_test\NLCD_2011_WV_90m_ArcMap_version2.tif "Rectangle 3 3 CELL" MAJORITY DATA
Start Time: Thu Jun 25 11:45:37 2015
Succeeded at Thu Jun 25 11:45:40 2015 (Elapsed Time: 2.46 seconds)

Input:

Output (noData in red):

Upon further inspection, I realized that the nodata cells only occurred where there was a tie for majority. I contacted an Esri support person who told me that this is the expected behavior.
I can understand why the algorithm was created like this, but I'd argue that it would be much more useful to have two options - one that would set ties to nodata and one that would, say, break ties randomly. I needed the output surface for species distribution modeling and a version with a bunch of nodata is useless to me. I'm also unsatisfied with Esri's other methods of resampling categorical data, such as the majority option in the Resample tool or resizing by snapping, as they drop a lot of information from the original raster by looking only at a limited window around the centroid of the new cell.
Fortunately, I was able to create the surfaces I needed using R - specifically, the aggregate tool in the raster library.


texts: How to fit an OpenLayers view to an extent that spans the map extent (wrap) or IDL?

texts:
How to fit the view to an extent that spans the map extent (wrap-point) or IDL?

Please excuse my semantics if they are incorrect, I'm new to GIS.

I've set up some code that does a zoom to extent based on these OpenLayers examples:

http://openlayers.org/en/latest/examples/center.html

http://openlayers.org/en/latest/examples/draw-features.html

The code (see below), works fine for extents that stay within the extents of the world map, e.g. Africa.

It does not work as I would like for my Pacific region, which spans the IDL and the point at which the map "wraps" in OpenLayers.

The actual behaviour is that the line between the 2 points goes "the long way", showing me almost the whole map, instead of the short way which gives only the Pacific. See screenshots after code sample.

var Pacific = {
            bottomLeft: [112.67578, -49.42705],
            topRight: [-125.85937, 45.78093]
        };

var Africa = {
            bottomLeft: [21.62109, -33.92741],
            topRight: [60.11719, 37.80327]
        };

var region = Pacific; //Africa;

myLib.zoomToExtent(region.bottomLeft, region.topRight, function (){;});

zoomToExtent: function (bottomLeftLonLat, topRightLonLat, onZoomCompleteCallback) {
        var map = _getMap(_domTarget);
        var view = map.getView();

        var size = /** @type {ol.Size} */ (map.getSize());

        var ptA = ol.proj.fromLonLat(bottomLeftLonLat, _defaultProjection);
        var ptB = ol.proj.fromLonLat(topRightLonLat, _defaultProjection);

        geometry = new ol.geom.Polygon(null);

        var start = ptA;
        var end = ptB;

        geometry.setCoordinates([
          [start, [start[0], end[1]], end, [end[0], start[1]], start]
        ]);

        view.fit(geometry, size, { padding: [1, 1, 1, 1], constrainResolution: false });

        _numTilesLoaded = 0;
        _setBaseLayerLoadedCallback(map, onZoomCompleteCallback);
}

Actual effect of the above code (not exact, for example):

Desired effect (not exact, for example):


texts:
New to the GIS, I want to upload these file formats to PostgreSQL server
.tfw
.tif
.tif.aux.xml
tif.vat.cpg
tif.vat.dbf
tif.xml
tif.ovr

I used the standard org2ogr way to upload files but it wouldn’t recognise any these formats, as none of the drivers listed in the error message support these file formats.
I’ve also noticed I don’t have raster2pgsql.exe in my Postgresql/13/bin/ directory. What might be the solutions?
(Edited to include gdalinfo as requested)
Driver: GTiff/GeoTIFF
Files: C:\filepath\filename.tif
       C:\filepath\filename.tif.ovr
Size is 43201, 15779
Coordinate System is:
GEOGCRS["WGS 84",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["geodetic latitude (Lat)",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic longitude (Lon)",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4326]]
Data axis to CRS axis mapping: 2,1
Origin = (-180.000000000000000,72.000059984999993)
Pixel Size = (0.008333333332238,-0.008333333332214)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (-180.0000000,  72.0000600) (180d 0' 0.00"W, 72d 0' 0.22"N)
Lower Left  (-180.0000000, -59.4916067) (180d 0' 0.00"W, 59d29'29.78"S)
Upper Right ( 180.0083333,  72.0000600) (180d 0'30.00"E, 72d 0' 0.22"N)
Lower Right ( 180.0083333, -59.4916067) (180d 0'30.00"E, 59d29'29.78"S)
Center      (   0.0041666,   6.2542267) (  0d 0'15.00"E,  6d15'15.22"N)
Band 1 Block=43201x1 Type=Byte, ColorInterp=Gray
  Min=1.000 Max=4.000
  Minimum=1.000, Maximum=4.000, Mean=nan, StdDev=nan
  NoData Value=255
  Overviews: 21601x7890, 10801x3945, 5401x1973, 2701x987, 1351x494, 676x247, 338x124
  Metadata:
    STATISTICS_MAXIMUM=4
    STATISTICS_MEAN=nan
    STATISTICS_MINIMUM=1
    STATISTICS_STDDEV=nan


texts: Handling "bad layers" with Esri File Geodatabases?

texts:
File Geodatabases are not connected to the same as other OGR data sources. In the Add Vector Layer dialog, the user has to select Directory as the Source Type, then further specify the type (as Open FileGDB or Esri FileGDB), then use the file browser to select the top-level folder of the FileGDB.

Example:

When one moves files around, data may not be where QGIS expects to find it. On opening the *.qgs file, the Handle Bad Layers dialog appears. However, this opens a file browser that appears to be the same one used to add a file-type OGR data source, such as Esri Shapefile, GeoJSON, KML, GeoPackage, etc. As this dialog is not used to add a FileGDB layer in the first place, it does not appear to work in the Handle Bad Layers dialog. When you navigate to the top-level folder and attempt to select it, the file browser just drills down into the FileGDB hierarchy.

Is there a workaround for this behavior? How can I reconnect data sources (what the Handle Bad Layers dialog does) for Esri File Geodatabases when files have been moved and paths have changed?


texts: How do I distinguish freeway onramps from offramps?

texts:
Freeway entrances and exits are tagged as highway=motorway_link. However, for some analysis I'm doing, I need to distinguish between entrances ("on ramps") and exits ("off ramps").

Is there an automated way to do this with existing OSM tags?


texts: arcpy field mapping for a spatial join: keep only specific columns

texts:
I want to perform a spatial join to find the join point (a node in the code below) that is closest to a target point (a consumer in the code below). In the result, I only need certain attributes from both input features:

node: "ID", "name"
consumer "ID", "a_value"

For doing this in arcpy it seems that I need a FieldMappings object to specify which columns I want to keep, but I can't figure out how. There are few examples in the documentation, but they appear to be more complex than what I need (merge rules an all that)

How do I perform a spatial join in arcpy where the result only contains specific columns?

fieldmappings = arcpy.FieldMappings()
fieldmap_consumer = arcpy.FieldMap()
fieldmap_node = arcpy.FieldMap()

# How to specify columns to be included in join result?

arcpy.SpatialJoin_analysis(target_features = paramConsumerFeatures,
  join_features = paramNodeFeatures,
  out_feature_class = paramWS + "\\" + paramOutFile + "_tmpjoin",
  join_operation = "JOIN_ONE_TO_ONE",
  join_type = "KEEP_ALL",
  field_mapping = fieldmappings,
  match_option="CLOSEST",
  distance_field_name="join_dist")


texts: Query to select attributes except specific fields

texts: What is the expression in ArcGIS to select values except for specific attributes?


texts: Prevent drawing of overlapping point symbols - ArcGIS Server/Javascript API

texts:
Is there any way in ArcGIS to prevent symbols from being drawn when they overlap each other?  I know this is possible for Labels using Placement Properties/Label Weight, but I can't seem to find corresponding functionality for symbology.  The resulting map document will be published as an ArcGIS Server Map Service and displayed using the ArcGIS Javascript API.  The data will be updated frequently, thus the service must be drawn dynamically (no caching).

Basically, I have a point feature class that can be very dense or very sparse depending on zoom scale.  When zoomed out, there are too many symbols to display simultaneously in a coherent way.  In an ideal world, I'd like to be able to prevent cluttering by specifying a buffer distance (in pixels) around each symbol where no other symbols can be drawn from the same layer.  I don't want to cluster the symbols together - I only want to show a single point to represent the area.

My symbology is also fairly complex, and involves graduated symbols based on data value (as well as value-based rotation) as well as several labels with specific placement:

Short of thinning/downsampling the dataset before symbolizing, what are my options (if any) for reducing symbol density?


texts:
In SQL Server 2012 I can do this:

-- Construct a new FullGlobe object (a WGS84 ellipsoid)
DECLARE @Earth geography = 'FULLGLOBE'

-- Calculate the area of the earth
SELECT PlanetArea = @Earth.STArea()

Is there a similar construct in Oracle Spatial?

Alternatively what would be the best way to represent the entire earth in an 8307/4326 SRID?

Is there anything in the WKT standard that speaks to this?


texts: Mapbender App Publishing (disable localhost)

texts:
This is surely an elementary question but I am new in the Mapbender world. Using OSGeoLive 13 on the Oracle VM Virtual Box, I have created a WMS (1.3.0) in GeoServer (2.15.1). In Mapbender (3.0.8.1), I have created a new application and did all of the necessary configuration, defined the layersets and the layouts.

The problem is I see the app just as localhost (also in GeoAerver) so I cannot call the app from another computer. My aim is to publish the app for the broad public and not just as localhost. How I can do this; which settings are necessary?

I don´t know if this has something to do with the virtual machine itself or if I need to do something in GeoServer/Mapbender. Or is this in OSGeoLive in general not possible and do I need to install GeoServer and Mapbender manually? I am not sure if the OSGeoLive contains automatically after installation the Apache or no (maybe this could be the issue). I know how to create the WMS layer incl. styles etc. or configure the app in Mapbender but I don´t know to publish is as an open web-service for everyone...


texts: QGIS: Terrain accessibilty in hours/days usind DEM/CSV data

texts:
I am an Archaeologist. I am interested in settlement dispersal and land use strategies in late antiquity and early middle ages.

Therefore, I want to analyze a certain terrain by:
- least cost polygons or
- terrain accessibility by time from a certain point.

I am using QGIS 2.18.6
ASTER DEM 30m (tiff)
.CSV data for my settlement points.

GRASS and SAGA don´t want to operate.
I calculated the slope as friction for r walk/ r drain.

I am surely not an expert and there might be some basic calculations missing.
Maybe you can help me along or have another strategy how to draw my polygons.


texts: Change encoding type for all files or whole project

texts:
How can we change encoding mode for the whole project rather than doing to single shapefiles in QGIS 3.4?

I wanted to change "UTF-8" to "Shift_JIS" to all project vector files.


texts:
This is a tricky one.

A little context first: Bentiu wiki page: https://en.wikipedia.org/wiki/Bentiu + OSM map of Bentiu refugee camp, South Sudan: http://www.openstreetmap.org/#map=16/9.3314/29.7920

My task is to come up with a map of the dwellings and the facilities in the camp, each with their corresponding address.

For this I have (1) a schematic grid of the address system (.pdf), (2) indications on how the address system works (start numbering the buildings from the NW corner to E then S), and (3) the .shp of the individual houses and buildings extracted from aerial imagery analysis.

Problems: (a) while the grid assumes that buildings are perfectly-aligned squares, the shp shows different shapes of dwellings that partially follow the spatial organisation set up for the address system; (b) the drawing and the numbering are currently being used for healthcare provision and cannot be modified to fit the structure of the dwellings.

So far I have tried a couple of tricks unsuccessfully, mainly because the drawing is only a representation therefore cannot be georeferenced. I tried to assign every house in the .shp a number with $rownum but this did not work either because the order was not respected.

I'm more than willing to discuss any ideas!


texts:
I'm running GDAL 1.11.2, released 2015/02/10.

The following displays all the attributes in my shapefile on stdout:

  /opt/gdal-custom/bin/ogr2ogr -f CSV /vsistdout/ foo.shp -sql "SELECT *  FROM foo "

However, if I add the WHERE clause  WHERE ST_Intersects(GEOMETRY,BuildCircleMbr(0, 0, 10000,  4326)) I get ERROR 1: Undefined function 'ST_Intersects' used.

What am I doing wrong?


texts:
I am using leaflet 1.0.3. Its a simple map having a geojson layer and no basemap at all. I want to change zoom level step from 1 to 0.1 or 0.25 so that when user scroll or press zoom in or zoom out, map zoom level should change e.g. from 6 to 6.1 instead of 7.
I tried zoom snap but it only change the value from integer to decimal but the step value remain the same (i.e. 1)


texts: How to get values from a raster attached to a line in ArcGIS 10.1

texts: In MapInfo, using Vertical Mapper, I can do a line inspection and get the Max, Min, Mean etc values for a polyline from one or more open grids which is then added as an attribute to the polyline. I'm trying to do a similar process in ArcGIS 10.1 but have been unable to work out which tool I would need to run. Any suggestions, I have access to most extensions including Spatial Analyst and have ArcInfo level licencing.


texts: How to import CSV layer to QGIS using PyQGIS in Windows?

texts:
I'm using QGIS 2.18.14 Las Palmas on Windows 8.1
I tried to import CSV layer to QGIS using this script:

uri = 'file://D:/26mei.csv?d?delimiter=%s&xField=%s&yField=%s' % (",", "Lon", "Lat")
layer_csv = QgsVectorLayer(uri, 'paramHilal', 'delimitedtext')
QgsMapLayerRegistry.instance().addMapLayer(layer_csv)

The last line is to add the CSV layer to layer panel in QGIS, but it didn't appear. It could be because the layer is not valid.
I knew because I have tried layer_csv.isValid() and it return False. But I don't know where my mistake is.
Does anyone knows why my script is not work?

PS: I have tried this on QGIS 2.8 Wien on Ubuntu 16.04 and it worked just fine.


texts: Saving Feature Class Symbology into Layer or Object to apply to another Feature Class?

texts:
I am creating a tool to automate a data standardization process. Part of this process is saving the feature class (variable "filename") as a layer file, projecting it into another coordinate system (variable "wmsde"), and then saving that one as a copy as well. I want the user to symbolize the feature class and then run the script/tool which would take the symbology they designated and apply it to the layers created. I realize it is incredibly easy to right click the layer and click Save As Layer, but I want as much as possible of this process to be automated. I am trying to use representations but a layer with the symbology needs to exist before the tool is used, which would be fine if I could find a way to save the symbology to a layer file. I have tried arcpy's Make Feature Layer and Save to Layer File (Management) to create a layer from the data but it does not save the symbology. I read on an answer to another question on here that this may be normal. I think the representation tool is working because  the second layer's symbology is the same color as the first's (and they are usually assigned random colors). However the Drop Representation tool is failing and I'm not sure why. I get error 000252. Here is the (edited) code I have for it. Other tools that create symbology for a layer such as Apply Symbology From Layer (Data Management) also requires a symbolized layer.

Does anybody have any ideas as to what I can do to save the user's symbology designations or even another route aside from representations?

arcpy.env.workspace = r'C:\Data\Data.gdb'
worksde = r"C:\Users\me\AppData\Roaming\ESRI\Desktop10.2\ArcCatalog\DatabaseName.sde"
lyrFolder = r"C:\Layers"
ws = arcpy.env.workspace
filename = r"PlayAndTest"
sdefilename = r"PlayAndTest"
infile = os.path.join(ws, filename)
sdefile = os.path.join(worksde, sdefilename)

# Create a DD layer file
print "Creating DD LYR"
DDlyr = os.path.join(lyrFolder, filename + "_DD.lyr")
arcpy.SaveToLayerFile_management(filename, DDlyr)

# Reproject as WM
print "Projecting to WM"
wmsde = sdefile + "_WM"
prjfile = os.path.join(arcpy.GetInstallInfo()["InstallDir"],
                   "Coordinate Systems\Projected Coordinate Systems\World\WGS 1984 Web Mercator (auxiliary sphere).prj")
spatialRef = arcpy.SpatialReference(prjfile)
arcpy.Project_management(infile, wmsde, spatialRef)

# Make the WM Layer
WMlyr = os.path.join(lyrFolder, sdefilename + "_WM.lyr")
arcpy.SaveToLayerFile_management(sdefilename, WMlyr)
arcpy.ApplySymbologyFromLayer_management(WMlyr,DDlyr)


texts: Using Calculate Value to get extent of feature to input into Generate Tessellation tool

texts:
Here is what my model looks like:

Takes Excel with lat/long and plots points
Splits by name field
Iterate
And buffer into separate .gdb
Calculate value then pulls from the buffer.gdb

This is where I get lost:

Inputs extent into Generate Tessellation

I found this code from a separate question on here, but I haven't figured out how to input into the Calculate Value and then add the resulting extent into the tessellation.
arcpy.Describe(fc).Extent

Here is the other post; Generate Tessellations looping through all shapefiles from folder in ArcGIS Pro ModelBuilder?


texts: network related issues arising when running an arcpy script

texts:
Essentially I wrote an arcpy script that is within a toolbox that I ported over to a different machine to test out after the script worked locally.

The script essentially utilizes the psycopg module to connect to PostGIS and subsequently geoserver which I ran locally as opposed to over a server.

The new machine that I am working from in which I am running geoserver over the server is giving me errors as far as running the script. Errors which i didn't experience running on a local instance of geoserver.

The local machine is running ArcGIS 10 service pack 2 and the other machine is running ArcGIS 10 service pack 3.

Here is the error I'm getting and I think its J2EE related:

  AssertionError: Tried to create PostGIS Layer but got 500: binding :
  binding : binding : binding

  ---- Debugging information ----

  message             : binding : binding

  cause-exception     :
  com.thoughtworks.xstream.mapper.CannotResolveClassException

  cause-message       : binding : binding

  class               : org.geoserver.catalog.impl.FeatureTypeInfoImpl

  required-type       : org.geoserver.catalog.impl.AttributeTypeInfoImpl

  line number         : 1

I've been able to narrow down that the issue is with me giving it feature types when I need to be giving attribute types as far as uploading to geoserver to display. However considering this is an arcpy script that I've written in Python, is there any way to salvage this without going to Java, as I know that the com.vividsolutions.jts.geom series of classes include attribute defining functions such as getAttributes etc.

The reason being the way that I originally retrieved all the data necessary in order to go from my PostGIS installation and post to geoserver was due to extracting the PostGIS data via queries and then posting them in the geoserver required format. So, I think I may be able to if I can figure out a Postgres series of queries that would get me the attributes instead of the features.


texts: Joining Google Docs table to feature class at ArcGIS Online to use in Operations Dashboard?

texts:
I am using ArcGIS operation Dashboard.

I get the Data Source from excel file that i have published using google docs as a dynamically data source.

I need to join the excel table to hosted feature service on Arcgis-Online and give it diffrent colors and diffrent labels by the joined attributes.

For example - Building Layer Joined to money collection table - green buildings already paid, red building didn't paid yet. The color of the layer is changing by the attribute of the table. And when I update the excel table the layer will automatically be updated.

How can I join this table (google docs) to a feature class in AGOL and use it on Operations Dashboard ?


texts: Popup for ArcGISDynamicMapServiceLayer point data

texts:
I was going through this example from ESRI Resource
popup.
I was able to integrate it to my polygon layer but it is failing for point and line layer.I am not sure but

query.spatialRelationship = esri.tasks.Query.SPATIAL_REL_INTERSECTS;

could be the reason.I tried the other combination of the task but was not successful if anyone has another example or if anyone could suggest something in existing example then it will be very helpful.I am using ArcGIS server 10.2,Javascript API 3.1.


texts: Seeking free nautical route planner Data?

texts:
I don't really need an application, but rather an API. I plan to develop my own application and don't need someone else's GUI, just their information.

I want to know the distance a vessel at sea must travel between two points. Obviously this isn't just a straight line thing, if there is land on that straight line.

In my previous question, Geoist asked "You want something that takes into account most efficient routes given currents, tides, etc?".

My answer is, yes, please, if that is possible. If it costs extra, then I would prefer to start with a free solution and add that later.

Now, I can't expect everyone to be an application developer like me, so I can't expect everyone to understand the difference between an application and an API.

Briefly, an application is a whole, complete, stand-alone program, something like a Windows .EXE file or even a browser based program.

Whereas an API is an Application Programming Interface; it provides and interface which allows a different application to get at some of its inner workings (like a Windows DLL file, or just accessing someone's database, but not using their GUI (which is exactly what I want to do)).

In this case, I am looking for an interface which will allow me to query the distance between two GPS points, possibly with route planning, etc.

I agree with @chriserik route plotting at sea is not as simple as taking the distance "as the crow flies", because the crow in a boat which is subject to tides, winds, currents, etc

Note that I could just run OpenCpn and monitor the HTTP traffic from the program in order to reverse engineer the interface and then writ emy own program. But I won't do that on ethical grounds.


texts: name 'AddressLocator' is not defined, Python created Address Locator

texts:
When I run this tool, I get the error 'name 'locator_AddressLocator_n_' is not defined'.

My guess is that this is an issue with the path where the locator is trying to save to?  So, I am assuming I have not set the workspace or the scratch workspace?  These are just guesses though.

I will need to save each of the geocoding outputs, so should I just have the code create a new GDB and push the outputs to that location?  How could this be accomplished?

import os
import arcpy
from arcpy import env

env.overwriteOutput = 1

Here is the geocoder code:

arcpy.CreateAddressLocator_geocoding("General - Single Field", inputFC+"'Primary Table'", "'Feature ID' '' VISIBLE NONE;*KeyField GeoLabel VISIBLE NONE;'Additional Field' '' VISIBLE NONE;'Altname JoinID' '' VISIBLE NONE", locator_AddressLocator_n_, "")

UPDATE:

Created a GDB from paramter and set the workspace to that location, but I am still getting the same error when running the script.  name 'AddressLocator' is not defined

inputGDBname = arcpy.GetParameterAsText(3)
inputGDBloc = arcpy.GetParameterAsText(4)

#Create GDB for address locators
arcpy.CreateFileGDB_management(inputGDBloc, inputGDBname)

env.workspace = inputGDBloc+"/"+inputGDBname

UPDATE

New Error on the last line...

RuntimeError: Object: Error in executing tool

I have included all of my code this time, just so everything is visible to see what could be throwing the error.

import os
import arcpy
from arcpy import env

env.overwriteOutput = 1

#Create GDB for address locators
inputGDBname = arcpy.GetParameterAsText(0)
inputGDBloc = arcpy.GetParameterAsText(1)
arcpy.CreateFileGDB_management(inputGDBloc, inputGDBname)
env.workspace = os.path.join(inputGDBloc, inputGDBname)

#Get FC and Fields
inputFC = arcpy.GetParameterAsText(2)
inputField = arcpy.GetParameterAsText(3)
inputFieldLabel = "GeoLabel"
inputFieldName = arcpy.GetParameterAsText(4)

#Add the GeoLabel field to the Selected Feature Class
arcpy.AddField_management(inputFC, inputFieldLabel, "TEXT")

#Calculate the field values for the selected field and populate it with the inputFieldLabel + InputField
rows = arcpy.UpdateCursor(inputFC)
for row in rows:
    idt = row.getValue(inputField)
    exp = inputFieldName + ': ' + idt
    row.setValue('GeoLabel', exp)
    rows.updateRow(row)
del row, rows

Reference_Data = inputFC + ' Primary Table'
RefData = Reference_Data

arcpy.CreateAddressLocator_geocoding("General - Single Field", RefData, "'Feature ID' <None> VISIBLE NONE;*KeyField GeoLabel VISIBLE NONE;'Additional Field' <None> VISIBLE NONE;'Altname JoinID' <None> VISIBLE NONE", "Address_Loc_n_", "")


texts: using r.volume model GRASS plugin in QGIS

texts:
I am calculating volume of bathymetric data in the form of a clipped raster in GRASS with a Custom Lambert projection centered around 70degrees North, using the module r.volume in GRASS.

Side note:
I have concluded that my original raster appears have both depth of the ocean (in negative figures) and height above see level (in positive figures). However I have now clipped it to represent only information covering the water in the NE Atlantic. Please see at the bottom a copy of all the info in the metadata tab for this bathymetric file and do not hesitate to correct me if you notice something strange as this information has changed a lot since I imported it into GRASS. If you require the metadata of my original bathymetric file I can provide it.

My question is:

How do I know that the r.volume module in GRASS has accurately calculated the volume of the ocean for this area? Please see screen shots below of my data and of the r.volume results generated.

Additionally can someone please confirm that this output is in units m3 seeing as the LAMBERT projection is in meters? Is this correct?

Driver
GRASS raster provider
GISDBASE: /Users/aninahenggeler/Documents/Humboldt/Thesis /High seas proposal/GIS & Google Earth/GRASSDATA
LOCATION: allnewicesareasgeotiffLAMBERTLocation
MAPSET: allnewicesareasLAMBERTclippedgeotiffraster
MAP: LAMBERTCLIPPEDwithnewICESdensifiedareasGRASS
MAX_VALUE : 1.74800000000000000e+03
TYPE : 0
MIN_VALUE : -6.37400000000000000e+03
No Data Value
-2.14748e+09
Data Type
Int32 - Thirty two bit signed integer
Pyramid overviews
Layer Spatial Reference System
+proj=laea +lat_0=70 +lon_0=0 +x_0=0 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
Layer Extent (layer original source projection)
-3685610.5270899999886751,-3707057.7249250002205372 : 2075705.6314520000014454,2218572.3506209999322891
Band
Band 1
Band No
1
No Stats
No stats collected yet


texts: Measuring distance between two points considering elevation raster?

texts: Is it possible to compute the distance between 2 coordinates (x1,y1,z1) - (x2,y2,z2) considering an elevation raster (cells containing elevation)? On the image, there are 2 coordinates (red crosses). I would like the green distance. I know how to get the direct distance taking the elevation into consideration (red dashed line on the picture), but this does not take that there can be a mountain or a valley between the 2 points into account. I am using R but if there is a nice way to do it in QGIS that would be great too.


texts:
I have been downloading Sentinel2 and Sentinel3 images for Convolutional nerual networks but I keep getting black TIFF files.

Can you check my code?
I have been using the GEETOOLS by fitoprincipe, can be viewed here: link.

The code:

var start = ee.Date("2018-04-01")
var end = ee.Date("2018-12-30")

var region = ee.FeatureCollection(allRegions.filterMetadata("system:index","equals","00000000000000000000"));

var filter = ee.Filter.and(
  ee.Filter.date(start, end),
  ee.Filter.bounds(region))

var batch = require('users/fitoprincipe/geetools:batch');

var s2_dates = ee.ImageCollection('COPERNICUS/S2_SR')
  .filterDate(start,end)
  .filterBounds(region)

var dates = s2_dates
    .map(function(image) {
      return ee.Feature(null, {'date': image.date().format('YYYY-MM-dd')})
    })
    .distinct('date')
    .aggregate_array('date')

var s2 = ee.ImageCollection('COPERNICUS/S2_SR')
  .filterDate(start,end)
  .filterBounds(region)
  .map(function (image) {
    return image
      .normalizedDifference(['B8', 'B4'])
      .updateMask(image.select('QA60').not())
      .rename('S2')
      .clip(region)
      .set('date', image.date().format('yyyy-MM-dd'))

  })

var s3 = ee.ImageCollection('COPERNICUS/S3/OLCI')
  .filterDate(start,end)
  .filterBounds(region)
  .map(function (image) {
    return image
      .select(['Oa17_radiance', 'Oa08_radiance'])
      .divide(ee.Image([0.00493004, 0.00876539]))
      .normalizedDifference(['Oa17_radiance', 'Oa08_radiance'])
      .rename('S3')
      .clip(region)
      .set('date', image.date().format('yyyy-MM-dd'))
  })

var s3_dis = s3.distinct('date')

print("S2",s2)
print("S3",s3_dis)

batch.Download.ImageCollection.toDrive(s2, 's2_images', {name: '{system_date}'})
batch.Download.ImageCollection.toDrive(s3_dis, 's3_images', {name: '{system_date}'})

The reason I am using this tool is to prevent manually downloading every single image.


texts: Assigning attribute by location for multiple layers

texts:
I am trying to assign values by location of one of my polygon layer.
It is overlayed by multiple polygon layers and I want to assign the value of the layer that covers it the most.
Is there a tool or a way you can think of in QGIS?

To make it more clear, I have data as single lines and created buffers of the lines and then areas. The buffers overlap with the areas and I want to assign the yellow (garden) attribute to the areas mostly covered by that, also the same for other values.
So in the example in the image above the area marked with G will get the yellow value as it is mostly covered by that layer and not the pink one.


texts: Obtain virtual output directory path in geoprocessing service

texts:
I am publishing a geoprocessing service and I want to be able to access the virtual output directory path configured for the service.

Now I know that gp services will magically return virtual paths for output parameters of data type file, but is there a way to access the virtual output directory path string directly?

right now I'm doing something very hackish:

def get_results_virtual_path(resultsFilePath):
    file_url = urlparse.urljoin('file:', urllib.pathname2url(resultsFilePath))
    if 'directories' in file_url:
        return SERVER_VIRTUAL_DIRECTORIES + file_url.split(r'directories')[1]
    else:
        return file_url

Its basically hard-coding the virtual path.  Now I know the virtual path can be complicated by load-balancers and subdomains and such, but all I want to do is read the environment variable from the service properties.


texts: Using Google Maps API to dynamically display delivery-time zones?

texts:
I'm looking to replicate the functionality from the UPS Ground Map which will allow users to enter a zipcode, then dynamically show delivery times by region. I know I need travel times by state data to do this. But lets say I did, could this be possible to do over the Google Maps API if I had this data?

If yes, could anyone offer any advice on how to accomplish this? I'm not looking for code handouts, but any references/pointing me to useful documentation would be great.


texts: How to create a feature mask with a polygon?

texts: I need to create the SAME polygon shape shown in the image below, just 2 miles larger.  The buffer wizard tool is only creating circular buffers around the irregular shape below.  Any ideas on how to preserve the shape just 2 miles larger?


texts:
I am trying to zoom into a location on a single click. How can I zoom towards the point on single click?

I have tried using OpenLayers.Control.ZoomIn and adding a click handler, but for some reason its not registering the "click".

Here is a run-down of my code:

map = ...; // map initialization
... // added layers

// list of tools for map/feature manipulation
tools =
{
  // ... tools
  // Tool In question :
  zoomIn : new OpenLayers.Control.ZoomIn({
                // zoom_in needs a click handler, otherwise, its just button on the map
                handler : new OpenLayers.Handler.Click(
                    this, // control
                    {'click' : this.trigger}, // callbacks
                    {'single' : true, 'double' : false, 'pixelTolerance' : 0} // options
                ),
                // overriding ZoomIn.trigger so I can call console.log and alert
                trigger : function(e) {
                    alert("Zoom in, please!");
                    console.log(e);
                    this.map.zoomIn();
                }
          }) // end zoomIn
} // end tools

for(key in tools) {
   map.addControl(tools[key]);
}

// I only made zoomIn active for testing...
tools.zoomIn.activate();

So, the alert doesn't happen and neither does the console.log(e)... How come the callback function is not getting called when I click?


texts: In CartoDB, when I switch to Map View on any of my tables, I have the opportunity to change the basemap, including using "Yours".  In my attempts to connect to our MassGIS XYZ basemaps, I've filled up all 8 slots and I can't figure out how to clear these out?


texts:
I'm creating a new feature layer like that:

createNewFeature: function (layer, geometry, layerInfo) {
    var loJsonFS = {
      "displayFieldName": layer.displayFieldName,
      "fieldAliases": layer.fieldAliases,
      "geometryType": layer.geometryType,
      "spatialReference": layer.spatialReference,
      "fields": layer.fields,
      "features": [{
        "attributes": loAttributes,
        "geometry": loGeometry,

      }]
    };
    var loFS = new FeatureSet(loJsonFS);

    var featureCollection = {
      layerDefinition: {
        "geometryType": layer.geometryType,
        "fields": layer.fields
      },
      featureSet: loFS
    };

    var loFL = new FeatureLayer(featureCollection, {
      mode: FeatureLayer.MODE_ONDEMAND
    });

    this.map.addLayer(loFL);

}

But it does'nt show in the map. The attributes that receive are empty, and the geometry is based in a draw that user can do in the map.
What I need to show it in the map?


texts: Calculate the nearest distance between each pixel inside my region of interest to the closest streamline (a FeatureCollection) in Google Earth Engine

texts:
I need to find the nearest distance between every pixel in my area of interest to the streamline closest to it. I'm using the following dataset for the streamlines ee.FeatureCollection("WWF/HydroSHEDS/v1/FreeFlowingRivers").
Initially, I thought image.cumulativeCost() would do the trick, but it seems like increasing the maxDistance only increases the buffer around the streamlines. The resulting value of cumulative_cost is 0 when the area is within the 1km buffer of the streamlines, and 'masked' otherwise.
What are the values of the 'masked' pixels? Is it even feasible to achieve what I need using image.cumulativeCost()?
var streams = ee.FeatureCollection("WWF/HydroSHEDS/v1/FreeFlowingRivers");
var roi = ee.Geometry.Polygon(
        [[[31.79857587268242, 40.34237834590226],
          [31.79857587268242, 39.1087726665225],
          [35.77013348986992, 39.1087726665225],
          [35.77013348986992, 40.34237834590226]]], null, false);

//polylines representing river networks
var filteredstreams = streams.filterBounds(roi);

var streamImg = filteredstreams
  .reduceToImage({
    properties: ['LENGTH_KM'],    //picked a random property here just to convert the FeatureCollection
    reducer: ee.Reducer.first()   //into an image
});

//create a source image where the geometry is 1, everything else is 0
var streamsources = ee.Image().toByte().paint(roi, 1);
streamsources = streamsources.updateMask(streamImg);
var emptysource = ee.Image().toByte().paint(roi, 0);

//define distance
var maximumdist = 1000; //1km

//cost image
var thedist = emptysource.cumulativeCost({
  source: streamsources ,
  maxDistance: maximumdist,
});
print(thedist);
Map.addLayer(thedist);


texts:
I have a large (3,000+) set of drainage polygons that I derived from the ArcHydro Stormwater Processing workflow. I would like to aggregate them based on the geometric network elements they contribute to. I essentially need to chose a pour point, place a flag, trace upstream, then select polygon features by location (see figure).

Is there a way to automate this process given a pre-determined set of point features (pour points)? I could do it manually, but there are ~500.


texts: Selection subquery using another table in database ArcGIS

texts:
I need to 'extract' a list of variables based on the IDs from another table. Using a subquery seems to be the best method. I have a large list of variables and have outside calculated a second list. This second list does not have coordinate data attached to it, so I need to select the same values in the original list to pull the data out. The shapefile and table are both in the catalog.

I try to do a simple select by attribute. It gives me the

SELECT * FROM US_Maize_All WHERE: "pointid" = (SELECT "pointid" FROM Maize_Table)

but this only gives me a SQL expression error. The "pointid" parameters are selected from the option box. As far as I can tell this is the exact format as any other info I can find online. Has to be something I'm overlooking why this will not work.


texts:
I'm working with shapefile files, but I need to convert to ASCII grid, and I was trying to use QGIS, but the following error occurred:

File "/usr/share/qgis/python/plugins/GdalTools/tools/dialogBase.py", line 208, in onFinished
            msg = unicode(self.process.readAllStandardError())
        UnicodeDecodeError: 'ascii' codec can't decode byte 0xc3 in position 24: ordinal not in range(128)

I am new to the area, when I am doing the conversion, only my shapefile file is enough? Can I perform this procedure in GDAL?

I know there are related questions, but none ends up solving the problem.


texts: How to project from Greek Grid to ETRS LAEA

texts:
How to make a project transformation between Greek Grid (EPSG:2100) to ETRS LAEA L52_M10?

I'm using ArcGis 9.3.1, but when I use the project tool, it don't have anything in "Geographic Transfornation" section.

How is it possible to make the transformation with proj4 library?


texts: Optimal ArcGIS PDF export compression settings for raster topographic basemaps

texts:
There are a number of different compression options for exporting pdf files from ArcGIS (eg RLE, LZW, Deflate etc) but which is best for exporting raster topographic basemaps?

Generally such maps will have quite a lot of white space, along with areas of even tone (eg lakes, forests). I'm looking for guidance on what settings work well for folks - clearly I'm trying to minimise file size whilst preserving quality.

I've had reasonable luck with the 'adaptive' compression method, but relying on trial and error is quite difficult because there are a few parameters to play with: resolution, output image quality, and compression-specific quality as well.


texts: Add or subtract days to date field using the field calculator

texts:
I am using ArcGIS 10.1

I'm not a programmer

I want to use the field calculator to add 4 days to the current date if the status is "delivery"

or

subtract 5 days if its status is "order"

this result need to have in the field of "New Date"

current date  13/11/2014

If the status is "delivery" is new date 16/11/2014

If the status is "order" new date is 09/11/2014


texts:
I want to use osgeo tools in a virtualenv with Django. I already followed the directions in the Django manual succesfully on different linux systems. But this time I got stuck while intalling GDAL.

My system is Lubuntu 15.10.

$ tar xzf gdal-1.9.2.tar.gz
$ cd gdal-1.9.2
$ ./configure
$ make

After a while 'make' hangs with this message:

libtool: compile:  gcc -g -O2 -Wall -Wdeclaration-after-statement -DOGR_ENABLED -I/home/user/tmp/gdal-1.9.2/port -I/usr/include -c json_object.c  -fPIC -DPIC -o ../../o/.libs/json_object.o
In file included from /usr/include/string.h:634:0,
                 from /home/user/tmp/gdal-1.9.2/port/cpl_port.h:132,
                 from /home/user/tmp/gdal-1.9.2/port/cpl_conv.h:34,
                 from json_object.c:12:
json_object.c:29:9: error: expected identifier or '(' before '__extension__'
   char* strndup(const char* str, size_t n);
         ^
GNUmakefile:23: recipe for target '../../o/json_object.lo' failed

I could trace the problem back to this command

$ gcc -g -O2 -Wall -Wdeclaration-after-statement -DOGR_ENABLED -I/home/user/tmp/gdal-1.9.2/port -c /home/user/tmp/gdal-1.9.2/ogr/ogrsf_frmts/geojson/jsonc/json_object.c

which results in:

In file included from /usr/include/string.h:634:0,
              from /home/user/tmp/gdal-1.9.2/port/cpl_port.h:132,
              from /home/user/tmp/gdal-1.9.2/port/cpl_conv.h:34,
              from /home/user/tmp/gdal-1.9.2/ogr/ogrsf_frmts/geojson/jsonc/json_object.c:12:
/home/user/tmp/gdal-1.9.2/ogr/ogrsf_frmts/geojson/jsonc/json_object.c:29:9: error: expected identifier or ‘(’ before ‘__extension__’
char* strndup(const char* str, size_t n);

relevant lines in json_object.c are:

12 #include "cpl_conv.h"
13
14 #include "config.h"
15
16 #include <stdio.h>
17 #include <stdlib.h>
18 #include <stddef.h>
19 #include <string.h>
20
21 #include "debug.h"
22 #include "printbuf.h"
23 #include "linkhash.h"
24 #include "arraylist.h"
25 #include "json_object.h"
26 #include "json_object_private.h"
27
28 #if !HAVE_STRNDUP
29   char* strndup(const char* str, size_t n);
30 #endif /* !HAVE_STRNDUP */

I realize I am not working with the most recent version. It looks like Gdal-1.9.2 has a bug which caused my problem as commented by @huckfinn but finally I used the recent version Gdal-1.11.2 which does not cause any trouble. So problem is 'solved'


texts: How can I show only some classes of a layer in UMN-mapserver?

texts: I visualize a shape-layer in the mapserver. Now I want to make some requests for only some classes of this layer. I want for example one request for only the water-classes and later for only the wood-classes. Can I realise this with the mapserver?


texts: How to estimate overlap of certain pixel colors in a raster image and a linestring shape added on top of it?

texts:
I am trying to calculate how much of the electric feeders exist underneath tree covers. I am using NAIP infrared satellite imagery and have a shp file for a feeder circuit. Like so:

I need to find out how much of my green lines actually exist underneath the vegetation cover, and do segmentation of the line shape file based on that information. I have looked at deep learning capabilities within ArcGIS.
I am wondering if there is a simpler way of doing this using ArcGIS?


texts:
I am looking at applying a direction indicator to land parcels and houses (in MasterMap).

Effectively, I am looking at applying a direction indicator to the land and the buildings associated to the land (NW-SE, N-S, E-W etc).  Is there an easy way to do this?

I guess you'd need to also find the direction of the road it was situated on for completeness?


texts: Is there any tool that will let me visualize a WKT polygon on a world map such as google maps?


texts: GeoTools correct use of ReprojectingFeatureCollection

texts:
I have trouble understanding the use of the function ReprojectingFeatureCollection.
Here is my code:

try {

    Double buffer = 500.0;
    String attributeBuffer = "";

    ReprojectingFeatureCollection clip_repojected = new ReprojectingFeatureCollection(clip, CRS.decode("EPSG:3857"));

    myBufferM buffered = new myBufferM();
    SimpleFeatureCollection bufferedFeatures = buffered.createBufferM(clip_repojected.collection(), buffer, attributeBuffer);

} catch (OperationNotFoundException e1) {
    // TODO Auto-generated catch block
    e1.printStackTrace();
} catch (FactoryRegistryException e1) {
    // TODO Auto-generated catch block
    e1.printStackTrace();
} catch (NoSuchAuthorityCodeException e1) {
    // TODO Auto-generated catch block
    e1.printStackTrace();
} catch (SchemaException e1) {
    // TODO Auto-generated catch block
    e1.printStackTrace();
} catch (FactoryException e1) {
    // TODO Auto-generated catch block
    e1.printStackTrace();
} catch (IOException e) {
    // TODO Auto-generated catch block
    e.printStackTrace();
}

So my function gets a SimpleFeatureCollection as parameter and now tries to Reproject it to EPSG:3857. I have modified the gs:BufferFeatureCollection from the list of WPS processes of GeoServer and now call it myBufferM. I am trying to Buffer a FeatureCollection in meter instead of degrees (data's default unit). Previously I had the transformation to meter build into my Class myBufferM which worked good (transforming every feature geometry itself and then buffering it and reverse transforming it).

Now I am trying to transform the geometries beforehand using the Buffer function as it was before I modified it.

The question is: am I using the ReprojectingFeatureCollection correctly? Because the results I get now are wrong (in degrees).

If you need more code I will update my question. I have reduced it to the important parts. The total code would be kind of huge.


texts: Plot a series of subsetted pixels from a raster stack in R

texts:
I have created two raster stacks and would like to plot a row or column of pixels through the stack, but am not sure how to approach this.

I know I can subset by layer, but how do I subset vertically through the stack at a specific pixel location?

# Subsetted stack by layer
> tr_sub <- subset(tr_stack, 1:10)


texts: Triggering events when loading and finished loading from WMS

texts:
I'm currently creating a webgis application using OpenLayers 3.0.0, and I'm trying to display a loading gizmo when loading layers from a WMS. I've found examples for OL 2.1.0 that registers 'loadstart' and 'loadend' events.

However I haven't had any success with registering those events in OL 3.0.

Anyone got any pointers on how to do this or know if this is even possible to do in OL 3.0?

Thanks!


texts: New QGIS user, 2.6 Brighton.  My project CRS is default to Nad83ft / N.Carolina (US).  I loaded various polygon and line (provided to me) that have various CRS's but they still show correctly in place over USGS / ortho-imagery.  Oddly, setting each file's CRS to NAD83ft / N. Carolina moves them out of view. Clearly I'm missing something about CRS basics.


texts: Retrieve pixel value with geographic coordinate as input with gdal

texts:
I have a list of X,Y coordinates in UTM called coords. I also have a .tif of a digital terrain model (DTM) also referenced in UTM. I would like to use the Python wrapper for gdal to return the pixel values (i.e. elevation values) for each coordinate in coords.

In my search for an answer so far I have only found this answer, which is related, but does not just use a simple list of coordinates and is seemingly a bit more complex than the simple task I am trying to achieve.

My concern is not so much the particular script that needs to be written, but the methods or functions within gdal that can perform this type of function. Is there a simple function within gdal that takes a geographic coordinate and returns a pixel value?


texts: Trouble inserting name of feature in TOC into attribute table using arcpy.calculateField

texts:
I currently have about fifty rivers and creeks in a map, and none of the features contain the name of the waterbody. However, the features are named correctly in the table of contents

Ideally, I would have a script that created a field named 'NameFromLYR', and then took the name of the feature from the TOC and inserted it into the new field in the attribute table.

This is the script so far:

import arcpy

mxd = arcpy.mapping.MapDocument (r'C:\Users\athom\Documents\ArcGIS\2014\ScriptTest.mxd')

lyrList = arcpy.mapping.ListLayers (mxd)
for lyr in lyrList:
    layerName = lyr.name
    NameVar = lyr.dataSource
    arcpy.AddField_management (NameVar, "NameFromLYR", "TEXT")

    expression = "reClass(!NameFromLYR!)"
    codeBlock = """def reClass(lyrName):
        return (layerName)"""

    arcpy.CalculateField_management (NameVar, "NameFromLYR", expression, "PYTHON", codeBlock)
    print layerName
    del layerName
    del codeBlock
    arcpy.RefreshTOC()

del mxd
del lyr
del lyrList

The script works, but there is one flaw - it uses the name of the first feature in the TOC and inserts it into all the features.

As you can see in the image below, 'AOI' has been added into NameFromLYR field in all the features, not just the feature called AOI.

Does anyone have any hints on what is wrong with this script? It appears as though the name of the first feature in the TOC does not get deleted, as it stays with the script for a few minutes after the features in the TOC have been shuffled, or ArcMap and PyScripter has been turned off and then on again.

I realize that using an update cursor is the smart way (easier and faster), but I have spent enough time on trying to use the arcpy.CalculateField_management that it would be nice to know what is wrong.

EDIT:
The script work with an update cursor if the expression, codeBlock and arcpy.CalculateField_management are replaced with this:

with arcpy.da.UpdateCursor (NameVar, "NameFromLYR") as cursor:
    for row in cursor:
        row[0] = layerName
        cursor.updateRow(row)


texts: Importing text only out of string field to new field for labelling in ArcMap?

texts:
I want to import data out of a string field in my shapefile to a new field for labeling.

The string field contains river names and codes and I just want to display the text i.e river names


texts: Updating master layer with data imported from GPS in QGIS?

texts:
All the data I import comes in from GPX files in dribs and drabs.  I load the file into a new layer using GPSTool or "Add Vector Layer" but neither of these give me any control over the fields created in the new layer.  This makes merging data into the master table difficult.

Ideally I want to be able specify what features get imported and what they are named in the new layer and also be able to populate fields like source, etc. with static data or alternately have a tool that takes the raw imported layer and adds the contents to a master layer while mapping feature names and adding static values.

It seems to me that this must be a common workflow in the GIS world and I want to make sure that I am not missing something before start coding something -- I am thinking of a more sophisticated version of the Merge Layer plugin from Michael Minn.


texts:
Could you tell me why bbox loading all data beginning?

'bbox=' + extent.join(','); },  strategy: ol.loadingstrategy.bbox

The map then loads very long. How can I make it calculate correctly?


texts: Can someone explain in simple terms the licence on Google “satellite” images?

texts: What can I and can I not do with the Google satellite images?


texts: 3D Multipatch feature class to 3D GeoJSON in ArcGIS Pro?

texts:
I have generated some building blocks in CityEngine and exported them as "Multipatch feature class" in GeoDataBase, then I imported them into ArcGIS Pro and tried to convert them into 3D GeoJSON format, but it turns out that the "Feature to JSON" tool available in Toolbox does not accept Multipatch feature class, only polygons, polylines, and points. So I first tried to convert the 3D Multipatch feature into Multipolygon, I have tried "Multipatch footprint", which turned out to be losing all the Z values of the building floors. And "Feature Vertices To Points" and "Multipatch to Raster" also does not suit my purpose (I think).

So how can I convert Multipatch to Multipolygon without losing Z-value, or convert Multipatch into GeoJSON format in general?


texts: For what ist the epsg-file in the "share/proj" subfolder of the QGIS - installation used for?

texts:
In the "share/proj" subfolder of the QGIS-installation there is an epsg-file with  definitions of the epsg-codes like this:

# DHDN / 3-degree Gauss-Kruger zone 3
<31467> +proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +datum=potsdam +units=m +no_defs  <>

For reprojecting layers (permanent and on the fly) the definitions from this file seem to be irrelevant as the definitions of the srs.db seem to be used.

Does it make sense to update this file if custom changes had been done to the srs.db file or is this file not important for QGIS as the srs.db is used?


texts:
Does anyone know where I can find the (very old) WorldWind Java source code for version 0.6?
For example: http://patmurris.blogspot.com/2009/03/worldwind-java-06-source-available.html


texts: Getting geometry of a feature for reduceRegion function in GEE

texts:
I would like to do an analysis for a few countries. I want to write a loop that calculates the loss for a given slope for each country.
I have a problem with geometry in reduceRegion function. It doesn't accept the geometry I provide.
I already did the same analysis for a country(with peru.geometry()), and it worked. But now, I am getting the following error when I call countries from a list by .get(),
List (Error)
Image.reduceRegion, argument 'geometry': Invalid type. Expected type: Geometry. Actual type: Feature.

How can I get the geometry of a feature in GEE? and How can I add this geometry to the code?
here is the link for my script:
https://code.earthengine.google.com/3a7d3053935e83b15e216dbef2056522


texts: "File not found at..." error when loading KML into Google Maps

texts:
I have a KML file on a public URL like http://host.com/dir/file.kml

If I put that KML file on another server, and put the new URL in the search box on Google Maps, it loads the KML into google maps.

If I put the above URL into google maps like so: http://maps.google.com/maps?q=http://host.com/dir/file.kml

It does not load, and I get the following message:

File not found at http://host.com/dir/f....

Suggestions:

    Make sure the URL is spelled correctly.
    Make sure the file exists.

Therefore, I assume this is a configuration issue on the server. They are not being very helpful with debugging/configuring. Can anyone here offer advice on how to fix this?


texts: Error in Universal Kriging using gstat package in R

texts:
I am trying to perform universal kriging in R using the gstat package and running into the following error:

krige(formula = CDD_LT_stC ~ Imperv1km + Elev_new, locations = sp.cal_data, newdata = sp.studygrid, model = modelfit, nmax = 200)

Error in eval(predvars, data, env) : object 'Imperv1km' not found

Any ideas why is this happening and how to fix this?

However, I was able to perform ordinary kriging with:

krige(formula = CDD_LT_stC ~ 1, locations = sp.cal_data, newdata = sp.studygrid, model = modelfit, nmax = 200)


texts: Creating planning unit shape to use with Marxan?

texts:
There used to be a tool in the old QMarxan plugin to create the planning unit shapefile but the new version QMarxan Toolbox no longer has this tool.

Any ideas on how to create the planning unit shapefile using QGIS?


texts: LANDSAT Image classification and open source solution(s)

texts:
I have some experience with ERDAS IMAGINE in doing some supervised/unsupervised LANDSAT imagery classifications and now find myself in the position of needing to do the same but without the access to ERDAS IMAGINE. I've been looking at the GRASS literature and found this here and see that it has similar functionality/tools such as i.class/cluster/gensig as well as PCA and Tassled Cap analysis, which also hold some interest to me. I plan on looking at SAGA as well when I get the chance to see what it can do in respect to classification. I also found this thread as well here but much of what is there dates back a few years.

My question is: are there other Open Source alternatives similar to GRASS or SAGA that I should look into as well that you may have used? Just doing my homework before I commit to something...and I am really only interested in the image classification abilities, or lack thereof. I will be using LANDSAT data for an AOI of 700 sq/mi and would like to utilize the ASTER spectral library. The class number would probably be fairly low as it is not a very diversified landscape. Operating system is Windows 7, 32bit(!) and no budget to speak of. Time frame will be mid-summer, full leaf on.


texts: Need to combine the lines created by Vector grid in QGIS into one polyline

texts:
I am having a fairly specific issue. I am starting off by creating 100m2 squares around specific points on a map by adding a buffer to the point and then using this buffer as the extent for a Vector grid. This works fine in that it creates the squares as needed but when I export the shp files to create gpx files for use in the field with a garmin etrex 20 I find that instead of the file representing one track as a square, it is instead 4 tracks each representing one side of the square.

I have traced the problem back to when the squares are created by Vector Grid. The attribute table shows 4 different lines being created rather than one polyline. How do I merge these 4 lines into 1 so that I can then continue with the rest of the workflow to create the squares as needed.

I realise I could manually edit the final gpx file in Notepad but this seems a long winded way of solving the problem. I have also tried outputting as polygons rather than lines in Vector Grid but for some reason is then gives me 200m2 squares instead of 100m.

The squares are needed for field research in precise locations and it is critical that the centre point is followed and the researcher can easily follow the Garmin to find the relevant grid when in the landscape.


texts: How to find neighbouring polygons with same attribute value?

texts:
I have certain amount of polygons with different classes (numbers). I would like select/find all these polygons that have a neighbor-polygon with the same class.

Is there a tool in  QGIS, ArcGIS (basic)? What would be the right approach?

Here an example of the data set:


texts: Converting DMS to decimal degrees, what precision?

texts:
I'm working on a piece of software that converts Lat & Long from DMS (35° 19′ 20″ N) to decimal degrees.

At the moment, the code creates number with 16 digits after the decimal point, which seems absurd.

Since we never have more than 2 digits for the arcseconds, would it be ok to round the decimal degrees we calculate? I'd like to round to the nearest 4th digit.


texts: How do I retrieve an Arc ASCII grid file from NED 1/3 second data using Python and OWSLib

texts:
I have my bounds, and have been able to connect to the server http://imselev.cr.usgs.gov:80/wmsconnector/com.esri.wms.Esrimap/USGS_EDC_Elev_NED_3? in Python.

What I don't know how to do is send a request to download the area in Arc ASCII grid format.  I would be grateful for any suggestions, examples, or pointers to where I might find examples.


texts: Minimum Bounding Box by Area using ArcObjects

texts:
I'm trying to calculate the minimum bounding geometry by area for a polygon in arcobjects. I know that the the MinimumBoundingGeometry class is exposed but I was hoping that this functionality was available at the IGeometry level compared to the feature class level.

This is what I'm trying to calculate; RECTANGLE_BY_AREA, for an IGeometry. I know you can quickly get the envelope as well as the convex hull in ArcObjects but I dont know how to get the minimum bounding box by area. Any thoughts or direction would be great, thanks!


texts:
I am trying to add an interactive map HERE to my ArcMap? I am using ArcGIS 10 and the path I have been following so far is through "Add WMS Server" (Catalog -> GIS Servers -> Add WMS Server). Unfortunately, this is not  working for me. The way I see it, this can have two reasons:

Did I do something wrong here ?
This is not the way to go.

In case of the former, what is that I am doing wrong?

I am using the AppID and AppCode from my HERE developer account. I have tried different things for the URL. The last one I tried was simply https://www.here.com/?map=51.54916,5.64706,10,normal (I was running out of ideas). Usually it tells me there is a connection error.

In case of option 2, that this isn't the way to go, what should I do to display HERE maps in ArcMap? It surely isn't offered through Add Basemap/ArcGIS Online.

In QGIS it is possible to add HERE maps in multiple ways. See this post: HERE background maps in QGIS possible?
I hope there is a way around to get it fix in ArcMap as well.

I think that there is a third option, and that it isn't possible at all, period. But I sure hope somebody can tell me how I can get it done.


texts: latest version of QGIS fails to edit attribute table of GeoJSON file

texts:
I am using QGIS 3.6.4 and working with GeoJSON files. I have come across a very annoying problem. Whenever I try to add an attribute and assign values to it using the field calculator, the values will appear to update and the field looks correctly populated. However as soon as I exit out of the attribute table (while saving changes) and reopen it the entire column is converted to NULL.

I have tried playing around with the precision and width of the output field to no avail. I don't think this is related to coordinate system units since this is not a geometric operation of type $area or $length. It is more like Acol1+Bcol2.

Has anyone had such an issue before? Please help... I would have thought basic functionality like this is already well tested and bug free in QGIS....

EDIT

Other symptoms include:

new field disappears if I remove the layer and add it again from file location, despite having saved layer changes....
a separate geojson file called result.geojson.tmp is created and contains the modified attribute table. However my result.geojson is unchanged. Adding an attribute field and populating it should not create a new file but modify the current one in place... why does this happen?
I have also attempted this with QGIS 3.6.1 - Same problem.

The problem appears to be under discussion in the QGIS issue tracker


texts: How to link my QGIS plugin main class to mouse events from a custom Tool?

texts:
I want to execute the merge(self,aLayer) method after clicking with the tool I've created, i.e. after the mouse click I want the control to come to the main.py module how can I do that using Python programming?

I have two files one is mouseevent.py and the other one main.py.

class MouseClick(QgsMapTool):
    canvasDoubleClicked = QtCore.pyqtSignal(object, object)

    def __init__(self, canvas,layer,Vertex_list,mergeDlg):
        QgsMapTool.__init__(self, canvas)
        self.canvas = canvas
        self.layer=layer
        self.Vertex_list=Vertex_list
        self.mergeDlg=mergeDlg
        self.rb = QgsRubberBand(self.canvas, QGis.Line)
        self.rb.setColor(QColor(255, 0, 0))
        self.rb.setWidth(2)
        self.rb.reset()
    def canvasPressEvent(self, event):

        if event.button() == 1:
            print "helo"
            x = event.pos().x()
            y = event.pos().y()
            Lpoint = self.canvas.getCoordinateTransform().toMapCoordinates(x, y)
            print 'Then Left points are',Lpoint

    def canvasMoveEvent(self, event):
        pass

    def canvasReleaseEvent(self, event):
        #Get the click
        if event.button() == 2:
            print "hello right"
            x = event.pos().x()
            y = event.pos().y()
            Rpoint = self.canvas.getCoordinateTransform().toMapCoordinates(x, y)
            print 'Then Right points are',Rpoint
            self.Vertex_list.append(Rpoint)
            print 'The list after right click is',self.Vertex_list

    def canvasSelectionChanged(self,layer):

        pass
    def activate(self):
        pass

    def deactivate(self):
        pass

    def isZoomTool(self):
        return False

    def isTransient(self):
        return False

    def isEditTool(self):
        return True

Main.py:

def merge(self,aLayer):
    canvas=qgis.utils.iface.mapCanvas()
    aLayer=canvas.currentLayer()
    aLayer.removeSelection()
    selection=aLayer.selectedFeatures()
    print 'length is',len(selection)

def run(self):
    """Run method that performs all the real work"""
    mergeDlg = MergeDialog()
    Vertex_list = []
    global Vertex_list
    MergeByFreeHand.readConfigfile(self)
    AreaList=[]
    Value=99
    canvas=qgis.utils.iface.mapCanvas()
    aLayer=canvas.currentLayer()
    aLayer.removeSelection()
    qgis.utils.iface.setActiveLayer(aLayer)
    count=0
    global count
    tool = MouseClick(canvas,aLayer,Vertex_list,mergeDlg,count)
    qgis.utils.iface.mapCanvas().setMapTool(tool)


texts: How to populate Unique ID field after sorting?

texts:
I am trying to create an new unique id field in a feature class. I already have one field called "SITE_ID_FD", but it is historical. The format of the unique value in that field isn't what our current format is, so I am creating a new field with the new format.

Old Format = M001, M002, K003, K004, S005, M006, etc

New format = 12001, 12002, 12003, 12004, 12005, 12006, etc

I wrote the following script:

fc = r"Z:\test.gdb\testfc"

x = 12001

cursor = arcpy.UpdateCursor(fc)

for row in cursor:
    row.setValue("SITE_ID", x)
    cursor.updateRow(row)
    x+= 1

This works fine, but it populates the new id field based on the default sorting of objectID. I need to sort 2 fields first and then populate the new id field based on that sorting (I want to sort by a field called 'SITE' and then by the old id field "SITE_ID_FD")

I tried running the script from the python window in ArcMap after manually sorting the 2 fields in hopes that the python would honor the sort, but it doesn't. I'm not sure how to do this in python. Can anyone suggest a method?


texts: Exporting map to PDF from template using PyQGIS?

texts:
I'm trying to update a QGIS plugin to PyQGIS 3 but I have some troubles to update the code for exporting layout to PDF. The layout derives from a template file. This is the "old" code:

canv = iface.mapCanvas()
rend = canv.mapRenderer()
rect = QgsRectangle(iface.mapCanvas().extent())

qpt = u'D:\\templates\\qvet_h_template.qpt')
ff = open(qpt, 'r')
tmplt = ff.read()
ff.close()
#zoom on the map
tmplt = tmplt.replace("iymin", str(rect.yMinimum()))
tmplt = tmplt.replace("iymax", str(rect.yMaximum()))
tmplt = tmplt.replace("ixmin", str(rect.xMinimum()))
tmplt = tmplt.replace("ixmax", str(rect.xMaximum()))

doc = QDomDocument()
doc.setContent(tmplt)

comp = QgsComposition(rend)
comp.loadFromTemplate(doc)

pdfpath = u"D:\\outputmap2.pdf"
out = comp.exportAsPDF(pdfpath)

As far as i know, it is necessary replace QgsComposition class and its methods (loadFromTemplate() and exportAsPDF()) as they are no longer present in pyqgis-3. I read to use QgsLayout to replace QgsComposition class, I tried to write some code but with failure. A draft of the code is as follows (the first part is the same as the above code):

doc = QDomDocument()
doc.setContent(tmplt)
ele = doc.documentElement()

project = QgsProject.instance()
lay = QgsPrintLayout(project)
lay.addItemsFromXml(ele, doc, QgsReadWriteContext())

exporter = QgsLayoutExporter(layout)
exporter.exportToPdf(u"D:\\TestLayout.pdf", QgsLayoutExporter.PdfExportSettings())