texts: Batch converting netCDF to Raster using ArcPy?

texts:
I'm having trouble using the Dimension Values (Optional) argument in arcpy.MakeNetCDFRasterLayer_md.

Whatever I seem to do the output raster is always the first layer which is the default value.

Below is my attempt to export all layers to raster in the time dimension.

Code:

def extractAllNetCDF():

    variable = "RRt_10m"
    x_dimension = "lon"
    y_dimension = "lat"
    band_dimension = ""
    dimension = "time"
    valueSelectionMethod = "BY_VALUE"

    outLoc = "E:/New Folder/"
    inNetCDF = "E:/netCDFFiles/RRt.nc"

    nc_FP = arcpy.NetCDFFileProperties(inNetCDF)
    nc_Dim = nc_FP.getDimensions()

    for dimension in nc_Dim:

        top = nc_FP.getDimensionSize(dimension)

        for i in range(0, top):

            if dimension == "time":
                dimension_values = nc_FP.getDimensionValue(dimension, i)
                nowFile = str(dimension_values)

                arcpy.MakeNetCDFRasterLayer_md(inNetCDF, variable, x_dimension, y_dimension, nowFile, band_dimension, dimension_values, valueSelectionMethod)
                arcpy.CopyRaster_management(nowFile, outLoc + nowFile + ".img", "", "", "", "NONE", "NONE", "")

            print dimension_values, i

The print method at the end will show the dates as they should be and the index i is also moving alon so there is no reason to think that there are other problems with the code other than the Dimension Value being incorrect and reverting to the default.

Does anyone have any idea how to get the subsequent layers to export?

Is there any code online that has a specific example of this using the Dimension Value argument other than empty quotes?


texts: Has anyone else been unable to locate and install the GRASS plugin in QGIS 2.0.1?

texts: After uploading and installing QGIS 2.0.1, I am unable to locate and install the GRASS plugin.  The plugin is not listed in the 'Plugin Manager'.  The 'Plugin_Installer' is also shown in red, indicating that there is a problem with the plugin.  Has anyone else experienced the same problem?  The look and feel of 2.0.1 is great though, and there are some nice new features.


texts: ogr2ogr giving "ERROR 1: Failed connect to 192.168.4.18:8181; Operation timed out" error?

texts:
This is my command, which was working OK yesterday:

ogr2ogr -f "GeoJSON" ../geojson/myfile.geojson ../data/gml/myfile.gml
  -s_srs EPSG:27700 -t_srs EPSG:4326

The output is this:

ERROR 1: Failed connect to 192.168.4.18:8181; Operation timed out

I'm baffled, why would ogr2ogr need to connect anywhere?


texts: GIS system for displaying GDSII(*.gds) format?

texts:
I want to display vector graphics on raster graphics and I have vector scheme in GDSII format.

There are some viewers for GDSII like KLayout for example, but they don't support raster graphics.

Is there any out-the-box solution? If no, if I know how to parse my vector graphics(for example I can get lines or polygons) can I somehow convert them to some format that some GIS system takes?


texts: ArcGIS Python API FeatureLayer query by buffered point

texts:
I am really running into some problems with the Python API, I am trying to query a feature layer and return all results that fall within a certain distance of a supplied point.

According to the docs here: http://esri.github.io/arcgis-python-api/apidoc/html/arcgis.features.toc.html#feature what I am doing should work but instead it is just returning 1,000 results each time. I know that it should only be returning 1.

pt =  Geometry({"x" : -96.95335386999966, "y" : 32.87209363900047, "spatialReference" : {"wkid" : 4326}})

print (pt.is_valid) # True
print (pt.type) # POINT
print (pt) # True

results = policies.query(geometry_filter=pt, distance=10, units="esriSRUnit_Foot")

print(len(results)) # 1000

If someone could point me in the right direction to get that to work that would be amazing.

The query works with the REST API so I know it is not a problem with the layer or anything like that. It's just something with the Python API


texts: Manually selecting colors for different polygons within layer in ArcMap?

texts: There is a polygon layer with several polygons. They are displayed with some colour which can be changed simultaneously for all of them. Is there a way to manually select several polygons (without referring to their attributes in the table) and to pick for them some other colour (that is, by manually selecting polygons and groups of polygons one by one to assign them some different colours to display)?


texts: What is the process by which ONS creates "Generalised" and "Ultra Generalised" Clipped Boundary shapefiles?

texts:
The ONS (Office for National Statistics) in the UK makes available a large number of shapefiles based on geocoded statistical data sets and admin/demographic boundary shapefiles. These are deposited at http://geoportal.statistics.gov.uk. I'm aware that there are a range of ways to generalise and simplify polygons within shapefiles, but what I'd like to know is how, in particular, this is being done by ONS. I've experimented with several similar processes: using GRASS v.generalise / v.simplify and ogr2ogr -simplify but can't seem to replicate anything quite the same as this. In particular, I'm curious how this can be done to generate such small shapefiles. So, for example, on the "Combined Authorities (April 2015) Full Extent Boundaries in England" dataset, the shapefile is about 1.4MB, where as for "Combined Authorities (April 2015) Ultra Generalised Clipped Boundaries in England" the size is around 29kb.

The closest thing I could find to an explanation was:

  Full resolution - extent of the realm (usually this is the Mean Low Water mark but in some cases boundaries extend beyond this to include off shore islands);
  Full resolution - clipped to the coastline (Mean High Water mark);
  Generalised (20m) - clipped to the coastline (Mean High Water mark);
  Super generalised (200m) - clipped to the coastline (Mean High Water mark) and
  Ultra generalised (500m) - clipped to the coastline (Mean High Water mark).

Perhaps there is an obvious mechanism for doing this that I'm missing because I haven't been searching in ArcGIS. I'd be glad for any suggestions as to the teechnical specifics of an industry "standard" approach to generalisation that might be in use here by ONS.


texts: How to add new base layer to QGIS OpenLayers from Freemap server

texts:
I have a problem with adding a new layer to QGIS 2.6 OpenLayers and I need help. It is a source of maps www.freemap.sk

Examples of working with the server are:

http://wiki.freemap.sk/LayerAllInOne
http://wiki.freemap.sk/TileDistributionTest

Now I created a file freemap.html in the directory:

~/.qgis2/python/plugins/openlayers_plugin/weblayers/html

&lthtml xmlns="http://www.w3.org/1999/xhtml">
  &lthead>
    &lttitle>OpenLayers FreeMap Layer
    &ltlink rel="stylesheet" href="qgis.css" type="text/css">
    &ltscript src="OpenLayers.js">
    &ltscript src="OlOverviewMarker.js">
    &ltscript type="text/javascript">
        var map;
        var loadEnd;
        var oloMarker; // OpenLayers Overview Marker
        function init() {
            map = new OpenLayers.Map('map', {
                theme: null,
                controls: [],
//                projection: new OpenLayers.Projection("EPSG:4326"), new OpenLayers.Projection("EPSG:3857"),
                units: "m",
                maxResolution: 156543.0339,
                maxExtent: new OpenLayers.Bounds(-20037508.34, -20037508.34, 20037508.34, 20037508.34)
            });

            loadEnd = false;
            function layerLoadStart(event)
            {
              loadEnd = false;
            }
            function layerLoadEnd(event)
            {
              loadEnd = true;
            }

            var osm = new OpenLayers.Layer.OSM(
              "FreeMapTuristika",
              "http://a.freemap.sk/T/${z}/${x}/${y}.png",
              {
                numZoomLevels: 16,
                eventListeners: {
                  "loadstart": layerLoadStart,
                  "loadend": layerLoadEnd
                }
              }
            );
            map.addLayer(osm);
            map.addControl(new OpenLayers.Control.Attribution());
            map.setCenter(new OpenLayers.LonLat(19.41,48.52).transform(
                                               new OpenLayers.Projection("EPSG:4326"), new OpenLayers.Projection("EPSG:3857")
                                                ), 11 );
//           map.setCenter(new OpenLayers.LonLat(4485,2841),13);
            oloMarker = new OlOverviewMarker(map, getPathUpper(document.URL) + '/x.png');
        }
    &lt/script>
  &lt/head>
  &ltbody onload="init()">
    &ltdiv id="map">&lt/div>
  &lt/body>
&lt/html>

The problem is that QGIS does not open map:

While the browser opens fine:

When I watched the network traffic, I found that QGIS was sending erroneous requests to the server:

GET /T/1/0/-1.png HTTP/1.1
User-Agent: Mozilla/5.0 (X11; Linux x86_64) AppleWebKit/537.21 (KHTML, like Gecko) QGIS2 Safari/537.21
Accept: */*
Connection: Keep-Alive
Accept-Encoding: gzip
Accept-Language: en-US,*
Host: a.freemap.sk

HTTP/1.1 302 Found
Server: nginx/1.6.2
Date: Wed, 05 Nov 2014 09:56:41 GMT
Content-Type: text/html; charset=iso-8859-1
Location: http://a.freemap.sk/core_allinone.php?/T/1/0/-1.png
Cache-Control: max-age=0
Expires: Wed, 05 Nov 2014 09:56:41 GMT
X-Varnish: 10693744
Age: 0
Via: 1.1 varnish-v4
Content-Length: 235
Connection: keep-alive

&lt!DOCTYPE HTML PUBLIC "-//IETF//DTD HTML 2.0//EN">
&lthtml>&lthead>
&lttitle>302 Found&lt/title>
&lt/head>&ltbody>
&lth1>Found&lt/h1>
&ltp>The document has moved &lta href="http://a.freemap.sk/core_allinone.php?/T/1/0/-1.png">here&lt/a>.&lt/p>
&lt/body>&lt/html>
GET /core_allinone.php?/T/1/0/-1.png HTTP/1.1
User-Agent: Mozilla/5.0 (X11; Linux x86_64) AppleWebKit/537.21 (KHTML, like Gecko) QGIS2 Safari/537.21
Accept: */*
Connection: Keep-Alive
Accept-Encoding: gzip
Accept-Language: en-US,*
Host: a.freemap.sk

HTTP/1.1 500 Internal Server Error
Server: nginx/1.6.2
Date: Wed, 05 Nov 2014 09:56:41 GMT
Content-Type: text/html; charset=UTF-8
Content-Length: 0
X-Powered-By: PHP/5.6.2
Cache-Control: max-age=172800, public
Vary: User-Agent
X-Varnish: 10693746
Age: 0
Via: 1.1 varnish-v4
Connection: keep-alive


texts: Adding RGB values from a .txt file to a .las file using laspy

texts:
I have a text file that has XYZ and RGB values in it. To convert it to .las format i use

lfile = laspy.file.File(r"D:\las\tolas.las", mode='w',header=LasHeader(x_scale=0.01,y_scale=0.01,z_scale=0.01))
header = laspy.header.Header(point_format=2)
lfile.x = x/0.01
lfile.y = y/0.01
lfile.z = z/0.01
lfile.set_red = r
lfile.set_green = g
lfile.set_blue = b
lfile.close()

But the RGB values are not being set. Where did I go wrong. The las format is 1.2


texts: How do I use a angular router link inside a Mapbox popupmarker?

texts:
I trying to use an angular routerlink with mapbox popup. The popup is displayed when a click on the map. The router link fails whenever I click it. It shows no errors it just fails to do anything.
I have tried to use a angular click on the a tag so that I use angular router to navigate to a specific route but still fails.
Here is my code:
 var markerPopup =`<img src="`+this.imageBaseUrl+location.imagem+`" height="120px" class="card-img-top rounded" alt="...">

 4500 (check ins)

`+location.name+`

 `+location.bairro+`, Luanda

     Visite Page

;`

 var popup = new mapboxgl.Popup().setHTML(markerPopup);
popup.addClassName('bg-warning');
    var marker = new mapboxgl.Marker({
      color:"#4CAF50",
      draggable: false
      })
      .setLngLat([location.lng,location.lat])
      .setPopup(popup)
      .addTo(this.map);


texts: Where is Request type drop down in GeoServer authentication settings page?

texts:
I am reading the GeoServer tutorial about configuring digest authentication. The tutorial says that on the authentication page, the Filter chains panel, there is a Request type drop down, which I cannot find anywhere in my geoserver instance (2.8.3 from OSGeoLive 10).

According to documentation, the Filters Chain panel should look like this:

But mine looks differently like this:

Am I looking at the wrong panel, is the OSGeo version GeoServer different, or did something else go wrong?

By the way, the address on the browser is: http://localhost/geoserver/web/?wicket:interface=:24:1:::


texts: Eliminate line sections in a merge to use Qchainage plug in

texts: in QGis I merged several layers, each containing sections of a road, to show a longer route on a map. Now when I want to use the qchainage plug in the plug in starts the count over again in each 'section' of the road before i merged them. How can I 'smooth' the route so that the plug in sees it as a single road?


texts: How to calculate exact area of a loaded jpg map in to arcmap?

texts:
I calculated the polygon area through arcmap.but i'm not sure that the area is not actual,may be scale is not right.

so how to calculate the exact area of a jpg image?


texts: Georeference .txt file to SpatialGridDataFrame

texts:
I have a matrix file in a .txt format and I want to georeference this matrix using latitude and longitude.
How could I do this in R, so that I can transform my .txt file into a SpatialGridDataFrame object and then import it in GRASS GIS?


texts: Classify image using segmentation analysis using QGIS and OrfeoToolbox

texts:
I am working on a project where we have made a RGB orthomosaic using a drone which I want to classify using QGIS. I am not a fan of unsupervised classification, and I started exploring the possibilities of object-based classification. I tried the trial of eCognition; just a wonderful workflow.

Sidenote:
the rasters pixels are 1.5 cm x 1.5 cm. The total area is around 3 hectares. So a lot of data. I exported the raster to a raster pixel width/height of 5.0 x 5.0 cm,

Target:
classify the image using object-based analysis of the (I suspect) the Orfeo toolbox

Workflow:
segmentation (meanshift) OTB:

processing.runalg("otb:segmentationmeanshift","G:/Veld1/veld1_5cm.tif",0,5,30,0.1,100,100,0,0,"G:/Veld1/veld1_5cm.tif",True,True,100,0.1,"layer","ID",512,1,"None","G:/Veld1/SEGMENTATION/segm_5cm_tile512,size1000.shp")

This works fine.
I want to classify the image using the polygons from the segmentation mean-shift; just like in eCognition.

Maybe make a sort of ruleset like an ROI in the semi-automatic classify plugin for QGIS? Steps:

add fields and rule polygons in certain classes
classifier

Does anyone know a good workflow? To give a feeling for the data; I uploaded an image with the image and the segmentation output.


texts:
The following PostGIS query involving ST_Relate gives FFFFFF102 as its result:
SELECT ST_Relate(
  ST_GeomFromText('POLYGON EMPTY'),
  ST_GeomFromText('MULTILINESTRING((0 0,0 1),(0 1,1 1,1 0,0 0))'),
  1);

I would have expected it to give FFFFFF1F2 instead. My reasoning for expecting that result:
I'm passing in 1 as the boundary node rule parameter, which corresponds to the OGC/MOD2 rule. My understanding of the MOD2 rule is that a node will be considered as part of the boundary of a collection geometry if and only if it appears in the boundary of child geometries an odd number of times. In this case, 0 0 is part of the boundary of both children in the MULTILINESTRING collection, as is 0 1 (i.e. they appear an even number of times, so shouldn't be part of the boundary).
The 8th digit in the intersection matrix code indicates the dimension of the intersection of the exterior of the first geometry and the interior of the second geometry. The first geometry is empty, so its exterior is the entire plane. And assuming my understanding above is correct, then the boundary of the second geometry should be empty. So therefore there should be no intersection (i.e. the 8th digit should be F).
This is using PostGIS version POSTGIS="3.1.0 5e2af69". I have seen the same behaviour in GEOS version 3.8.0 as well.
Can anyone help me understand why I'm getting FFFFFF102 rather than FFFFFF1F2?


texts:
I'm testing out QGIS for the first time, and I'm trying to find a way to highlight index contours.

I'd like to be able to do so via the layer Symbology, or if necessary by separating out the index contours I want into a separate layer.

I've tried building a query to select every 10th contour using the generate_series expression, but I'm not well versed in crafting these and thus far I haven't been successful.


texts:
I am looking for some help on proper syntax in ArcGIS Pro/ Python 3. I am trying to read data from a feature class and write to another feature class. The attribute fields are different so I am reading from source and assigning it to output. I get this error at insertCursor.insertRow(newRow):

sequence size must match size of the row

def LoadData(fSource, fTarget):

    # setup search cursor
    fields = ['*']

    fldShape = "SHAPE@"

    # target field mapping
    fldCountry = 'Country'
    fldCommodity = "Product"
    constructionYear = "CNSTR_YEAR"

    try:
        with arcpy.da.InsertCursor(fTarget, ['*']) as insertCursor:
            with arcpy.da.SearchCursor(fSource, "*") as cursor:
                for row in cursor:
                    country = row[cursor.fields.index(fldCountry)]
                    geometry = row[cursor.fields.index(fldShape)]
                    commodity = row[cursor.fields.index(fldCommodity)]
                    constructionYear = row[cursor.fields.index(constructionYear)]

                    newRow = [country, geometry, commodity, constructionYear]

                    insertCursor.insertRow(newRow)

    except Exception:
            print(e.args[0])


texts: Landsat L1GS product processing to change detection analysis ready state

texts:
At high latitude area of Svalbard Norway, one can find Landsat 8 data processed to L1TP level. With previous Landsat missions, this is not the case. Landsat 7 data are processed to L1Gt level and Landsat 1 to 5 data are processed to L1GS level.

https://www.usgs.gov/land-resources/nli/landsat/landsat-levels-processing

Target would be to get these products to state, where time series for change detection analysis could be composed. Spatial inaccuracies could be mitigated with automated co-registration of images. I thought this could be done with AROSICS software, using the L8 OLI L1TP -data as reference images.

https://www.mdpi.com/2072-4292/9/7/676

https://pypi.org/project/arosics/

When the imagery reliably aligns, I have understood that L1Gt data could be used for change detection analysis. This leaves the question of L1Gs data. Can one then employ topographic correction with a Digital Elevation Model (DEM) to correct for relief displacement in the new co-registered L1Gs images? Or does the previous geometric corrections done to L1Gs data obstruct this process?


texts: IDRISI .vec file translator to ESRI shapefile format?

texts:
I've recently got hold of almost 20 years old IDRISI .vec files which I'd like to convert to ESRI shapefile. I've been struggling to find any translator for this file format so far. Does anybody know of a way to convert this type of file? Any help much appreciated!

I'm going to post cross post on the http://forums.arcgis.com/


texts: Do you know where I can find a geoserver docker image that also contains the latest GDAL extension/version?


texts: Saving a QGIS project using the API loses mapcanvas (and other) XML nodes

texts:
I'm using the following python script to open a QGIS project, modify some things, and save it under a different name:

from PyQt4.QtCore import QFileInfo
from qgis.core import *

def main():
    QgsApplication([], True)
    QgsApplication.setPrefixPath('/usr/share/qgis', False) #If this is set to true qgis can't find resources..
    QgsApplication.initQgis()

    proj=QgsProject.instance()
    proj.read(QFileInfo("source.qgs"))

    #do other stuff to the project

    proj.write(QFileInfo("target.qgs"))
    QgsApplication.exitQgis()

if __name__ == "__main__":
    main()

While this seems to work fine at first, the XML nodes <mapcanvas>,<layer-tree-legend> and a few others are omitted in the now smaller target file. I want to serve the target file via QGIS Server, but QGIS Server doesn't pick up on all the available layers due to these missing nodes in the XML file. The GetProjectSettings request doesn't show any of the available layers even though they're clearly still in the project.

Opening the target file in QGIS Desktop and saving it again works fine, and then the missing nodes are added again to the qgs file.

Does anyone know why this information is lost when saving a project, and what I could do to preserve or regenerate it?

Edit: added setPrefixPath


texts: Extracting elements plotted on Google Maps API

texts:
I have the website about the planet occultations by Moon.
https://eco.mtk.nao.ac.jp/cgi-bin/koyomi/occulx_p_en.cgi
There are "Occultation maps", which are very nice:

I tried to use the Web Scraper for Chrome, but it looks like it grabs the stuff to Excel only.
I would like to extract this occultation map to a KML file or similar. Is it possible?


texts: Can't Add Mapbox hosted DigitalGlobe Basemap in QGIS 3.2

texts:
I'm having trouble viewing a Digital Globe base map in QGIS. I'm following Mapbox's instructions from their API
https://mapsapidocs.digitalglobe.com/v2.0/docs/tiles (bottom section that describes connecting with WMTS). In QGIS 3.2.0, I add new WMS/WMTS layer (“Layer” -> “Add Layer” -> “Add WMS/WMTS Layer…”), and I add the new layer with the URL with my own Mapbox public token (i.e. https://api.mapbox.com/styles/v1/digitalglobe/cinvynyut001db4m6xwd5cz1f/wmts?access_token=MY_TOKEN) and I hit connect, but I keep getting the error message:

Failed to download capabilities:
Download of capabilities failed: Error transferring https:
https://api.mapbox.com/styles/v1/digitalglobe/cinvynyut001db4m6xwd5cz1f/wmts?access_token=pk.ey...bunch of token characters...=WMS&REQUEST=GetCapabilities

-server replied: Not Found

My token starts with pk.ey and a bunch of random characters so I'm sure its the right one.  I'm using the mapbox Pay-as-you-go plan, and I haven't added any money to the account. Do I need to pay for this service? I thought Mapbox offered some limited views of Digital Globe images for free.  Or maybe it's a technical problem.


texts: Converting large value float value to integer using ArcGIS Desktop

texts:
I have a raster with 32 bit unsigned values that I want to convert to a polygon shapefile.  I ran the raster calculator in ArcGIS to remove the decimals by multiplying by a constant and then running raster calculator again but this time converting to integer.
The output integer values do not match what the values were before converting to an integer.  I'm not sure what is happening but by doing it this way will not work especially if its coming up with totally different values.
Clarification
I am running ArcGIS 10.4.1
For example - the raster holds values from 7881.530000 to 53613.300000
Multiplying the raster values by 1000000 because of the 6 decimal places
Converting the raster to integer results in completely different values
Real Example:
Pixel value: 39961.472656
Multiply by Constant: 39961473024.000000 (this number is already different)
Converted to Integer: 1306767360
If the pixel values are much smaller, the results are fine
Example:
Pixel Value: 0.087640
Multiply by Constant: 87640.140625
Converted to Integer: 87640
The end result is to convert to a polygon shapefile where I will add a field and calculate the value to be based on the Integer Value where I divide it by the constant to get back to the original pixel value.


texts: ogr2ogr bringing up "memory allocation error"

texts:
I've got a shapefile which contains, at the moment, roughly 1500 polygons. The polygons all display fine using whatever GIS I throw it in (ArcGIS, QGIS).

However, I'm trying to use ogr2ogr to convert it to a table on SQL Server and it gives the following error:

  ERROR 1: INSERT command for new feature failed. [Microsoft][SQL Server
  Native Client 11.0]Memory allocation failure ERROR 1: Unable to write
  feature 8 from layer PARCEL_MP. ERROR 1: Terminating translation
  prematurely after failed translation of layer PARCEL_MP (use
  -skipfailures to skip errors)

Of course, I can use -skipfailures but when I do, it manages to export roughly 120 of the 1500 records (it differs each time - sometimes it's 115, sometimes 120, a bit more or less).

Obviously what I'm shooting for is 100%. The server is SQL Server 2014.

My code (running this in OSGeo4W) is:

ogr2ogr -overwrite -f MSSQLSpatial "MSSQL:server=MY_SERVER;database=MY_DATABASE;trusted_connection=yes" "PARCEL_MP.shp"

EDIT - also, the first feature that throws an error changes each time. The shapes aren't complicated and are quite small, with a small number of vertices.

Also, just tried using ogr2ogr to convert it to another shapefile which has worked fine. Think it must be in my server. Code I used was:

ogr2ogr -overwrite -f "ESRI Shapefile" "PARCEL_MP_2.shp" "PARCEL_MP.shp" -skipfailures


texts: Problem understanding response from MapServer WCS

texts:
I can't understand what I am doing wrong with my WCS request.

When submitting this request from a browser the result is an output called  mapserv.

http://localhost:8080/cgi-bin/mapserv?map=/usr/local/mapserver_work/wcs1.map&SERVICE=WCS&VERSION=1.0.0&REQUEST=GetCoverage&coverage=jeo&CRS=EPSG:4326&BBOX=25.670,35.813,44.824,42.106&WIDTH=1600&HEIGHT=1066&FORMAT=GTiff&


texts: query / filter on attribute tables with external tables joined

texts:
i tried to do some queries on joiend table fields ...

i've got/found some answers concerning version 1.8 - which result in "use select + export to create a combined dataset"
so i just did some tests with Layer / query in qgis 2.2 valmiera and 2.3-master(build from 2014-05-02)) and ended up with the same problems:

while

"layer properties / style / column expression" tool is working fine
based on "external" fields - even if i build complex expressions
and "table - select faetures by expressionn" doesn't seem to have any problems with where-clauses

"layer / query" returnes strange results with my demo-data; i.e.:

no problems with expressions unsing fields from the layer's attribute table
but - as soon as i address "external" fields
but - as soon as i address "external" fields

the test feature returns an error:
a simple query like "zsp_experimente_xls_WFL01_AREA" >= 1.52"
(where 1.52 is valid value out of the sample-list) returns "OGR 3 error 1: type mismatch or improper type of arguments to >= operator"
but "zsp_experimente_xls_WOHNFL01_AREA" LIKE "%" returns all records when testing the query ("the where clause teturned 1364 rows")
the confusing thing about that: the the sample tool privides correct results ! ! !
((which means to me that the select part of the statement works fine))

so finally i am still argueing that one of the essential benefits of the join-concept is, that you do not have to mix up geometry and attribute data -> which was a really great idea !!!!!!!
to select data by expression and export geometry and "added" additional data really undermines this concept by duplicating data !!!

Does anybody have any idea how to deal with that problem ?


texts: Finding distance between multiple points along a line

texts:
I'm using Arcmap 10.4 full license, all extensions, including Network Analyst.

I have a GTFS (transport) dataset from which I am trying to calculate distances between each point (alighting bus stops) to every other point (end point bus stops) along a specific bus route (quite a few of them).

I am a bit of a newby with network analyst (only taught it to myself 2 days ago, so please excuse any silly questions).

I have tried various options including spatially joining the bus stops to the bus route, splitting the bus route by the bus stops (makes all sorts of additional segments, which I can't understand why, because I dissolved the route in the first place), and doing route analysis, closest neighbour and cost distance in network analyst.

I can't figure out how to do this on multiple routes, and when I calculate it on one route (by selecting the route, and using a selection layer, ditto for the points), it doesn't follow the bus route but appears to be a linear euclidean distance measure, which I don't want, I need the actual distance along the bus route.

I've added all the data in using the gtfs to network scripting, but suspect I have gone wrong somewhere.

Is there perhaps a way to batch or do this for the whole dataset, rather than individually by route and bus stop?


texts: Convert point shapefile into a polygon which may have holes in the middle using R

texts:
I have multiple point shapefiles with 20,000 - 100,000 points and I would like to create an outline/polygon from them using R. There may be holes in the middle of the polygon which I want outlined as well. I have posted pictures of how I want this to look. I have searched a lot for an answer to this using R and haven't found one that has worked. Many other examples I see don't follow the outline as closely as I would like. The polygon that I did get to work came from ArcGIS using the aggregate points (Cartography) tool with a 10 meter distance, but I want to use R.
    ye <- read.csv(ye1,header = F, stringsAsFactors = F)
    ye$gid <- 1
    new_names <- c("easting", "northing","lat","long","yield","moisture","swath","distance","flow","interval","agl","transect","gpstime","utmzone","rmcode","gid")
    names(ye) <- new_names

     ye_shp <- SpatialPointsDataFrame(ye[,1:2],data = ye,proj4string = CRS("+proj=utm +zone=14 +datum=NAD83 +units=m +no_defs +ellps=GRS80"))

This is an agricultural field with 70,000+ points.

This is what I want the polygon to look like.


texts: How to programmatically delete all inner rings of a polygon layer?

texts: I have many polygon layers, each layer containing many polygons with inner rings. I know that you can delete a ring manually with QGIS (Edit -> Delete ring). I need the python equivalent to automatically process all my polygon layers.


texts: Is there an equivalent of GeoDjango for Flask or other Python web frameworks?

texts: I have experimented with GeoDjango, but I think nobody would be surprised if I said there's quite a steep learning curve. I prefer Flask, which I also know much more than Django anyway, but I couldn't find anything like GeoDjango for Flask. Is there any such package around that I might have missed?


texts: Adding Nodes to Linestrings at Overlapping Points with Polygon Layer in Spatialite

texts:
I have two Spatialite tables.

A table called 'lines' which contains a PK column called OGC_FID and a LINESTRING column called GEOM.
A table called polys which contains a PK column called OGC_FID and a POLYGON column called GEOM.

How do I update the 'lines' table GEOM so that the intersection points with any POLYGONS in the polys table are added to the LINESTRING geometries?

(Preferably using SQL queries only.  Temp tables are fine)  (Postgis answer acceptable)


texts: Clip_management to raster extent outputs larger raster than raster extent

texts:
Using Python and ArcMap 10.0 SP5, I'm trying to reproject and clip a large input raster to exactly the same grid as a smaller target raster (same extent, pixel sizes and projection).  The problem is that the script produces an output raster that has 1 more row than the target raster (e.g. 1324 rows instead of 1323 -- columns match perfectly, pixels appear to be on the same grid).  There is an extra row of pixels that overhangs the target image - everything else appears to be okay.

The following is my code for projecting and clipping to a target raster.  What's needed to ensure that the output extent is exactly the same as in target raster?  I'm new to ArcMap so it could be something obvious.

import arcpy
import tempfile
import os
import shutil

sourceImage = "C:/data/WorldClim_alt/WorldClim_Alt.tif"
inputImage = "E:/subset/nCom_2000/nCom2000101to107.tif"
outputFilename = "E:/subset/test/nComTest_clipped.tif"

# Snap to source image and source image extent
arcpy.env.extent = sourceImage
arcpy.env.snapRaster = sourceImage

# Get cell size
cellx = arcpy.GetRasterProperties_management(sourceImage, "CELLSIZEX")

# Get temporary directory and filename
tempDir = tempfile.mkdtemp()
tempFilename = 'E:/subset/test/' + basename(outputFilename) +'_proj.tif'

# Reproject
arcpy.ProjectRaster_management(inputImage,tempFilename,sourceImage, "#", cellx)

#
# Clip to source file
#

# Get envelope for clipping
top = arcpy.GetRasterProperties_management(sourceImage, "TOP")
left = arcpy.GetRasterProperties_management(sourceImage, "LEFT")
bottom = arcpy.GetRasterProperties_management(sourceImage, "BOTTOM")
right = arcpy.GetRasterProperties_management(sourceImage, "RIGHT")
envelope = str(left) + ' ' + str(bottom) + ' '+ str(right) + ' ' + str(top)

# Clip
arcpy.Clip_management(tempFilename, envelope, outputFilename, sourceImage,"#","NONE")

Edit in Response to Paul's comment on Clip help:

  From the help on Clip, it mentions

    If the clip extent specified is not aligned with the input raster
    dataset, this tool makes sure that the proper alignment is used. This
    may cause the output to have a slightly different extent than
    specified in the tool.

I think that the input to clip (the output from ProjectRaster) may be aligned with the in_template_dataset.  The snapRaster environment variable was set to the in_template_dataset before ProjectRaster was called.  Also, the output cell size of ProjectRaster was set to the in_template_dataset cell size.  Visual inspection of the grids seems to indicate that they are aligned (except for an extra row and column in the ProjectRaster output).

The only effect of Clip seems to be the removal of a single column.


texts: Rasterize causes black background, need white to match BioClim data

texts:
I am using the rasterize tool to convert a polygon layer into a .asc or .bil raster layer. As it is a UK study, there is a lot of sea that when rasterized, turns black. I need it to be white to match the BioClim data.

I've tried to convert using the rasterize tool, which led to the initial problem. I have then used the reproject tool to attempt to convert the data with null values of 9999, rather than 0. However, still no difference. It isn't simply a case of inverting the colour ramp either, as the actual file needs to store the data with a white background.


texts: Should I start a new project using the Silverlight API?

texts:
I'm about to start a new mapping project for a client, and am considering using the ArcGIS Server Silverlight API, since it gives a lot of the required functionality out-of-the-box (especially with the Silverlight Viewer).

A colleague suggested this is a bad idea, since "Silverlight is dead technology". He says that most new internal Microsoft projects are moving to HTML5 instead of Silverlight.

This is a small project so I don't have the luxury of developing everything from scratch, which the JavaScript API entails. Using the Silverlight Viewer I can have the core functionality running almost immediately.

Does anyone have any opinions/inside information on whether Silverlight is a good idea for a new ArcGIS Server website?


texts: Sharing GeoPackage for collaborative work (web, ftp, cloud...)

texts:
I would like to share a GeoPackage as part of a GIS project (QGIS) with one of my customer. Only the customer will edit the GPKG.

I could use a cloud service like Google, Dropbox or pCloud. But that assumes that the other one has got a cloud account too. Could I also use an HTTP or FTP-Server? How?

I prefer to use GeoPackage instead of PostGIS because it's easier to share and backup.


texts: Why do valid polygons repeat the same start and end point?

texts:
In the GIS world, for most modern formats and Standards (e.g., Shapefiles, WKB/WKT, GML, KML, etc.), valid polygons need to have closed linear rings, which is a linestring of coordinates where the first point is a repeat of the last point. For example, a triangle requires four points (not three).

POLYGON ((10 20, 30 60, 50 20, 10 20))

Who started this convention, and why? Is it some legacy of pre-Shapefile storage? (like how MS Windows still uses 2-byte CR+LF newlines?) Other non-GIS standards (e.g. SVG) don't require this repetition to encode polygons.


texts: Stop pending tiles requests in Openlayers 3!

texts:
I have created a Main-map in OpenLayers 3 with a single layer.

I use 'setUrl' to change the layers sourceUrl, when I want to change the map:

 this.mainMapLayer = new ol.layer.Tile({ source: new ol.source.XYZ({ url: this.url })});
 this.addLayer(this.mainMapLayer);

changeMainLayerSource(id): void {
    this.url = "http://MyTileServer" + id + "/{z}/{x}/{y}.jpg";
    this.mainMapLayer.getSource().setUrl(this.url);
}

It works, but it doesn't abort/stop pending tiles when the url is changed.
If I change the sourceURL, while some tiles are still pending, the map disappear, until all the tiles are done loading.

How do I avoid this?

Is there any way, to stop/abort pending tiles?


texts:
I converted the QGIS polygon to QGIS line by geometry tools which originally provided in the vector tools tab.
the result is ok, but when I export it to DXF and open it with the AutoCAD, the lines are somehow hidden or something else that results in some errors.
so I could not use the file...
What should I do to solve it?


texts: Parallelize Extracting Large Raster Values (one band) to Points

texts:
I have sf data with single point geometries (points) and a 23,000 x 18,000 elevation RasterLayer with a single band, loaded from a geotiff (elev). They have the same CRS. I want to extract the elevation values to the points.

When I try running extract(elev, points), since my raster is very large I get:

  Error: cannot allocate vector of size 3.1 Gb

I followed advice here R - multicore approach to extract raster values using spatial points and tried to parallelize, but I only have a single raster, so chances are the work won't get divided. Still I try the following:

library(snowfall)
sfInit(parallel=TRUE, cpus=parallel:::detectCores()-1)
sfLibrary(raster)
sfLibrary(sp)

# Run parallelized 'extract' function and stop cluster
e.df <- sfSapply(elev, raster::extract, y=points)
sfStop()

I get the same result:

  Error: cannot allocate vector of size 3.1 Gb

All the examples I've found so far, such as the above, store multiple rasters in a list. How do I parallelize extraction from a single raster like this? Or is there a way to make use of sparse matrices and memory management as noted here?

I'm on Windows 10 running RStudio 1.1.419 (64 bit). 8GB RAM. Intel Core i5-6200U CPU @ 2.3GHz 2.4GHz


texts: Finding difference between two vector layers using QGIS

texts:
I have two vector layer land cover for 2016 and 2018. I want to get  result layer that will show changes in land cover e.g. the forest has changed into a meadow.
What option in QGIS can I use to do that?


texts: Creating NDVI from multispectral imagery?

texts: I have LISS IV image data in three bands (band4, band3 and band2). I wanted to know the required steps to obtain the NDVI image from the given data using ERDAS IMAGINE 9.2.


texts: Intersection between two 3D planar rectangles

texts:
I am doing research in roof reconstruction from different data sources.

After previous steps I receive a (plane) 3D rectangle by roof-part, stored with a 2D polygon, a 3D centroid point and a normal vector (to the plane) in PostgreSQL.
As the next step I'd like to cut the (overlapping) planes and calculate the intersection lines (e.g. ridge lines) so they fit together. It would be possible to calculate a 3D polygon / 3D surface from the stored roof-part data.

How should I calculate this? ST_3DIntersection from PostGIS seems propably not to be the right way... Maybe the Intersection3D function from SFCGAL? Unfortunately I've never used this lib.


texts: Is there a CASE Tool for modelling spatial databases?

texts: Does anyone knows if there is a CASE Tool that would be possible for modelling spatial tables? I need to create a model and export it to SQL or a database directly. I'm using PostGIS.


texts: Query about fid in QField causing difficulty digitising

texts:
I have been using QField for a while and  recently got a  message under fid saying NOT NULL, Unique and I can no longer digitise anything on QField. I read that you can make a change under layer properties but the fid field is not showing in any of my layers in QField to make that change.
I am not that experienced so simple directions on how to fix it would be great.


texts: Whats the status for OSM in print composer QGIS 2.0?

texts:
I have a Stamen Toner OSM background on a report am doing. It seems impossible to use "Print" or "Save as image" functions since the OSM is being projected wrong with 10 km or so.. Right now I'm taking screenshots from there as a work around but it feels kind of wrong?
Screenshot with correct OSM projection:

Save from print composer with wrong OSM projection: (Buffers should be located at the lower right corner


texts: Explaining laser scanning distortions on water surfaces?

texts:
I have received tif Laserscan images interpolated from Ascii data. I did not process the data so cannot say which processes or tools were used to interpolate the data. There are strange distortions on water surfaces which I cannot understand and seek an explanation.

Image 1 shows the Surface model. There are smooth elevations around vegetation. (Maybe due to the Spline interpolation technique?)

Image 2 shows the same areas but in the Terrain model. Here it seems the interpolation technique was different as the smoothing is missing.

Image 3 shows the las data.

Are the raised areas due to reflections of the trees on the water?

If so, why are reflections registered and not somehow eliminated?


texts: Why Google Maps is not popular with tech stack like Leaflet, OpenLayers, nebula.gl?

texts:
Google Maps is really good when it comes to places data and navigation and still is one of the most popular options for normal uses but whenever I try to look up a GIS tech stack they uses Mapbox, OSM or something else for mapping but not Google Maps

I want to know why Google Maps is not popular with these frameworks?

Few to name

Leafletjs
Openlayers
Nebula.gl
Kepler.gl


texts: How to get a slippy map from OziExplorer .map and .png tiles?

texts:
I have a bunch of topographic maps in the format of .map files and 4000x4000px png files.I would like to view them as a slippy map and be able to select an area and export that as one big image. I am wondering what my options are. QGis would be nice, or Mobile Atlas Creator, which supports local tiles as source for the following formats: OSMTracker, AndNav, BackCountry Navigator, Maverick and OSMAND.

It would be nice if at some time I could convert them to something more standard, like 256px tiles under a coordinate dir structure, or mbtiles so I could view it with tilemill or their server or map on a stick software.

So what are my options to getting these files to open in some slippy map way and being able to export specific areas as images.

I don't care what formats I have to convert to or which software/servers I have to install.


texts: Exporting from Layer to Shapefile or feature class?

texts:
I am aware that a layer contains style information and a shapefile has vector information, but I am talking about a layer which can be interpreted as shape information.

From my understanding there is some option that lyr file can be some shortcut of other shp and when I upload the lyr it show all the shape information that has in the original one, but now I want to export this lyr so I can use it without harming it, but when I do "right click" then "export data" it's make the shp with field problem like when the field type is a f-code it's make it text and show just the code and not what he present and I want it to be just like the original.


texts: How to fill composer map with canvas data and export composition as png using Python script in QGIS

texts:
I'd like to create a Python script in QGIS that will take some data (shp + tif, in my example), use tamplate Map Composer (from a file) and export created composition to a png image.

With almost no experience in programming (less than basic Python knowledge) I google some code snippets and try to make them work together. I took the Map Composer code from some previously answered question:
Save Print/Map QGIS composer view as PNG/PDF using Python (without changing anything in visible layout)?

I managed to load the data and Map Composer template (with a single map and legend items defined) but my export png image has an empty map frame (no vector/raster data within the frame). The legend item looks fine however.

Any help with making this code work?

from qgis.core import *
import qgis.utils

from PyQt4 import QtCore, QtGui
from qgis import core, gui

# ADD VECTOR LAYER

data_folder = "D:/QGIS/dane/"
granica = "granica_SZ VI_UTM34.shp"
granica_name = granica[0:-4]
granica = data_folder + granica
granica_style = "granica_style.qml"
granica_style = data_folder + granica_style

granica = iface.addVectorLayer(granica, granica_name, "ogr")
granica.loadNamedStyle(granica_style)
if not granica.isValid():
  print "Granica failed to load or already loaded!"

qgis.utils.iface.legendInterface().setLayerVisible(granica, True)

# ADD RASTER LAYER

landsat = "SZ_VI_LC8-190-022_2015-111_LGN00_OLI_TIRS_atm.TIF"
landsat_name = landsat[0:-4]
landsat = data_folder + landsat
landsat_style = "landsat_style.qml"
landsat_style = data_folder + landsat_style

landsat = iface.addRasterLayer(landsat, landsat_name)
landsat.loadNamedStyle(landsat_style)
qgis.utils.iface.legendInterface().setLayerVisible(landsat, True)

iface.zoomFull()

# MOVE RASTER DOWN (change order)

canvas = qgis.utils.iface.mapCanvas()
layers = canvas.layers()

root = QgsProject.instance().layerTreeRoot()

landsat_old = root.findLayer(landsat.id())
landsat_move = landsat_old.clone()
parent = landsat_old.parent()
parent.insertChildNode(2, landsat_move)
parent.removeChildNode(landsat_old)

# USE MAP COMPOSER TEMPLATE TO EXPORT IMAGE

from qgis.gui import QgsMapCanvas, QgsLayerTreeMapCanvasBridge
from PyQt4.QtXml import QDomDocument
from PyQt4.QtGui import QImage
from PyQt4.QtGui import QPainter
from PyQt4.QtCore import QSize

template_path = data_folder + 'template.qpt'
template_file = file(template_path)

# Set output DPI
dpi = 300

canvas = QgsMapCanvas()

template_file = file(template_path)
template_content = template_file.read()
template_file.close()
document = QDomDocument()
document.setContent(template_content)
ms = canvas.mapSettings()
composition = QgsComposition(ms)
composition.loadFromTemplate(document, {})

# You must set the id in the template
map_item = composition.getComposerItemById('map')
map_item.setMapCanvas(canvas)
map_item.zoomToExtent(canvas.extent())
# You must set the id in the template
legend_item = composition.getComposerItemById('legend')
legend_item.updateLegend()
composition.refreshItems()

dpmm = dpi / 25.4
width = int(dpmm * composition.paperWidth())
height = int(dpmm * composition.paperHeight())

# create output image and initialize it
image = QImage(QSize(width, height), QImage.Format_ARGB32)
image.setDotsPerMeterX(dpmm * 1000)
image.setDotsPerMeterY(dpmm * 1000)
image.fill(0)

# render the composition
imagePainter = QPainter(image)
composition.renderPage(imagePainter, 0)
imagePainter.end()

image.save(data_folder + "out3.png", "png")

QgsApplication.exitQgis()


texts: Determining if trees within forest gaps are clustered using R?

texts:
The attached dataset shows approximately 6000 saplings in approximately 50 variable sized forest gaps.  I am interested in learning how these saplings are growing within their respective gaps (i.e. clustered, random, dispersed).  As you know, a traditional approach would be to run Global Moran's I.  However, aggregations of trees within aggregations of gaps seems to be an inappropriate use of Moran's I.  I ran some test statistics with Moran's I using a threshold distance of 50 meters, which produced nonsensical results (i.e. p-value = 0.0000000...).  The interaction among the gap aggregations are likely producing these results.  I have considered creating a script to loop through individual canopy gaps and determine the clustering within each gap, although displaying these results to the public would be problematic.

What is the best approach to quantify clustering within clusters?


texts: Getting the depth of properties in ArcGIS

texts: Is there a way to get the depth of a parcel in ArcGIS Basic (or QGIS)? Depth meaning the distance from the road (frontage) to the rear of the property. I should also mention that not all parcels are rectangular or square, many parcels are irregular shapes.


texts: Using grid_terrain on a LAScatalog, when the first tile (chunk) has no ground points

texts: I am getting an error of "r: No ground points found. Impossible to compute a DTM." when running grid_terrain on a catalog. I think it is because my first tile is all water. Is there an option in LAScatalog to skip chunks based on class?


texts: Setting template or drawing layout to print map ol5

texts:
I have created a portal which export map and then print. I just wanted to add template or draw layout in canvas for that pdf as shown in the picture below.

Here is my code:

//print map
var dims = {
  a0: [1189, 841],
  a1: [841, 594],
  a2: [594, 420],
  a3: [420, 297],
  a4: [297, 210],
  a5: [210, 148]
};

var exportElement = document.getElementById("export-pdf");
exportElement.addEventListener(
  "click",
  function(e) {
    exportElement.disabled = true;
    document.body.style.cursor = "progress";

    var format = document.getElementById("format").value;
    var resolution = document.getElementById("resolution").value;
    var buttonLabelElement = document.getElementById("button-label");
    var label = buttonLabelElement.innerText;
    var scale = 4000;

    var mapView = map.getView();
    var mapProjection = mapView.getProjection();
    var mapResolutionAtEquator = mapView.getResolution();
    var viewCenter = mapView.getCenter();
    var mapPointResolution = ol.proj.getPointResolution(
      mapProjection,
      mapResolutionAtEquator,
      viewCenter
    );

    var dim = dims[format];
    var width = Math.round((dim[0] * resolution) / 25.4);
    var height = Math.round((dim[1] * resolution) / 25.4);
    var size = map.getSize();
    var extent = map.getView().calculateExtent(size);

    var mapView = map.getView();
    var currZoom = mapView.getZoom();
    var mapProjection = mapView.getProjection();
    var mapResolutionAtEquator = mapView.getResolution();
    var viewCenter = mapView.getCenter();
    var mapResolutionAtEquator = mapView.getResolution();
    var mapPointResolution = ol.proj.getPointResolution(
      mapProjection,
      mapResolutionAtEquator,
      viewCenter
    );
    var mapResolutionFactor = mapResolutionAtEquator / mapPointResolution;
    var pdf_name = dist + "_" + tahsil + "_" + village + ".pdf";

    map.once("rendercomplete", function(event) {
      var canvas = event.context.canvas;
      showText(event.context, "left-top", dist, tahsil, village);
      drawLegend(event.context, "right-bottom");
      drawNorthArrow(event.context, "right-top");

      WriteScaletoCanvas(event);
      var data = canvas.toDataURL("image/jpeg");
      var pdf = new jsPDF("landscape", undefined, format);
      pdf.addImage(data, "JPEG", 0, 0, dim[0], dim[1]);
      pdf.save(pdf_name);
      map.setSize(size);
      // map.getView().fit(extent, size);
      mapView.setZoom(currZoom);
      exportElement.disabled = false;
      document.body.style.cursor = "auto";
    });

    var printPointResolution = (scale * 25.4) / (resolution * 1000); // edit1: corrected
    var printResolutionAtEquator = mapResolutionFactor * printPointResolution;
    var printZoom = mapView.getZoomForResolution(printResolutionAtEquator);

    map.setSize([width, height]);
    // mapView.fit(extent, (map.getSize()));
    mapView.setZoom(printZoom);
  },
  false
);

For example:

How to create this layout pdf?


texts: Error type 000582 when trying to create buffer

texts:
I am trying to create buffers around a point FC that has an acreage field. The point FC has some rows that are null and some that report an acreage of 0.
I'm running into an issue when it comes to the arcpy.sf.CreateBuffers line where I get ERROR 000582. It looks like there is a null value somewhere that is causing a problem, but I can't figure out where. I put part of the error message and the code below. I tried making the cursor code set up a field that has no null or 0 values to avoid any issues. The field going into the buffer is type Float.
Does anyone have any idea what could be causing the issue? I tried running the process manually in Pro and got the same error.
  File "c:\program files\arcgis\pro\Resources\ArcToolbox\scripts\createbuffers_sa.py", line 31, in <module>
    units=get_value(4),

  File "c:\program files\arcgis\pro\Resources\ArcToolbox\scripts\gautils\utilities.py", line 19, in __call__
    return self.func(*args, **kwargs)

  File "c:\program files\arcgis\pro\Resources\ArcToolbox\scripts\geoanalyticssoap.py", line 32, in __init__
    self.token = auth_params.get("token")

AttributeError: 'NoneType' object has no attribute 'get'

The code that is causing the above error is as follows:
curFields = ["ApproxAcres", newField]
with arcpy.da.UpdateCursor(memBF, curFields) as cursor:
    for row in cursor:
        if row[0] is not None:
            if row[0] != 0:
                sqMeter = row[0] * 4046.86
                row[1] = (math.sqrt((sqMeter / math.pi)))
                print("Result is {0}".format(row[1]))
            elif row[0] == 0:
                row[1] = 8
            cursor.updateRow(row)
        if row[0] is None:
            cursor.deleteRow()
del cursor

# Buffer options
bufName = "Brownfields_Buf"
outBufLoc = os.path.join(wsGDB,bufName)
outBuf = arcpy.sfa.CreateBuffers(memBF, outBufLoc, "", newField, "METERS", "DISSOLVE")


texts: What network protocol does QGIS use for WMTS tile retrieval?

texts:
My unpublished android app happily maps Land Information New Zealand vector property boundaries above a backdrop of aerial photographs.

It directly obtains the individual tiles from LINZ using HTTP GET - for example:

http://tiles-a.data-cdn.linz.govt.nz/services;key=eb4828b0ef9144ada0bf7c4518182c33/tiles/v4/set=2/EPSG:3857/10/1006/621.png

Vector map data may be exported to shape files. To give my future users a desktop experience I recommend they load their shape files into QGIS then add aerial photographs as a WMTS layer using the following URL:

https://data.linz.govt.nz/services;key=eb4828b0ef9144ada0bf7c4518182c33/wmts/1.0.0/set/2/WMTSCapabilities.xml

This works well. The problem is that it works too well! [It's assembling the raster backdrop more quickly than my app.]

I investigated further:
(1) Start a WireShark ip trace
(2) Access the above 'Z/X/Y.png' URL from my browser [to prove the trace]
(3) Clear the QGIS WMTS cache.
(4) Open the QGIS project containing the above WMTS layer.
(5) Snap the trace result (image is below)

The only HTTP GET request I see is the one manually issued by my browser.

How is QGIS retrieving the raster backdrop?

Is it inhaling the tiles individually via HTTP GET protocol (in which case Wireshark misses them) or by some other TCP/IP protocol that I have no experience of?


texts: Cell shift when converting raster to array

texts:
I have recently been having some issues converting rasters to numpy arrays. The RasterToNumPyArray tool seems to be adding or subtracting rows or columns depending on which dataset I'm working on.

For instance I have a raster that contains 114 rows, and 160 columns. After converting it to an array it has 115 rows, and 160 column. If I convert this array back to a raster the data appears to be shifted up.

I have tried specifying the number of rows and columns to process in the RasterToNumPyArray tool. The row and column count issue is fixed but the output raster is still shifted up one row.

This issue happens whether the data is in a geographic or projected coordinate system.

Any idea why this would be happening? I am including a python script where I was seeing this issue.

import arcpy
import numpy as np

from arcpy import env
from os import path

env.overwriteOutput = True
env.outputCoordinateSystem = arcpy.SpatialReference(4326)

input_ras = r'PATH_TO_INPUT_RASTER'
input_ras_name = path.basename(input_ras)

ras_desc = arcpy.Describe(input_ras)
rows = ras_desc.height
cols = ras_desc.width
x_cell_size = ras_desc.meanCellWidth
y_cell_size = ras_desc.meanCellHeight

lower_left_corner = arcpy.Point(ras_desc.extent.XMin, ras_desc.extent.YMin)

print input_ras_name
print 'rows:', rows
print 'cols:', cols
print lower_left_corner

arr = arcpy.RasterToNumPyArray(input_ras, lower_left_corner, cols, rows, -9999)
print '\narry based on input raster'
print 'rows:', np.shape(arr)[0]
print 'cols:', np.shape(arr)[1]

print '\nconverting array back to a raster'

output_ras = arcpy.NumPyArrayToRaster(arr, lower_left_corner, input_ras, input_ras)
output_ras.save(r'PATH_TO_OUTPUT_RASTER')

edit: I explored the row and column count a little more. The following screen shot shows the row and column count three different ways - one of which is different. Could this mean anything?

edit 2: The problem does not seem to occur after converting my tif to a GRID and then processing the data. Ideally thought I would not want to perform this intermediate step.


texts: Relative Radiometric Normalization of multisensor Landsat Level 2 images using R

texts:
For my time series analysis, i need to use images of the three Landsat sensors (TM, ETM+ and OLI). I chose the already atmospheric corrected surface reflectance/Level-2 products. Due to the radiometric differences between the sensors, i still have to normalize them (i think) as a pre-processing step.

When i wanted to do this in R using a regression method via the relnorm-Function (landsat-package), i got the following Error Code despite the images having exactly the same values for dimensions and extent:

Error in array(x, c(length(x), 1L), if (!is.null(names(x))) list(names(x),  :
  length of 'dimnames' [1] not equal to array extent

I use the following code to clip and normalize identical bands of the TM-product of different dates. The images don't contain any clouds, thus no mask is used:

library(landsat)
library(lmodel2)
library(RStoolbox)
library(raster)

ma <- raster("D:/Working_Directory/LT05_L1TP_144051_20041216_20161127_01_T1_sr_band1.tif")

to <- raster("D:/Working_Directory/LE07_L1TP_144051_19991109_20170216_01_T1_sr_band1.tif")

shp <- readOGR("D:/Working_Directory/transect.shp")

#Clipping both bands:
clipped_ma <- crop(ma, extent(shp))
clipped_ma <- mask(clipped_ma, shp)
clipped_to <- crop(to, extent(shp))
clipped_to <- mask(clipped_to, shp)

#Trying to use the relnorm-function (no mask, OLS-method):
normalized <- relnorm(master=clipped_ma, tofix=clipped_to, method = "OLS", nperm = 0)

What is wrong with the code or the data(-format), that is causing the aforementioned Error?


texts: How to Configure Geoserver WFS to work with ArcGIS for Desktop?

texts: I am trying to figure out how to get a WFS service working with ArcMap.  I have Geoserver up and running and I am able to add a WFS layer in QGIS but not ArcMap.  I have tried to use the Data Interoperability Tool in ArcMap to connect to my server but there are no layers displayed in select Feature Types dialog.  I am able to connect to a WFS server so I'm assuming that I have overlooked some detail about how ArcMap wants Geoserver to be configured.


texts:
Being in a bilingual (English/French) organisation, I want to add dynamic text showing the date in both languages. I have a line that shows English month/year but since my ArcGIS is in English, by default, the month is written in English.

Is there a way or a tag that sets language?


texts: Redefine ArcGIS Krovak projection for MapInfo

texts:
I have an ArcGIS datum/projection (as below) and trying to write as a line into the Mapinfow.prj. Problem is that the MapInfo documentation for Krovak does not match the ArcGIS description - either in number or name of parameters (lesser number of parameters in MapInfo)
Was hoping for a line for Mapinfow.prj which incorporates the following…..

GEOGCS["GCS_GDA_1994",
DATUM["D_GDA_1994",
SPHEROID["GRS_1980",6378137.0,298.257222101]],
PRIMEM["Greenwich",0.0],
UNIT["Degree",0.0174532925199433]],
PROJECTION["Krovak"],
PARAMETER["False_Easting",253276.6182],
PARAMETER["False_Northing",1295660.263],
PARAMETER["Pseudo_Standard_Parallel_1",78.5],
PARAMETER["Scale_Factor",1.0003],
PARAMETER["Azimuth",105.74038],
PARAMETER["Longitude_Of_Center",120.531377],
PARAMETER["Latitude_Of_Center",-27.24038],
PARAMETER["X_Scale",1.0],
PARAMETER["Y_Scale",-1.0],
PARAMETER["XY_Plane_Rotation",259.05],
UNIT["Meter",1.0]]


texts: Use ST_startpoint and linear referencing on a self intersecting linestring that includes a road interchange loop ramp

texts:
I have a linestring (representing a route along roads) that has a self-intersecting loop.  As it has a self-intersecting loop it is not valid as a linestring. But I need it to be a linestring to determine its startpoint and endpoint in PostGIS and to carry out linear referencing for stops along the route.

Any ideas - must I convert to 3D?


texts: Can't open Google Earth Engine scripts from journal articles

texts:
I am trying to launch GEE scripts from a few journal articles, but once opening the code editor, I can't access the assets other authors worked with to see the full analysis. Here is an example, this code is from an article published in Remote Sensing : https://code.earthengine.google.com/e3a2d589395e4118d97bae3e85d09106
The error messages state that the collection assets ('ft:1POjOubAf8ZSvmMW3Ml3Ka3EUBt-5wY1eUg_J-04E' and 'ft:1goDR2acd9mOOaY3ha092MU_NOOMRvTh0sIApn_od') are not found.
I believe the code editor doesn't load the assets used by the author, does that mean that the author didn't allow the assets to be public?


texts: Set ESRI JavaScript Api's starting point on page load

texts:
I have the following block of code in my .js file:

require([
  "esri/Map",
  "esri/views/MapView",
  "dojo/domReady!"
], function(Map, MapView) {
  var map = new Map({
      center: [-46.807, 32.553],
      basemap: "streets"
  });

  var view = new MapView({
    container: "viewDiv",  // Reference to the DOM node that will contain the view
    map: map               // References the map object created in step 3
  });
});

The results is, as expected, a map of the whole world. How can I set the starting extent of the map to be my city, so clients will know my company's location. I've tried using the extent option, but it doesn't work. Any suggestions?


texts: PostGIS extension on Debian not found: could not open extension control file

texts:
I'm trying to set up the PostGIS extensions on a Postgres 9.1 database on Debian Wheezy, but getting an error:

=# create extension postgis;
ERROR:  could not open extension control file "/usr/share/postgresql/9.1/extension/postgis.control":
No such file or directory

I am not really sure what packages I need to install on Debian to get this working. I've done:

$ sudo apt-get install python-gdal
$ sudo apt-get install postgis
$ sudo apt-get install postgresql-9.1-postgis

all of which appeared to install fine.

What am I missing?

I've tried to install related packages that I've seen on here, but with no success:

$ apt-get install postgresql-9.1-postgis-scripts
Reading package lists... Done
Building dependency tree
Reading state information... Done
E: Unable to locate package postgresql-9.1-postgis-scripts
E: Couldn't find any package by regex 'postgresql-9.1-postgis-scripts'

Using autocomplete, apt-get can only see postgresql-9.1-postgis. I've even tried sudo apt-get install postgis* but it hasn't made a difference.


texts: How do I prevent 3rd party application from loading custom ".style" files?

texts:
Backround: I run Geosoft Target on my machine which has some compatibility functionality with ArcGIS. (mxd conversion, export to, etc.)  I use the Geosoft Target application "stand-alone" and never use the ArcGIS funcions.

Issue: Target installs a ".style" file in ints instalation directory which appears in the Styel manager.  If i uncheck the GEosoft Target ".style" file do clear it from the Style Manager it succeeds only in the current local session - the next time i open any other ArcMap doceument including the one in which i disabled it, the style file is automatically re-loaded.

The one possible patchwork solution would be to delete the Geosoft ".style" file but I am more interested in doing it properly in ArcGIS - is there a global setting that prevents 3rd party syle files from loading?


texts: Generalising lines of varying detail in QGIS

texts:
I am trying to generalise a line in QGIS so that I can export it as a vector without it having an unnecessarily huge file size. The line shows the GPS track of a road trip, and as such varies between having straight lines connecting a few points over a relatively large distance (eg where I was driving at 75mph along a straight road) and many points very close to each other (ie where I've stayed in the same building, but the GPS has wandered slightly and linked up several points close to each other).

The images below show the output from the Vector --> Geometry Tools --> Simplify tool (with a basemap for scale context): the original GPS track is red and the simplified version is black. As you can see, although the tool straightens out journeys of hundreds of km (more than I would want it to), if you zoom into the places where points have been more densely recorded (ie where I've stayed in the same place for an extended period of time), much of the unwanted noise is still there (albeit simplified).

Is there a tool for QGIS that will generalise in a way that doesn't simply remove a certain number of vertices without consideration of the distance between them, but will remove the noise of the type seen in the second screenshot? I am ultimately looking to map the whole route (~5,000 miles), so don't need to see where within a particular town I walked, but would like to still preserve the routes driven along roads (ie the red in the first image) to a reasonable degree.

I have tried various menu option configurations and an additional generalisation plugin, but all have yielded similar results.


texts: Getting SAGA-GIS Python Bindings without compiling from source

texts:
I am very new to SAGA-GIS, but I recently developed a workflow using SAGA for processing LiDAR data.  Now I need to batch the process using Python.  I installed SAGA from the Windows installer and I already have GDAL and other open source packages installed as well.

When I create a python script for a module using the SAGA GUI (right-click on module), my Python installation (Python 2.6 installed with ArcGIS 10.0) does not recognize the module called saga_api.  I found saga_api.dll in my SAGA install directory, but not saga_api.py.

All the instructions I have found for installing the SAGA Python Bindings require compiling SAGA from source, but I would rather not mess with that since I already have it installed an working.  Not to mention, I worry that compiling all the dependencies, like gdal, will break the installations of those packages I already have installed.

Can anyone offer any insight?


texts: TileStache + OpenLayers TMS - Tiles are rendered disorderly

texts:
I newly turned to TileStache to render Mapnik tiles because of the caching option. Before, I was using Mapnik python bindings only and it was working perfectly. Now that I use TileStache, the tiles are rendered randomly at wrong position in wrong order, just like pieces of a puzzle. Only at zoom 0 the features are on the right spot on the map.

This make me think that the problem could be with the row (y) and column (x) used to generate the tiles with coord = ModestMaps.Core.Coordinate(y, x, zoom) or with the projection wich is EPSG 900913.

I use OpenLayers.Layer.TMS for the mapping.

Here some codes:

Tiling script:

def tile(request, version, shapefile, zoom, x, y):
        config = {
            "cache": {
            "name": "Disk",
             "path": "/mygeosite/stache",
             "umask": "0000"},
          "layers": {
            shapefile.filename: {
                "provider": {"name": "mapnik", "mapfile": "c:\\mygeosite\\geodjango\\test.xml"},
                "projection": "spherical mercator"
            }
          }
        }

        # like http://tile.openstreetmap.org/1/0/0.png
        # x = row , y = column , zoom = zoom level
        coord = ModestMaps.Core.Coordinate(y, x, zoom)
        config = TileStache.Config.buildConfiguration(config)
        bytes = TileStache.getTile(config.layers[shapefile.filename], coord, 'png')
        return HttpResponse(bytes[2], mimetype="image/png")

mapnik.xml

<?xml version="1.0" encoding="utf-8"?>
<Map srs="+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs +over" background-color="rgb(115,145,173)">
    <Style name="featureLayerStyle">
        <Rule>
            <PolygonSymbolizer fill="rgb(247,237,238)"/>
            <LineSymbolizer stroke-width="0.5"/>
        </Rule>
    </Style>
    <Layer name="featureLayer" srs="+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs +over">
        <StyleName>featureLayerStyle</StyleName>
        <Datasource>
            <Parameter name="dbname">ezmapping</Parameter>
            <Parameter name="geometry_field">geom_multipolygon</Parameter>
            <Parameter name="geometry_table">&quot;shapefile_feature&quot;</Parameter>
            <Parameter name="password">ogrant22</Parameter>
            <Parameter name="simplify_geometries">1</Parameter>
            <Parameter name="srid">900913</Parameter>
            <Parameter name="table">(select geom_multipolygon, id from &quot;shapefile_feature&quot; where shapefile_id=68) as geom</Parameter>
            <Parameter name="type">postgis</Parameter>
            <Parameter name="user">postgres</Parameter>
        </Datasource>
    </Layer>
</Map>

OpenLayers.js:

map = new OpenLayers.Map('map',
          {maxResolution: 156543.0399,
            projection: new OpenLayers.Projection("EPSG:3857"),
            units: 'm',
            maxExtent: new OpenLayers.Bounds(-20037508.34,-20037508.34,20037508.34,20037508.34),
            numZoomLevels: 21});
tiledLayer = new OpenLayers.Layer.TMS('TMS',
                "{{ tmsURL }}",
                {serviceVersion: "1.0",
                layername: "{{ shapefile.id }}",
                type: 'png'});
map.addLayer(tiledLayer);


texts: How to serve an image layer with tilestache

texts:
I have some .png images with associated .wld and .tfw files.

I'd like to display this in OpenLayers, and I already have a Tilestache server running, so it would be good to leverage this.

OpenLayers supports drawing an image with ol.layer.Image, but this requires you to provide the geographic boundaries of the image in the Javascript.  What I would like to do is have this done on the server side, through Tilestache if possible.

Does Tilestache provide some way of doing this?


texts: Looking for API providing accurate timezone data

texts:
I'm looking for an API providing accurate timezone data.

Thus far managed to find earthtools api but daylight saving time info isn't very accurate.

I need to lookup by latitude/longitude or as an alternative (either one or the other is fine) based on <city,country>

Any ideas?


texts: Using Contour plugin, progress bar sticks at zero

texts:
I am using Quantum GIS 1.8.0. I have a base map sourced from OSGB Opendata VectorMap District (vector format), and that displays fine. I also have created a spreadsheet showing journey times to a chosen destination, with the grid references of 4623 postcodes within 10km of that destination. The first few lines of this spreadsheet are:
        from    easting northing  time
        DH1 1AB 430317  543630     527
        DH1 1AD 430239  543638     534
        DH1 1AE 430220  543685     545
        DH1 1AF 430340  544013     555
        DH1 1AG 430323  543876     556
I have imported this as a delimited text layer (is that the right terminology) and my map shows the locations of all these postcodes.
What I want to do is to plot contour lines showing the areas that are within 5 / 10 /15 / 20 minutes of the destination and shade them. (Is this called an isochron?) When I try to use the contour plugin to do this, it goes to 99% cpu and the progress bar does not move off zero. I should perhaps add that postcodes are packed quite densely in settlements and quite far apart in open countryside.
There must be a better way to do this, but as a newby to this software I need some advice please.


texts: How to create WMS using mapproxy with Mapnik XML file as input?

texts:
I have Mapnik generated XML file for mapproxy input, but I didn't understand how to manage this into mapproxy.
I want to create WMS using mapproxy with Mapnik XML file input.
My XML file: myxml_mapnik.xml is:
<?xml version="1.0" encoding="utf-8"?>
<Map srs="+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs" background-color="rgb(255,255,255)">
  <Style name="mapStyle">
    <Rule>
      <PolygonSymbolizer fill="rgb(0,100,0)"/>
    </Rule>
  </Style>
  <Layer name="mapLayer" srs="+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs">
    <StyleName>mapStyle</StyleName>
    <Datasource>
      <Parameter name="dbname">gis</Parameter>
      <Parameter name="host">localhost</Parameter>
      <Parameter name="password"/>
      <Parameter name="port"/>
      <Parameter name="table">polygons</Parameter>
      <Parameter name="type">postgis</Parameter>
      <Parameter name="user">apple</Parameter>
    </Datasource>
  </Layer>
</Map>

Now I want to take this as a input in mapproxy
How can I do this?


texts: When selecting features which are duplicated it does not appear right in QGIS

texts:
I have a road network data. In this data, I have segments which has a twin segment in the opposite way in the same layer. It is being considered as a single carriageway. So their just start points and end points crosses (which will lead me to determine the direction). Thus I have two issues,

When I try to select one of them, its twin segment could not be seen clearly. For instance, the segment with id 3823233 can be seen clear and 38223228 cannot. See below images,

Moreover, when I try to zoom out a little bit the selections just disappear. I think it is also because of the twin geometry issue.

I want to make the selected one appear on the top all the time. Any suggestions to deal with ? (Note that I am going to use my own selection tool in a project. So, I am open to any PyQGIS solution.)


texts: Custom forms feature in QGIS: QMessageBox.warning closes QDialog as well

texts:
I am attempting to work with custom forms in QGIS. I have included a validate function in my code that will pop a QMessageBox() if the validation fails. My problem is that when I close the QMessageBox the QDialog also closes without giving me the chance to modify its contents. Here is my code:

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

form = None
name_field = None

def form_open(my_dialog, layer_id, feature_id):
    global form
    global name_field

    global name_field
    form = my_dialog
    name_field = form.findChild(QLineEdit, "field")
    button_box = form.findChild(QDialogButtonBox, "buttonBox")
    button_box.accepted.connect(validate)
    button_box.rejected.connect(form.reject)

def validate():
    #This function validates that all of the required fields have been completed correctly.
    results = list()
    results.append(validate_field(name_field, "str"))

    if False in results:
        msg_box = QMessageBox()
        msg_box.setIcon(QMessageBox.Warning)
        msg_box.setWindowTitle("Field")
        msg_box.setText("Please complete the Field")
        msg_box.exec_()
    else:
        form.accept()

I won't to be able to close the QMessageBox and go back to my initial QDialog window. Anyone kind enough to take a look?


texts: Conversion of GeoTIFF to NetCDF file doesn't create a NetCDF file but my database thinks it does

texts:
I am using miniconda and the Open Data Cube project. I am trying to ingest some Landsat_5 data into a datacube, where the geotiff are converted to netCDF files. Using the datacube code I used the ingest tool:

datacube -v ingest --executor multiproc 12 -c C:\envs\DataCube_Files\LS5_USGS_UTM17N.yaml

It creates NetCDF files and they look correct, however when running through my jupyter notebook code where I load the landsat dataset with:

landsat_dataset = dc.load(x=x_extents,
                      y=y_extents,
                      crs='EPSG:32617',
                      platform = platform,
                      product = product,
                      measurements = ['red', 'green', 'blue', 'nir', 'swir1', 'swir2', 'pixel_qa'])

I encounter an error where its looking for a LS5...6_48...nc file that doesn't exists. I have the LS5...6_49...nc , 6_50, and onwards files but not the 6_48. It occurs for multiple files and not just one. I think maybe its an error when converting the geotiff to NetCDF but not sure how to fix it. I've tried making sure my libraries are updated and it looks all good. Error is as follow:

RasterioIOError: 'NetCDF:D:\Cameron\CMAs\Toronto\NetCDF_Tiles\LS5_TM_USGS_SR\LS5_TM_USGS_SR_32617_6_48_19890111152725138088.nc:red' does not exist in the file system, and is not recognized as a supported dataset name.

Any suggestions? I've reduced my images to 15 from 8000 and have completely dropped the datacube and made a new one with only those 15 and still encounter the same problem.


texts: Project a point onto a line feature at a given angle in PostGIS - ST_Project only works for geography?

texts:
I have a point, a line and a known angle. I want to project the point onto the line at that angle. In the below picture, I want to find the red point:

I think this may be possible by somehow using ST_Project, but it seems that it only works with geoography, not geometry. Even the example on the doc page says:

SELECT ST_AsText(ST_Project('POINT(0 0)'::geography, 100000, radians(45.0)));

                 st_astext
--------------------------------------------
 POINT(0.635231029125537 0.639472334729198)
(1 row)

So, projecting a point by 100000 meters puts it at (0.635231029125537 0.639472334729198), which is definitely not right in a cartesian system.

Is there any way to wrangle it into use for geometry instead?


texts: How to clear the cache of QGIS 3.10 with Python?

texts:
I want to add new attributes to a new and empty shapefile layer.
The problem is, that QGIS 3.10 takes attributes from the shapefile layer which was used before when loading it to the GUI. So there are attributes in the table where no attributes are not yet added.
The next problem is, that when I want to edit another layer, for example a point layer and I edited a polygon layer before, QGIS loads the new point layer as polygon.

Is there a possibilty to clear the cache of QGIS before editing a new layer?

Please find below my code.

processing.run("native:reprojectlayer", {'INPUT':'C:/shp_P.shp',
    'TARGET_CRS':QgsCoordinateReferenceSystem('EPSG:32633'),
    'OUTPUT':'C:/edit_N.shp'})

#load layer
P_n =  QgsVectorLayer('D:/SFT_editierbare_Layer','edit_N','ogr')
P_n.setCrs(QgsCoordinateReferenceSystem(32633))
QgsProject.instance().addMapLayer(P_n, True)

P_n.startEditing()
P_n.dataProvider().addAttributes([
    QgsField ('obj_id',QVariant.Int, 'integer', 10),
    QgsField ('DKM_bh', QVariant.String, 'character', 200),
    QgsField ('DKM_kg', QVariant.Int, 'integer', 10),
    QgsField ('DKM_gnr',QVariant.String, 'character', 20),
    QgsField ('FK_rev_nr',QVariant.String, 'character', 200),
    QgsField ('FK_abt_nr',QVariant.Int, 'integer', 20),
    QgsField ('FK_uab_tx',QVariant.String, 'character', 20),
    QgsField ('ALS_seeh', QVariant.Int, 'integer', 8),
    QgsField ('ALS_neigpr', QVariant.Int, 'integer', 8),
    QgsField ('ALS_exptxt',QVariant.String, 'character', 20),
    #QgsField ('ALS_exint',QVariant.Int, 'integer', 8),
    QgsField ('oberh', QVariant.Double, 'double', 3, 1),
    #QgsField ('FK_alter', QVariant.Int, 'integer', 3),
    #QgsField ('nhpro', QVariant.Int, 'integer', 3),
    #QgsField ('lhpro', QVariant.Int, 'integer', 3),
    #QgsField ('vha_nh', QVariant.Double, 'double', 10,3),
    #QgsField ('vha_lh', QVariant.Double, 'double', 10,3),
    #QgsField ('vha_sum', QVariant.Double, 'double', 10,3),
    QgsField ('x_coord', QVariant.Double, 'double', 10,2),
    QgsField ('y_coord', QVariant.Double, 'double', 10,2),
    QgsField ('baumart',QVariant.String, 'character', 254),
    QgsField ('anz_baeume',QVariant.Int, 'integer', 10),
    QgsField ('efm_schatz', QVariant.Double, 'double', 10,3),
    QgsField ('bemerkung',QVariant.String, 'character', 254),
    QgsField ('data_img',QVariant.String, 'character', 254),
    QgsField ('data_doc',QVariant.String, 'character', 254),
    #QgsField ('row_privat', QVariant.String, 'character', 4),
    QgsField ('user_ins',QVariant.String, 'character', 20),
    QgsField ('user_mod',QVariant.String, 'character', 20),
    QgsField ('date_ins',QVariant.Date),
    QgsField ('date_mod',QVariant.Date)])
P_n.updateFields()
P_n.commitChanges()


texts: Working as a Python Development Project Team with ArcGIS?

texts:
We have a development project in Python (ArcGIS 10).  This project involves a mixture of toolboxes, map templates, layer files, file geodatabase templates (acting as templates being imported into a map by scripts) and various other things.

We use Eclipse as our source editor and SVN as our Source Code Repository.

Though we have a problem with keeping all the files (that aren't py files) in a synchronised project by everyone.  The toolbox routinely gets messed up by multiple people editing the toolbox and then template files get adjusted and then not updated for other people as they are not checked back in.

How are people in organisations with more than one python developer on a company toolbox project ensuring that the project and all the different files are get versioned and managed correctly?  Or is it a case through everything goes into Eclipse (including template layers and GDB's used by the scripts) into the project and hope that people check out the files correctly?


texts: Integrate every single street into ArcGIS network analysis?

texts:
This is my first time working with a network analysis (ArcMap 10) and there are some problems which haven't been solved until now.

I have bus stops and Streets (with length, slope, time from-to for a pedestrian (2.5 km/h)).

I want to figure out from which point a person in a village can reach the bus stop within 5, 10 and 15 min. ( I did an ESRI tutorial about this topic, and actually it shouldn't be so difficult, but the tutorial data was well prepared)
At first I don't want to add slope restrictions. I just want to generate polygons around the bus stops with 5,10 and 15 min.

The result is quite weird. It seems that some streets have been analysed but others haven't been. (see screen)

All roads are connected and they have no further restrictions so i wonder why they are not displayed in analysis.


texts: How to check what application uses my services - arcgis server 9.3

texts: I serve a lot of map services via my arcgis server 9.3. Can I determine somehow which webapplication or desktop clients are using or were using my services? I know that I can view apache access log and by this way execute this task, but I think that it is not most elegant solution.


texts:
I'm stuck on a problem with photo2shape in QGIS 2.1x (I already tried several versions).

I installed Exif Reader and photo2shape worked fine for me. Since last week QGIS crashes as soon I start creating the photo point layer. I tried several QGIS versions (2.14, 2.16, 2.18) and reinstalled the exifread plugin too.

No idea how to proceed...


texts:
sorry for the changes but this is the overview of my problem

I'm using 10.2.2 arcgis

I have a geodatabase that every day enter new feature class
I want to create a tool that lists all the feature class I and through a checkbox feature class I decide to be added to my dataframe

Currently I have the following script

import arcpy
import os
arcpy.env.workspace = "C:/Users/Usuario/Documents/Modelo Geomarketing/base/Infomacion.gdb"
datasets = arcpy.ListDatasets(feature_type='feature')
datasets = [''] + datasets if datasets is not None else []
for ds in datasets:
    for fc in arcpy.ListFeatureClasses(feature_dataset=ds):
        path = os.path.join(arcpy.env.workspace, ds, fc)
        print path

results

C:/Users/Usuario/Documents/Modelo Geomarketing/base/Infomacion.gdb\Datos\Buses
C:/Users/Usuario/Documents/Modelo Geomarketing/base/Infomacion.gdb\Datos\Gasolineras
C:/Users/Usuario/Documents/Modelo Geomarketing/base/Infomacion.gdb\Datos\Bancos

I want to know what changes to make to add to my data frame by Means a checkbox
feature class result


texts: Prepare a QGIS plugin installation package (for windows)

texts:
I'm developing a plugin for QGIS that depends on python libs that are not default installed, hence could easily be installed by the network installation.

When the plugin is a little bit more mature I will upload it to github and try to get it included into the available plugins (since I'm not a programer I'll expect it to take some time before it is included).

In the meantime, is it possible to make some kind of .bat file (or similar) that install all python libs that my plugin depend on? The users of my plugin is (hopefully) farmers that (probably) does not have a huge experience of installing complex programs so I want to make it as easy as possible for them.


texts: QGIS mean error - RMSE - QGIS georeferencer

texts:
I would like to know if the mean error computed by QGIS georeferencer is referring to the RMSE (root-mean-square error). I calculated myself with the residual pixels X and Y but I don't have the same results.
I find 4.48 with my calculation and QGIS gives me 5.79.


texts: Deleting polygons with certain height and width in Python script with QGIS?

texts:
I have a polygon layer with bounding boxes (boundingLayer). This layer has attribute columns "HEIGHT" and "WIDTH". I want to delete all bounding boxes which have HEIGHT and WIDTH less than 200m. How should I proceed?

boundingLayer = QgsVectorLayer("/home/usr/Desktop/BoundingBoxes.shp", "BoundingBoxes", "ogr")
QgsMapLayerRegistry.instance().addMapLayer(boundingLayer)
provider = boundingLayer.dataProvider()
updateMap = {}
# To find the right id of the fields I need:
idHeight = provider.fields().indexFromName( 'HEIGHT' )
idWidth = provider.fields().indexFromName( 'WIDTH' )
iter = provider.getFeatures()
for feature in iter:
    Height = feature.attributes()[idHeight]
    Width = feature.attributes()[idWidth]
    if Height < 200:
    ..............


texts: Setting up drawing order for point shapefile

texts:
I am having a difficulty with creating polyline feature from a point dataset that represents a boundary. This dataset consists of 5084 points which are randomly ordered. I have no Network Analyst extension available. I tried to use X Tools Pro and 'Make Polylines from Points' tool, but not having the right drawing order set up the output is just a random cross of polylines.

I have also tried 'Tracing' in the edit session but I've found that it wasn't possible to trace points (only lines or polygons) - please correct me if I'm wrong.


texts: I am using QGIS 2.8. For the documentation of QgsExpression I see the evaluate method on feature as deprecated. What am i supposed to use instead of this methods ? I have some expression and need to check if expression evaluate to True or False on given feature.


texts: digitizing (tolerence) error in line data

texts:
I have designed a large number of boundaries. For accuracy, I took lines in patches from different data sets and merged them into one. Now I want to turn them into polygons from lines. I every second boundary I have to correct digitizing error manually. Either boundaries have gap between them, or they are overlapped or two nodes cross and a gap is left behind. Is there some way that the boundaries may be adjusted that polygons - out of lines - are made. I don't want to correct each boundary manually as they are a lot. Plus I am not a programmer. anything that can be done via arcgis 9.3 or quantumGIS, which i recently installed... :D

Regards,


texts: Time Manager runs with different speed in QGIS 3.10

texts:
I created a map to show my data with time manager. The feature type is polygon.
Changes in "Show frame for x milliseconds" affect the time the frames are shown. But the first frame is shown for x milliseconds and the next frames are shown for much less milliseconds. It almost seems it needs to much time to render and then it runs through the frames.
When I export and generate a GIF. That works fine.
Are there any settings I can change to improve this?


texts: Adding a WMS connection from a paid service

texts:
I'm trying to add a WMS connection to QGIS from a paid service which I subscribed already. Link sample as follows.

http://data.mar-nets.com/08d6d63b1d3b1e3e8541e9805852162a?LAYERS=eslsonar&FORMAT=image%2Fpng&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap&STYLES=&SRS=EPSG%3A900913&BBOX=15028131.255,4070118.8815625,15184674.288906,4226661.9154687&WIDTH=256&HEIGHT=256

I added
"http://data.mar-nets.com/08d6d63b1d3b1e3e8541e9805852162a?LAYERS=eslsonar&FORMAT=image%2Fpng&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap" in the URL section and entered my username and password as authentication. 4 layers will appear but all empty. How should I fix it? is my URL correct?


texts: Downloading kmz/shapefile of differences between 2 compared image dates from planet.com?

texts: I'm interested in downloading a KMZ or shapefile of the difference feature when you compare days on planet.com.  The images are huge composites so I don't think downloading 7,000 items for each date and then using difference function in ArcMap is feasible.


texts: QGIS define hexagon area around mapped points

texts:
I am using QGIS to map a number of radio antennas (for mobile network).
Once I have the layer of the sites (done), I would like to assign to each site a coverage area (envelope). I could achieve that using the "fixed distance buffer" function of fTools

However, I couldn't make the envelope have a Hexagon shape.

When I specify the number of segments:
1: I obtain a square
2: I obtain an octagon ...

How can I do to have a Hexagon?

Thank you


texts: How to set fixed distance labels from data defined centroid in QGIS (3.4)

texts:
I have a labelling issue in QGIS (using 3.4, but haven't tried in other versions). I have polygons for which I show an id number as well as some data information in the labels. I've used rule-based labelling for this: I set two labels and gave them different offsets to the centroid. This results in a two-line-label with on the top part the information you want to display and on the bottom part the id number of the polygon.
This is when not using columns to define the X and Y position.

Sometimes the centroid is not where you want it to be, and so you can use X and Y columns to change the position of the label. Problem is that using these columns negates the offsets I set in the labels, so that they are renderend on top of each other. It doesn't help to make the offsets bigger: they just have no effect anymore, I assume because these settings do not affect the data defined label position.

Is there any way to define offsets to a data defined centroid for Polygons? Or does anyone have a different way to achieve the labels that I want without using offsets?

Example: I changed the position of the label of polygon 20 but now it's broken:

I have two rule-based labels applying to all polygons:

I have settings for offsets but these don't apply to the data-defined position and I seem to have no options to do that in that part of the settings:


texts: Clip an interpolated raster using points or create an extent polygon from points using ArcGIS for Desktop

texts: I have a point feature class that I have imported from lat longs in an Excel spreadsheet. I have created a raster image from these points however areas outside of the point set extent are interpolated. I would either like to find a way to only interpolate within the points extent or clip the raster using the points / create an extent polygon from the points and use that to clip the raster.


texts: How to clip MrSID Files in R and export to GeoTiff?

texts:
I would like to use R to:

Read in an MrSID raster (.sid format)
Clip the extent using a shapefile
Save the new clipped raster as a geotiff

I am fine with steps #2 and 3, but I can't find any information about reading an MrSID file to R. Here is an example of what my workflow would look like if the raster was a .tif instead of a .sid:

# load packages
require(raster)
require(rgdal)

# Step 1: read in raster from .tif file
r <- raster("C:/.../raster.tif")
# Q: how do you do this with a .sid file?

# read in shapefile used for clip
shp <- readOGR(dsn="C:/.../folder", layer="shapefile")

# Step 2: clip raster to shapefile
r.clip <- crop(r, shp)

# Step 3: save clipped raster as tiff
writeRaster(r.clip, "C:/.../raster_clip.tif", format="GTiff")

I know this could be done relatively easily in something like QGIS, but I have a large number of shapefiles I would like to use as clipping extents and then do some follow-up data analysis, so I would like to use R to automate it.

I've found other questions regarding GDAL and MrSID (e.g. here), but none related to R, and I'm kind of a GDAL newbie so I haven't been able to successfully translate the knowledge.


texts: QGIS 3.4.6 Cannot create fields in new shapefiles

texts:
Running 3.4.6 on Windows 10.  Until today I have had no problems creating new shapefiles and fields. Now when I create a new shapefile everything in the new dialog box works as usual when I create fields, then when done I click OK to save and open the attribute table the fields are not there and the visible fields are "Actions" and "_1", "_2" etc. depending on the number of fields originally created. I get the same result when I try to add fields. I have tried this with multiple new shapefiles. How do I get the normal behavior for adding new shapefiles and fields?

I am getting an SIP module error disabling Python on startup. This is probably the source of the problem but I can't seem to fix that. I have posted question and photo of error messages at Couldn't load SIP Module but if anyone here can help please let me know!

]


texts: How to create a googlemaps marker from kml file information

texts:
To clearify:

I tracked myself via my google account...
exported the location history to an xml file.

Now I try to extract the coords from a particular time and location but the format confuses me (

<when>2010-05-28T16:02:44Z</when><gx:coord>0.265749 44.3140135 93</gx:coord>

)

if I insert the value on a googlemaps site that accepts lat and long I end up in Africa, but It is supposed to be France?

How can this information be used to get the correct lat/longtitude values?


texts: Does anyone know a way to use arcpy.describe to see if there is a colormap on a raster? Can't find anything in the help.


texts: Google Earth Pro - Automatically load a secondary database

texts: Is it possible to configure Google Earth Pro so that it automatically loads a secondary database (e.g QLD Globe) when starting up?


texts: Unable to get output for WMSGetFeatureInfo in OpenLayers

texts:
Can anyone please help me on using the WMSGetFeatureInfo in OpenLayers. I am able to add WMS Layers on top of OpenStreetMap. But while using the FeatureInfo functionality I am not getting the results. I am stuck here since long.
I am using the below code in my Application:
var laers = ['Network','HTNetwork','OtherObjects','LandbaseObjects']
              var  info = new OpenLayers.Control.WMSGetFeatureInfo({
                url: 'http://nv-qa:8080/gssogc/ogc?',
                title: 'Identify features by clicking',
                queryVisible: true,
                bbox:bbb,
                i:ee.xy.x,
                j:ee.xy.y,
                info_format:'application/json',
                query_layers:laers,
                format:'image/png',
                crs:'EPSG:4326',
                layers:laers,
                width:map.getSize().w,
                height:map.getSize().h,
                eventListeners: {
                    getfeatureinfo: function(ee) {

                        map.addPopup(new OpenLayers.Popup.FramedCloud(
                            "chicken",
                            mapb,
                            null,
                            ee.text,
                            null,
                            true
                        ));

                    }

                }
            });

            map.addControl(info);

            info.activate();
    });


texts: How to snap and move circles by the center?

texts:
I have a point layer, and a polygon layer with buffers of a given radius around each point. I need to move points, and check that each buffer is not above forbidden land/plots.

In autocad, you can move a circle snapping its center, so you can see immediately the circle limits, and make sure it is not above forbidden land.

Here in qgis, I can only :

move my point
buffer my point
check the buffer doesn't come over forbidden plots.

and if it does, do the process all over again. A pain !

So I would like to :

move at the same time point and buffer by the centre, although it is probably impossible as it's on two different layers,
move a circle by its center, so I can directly snap it to my point and check the circle position.

Any ideas ?

I guess, it could be a plugin doing polygon centroid + move item snapping by its centroid. I didn't find anything of the kind.


texts: Mapping tabular data based on city/province name?

texts:
I am very new to GIS and ArcGIS.

I have an Excel file containing data on  cities (only names) and provinces (of those cities) for 158 countries around the world. The data also contains some economic (e.g. employment level) and demographic data (e.g. population).

I would like to visually analyse (map) the data via ArcGIS.

Do I need to add some geographic data (eg. latitudes and longitudes) into my Excel spreadsheet so that ArcGIS read it and, if yes, where do I get it from?