texts: Adding vertex to line at point location using ArcPy?

texts:
I have a line layer and point layer on my ArcGIS project. The points are on the lines. I want to create vertex on lines at this point coordinates.

import arcpy

point_rows = arcpy.da.SearchCursor("D:\Test.gdb\MyPoints", ['OID@', 'SHAPE'])

rows = arcpy.UpdateCursor("D:\Test.gdb\MyLines")

for row in rows:
    geom = row.SHAPE
    for part in geom.getPart():
        for pnt in point_rows:
                row.X = pnt.X
                row.Y = pnt.Y
                rows.updateRow(row)

I could not add vertex on points, if point is on the line. How can I do this using arcpy?


texts: OpenLayers displaying Polygon with Icon style

texts:
I was working on a React project with an OpenLayers map and I was having problems with showing Polygons on the map. When the map is too zoomed out, the Polygons would be to small to see them. So then I tried to display them as pin Icons. By adding a Icon style to the Vector style and only displaying that Vector at a certain zoom level. But I can't get the pin Icons to show on the map.

This is the function where I add the new layer to the map. Note that I have disabled the maxZoom property in the Vector just for testing purposes. So that the Icons should be visible at all zoom levels. But when I get everything working, this would be undone.

import PinIcon from "../../assets/icon.png";

class Map{

...

addUsersPlotBoundriesLayer(newFeatures) {
    if (newFeatures.length > 0) {
      let vectorSource = new VectorSource({
        format: new GeoJSON(),
        minScale: 15000000,
        strategy: bboxStrategy,
      });
      newFeatures.forEach((newFeature) => {
        let feature = new Feature({
          geometry: new Polygon([newFeature.geometry]),
        });
        vectorSource.addFeature(feature);
      });
      let iconStyle = new Style({
        image: new Icon({
          anchor: [0.5, 46],
          anchorXUnits: "fraction",
          anchorYUnits: "pixels",
          src: PinIcon,
        }),
      });
      let vector = new Vector({
        //maxZoom: 13,
        source: vectorSource,
        style: iconStyle,
      });
      vector.setZIndex(10);
      vector.set("name", "plotUserBoundriesLayer");
      this.map.addLayer(vector);
    }
  }

...
}


texts:
I am trying to overlay a jpg image over an OpenLayers map. I have been to the example on their site and looked all over for a solution with no luck. So I am guessing that their are rules for setting an openlayers.layer.image as an overlay layer. I have tried jpg, png, gif. They all work fine as a base layer with the following code:

var image_layer = new OpenLayers.Layer.Image( 'Wallpaper', 'images.png', new
OpenLayers.Bounds(-90,-45,90,45), new OpenLayers.Size(318,159), {numZoomLevels:7,
maxResolution:.625} );

But once I include:

var image_layer = new OpenLayers.Layer.Image( 'Wallpaper', 'images.png', new
OpenLayers.Bounds(-90,-45,90,45), new OpenLayers.Size(318,159), {numZoomLevels:7,
maxResolution:.625, **isBaseLayer: false**} );

the layer is grey out within the layer switcher

I did figure out that if I include two different image layers - I can specify the second as an overlay. However, once I add a WMS, that overlay greys out!

What am I doing wrong - or better yet - what is OpenLayers doing wrong?


texts: PostGIS ERROR: could not load library "/usr/lib/postgresql91/lib64/postgis-2.0.so"

texts:
OS is openSUSE 12.1, PostgreSQL - 9.1, PostGIS - 2.0. Postgre/PostGIS are installed from Application:Geo repository.

After recent updates when I try to connect to my database I get the following error:

    ERROR:  could not load library "/usr/lib/postgresql91/lib64/postgis-2.0.so":
 /usr/lib/postgresql91/lib64/postgis-2.0.so: undefined symbol: pj_get_spheroid_defn

I searched the web for this error (hint-1, hint-2) and seems that it may have something to do with ld.so.config (or ldconfig command). But I can't figure out what I have to do in my case.


texts: ArcGIS field calculator to calculate most recent date of 4 date fields?

texts:
Currently I am trying to use the field calculator to generate the most recent date values out of four fields into one field. I attempted to use the max function to calculate the most recent dates but the calculator doesn't seem to be reading the values as dates. For an example, if I have a date that reads 12/9/2016 and one that reads 7/20/2017, the field calculator chooses the first date because it reads the 9 (in the 12/9/2016 date) over the 20 (in the 7/20/2017 date) as the max. The case seems to be that the field calculator doesn't read the values as dates.

Any possible python codes for the codeblock to solve this issue? I'm thinking the datetime function is necessary in this case?


texts: Can't see map after renaming fields using Table Manager plugin

texts:
I'm quite the beginner at using QGIS and I'm not sure what I did.

I tried renaming two of my columns/fields using the plugin Table Manager, and after changing the names, my map no longer appears.  I feel that it should be an easy fix but can't figure out what to do at all.


texts: Extract training data with correspond to a reference image

texts:
While extracting the pixel values from a given image, the following technique works, with given geometries (feature collection fc).
var training = image.sampleRegions({
  collection: fc,
  properties: ['land cover'],
  scale: 100});

However, my requirement is that, I want to extract pixel values, using a corresponding image, instead of fc
How can I achieve it?
I'm dealing with two classes (0 and 1) for classification.
image is a multi-channel image.
I need to extract all pixel values with correspond to reference image's 0 pixels as 0 feature value and 1 as feature value 1.
The reference image has 4 values, [0,1,2,3].


texts:
I search a possibility to save querys, that other users can use them.

For example:
I have a polygon layer with different columns:area, id and owner.
Now I want to search in my layer a specific owner or an id.
Normaly I use to query tool for that, but I want to pass the project to other users of my company.
These users are not familiar in handling the query tool. Therefore I search an easier way for query, maybe with the graphical modeler?
But I´m new into the modeler and don´t know how to create a model, in which I can change my search parameters at the time I start the model.

Or is there any easier way for querys in QGIS? Can I save my querys in a project and hand it to other users?


texts: I want to install a dt-classifier plugin I have downloaded from a web page (http://gis-lab.info/qa/dtclassifier-eng.html). I unzipped it to the qgis plugin directory but it won`t show up in my plugin manager in Q-GIS.What can I do?


texts: plotKML RasterBrickTimeSeries (timeseries charts) in parallel and regions

texts:
I am trying to plot a very large timeseries charts associated with RasterBrickTimeSeries-class in plotKML. The idea is to divide the RasterBrickTimeSeries into regions, and to use snowfall multicore as well.

A reproducible example is below:

library(plotKML)
library(raster)
library(sp)
library(snowfall)
data(LST)
setwd("/Users/x/test")
#create a sample point for each cell
pnts <- LST[,47:48] # lon and lat columns
pnts <- cbind(1:nrow(pnts), pnts)
coordinates(pnts) <- ~lon + lat
proj4string(pnts) <- CRS("+proj=longlat +datum=WGS84")

gridded(LST) <- ~lon+lat
proj4string(LST) <- CRS("+proj=longlat +datum=WGS84")
dates <- sapply(strsplit(names(LST), "LST"), function(x){x[[2]]})
datesf <- format(as.Date(dates, "%Y_%m_%d"), "%Y-%m-%dT%H:%M:%SZ")
## begin / end dates +/- 4 days:
TimeSpan.begin = as.POSIXct(unclass(as.POSIXct(datesf))-4*24*60*60, origin="1970-01-01")
TimeSpan.end = as.POSIXct(unclass(as.POSIXct(datesf))+4*24*60*60, origin="1970-01-01")
## get the dates from the file names:
LST_ll <- brick(LST[1:5])
LST_ll@title = "Land Surface Temperature images"
LST.ts <- new("RasterBrickTimeSeries", variable = "LST", sampled = pnts,
              rasters = LST_ll, TimeSpan.begin = TimeSpan.begin[1:5],
              TimeSpan.end = TimeSpan.end[1:5])
data(SAGA_pal)
## plot images in Google Earth:
plotKML(LST.ts, colour_scale=SAGA_pal[[1]])

Till here the code is fine, a chart is plotted for each point

Trying the kml.tiles on pnts works:

tiles.LST.ts <- kml.tiles(pnts, block.x=0.05,
                     size=0.4)

Tiles are produced and work fine, but of course it plots only the point locations  not the actual charts.

I am not much interested in plotting the raster themselves, rather the timeseries charts only. Hence, plotKML.GDALobj might not be the best option

Is there a way to use kml.tiles on the LST.ts (i.e. RasterBrickTimeSeries) to plot the timeseries charts in parallels and regions?


texts: I am the GIS administrator for a 911 call center. Currently we use ArcGIS for our mapping, and my director recently told me that I have a budget every year to purchase GIS/GPS equipment or training. I am looking to for useful equipment that would assist in collecting street or address data, or really anything that would assist mapping for a 911 center. Would anyone on here have any kind of experience with something like that?


texts:
I am unable to solve the following python error in the General tab of the Log Messages when I startup QGIS:
2021-03-29T09:19:03     WARNING    Traceback (most recent call last):
    File "C:/PROGRA~1/QGIS3~1.18/apps/qgis/./python/plugins\processing\script\ScriptAlgorithmProvider.py", line 121, in loadAlgorithms
    alg = ScriptUtils.loadAlgorithm(moduleName, filePath)
    File "C:/PROGRA~1/QGIS3~1.18/apps/qgis/./python/plugins\processing\script\ScriptUtils.py", line 66, in loadAlgorithm
    spec.loader.exec_module(module)
    File "", line 724, in exec_module
    File "", line 860, in get_code
    File "", line 791, in source_to_code
    File "", line 219, in _call_with_frames_removed
    File "C:\Users\Matthew\AppData\Roaming\QGIS\QGIS3\profiles\default\processing\scripts\plugin_upload.py", line 20
    print "Connecting to: %s" % hidepassword(address)
    ^
    SyntaxError: invalid syntax

In the processing tab of the log messages I get the following as well:
2021-03-29T09:19:03     CRITICAL    Could not import script algorithm 'export_to_CSV' from 'C:\Users\Matthew\AppData\Roaming\QGIS\QGIS3\profiles\default\processing\scripts\export_to_CSV.py'
    module 'processing' has no attribute 'getObject'


texts: Deleting duplicates from multiple tables using MapInfo

texts: Using the delete duplicates tool is there a way in Mapbasic to have it run through a large number of tables, or do you have to manually select each table and run the tool?


texts: Function to report all polygons and their data within another polygons?

texts:
I have a layer with different types of vegetation and another layer with estates. Is there a function in ArcGIS 10.3 to create a report of each estate where it lists what It finds on a selected estate?

It is worth mentioning that I have already used union to cookie cut the "nature polygons" to the estate, which means there can be cases where there are 3-4 different nature types on a estate which a lot of information assosiated to each of them


texts: SQL Server Geography Query works but Geometry Query Does Not - Different Tables, Similar Schemas

texts:
I am trying to get a Geometry query to work. A similar Geography query works fine but I must work with a table that uses a Geometry type. Although the Geography version returns lots of records as expected, I cannot get the Geometry version to return any records. Both tables have exactly the same latitude and longitude records .

This Geography query works fine:

DECLARE @home GEOGRAPHY
SET @home = GEOGRAPHY::STPointFromText('POINT(-0.7799193 51.3083162 )', 4326);

SELECT OutwardCode, InwardCode, Latitude, Longitude
FROM dbo.PostCodeData
WHERE GeoLocation.STDistance(@home) <= (5 * 1609) -- 1609 = approx metres in 1 mile

The table schema is:

+-------------+--------------+
|    Field    |     Type     |
+-------------+--------------+
| OutwardCode | Varchar(4)   |
| InwardCode  | Varchar(3)   |
| Latitude    | Decimal(9,6) |
| Longitude   | Decimal(9,6) |
| GeoLocation | Geography    |
+-------------+--------------+

Example Table Data:

+-------------+------------+------------+----------+------------------------------------------------+
| OutwardCode | InwardCode | Longitude  | Latitude |                  GeoLocation                   |
+-------------+------------+------------+----------+------------------------------------------------+
| GU14        | 9HL        | -0.7803759 | 51.30818 | 0xE6100000010C01A4367172A7494027C522E1D6F8E8BF |
+-------------+------------+------------+----------+------------------------------------------------+

This Geometry query returns no records (I have exactly the same latitude and longitude records in the database but have Geometry as a centre point for the street and Postcode is a joined version of OutwardCode and InwardCode):

DECLARE @home GEOMETRY
SET @home = GEOMETRY::STPointFromText('POINT(51.3083162 -0.7799193)', 0);

SELECT Postcode, Latitude, Longitude
FROM dbo.OS_Locator
WHERE Centre.STDistance(@home) <= (5 * 1609) -- 1609 = approx metres in 1 mile

The table schema is:

+-----------+--------------+
|   Field   |     Type     |
+-----------+--------------+
| Postcode  | nvarchar(10) |
| Latitude  | Decimal(9,6) |
| Longitude | Decimal(9,6) |
| Centre    | Geometry     |
+-----------+--------------+

Example Table Data:

+----------+-----------+-----------+------------------------------------------------+
| Postcode | Latitude  | Longitude |                     Centre                     |
+----------+-----------+-----------+------------------------------------------------+
| GU14 9HL | 51.308304 | -0.779928 | 0x346C0000010C00000000549C1D410000000018330341 |
+----------+-----------+-----------+------------------------------------------------+

Where am I going wrong?


texts: Offline editing slides the attribute data over one column and ads an "fid" column

texts:
My QGIS 3.4 Madeira project is attached to a PostGIS database (Azure) and looks good. It contains two layers, one with points and one with polylines.   Each item in the Pointz layer has about ten attributes.  When I use the standard Offline Editing plugin, it slides the attribute data over one column and ads an "fid" column - in the pointz layer only.   The data is all there but in the wrong column.

How do I stop this or is this a bug in the Offline Editing plugin?

This is how the attribute columns are supposed to be -

This is what happens when I use the Offline Editing Plug-in.

I tried to add the fid attribute column as described below just to see what would happen.  It created the fid_1 column as predicted.

It must be something in the creation of the offline file that the mapping gets scrambled but how does it figure out how to unscramble them when Offline Editing Synchronizes?

I have done this on two different databases now and it is the same.

I have created a much smaller shapefile from scratch with exactly the same attributes and it works correctly.

It may be a size issue?

I do not know what is happening exactly.

Does anyone have an Idea on how to troubleshoot?

I do not think there is anything wrong with my shapefile.  When I look at the attribute table it looks OK.


texts: Reproject and set central meridian and standard parallel in Google Earth Engine

texts:
I am creating a time-lapse video of Sentinel-1 imagery. In order to compare it well with some other images I have, I want the projection to be in Antarctic Stereogrpahic (EPSG:3031) and I want the standard parallel set at -75 and the central meridian at 169.
Could someone please advise how to do this. Without any re-projection lines, my code successfully produces the video in the default projection (https://code.earthengine.google.com/6ad5f05b2dc59e8d66e164c89843e888)
Based on this topic I tried adding:
var proj = ee.Projection('EPSG:3031');
var newcol = newcol.map(function(image){
  return image.reproject(proj,null, 40);// , null, 10
})

In that case the video was still produced but it took way longer to process while failing to reproject.


texts: OGC API Processes returns error 500 or 400

texts:
I am testing the following OGC API endpoint on demo pygeoapi server
https://demo.pygeoapi.io/master/processes/hello-world/jobs
following the documentation page 32
I post the following JSON body:
{
    "id": "dzed5zef45zfze",
    "inputs": {
        "literalInputId":
            {
                "dataType": {
                    "name": "string"
                },
                "value": "Marsel"
            }

    },
    "outputs": {
        "literalOutputId": {
            "transmissionMode": "value"
        }
    },
    "response": "document",
    "mode": "sync"
}

with the following python command:
>>> r = requests.post('https://demo.pygeoapi.io/master/processes/hello-world/jobs', data=json.dumps(body))
>>> r
<Response [500]>

or
>>> r = requests.post('https://demo.pygeoapi.io/master/processes/hello-world/jobs', data=body)
>>> r
<Response [400]>


texts: QGIS difference of a layer with itself is not empty. Why?

texts:
I have a Multiline layer A with about 900 Features loaded (only a few of those are actually multi part geometries)  in QGIS 2.14.
I try to do Vector  -> Geoprocessing Tools -> Difference with Layer A as the input vector layer and Layer A as the difference layer, expecting to get as a result an empty layer.  But, the Output Shapefile has ~150 Features in it.  Topology Checker says, there are no invalid geometries, but there are about 90 pseudo nodes in the layer.  I started snooping around and found other line shapes that behave this way.  Have I fundamentally misunderstood what a difference is, or is there some specific geometry problem that causes difference to behave oddly?  The use case is silly, but I started with comparing an exported shape to its source, got weird results, and then tried to compare the source to itself as a sanity check.  Now I feel I might have lost sanity.  Help?

I get the same result in 2.8, so this is not a 2.14 bug.


texts: Creating shared, custom border with unique colors for each polygon using ArcMap

texts:
In ArcGIS Desktop 10.7 I am attempting to make a custom border for each polygon where there is a:

shared, custom border
unique border color
transparent fill

Here is what I tried: In Properties under the Symbology tab I selected Categories: Unique Values. The default settings shows <all other values> and I edited the symbol to create the custom border. Then I 'Add All Values' and the features populate with the custom border. However, the selected Color Ramp does not apply to the border but instead to the fill. Each border is black with a unique fill.


texts: Interpolation for unevenly spaced data points

texts:
I am currently having some trouble finding a good interpolation technique for unevenly distributed spatial data points. Let's give a quick example :

I have a 1000 points dataset. As several surveys were made in the past, several spatial spacing were used. Some of them use a 200 meter grid (red bounding box), whereas other used a 1 km grid ( blue bounding box, see image below for a quick drawing). I merged all the datasets and now i would like to interpolate the data in an even grid. The data here just have spatial components (x and y in projected system) associated with the variable of interest, there is not time component to this data.

I usually use ordinary kriging or IDW, but they do not feel quite appropriate in this case, as they will be influenced by the proximity change  of neighbouring data points.
My search on the internet did not lead me to a solution.

Edit : Those methods seem inappropriate for me at the moment  as a search algorithm based on a diameter will use more samples in areas of high density compared to low density regions, that may bias interpolation, in my opinion.

Did you already encounter such a problem with spatial data ? How did you manage to handle it ?


texts: Difference between overlapping polygons using OGR

texts:
I’m trying to do a simple difference between two polygons using ogr, but I can’t make it work. The input layers have the same CRS. The input polygon is a multipolygon, and the overlay polygon is the result of merged polygons.
I’m using this code:
ogr2ogr diff8.shp project_landuse.shp -dialect sqlite -sql "SELECT ST_Difference(project_landuse.Geometry, 'merge_terrain_dissolve.Geometry) AS Geometry FROM project_landuse a LEFT JOIN 'merge_terrain_dissolve.shp'.'merge_terrain_dissolve' b

The result looks like this:

And should look like this:


texts: Varying layout map scale based on layer visibility in QGIS?

texts:
I am looking for a method to vary the map scale of a QGIS layout according to the visibility of a layer.

The output should show the minimum scale at which the layer can be seen in the layout from a range of predefined scales (e.g. 1:1000, 1:2500 etc.).

I would like the method to work with the QGIS Atlas so that the output can be exported efficiently.

As an example, the images below show the layout at 1:1000 (site boundary only), 1:2500 (site boundary and surface water) and 1:5000 (site boundary and surface water). I would like the layout scale to automatically adjust depending on the minimum scale at which surface water can be seen from a range of predefined scales. So in this case, produce a layout at 1:2500 and not 1:1000 or 1:5000.


texts:
I am developing a SOE for ArcGIS Server 10 SP 5 in Java.  The SOE makes use of some GeoTools functionality - specifically the HeatmapProcess.

I have written a small api (with all GeoTools, no ArcObjects) and packaged it as a standalone jar.  The standalone jar works perfect in projects outside of ArcGIS Server.  However, when used inside the SOE, I am eventually get hit with a "Provider com.sun.media.jai.imageioimpl.ImageReadWriteSpi not found" error message.

I am traditionally a .NET guy, so Java's class loading mechanisms are a little alien to me.  I am fairly confident I have the class-path of my SOE set properly - as all other GeoTools classes instantiate just fine.  Looking at the ESRI install, I can see that some version of JAI is included in the java/lib/shared directory - but it is an older version than what I have included in the class-path of my SOE.

Again admitting my lack of knowledge about all things Java to include the JVM, I am assuming that the version of JAI shipped with ArcGIS Server is taking precedence over the version I have included in the class-path of my SOE - through how ArcGIS Server starts the JVM.

So I have a couple of questions for you Java guru's out there:

Is there some way I can update the version of JAI loaded by ArcGIS Server?
Is there some way I can make ArcGIS Server start the JVM giving precedence to the class-path I define?


texts: Access specific feature property in local GeoJSON file?

texts:
I am working with a GeoJSON feature collection file. It has >100 features, each feature has a set of coordinates labelled "geometry" and a set of characteristics labelled "properties". I have been able to display the geometry (polygons) on my map with the following code:
var geojsonLayer = new ol.layer.Vector({
    source: new ol.source.Vector({
        format: new ol.format.GeoJSON(),
        url: 'http://myserver.com/mygeojson.js'
    }),
    style: new ol.style.Style({
    stroke: new ol.style.Stroke({
      color: "red",
      lineDash: [4],
      width: 3
    }),
    fill: new ol.style.Fill({
       color: 'rgba(0, 0, 255, 0.1)'
    })
  })
});

geojsonLayer.setZIndex(100);
map_obj.addLayer(geojsonLayer);    // map_obj is the map defined elsewhere.

There are some more coordinates stored in the properties of the GeoJSON which I want to access and display, I'm also interested in a couple of other parameters. How do I go about accessing these features in OpenLayers? Or maybe this is not possible with OL and I need to parse it using fetch or something similar?
The file is large, so below I have indcluded just one feature:
{
  "type": "FeatureCollection",
  "features": [
    {
      "geometry": {
        "type": "Polygon",
        "coordinates": [
          [
            [
              5.785,
              60.497
            ],
            [
              6.026,
              60.502
            ],
            [
              6.103,
              60.584
            ],
            [
              6.137,
              60.748
            ],
            [
              6.093,
              60.83
            ],
            [
              5.971,
              60.827
            ],
            [
              5.833,
              60.742
            ],
            [
              5.785,
              60.497
            ]
          ]
        ]
      },
      "type": "Feature",
      "properties": {
        "NumIdBirth": 2,
        "CType": 6,
        "LonTrajCellCG": [
          5.9569,
          5.97
        ],
        "CRainRate": 0,
        "ConvTypeMethod": 1,
        "ConvType": 1,
        "SeverityIntensity": 0,
        "ConvTypeQuality": 3,
        "MvtDirection": 299,
        "CTPressure": 60600,
        "NumIdCell": 1,
        "CTReff": 0,
        "LonG": 5.9569,
        "LatTrajCellCG": [
          60.6425,
          60.62
        ],
        "ExpansionRate": 57.384,
        "CTPhase": 1,
        "LatG": 60.6425,
        "BTmin": 268.25,
        "BTmoy": 273.98,
        "CTCot": 23.57,
        "CTCwp": 0.237,
        "NbPosLightning": 0,
        "SeverityType": "",
        "PhaseLife": 2,
        "MvtSpeed": 0.645,
        "Surface": 240,
        "Duration": 900,
        "CoolingRate": -11.2104
      }
    }


texts: Get regions containing a coordinate in Overpass-API

texts:
I am using the Overpass-API to access OSM data. I want to retrieve all relations or ways (don't really know which of them, since I'm very new to this...) containing a given coordinate. I thought of using the way(around:10, lat, lon) selector of the API but this seems to only give ways where one defining node is inside a circle with radius 10 around the position.

Is there any way to come around this?


texts: Selecting first distinct values and associated geometry using python

texts:
I have a layer I have created and assigned unique IDs. I only have a standard Arc license and require only one of each of the rows with a unique ID. Is there any way to select these using python? Or a work around where I can label them and select from there? These items do not always have duplicate geometry.

An example of what I hope to achieve is below:

(I have added the ORDER_ column and populated manually purely as an example)


texts: ArcMap Add-in with app.settings not recognizing app.config changes?

texts:
I have developed an ArcMap Add-in which requires a configuration file. After spending a while trying to read config values from a single App.Config file (and always getting null) I believe that the Add-in cannot read values from here as it is a class library, and looks for the calling application (ArcMap)'s config file when I ask for a key's value (hence the null).

To get around this I used an App.Settings file, which the application can read fine. Creating this also introduces an App.Config file into the environment and Visual Studio seems to keep the two files in-sync during development.

Now that the Add-in is being deployed I need to be able to change configuration values (e.g. log file location). I have tried opening / extracting the .esriaddin file and updating the App.Config file in there but the Add-in retains the same configuration values it had when compiled. I know that the new App.Config values are being persisted in the .esriaddin file because I can view them again after closing the archive.

Does anyone know a reliable way to configure an Add-in and permit this configuration to be updateable once deployed? Any suggestions very welcome as it seems ridiculous that I should need a custom config file for this.

App.Settings values are at the application level, and currently both App.Settings and App.Config have build action: none / do not copy.


texts: QGIS Compiling 3.6.0 Troubleshooting about PyQt5, SIP, Qscintilla on Ubuntu 18.04 LTS

texts:
I work on UBUNTU 18.04 LTS/bionic, I am trying to compile QGIS 3.6.0, for an upgrade from QGIS 3.2. 0. I am facing the same trouble as several people the release  compatibility between SIP and PyQT5 an more maybe.

For QGIS 3.2.0 compilation I worked with SIP 4.19.7 and PyQt5 5.11.3, I succeed the compilation.

For QGIS 3.6.0 I try PyQt5 5.12 and SIP 5.19.15 (GPL release).
I download source tarball from RIVERBANK repository.
I compile each of them, PyQt5 5.12 and SIP 5.19.15,  with configure.py and PyThon 3.6.7. Everything seems OK installation are done. I remove any garbage (I think) about other SIP installation. I installed also QSCIintilla, release 2.11.1. The Qt release is 5.9.5

When, for QGIS 3.6.0,  I prepare the Makefile with  CMake, I have the following messages :

Traceback (most recent call last):
  File "/home/mbariou/WKSCOMPILE/QGIS/QGIS-final-3_6_0/cmake/FindPyQt5.py", line 34, in <module>
    import PyQt5.pyqtconfig
ModuleNotFoundError: No module named 'PyQt5.pyqtconfig'

In fact I have not a such module but I have:

/usr/local/py356/lib/python3.5/site-packages/pyqtconfig/config.py

More, also this error message from CMake :

Found SIP version: 4.19.15
Traceback (most recent call last):
File "/home/mbariou/WKSCOMPILE/QGIS/QGIS-final-3_6_0/cmake/FindQsci.py", line 45, in <module>
from PyQt5.Qsci import QSCINTILLA_VERSION_STR
RuntimeError: the sip module implements API v12.0 to v12.5 but the PyQt5.QtCore module requires API v12.6

Qsci sip file not found - disabling bindings for derived classes

I met the same difficulty when I compiled QGIS 3.2.0 the conflict was the same with other release, bit I succeeded in a random way it seems, in fact I downgrade to the Ubuntu 18.04/bionic distribution about PyQt5 and SIP.

But the following message is recurrent :

RuntimeError: the sip module implements API v12.0 to v12.5 but the PyQt5.QtCore module requires API v12.6

Can somebody deciphers this conflict and what are the release rules to respect between :

PyQt5,
sip,
Qt,
QGIS ?


texts: Where is my mistake in extracting percentages for Urban Atlas 2018 in QGIS?

texts:
When extracting percentages for Urban atlas and adding up all the percentages for one address I either end up at ~60% (80 percent of the times) and also sometimes at 120% (the remaining 20%). Do you know what I am doing wrong? Here is what I did:
I have the x/y coordinates of addresses in GB, Germany and Ireland and downloaded the 2018 versions of Urban Atlas (all cities for Ireland, Germany and GB). I included the addresses by the step "add delimited text layer" and made buffers with the radius of 200m, 500m, 1km and 2km (using MMQGIS). I merged all Urban atlas pieces to one big map. Then I reprojected both layers to EPSG:3035. I made a new variable in the buffer shapefiles called "catchArea" that included the $area of the current buffer. Then I intersected the buffer-shapefile with the all-in-one UA-map. After that step I made a new variable in the intersected file with $area (called classArea) and then made another variable with the formula: (classArea/catchArea)*100.
Now adding up all the percentages I don't end up at 100% as it should be. This error does not exist if I have a buffer with a radius of 10m but I do not get less percentages with the 2km buffer. (So it is 60% in 200m and 60% in 2km most of the times).
I now started to make an extra step with dissolving the all-in-one map first and now intersect it. Dissolving worked fine but intersection now runs for 5 days now and there is no end in sight. Even if this is my solution it takes too long to be a solution and I am hoping for your swarm intelligence to tell me where I have made a mistake.


texts:
I have a table of lines with their geoinformation in PostGIS. All of them are part of one roadline.

In QGIS, I want to query all geo coordinates that are included between two random points on the road.

Is it possible?


texts: When iterating through TOC and erasing using ArcPy parameters are not valid?

texts:
I am trying to Iterate through my TOC and "Erase" part of the feature layers by using a predefined mask. However I keep getting an error that tells me the "Mask" is either not supported, or does not exist. I can say that the feature exists and that the format is supported, as it does work when I use "Erase" within ArcMap itself. Hence it must be something in my code.

I am a beginner at best when it comes to Arcpy/Python.

I have checked the path and this is correct.

Code:

import arcpy
import arcpy.mapping

Project_number = arcpy.GetParameterAsText(0)
Export = r"G:\\zn\\NM\\GGB\\Applicaties-Admin\\Kerngis\\ExpImp\\{}\\Export_{}\\Mask".format(Project_number, Project_number)
Import = r"G:\zn\NM\GGB\Applicaties-Admin\Kerngis\ExpImp\{}\Import_{}\Erased".format(Project_number, Project_number)

mxd = arcpy.mapping.MapDocument("CURRENT")  # Uses your currently open MXD
for df in arcpy.mapping.ListDataFrames(mxd, ''): # Loop through dataframes
    for lyr in arcpy.mapping.ListLayers(mxd, '', df): # Loop through layers

    arcpy.Erase_analysis (lyr, Export, Import)


texts: Calculate the shortest distance between polygons QGIS 3.12

texts:
I need to calculate the shortest distance from the boundary of one polygon to the boundarylines of a set of polygons as shown in the image below

I have followed the instructions according to this question Calculating shortest distance between multiple Polygon borders to a common Polygon border?.
It seemed to provide a good solution using the field calculator to generate a new feature in the attribute table (eg the attribute table of the blue polygons) using the following formula
distance($geometry,geom_from_wkt( 'Polygon ((........))' ))

as discribed in the above link. However, when I randomly checked some of the distances using the measuring tool, while many of my tests were consistent, I found several that were incorect e.g 196.11 km vs 216.9 km.
All the shapefiles and my project are in the same projection system.
I cannot figure out why there are the errors or if this is the best way to calculate these distances.
Ideally I would like tracklines generated to show the shortest distance, however this may be over complicated.


texts: Creating spatial-join which will calculate sum of values in certain field behind attribute in ArcMap?

texts:
When using the "Update Column" tool within MapInfo Pro, you have the ability to update a column in one table with the sum of values from a column in another table; joining them only where objects in the second table lie within objects from the first table. I've attempted to create an example to help me explain this...

Table 1 shows the total number of various sewer incidents reported, per day e.g. there were 3 blockages and 4 odours reported on 01/02/03.

The Geoplan displays this information visually and shows us the locations of these incidents within each of the three cities.

Table 2 is what I'm aiming to achieve. The sum of the total number of incidents reported for each separate category - based on their geographical location.

I know this is possible using MapInfo Pro - however, I'm finding more and more companies are beginning to use ArcGIS. I'm a lot more familiar with the ArcGIS software than I am with MapInfo, and am aware of the "Join Data" tool in ArcMap. However, from what I can see you can only join data either: based upon a common field within the two tables you're joining, or based upon the spatial location of data, joining where points fall within the boundary of a polygon.

The second method is more relevant for me, and it gives an option for the sum of the points within each polygon to be calculated. However, this will calculate the total number of attributes... not the sum of the values within a certain field behind those attributes e.g. so it would calculate 2 flooding incidents in City 1, when there are in fact 4.

Any ideas on how this would be possible?

I'm using ArcGIS 10.3.1.


texts: How can one compute the proportion of a country's area that lies within a given distance from the country's border?

texts:
A lawyer friend of mine is preparing a legal argument that concerns the border region of Germany: the part of Germany that is within a given distance (e.g.,  30 km) of the country's border. He needs to know as precisely as feasible what proportion of Germany's total area lies within that border region.
As a person unfamiliar with GIS, I have been unable to find even a hint as to how this proportion might be computed.

I seek any hints, even if they are high-level like "Use function X in software Y" or "Ask company Z to compute it for you".


texts: Splitting strings into new rows, retaining geometry

texts:
UPDATED BELOW

I have 9 shapefiles, each shapefile has 700,000+ rows and I have roughly 12 attribute fields in each shapefile.

There are three attribute fields, common to every shapefile, that I am interested in - Type, Code, Area. In the majority of records, there is one entry for each field (| B | X2 | 3.45 |) but in some cases you have more than one entry in a field (up to 8) separated by comma;

| A, B, B | X2, ZZ9, GGR | 1.2, 3.1, 0.23 |

It is nearly always the case that if there is 3 records in one field, the other two fields concerned will also have 3 records. If this does not apply, then the data is incomplete and can be ignored. It is also the case that the 1st entry in one field is associated with the 1st entry in another field.

I have found how to Splitting string in Python parser of ArcGIS Field Calculator? and ArcMap 10.1 Python split inconsistent string but that results in deconcatenating the string by the comma delimiter and creating new attribute fields for every comma separated entry.

The thing i need is for the deconcatenated data to create new records NOT new fields. In the instance above;

 1. | A, B, B | X2, ZZ9, GGR | 1.2, 3.1, 0.23 |

would become

 1. | A | X2 | 1.2 |
 2. | B | ZZ9 | 3.1 |
 3. | B | GGR | 0.23 |

the top one is the original record with the superfluous data removed into two new records below it. I need the new records to retain all the other information from the other 9 fields that don't need to be touched and also the geometry as they refer to the same polygon space.

thanks very much, I'm stuck!

(Answers for ArcGIS 10.1, FME, R, Python all gladly accepted)

EDIT 9th APR:

i used the following code to split concatenated entries (by a comma) in fields 'TYPE' and 'AREA' into new rows (there are between 1 and 8 concatenated entries in TYPE and AREA):

import arcpy
searchC =  arcpy.da.SearchCursor("FeatureClass", ("TYPE","AREA", "SHAPE@"))

insertC = arcpy.da.InsertCursor("newFeatureClass", ("TYPE","AREA", "SHAPE@"))

for row in searchC:
   A = [x.split(",") for x in row[:-1] ]

   for i in range(min(len(A[0]),len(A[1]))):
       newRow = [A[0][i] , A[1][i], row[2]]
       insertC.insertRow(newRow)
del insertC

this produces the correct amount of new rows based upon the content of the original concatenated rows.

but when doing this, i only get the 2 fields TYPE and AREA next to a new FID. each entry also has a uniqueID that has to come across as well, so i tried (amateurishly) this;

import arcpy
searchC =  arcpy.da.SearchCursor("FeatureClass", ("UNIQUEID", "TYPE", "AREA", "SHAPE@"))

insertC = arcpy.da.InsertCursor("newFeatureClass", ("UNIQUEID", "TYPE","AREA", "SHAPE@"))

for row in searchC:
    A = [x.split(",") for x in row[:-1] ]

    for i in range(min(len(A[0]),len(A[1]),len(A[2]))):
        newRow = [A[0][i] , A[1][i], A[2][i], row[3]]
        insertC.insertRow(newRow)
del insertC

However, this now leaves the original concatenated row plus one new row that takes only the first entry from Type and AREA. UNIQUEID does populate both rows.

Any help? i want to split out the concatenated entries but make sure the UNIQUEID that applies to all is carried across, thanks very much


texts:
I have a problem using ArcGIS. I have generate a DEM raster with a 3D online software and exported it in ArcGIS. So I have a raster of my surface. After that, I have some GPS point (x,y,z) that we have take on the zone. To correct my DEM I want to replace my dem value by the point value where DEM=point. I have tried to use the raster calculator like this: con(DEM=points,points,DEM) but it gives me only the value of the points and not the DEM layer.

How can I integrate the points value to my DEM (I have the two in raster form)?


texts: Digitizing single land use (agriculture) from raster in ArcMap?

texts:
I am working on digitizing the Agriculture land cover in a certain region. I did something very similar before and couldn't find a good alternative method beside manually digitizing, which is utterly time-consuming.

I have tried several methods. but none seems to work perfectly.

Image Classification
This is part of the original map

What I got from the Image Classfication Result:

And here is the legend:

The big problem here is that the tool seems to confuse some of the personal yards with Ag Land, which is totally understandable as some of the Ag land is like this:

But I think this can still be a foundation to base my further work on.

I tried to extract certain classes from a landcover raster and then erase tree canopy off it. The results are not quite optimal as you can see here:

Too many "noises" and the boundary are too coarse and edgy.

Is there any solution?


texts: Calculating mean and standard deviation for raster values using whole pixels within a polygon

texts:
I have a raster image and a corresponding shapefile. I would like to calculate the statistics of the values within polygons which have a common attribute but only using pixels which are fully contained in the polygons.

I would like to do this in R, however I have so far only managed to load the image and shapefile using readGDAl and readOGR respectively!


texts: Installing Postgis on Ubuntu 12.04: unmet dependencies libgdal1

texts:
After sudo apt-get install postgresql-9.1-postgis I got that massage:

  The following packages have unmet dependencies: postgresql-9.1-postgis
  : Depends: libgdal1 (>= 1.9.0) but it is not going to be installed E:
  Unable to correct problems, you have held broken packages.

I tryed the ubuntugis-unstable and ubuntugis-stable repositories. I have QGIS 2.2 installed.


texts:
I have a shapefile with points, like this:
x, y, id, data
-1.3345,-0.4760, 1, data.csv

The data.csv contains:
date, data
2020-01-01, 1
2020-01-02, 2
2020-01-03, 3
2020-01-04, 4

How can I plot the data from the CSV file in QGIS?


texts: Measuring roundness using buffers around centroids

texts:
I'm trying to measure the compactness (roundness) of some 800 polygons of different shape and size. I created the centroids and what I want to do now is to buffer those centroids with the exact same area of their polygons.
Sorry if this is too simple but I just can't find the way to do it.


texts: How to convert from a webmap on ArcGIS to downloading shapefile format?

texts:
I have a ArcGIS item in the format of a Web Map and would like to download it as a shapefile (or similar format, such as KMZ or CSV, would also be acceptable)

Is this possible, and if not, why?

By the way, here's the menu containing options that I see:


texts: QGIS 3.2 Layout Map render/cache/rectangle no longer available

texts: In previous version of QGIS when using the print composer there was an option to have the map displayed as a rectangle, cached or render. This was particularly useful when the map showed lots of information and took a long time to load/redraw. I am setting up a template to auto generate a lot of maps using an atlas but can no longer find this option and it takes quite a while to render the map, even though i don't need to see it, just want to export it. Does anyone know if this feature has been removed completely, hidden somewhere else or will be implemented back in further versions?


texts:
Android defineds a GPS device capability I couldn't find reference to
It's defined as GPS_CAPABILITY_ON_DEMAND_TIME in gps.h
The description in the comment explains:
GPS on demand time injection.
What does it mean ?


texts: Using Spatial Join to copy a single attribute

texts:
I am attempting to use Spatial Join to combine tax parcels and zoning information, such that the parcels contain the designation of the underlying zone. However, I only want to take the "esmZONE" attribute from the zoning feature class, without any of the other information. The following script will copy the info from "esmZONE" into "TEST" in the joined feature class, but it also transfers the rest of the attributes from the zoning FC. Can anyone provide some direction as to how I'd remedy this?

import arcpy
import os

base_gdb = r"\\Esm8\gis\python_practice\bellevueBaseData.gdb"
arcpy.env.workspace = base_gdb
target = os.path.join(base_gdb, "inputData", "bellevueParcelsBase")
join = os.path.join(base_gdb, "inputData", "bellevueZoning")
outfc = os.path.join(base_gdb, "outputData", "bellevueTest1")

##good_fields = ["esmZONE"]
##
##field_names = [f.name for f in arcpy.Describe(zoning).Fields if not
##               (f.type in ["OID","Geometry"] or
##               f.name in ["Shape_Length","Shape_Area"] or
##               f.name in good_fields)]
##
##if field_names:
##    arcpy.DeleteField_management(zoning, field_names)

#Create FieldMap and FieldMappings objects
fieldmappings = arcpy.FieldMappings()
fieldmappings.addTable(join)
fieldmappings.addTable(target)

zonefield = fieldmappings.findFieldMapIndex("esmZONE")
fieldmap = fieldmappings.getFieldMap(zonefield)

parcel_field = fieldmap.outputField

parcel_field.name = "TEST"
parcel_field.aliasName = "TEST"
fieldmap.outputField = parcel_field

fieldmappings.mergeRule = "Join"
fieldmappings.replaceFieldMap(zonefield, fieldmap)

arcpy.SpatialJoin_analysis(target, join, outfc, "", "KEEP_COMMON", fieldmappings, "WITHIN")


texts: Using python to batch process lasclip - LAStools/QGIS

texts:
I am trying to use the python console to batch quite a few processes from lastools. I'm having trouble with the arguments to lasclip.

clipped_las = processing.runalg('lidartools:lasclip', True, flight_las, clipper, False, 0, 0, outputname)

is complaining:

Wrong number of parameters
. . .
OPERATION(what to do with points)
     0 - clip
     1 - classify

But I told it to "0" (in the 5th argument) i.e., clip.

I guess I don't understand what ParameterSelection, the type (?) associated with OPERATION is. I don't want to make this a tool with parameters, just want to run it once with the parameters specified in the code.

What ought to go in that 5th spot?


texts: Iterating through every raster cell and calculating distance to other cells using Python?

texts:
I have a binary raster (.tif) in ArcGIS that shows builtup areas (value=1) and non builtup areas (value=0) like this:

In order to make a proximity analysis, I need to iterate through the raster with arcpy/python and for each cell (i) get the distance to EVERY other cell (j) with value=1 within in a given radius around the current cell (i) and write all those distance values to a list/array.

Is there a performant way to do this using Python?

I know there is the EucDistance tool in ArcGIS:

but it may not be the smartest, to create several thousand new rasters.

My plan was to use a RasterToNumpyArray and then calculate on the numpy array

So far I tried the following Python code:

import numpy
from scipy.spatial import distance

def dist_euclid(cell_a, cell_b):
    dist = distance.euclidean(cell_a, cell_b)
    return dist

'''Masking is important to speed up processes'''
def cmask(coords, radius, array):
    a, b = coords[0], coords[1]
    nx, ny = array.shape
    y, x = numpy.ogrid[-a:nx - a, -b:ny - b]
    # create mask where each value in radius around coords is true and everything else is false
    mask = x * x + y * y <= radius * radius

    # reshape mask to rectangular shape
    i, j = numpy.where(mask)
    indices = numpy.meshgrid(numpy.arange(min(i), max(i) + 1),
                             numpy.arange(min(j), max(j) + 1),
                             indexing='ij')
    masked_array = array[tuple(indices)]
    return masked_array

def calcCellDistancesForEachCellInRaster(in_array, buffer_size):
    '''For each cell make slice sub-array around the cell and calculate distances to every other cell (value=1)'''
    dist_arr = []
    for (source_x, source_y), source_value in numpy.ndenumerate(in_array):
        if source_value == 1:
            masked_array = cmask([source_x, source_y], buffer_size, in_array)
            for (target_x, target_y), target_value in numpy.ndenumerate(masked_array):
                if target_value == 1 and ([source_x, source_y] != [target_x, target_y]):
                    dist = dist_euclid([source_x, source_y], [target_x, target_y])
                    if dist <= buffer_size:
                        dist_arr.append(dist)
                # print("dist from %s to %s = %s" % ([source_x, source_y], [target_x, target_y], dist))
                else:
                    pass
        else:
            pass
    return dist_arr

# Create a simple array from scratch using random values
myArray = numpy.random.randint(2, size=(100, 100))
myArray.shape = (100, 100)

dist_arr = calcCellDistancesForEachCellInRaster(myArray, 2)


texts: How to make a certain attribute Toggelable?

texts:
I have imported data from an excel spreadsheet containing lat,long coordinates associated with other attributes such as current operator of each point (well).  I imported as a tab delimited file, then saved this as a shapefile.

I want to be able to organize the files by operator so that they can be toggled on and off (by operator).  I have been doing this by opening the attribute table, making a selection and then saving the selection as a new shapefile for each operator, which works, but I have reached the limit of shapefiles for a single project, and am wondering if there is another way to accomplish this.


texts: Shapefile no map and no attributes in QGIS

texts: I'm new to QGIS and I'm using QGIS 3.14. I was following the tutorials, and I got stuck in adding my very first new vector layer. Every time I tried to do the option new Vector layer, no map appears and the the attribute table is empty.  I tried using other shapefiles, the same thing happens, no map and no attributes. I have also tried using older versions but none works. What do I do?


texts:
I'm looking for an option to import nautical maps (mostly for depth data) into QGIS. So far, I managed to import NOAA .KAP files, but coverage and resolution of some areas are not very good. Is there a chance to import data from Garmin or Navionics into QGIS (I know that this is an expensive option). Or does anyone know a free alternative?

I'm using QGIS 3.4.10 and am a beginner, complex operations to import data are not an option for me right now (except with very detailed instructions for QGIS 3.4.10).


texts:
I am facing a problem when I am trying to execute few queries in PostGIS. It returns the following messages.

    ERROR:  GetProj4StringSPI: Cannot find SRID (4326) in spatial_ref_sys

********** Error **********

ERROR: GetProj4StringSPI: Cannot find SRID (4326) in spatial_ref_sys
SQL state: XX000

What I am missing and what should I do?


texts: Importing toolbox and using one of its tools in ArcGIS Pro

texts:
I have a python script where I need to import a toolbox which contains a number of tools so that I can use one of them.
.............
Name of the toolbox: Generic_Tools.tbx
Alias of the toolbox: GenericToolbox
Name of the module I need to use: MyTool
.............
And these are the python lines I am trying with to import the toolbox and use the tool:
#Importing
arcpy.ImportToolbox(r'..............\Generic_Tools.tbx')

#Using the tool
arcpy.GenericToolbox.MyTool("TempLayer", buffer_fc, buffer_dist, "Meters", OutputWorkspace, "false")

When I run the script it doesn't throw any error, but the tool MyTool is not run at all. Basically, the line where I am calling the module is just skipped.
In ArcMap this is all working fine. I use the following line to call the tool:
arcpy.MyTool("TempLayer", buffer_fc, buffer_dist, "Meters", OutputWorkspace, "false")

But if I try with that same line in ArcGIS Pro, I get the following error:
AttributeError: module 'arcpy' has no attribute 'MyTool'


texts: What does "PostgreSQL layer has no primary key." when loading layer from database manager mean?

texts:
I have a complex join that I have working in the database manager but when I try and load it nothing appears and I found this in the logs:
    PostgreSQL layer has no primary key.
What is actually wrong?


texts: Adding values in X column to list using ArcPy?

texts:
My apologies if this question has been asked before.  I've spent some time searching for a similar question with no results.

Is it possible to add each value under X column to a list using ArcPy?  I think I could brute force it by using cursor and writing each to a list but I'm hoping there is an easier way.

I'll admit I'm fairly new to Python (and ArcPy) but I'm trying to explore the possibilities with several ideas I have in mind.


texts: How to get labels to diplay in Google Earth after saving to KML?

texts: When saving a QGIS layer as kml, is it possible to have labels be displayed when opened in Google Earth?


texts:
from qgis.core import *
from qgis.gui import *
from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry
from math import *
from PyQt4.QtCore import QFileInfo
import os
from osgeo import gdal
from array import *
import ctypes

entries = []
fileName = "E:/WRA/pet/mothmean81j82m.tif"
fileInfo = QFileInfo(fileName)
baseName = fileInfo.baseName()
rlayer = QgsRasterLayer(fileName,baseName)

boh1 = QgsRasterCalculatorEntry()
boh1.ref = 'boh@1'
boh1.raster = rlayer
boh1.bandNumber = 1
entries.append( boh1 )

if 'mjune@1' > 99.6:
    print "not calculate"
else:
    juneop = '( 16*mjune@1)'

print juneop

calc = QgsRasterCalculator(juneop, 'E:/WRA/pet/juneop11.tif', 'GTiff', rlayer.extent(), rlayer.width(), rlayer.height(), entries )
calc.processCalculation()

The if condition is not properly working exactly if condition is false but output is true .
what is the exact expression of greater than to  bands value ('mjune@1'> 99.6:).


texts: Register custom python function in QGIS Server 3.10

texts: I have some custom python functions I have defined in QGIS which I use in print composer. I tried copying the python folder into the QGIS Server path but they do not get propagated. I have seen the following http://osgeo-org.1560.x6.nabble.com/QGIS-Developer-Custom-python-expressions-in-QGIS-Server-td5362875.html but the answer seems vague to me. Any hints / documentation or  that I could read on this anyone out there with an example that I could look at.


texts: Print Composer: Finder window disappearing when pressing save button?

texts:
I recently updated to 2.6.1 (using mac OS 10.10.1). When I'm in print composer, and I try to save my map (or open a new file, or do anything that brings up a finder window), the window disappears before I can select anything. Which means I can't save any files. I've checked to make sure it's not just popping up behind other windows, and as far as I can tell, its' not.

I've tried redownloading and installing GDAL and Matplotlib just in case. I've restarted (both QGIS and the computer).

Any ideas?


texts: Google Earth Engine pull least cloudy Landsat scene by location using conditional

texts:
My goal is to write code that will pull all Landsat 8 scenes available for a given month in a given year for the US, add a cloud score band to each image, and if there are multiple images for the same location (which I plan to define using the WRS_PATH and WRS_ROW properties), pull the least cloudy image into a new image collection.

So far, I have written/pulled code to filter the Landsat scenes and add a cloud score band to them. I am stuck with the conditional.

See my code below - any ideas??

// Load a region representing the United States.
var US = ee.FeatureCollection('ft:1tdSwUL7MVpOauSgRzqVTOwdfy17KDbw-1d9omPw')
.filter(ee.Filter.eq('Country', 'United States'));

// SimpleCloudScore, an example of computing a cloud-free composite with L8
// by selecting the least-cloudy pixel from the collection.

// A mapping from a common name to the sensor-specific bands.
var LC8_BANDS = ['B2',   'B3',    'B4',  'B5',  'B6',    'B7',    'B10'];
var STD_NAMES = ['blue', 'green', 'red', 'nir', 'swir1', 'swir2', 'temp'];

// Compute a cloud score.  This expects the input image to have the common
// band names: ["red", "blue", etc], so it can work across sensors.
var cloudScore = function(img) {
  // A helper to apply an expression and linearly rescale the output.
  var rescale = function(img, exp, thresholds) {
    return img.expression(exp, {img: img})
        .subtract(thresholds[0]).divide(thresholds[1] - thresholds[0]);
  };

  // Compute several indicators of cloudyness and take the minimum of them.
  var score = ee.Image(1.0);
  // Clouds are reasonably bright in the blue band.
  score = score.min(rescale(img, 'img.blue', [0.1, 0.3]));

  // Clouds are reasonably bright in all visible bands.
  score = score.min(rescale(img, 'img.red + img.green + img.blue', [0.2, 0.8]));

  // Clouds are reasonably bright in all infrared bands.
  score = score.min(
      rescale(img, 'img.nir + img.swir1 + img.swir2', [0.3, 0.8]));

  // Clouds are reasonably cool in temperature.
  score = score.min(rescale(img, 'img.temp', [300, 290]));

  // However, clouds are not snow.
  var ndsi = img.normalizedDifference(['green', 'swir1']);
  return score.min(rescale(ndsi, 'img', [0.8, 0.6]));
};

// Filter the TOA collection to a time-range and add the cloudscore band.
var collection = ee.ImageCollection('LC8_L1T_TOA')
    .filterDate('2016-01-01', '2016-01-31')
    .filterBounds(US)
    .map(function(img) {
      // Invert the cloudscore so 1 is least cloudy, and rename the band.
      var score = cloudScore(img.select(LC8_BANDS, STD_NAMES));
      score = ee.Image(1).subtract(score).select([0], ['cloudscore']);
      return img.addBands(score);
    });

// Check: How many landsat scenes are in this collection?
var count = collection.size();
print('Number of Landsat scenes in collection:', count);

// Create a function that says if two Landsat scenes have been pulled for the same location, selected the scene with the lower cloud score. If only one scene has been pulled, select that. Add all selected scenes to a new ImageCollection.
var conditional = function(image) {
  return ee.Algorithms.If('WRS_PATH'

  )
}

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


texts: How to change the coordinates of Shapefiles from another GIS platform to QGIS

texts:
I have over 800 shapefiles which when paired with google maps/bing maps, the point is off by at least 1000 km. Been scratching my head to make this work by changing the default CRS. (Previously was in Atlantic Ocean)
Is it possible to change the longitude and latitude manually for all 800 shapefiles? (its a long-winded process) or is there any other method I could employ for this?
Currently using CRS Layer : GDM2000 / Johor Grid


texts:
I am working with Modis Data (MOD02, MOD11, MOD03, MOD35).
In particular I am looking for a solution to reproject the subdatasets below, into albers equal. Because I am using NLCD dataset which uses albers.(PROJ4: +proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23 +lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs)  area in 1000 m resolution. I will run a model for whole USA so I thats why I chose Albers.

Can somebody help me? By the way I am aware of the MRT Swath but I want to do this in Python environment.

So far I tried this but there was shifts in the pixels. (Alignment problems).

gdalwarp -tr 1000 1000 -of GTIFF -tps -t_srs "+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23  +lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs” -r cubic 'HDF4_EOS:EOS_SWATH:
"MOD021KM.A2007181.1750.005.2010199133141.hdf”:MODIS_SWATH_Type_L1B:EV_250_Aggr1km_RefSB’ x.tif

This one works for the MOD021KM but I couldn't work it out for the rest.
I am confused.

These are the subdatasets that I want to extract as geotiffs.

HDF4_EOS:EOS_SWATH:"MOD021KM.A2007157.1700.005.2010197153031.hdf":MODIS_SWATH_Type_L1B:EV_250_Aggr1km_RefSB
HDF4_EOS:EOS_SWATH:"MOD021KM.A2007157.1700.005.2010197153031.hdf":MODIS_SWATH_Type_L1B:EV_500_Aggr1km_RefSB
HDF4_EOS:EOS_SWATH:"MOD03.A2007157.1700.005.2010196183307.hdf":MODIS_Swath_Type_GEO:SensorZenith
HDF4_EOS:EOS_SWATH:"MOD03.A2007157.1700.005.2010196183307.hdf":MODIS_Swath_Type_GEO:SolarZenith
HDF4_EOS:EOS_SWATH:"MOD11_L2.A2007157.1700.005.2007158222249.hdf":MOD_Swath_LST:LST
HDF4_EOS:EOS_SWATH:"MOD11_L2.A2007157.1700.005.2007158222249.hdf":MOD_Swath_LST:View_angle
HDF4_EOS:EOS_SWATH:"MOD35_L2.A2007157.1700.005.2010197153939.hdf":mod35:Cloud_Mask


texts: Getting black square after merging two rasters in QGIS

texts:
After merging two SRTM TIFF files all I get is a black rectangle.

Merged raster has zeros next to both black and white squares under layer name.

I have no clue what I am doing wrong. Any suggestions?

I downloaded 2 tiles from here: https://dds.cr.usgs.gov/srtm/version2_1/SRTM3/Islands/

Tiles: S44E170 and S44E171

I opened both tiles in QGIS Desktop 3.4.3

And after merging them I get this:


texts:
I have a SpatialPointsDataFrame and a raster for a region. I would like to retain only those points that occur within the elevational limits of the raster. Unfortunately, no points are being filtered. I still see points that occur outside of the 'masked extent' of the raster, although they occur within the overall extent of the raster.

Raster below:

> newelev
 class       : RasterLayer
 dimensions  : 642, 382, 245244  (nrow, ncol, ncell)
 resolution  : 1000, 1000  (x, y)
 extent      : 461951, 843951, 892583.3, 1534583  (xmin, xmax, ymin, ymax)
 coord. ref. : +proj=utm +zone=43 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
 data source : in memory
 names       : alt
 values      : 1400, 2585.524  (min, max)

Below is my SpatialPointsDataFrame

> pr
class       : SpatialPointsDataFrame
features    : 262
extent      : 529866.4, 798884.4, 945186.2, 1498039  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=43 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
variables   : 2
names       :   X, pr
min values  :   1,         1
max values  : 899,         1

I use raster::intersect to obtain only points that occur within the elevational limits of the above raster.

> trial <- raster::intersect(pr,newelev)
> trial
class       : SpatialPointsDataFrame
features    : 262
extent      : 529866.4, 798884.4, 945186.2, 1498039  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=43 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
variables   : 2
names       :   X, pr
min values  :   1,         1
max values  : 899,         1

Below, you can see the plot, where points occur outside the masked extent of the elevation raster. Any suggestions?


texts: can't generate/display signed float values using python/gdal

texts:
I'm trying to create a raster with signed float values as the pixel values and display the results in MapWindow on windows.

What is the proper gdal band type to use?  I've tried both gdal.GDT_Byte and gdal.GDT_Float32 and using both methods I can't seem to get the pixels displayed.  If I change the pixel value to the abs() value they will display.

dataset = driver.Create(filename, cols, rows, number_of_bands, gdal.GDT_Byte)

gdal.GDT_CFloat32 allows for negative numbers but doesn't properly display the floating point value past the point.


texts: Updating field of point which fits in region using ArcMap?

texts:
I have a problem when I was working with ArcMap:
I want to get value from a field of a region and update for points which fit in that region,
I don't known what tools can do it.

What tool should I use?


texts: Convert list of Java objects with JTS geom to GeoJson

texts:
I have a list of Java POJOs that all have a JTS geometry associated with them (Point) - a list of POIs, basically.

I want to output them as GeoJson (FeatureCollection, Feature etc.)

jackson-datatype-jts (https://github.com/bedatadriven/jackson-datatype-jts) will encode JTS Points into GeoJSON, but not a list of objects.

geojson-jackson (https://github.com/opendatalab-de/geojson-jackson) will let me make FeatureCollections and Features with geometries, but uses it's own classes for Geometries so I'd have to convert all of JTS to theirs.

That's not that complicated, but is there any Java library that will let me create proper GeoJson documents out of Java objects with JTS geometries? I started looking at GeoTools but can't really find any example that does that.


texts: Google / Mapbox static maps displaying custom image with certain bounds

texts:
I have a process which generates a simple Google Static Map (can also be a MapBox static image as the code is pretty much interchangeable).

It's recycled code from a project I found years ago.

The final map url creation string is something like this:

var mapurl = "https://maps.googleapis.com/maps/api/staticmap?center=" +
            center[1].toFixed(6) + "," + center[0].toFixed(6) +
            "&zoom=" + zoominfo.zoom + "&size=" + mapsize.width + "x" + mapsize.height + "&scale=4&format=png&key=APIKEY" + mapStyle

Lets say a URL string is this:

https://maps.googleapis.com/maps/api/staticmap?center=42.803331,4.744135&zoom=7&size=640x640&scale=4&format=png&key=APIKEY&maptype=hybrid

And produces this map:

What I'd like to do is set up a bounding box for a specific set of lat/lon co-oridnates, and if the the url generated comes within those bounds, a seperate "customised" image is instead shown as an overlay.  Think of this as essentially a custom map tile, nothing to do with google or mapbox static maps.

This map tile represents a particular region on earth, aligned to specific co-ordinates.

The end result would be something like this:

The closest thing I've seen is an implementaion of Google Maps GL JS where an overlay can be added: https://docs.mapbox.com/mapbox-gl-js/example/image-on-a-map/

This is essentially what I'm trying to do, albeit with Static Maps.  Static Maps does not allow you to add custom raster overlays.


texts: Poor fit of variograms in spatial GLS models

texts:
I run into some problems when investigating correlation structures in spatial GLS models, precisely that fitted variogram does not fit residuals.

Template data can be downloaded here: http://harmackova.wz.cz/other/temp_data.csv

R code:

temp.data <- read.csv('temp_data.csv')[-1] # read data

model.gls <- gls(response ~ predictor.1 + predictor.2 + predictor.3 + predictor.4, data = temp.data) # simple gls
model.gls.spat <- gls(response ~ predictor.1 + predictor.2 + predictor.3 + predictor.4, data = temp.data, corr = corExp(form = ~x+y)) # gls with correlation structure (Exponential model has the lowest AIC)

plot(Variogram(model.gls, form = ~x+y)) # plot variograms
plot(Variogram(model.gls.spat, form = ~x+y))

The resulting variogram looks like this:

It does not fit at all.
When inspecting other correlation structures they also show poor fits.

Then I took a look at Moran's I and how well did this model dealt with spatial autocorrelation:

corr <- data.frame(no = 1:20)
corr$gls <- correlog(x = temp.data$x, y = temp.data$y, z = model.gls$residuals, increment = 2, resamp = 100, quiet = T)$correlation
corr$gls.spat <- correlog(x = temp.data$x, y = temp.data$y, z = model.gls.spat$residuals, increment = 2, resamp = 100, quiet = T)$correlation
plot(corr$gls ~ corr$no, pch = 1, xlab = 'Distance', ylab = 'Moran I')
points(corr$gls.spat ~ corr$no, pch = 16)
legend('bottomleft', c('GLS', 'spatial GLS'), pch = c(1,16), cex = 0.7)

Residual spatial autocorrelation even increased in the spatial model.

I also tried using LME instead of GLS, but results were similar:

model.lme <- lme(response ~ predictor.1 + predictor.2 + predictor.3 + predictor.4, data = temp.data, random = ~ 1|region)
model.lme.spat <- lme(response ~ predictor.1 + predictor.2 + predictor.3 + predictor.4, data = temp.data, random = ~ 1|region, corr = corExp(form = ~x+y))

plot(Variogram(model.lme.spat, form = ~x+y))

corr$lme <- correlog.nc(x = temp.data$x, y = temp.data$y, z = model.lme$residuals, increment = 2, resamp = 100, quiet = T)$correlation
corr$lme.spat <- correlog.nc(x = temp.data$x, y = temp.data$y, z = model.lme.spat$residuals, increment = 2, resamp = 100, quiet = T)$correlation
plot(corr$lme ~ corr$no, pch = 1)
points(corr$lme.spat ~ corr$no, pch = 16)

I would like to know what causes this behaviour and how can I better fit some correlation structure to actually match my data.


texts: I have a raster layer with travel times from a specified point in my city's downtown. What I'd like to do is select all of the census blocks that lie at least partially within 45 minutes of downtown (i.e. I'd like to select all of the vector features that contain at least one raster pixel with a value of less than 45). I'd like to do this in QGIS if at all possible; suggestions for other FOSS software are also welcome.


texts: QGIS - how to symbolise trajectory segments by their start time

texts:
I have a trajectory layer in PostGIS of type LINESTRINGM. I want to visualize it in QGIS and symoblize such that each segment of the trajectory is classified by its start time. the layer field are as follows:
trajectory geometry,
tripid text
min_t timestamp without time zone,  -- start time
max_t timestamp without time zone,  -- end time
pt_cnt numeric  -- no.of points in trajectory

What I have done:
I first split the trajectory into its segments by using geometry generator function segments_to_line in QGIS.
In the second part I need to style the segments with an expression such that it is symbolized/classified by hour interval on the field min_t. Currently I am doing interpolation using scale_linear function as below.
ramp_color(
  'RdYlBu',
  scale_linear(

      hour(min_t)
  ,
    0,24,
    0,1
  )
)

But I want to actually define classes as given below. How can i define classes with an expression?
0h - 6h (night)
6h - 12h (morning)
12h to 18h (afternoon)
18h to 24h (evening)
I am totally new to QGIS expression and need help here.


texts: Calculating euclidean distances of forest areas to road infrastructure in ArcGIS for Desktop?

texts:
I have prepared a forest cover map at 100 m resolution and have also created a euclidean distance layer of road infrastructure for that area at 100 m resolution. I would like to know how I can calculate the distance between each forest cell to the nearest road infrastructure. Should I overlay these two layers to find out which percentage of forest falls within a certain distance from the road?

Preferably in ArcGIS for Desktop.


texts: resize divIcons (svgs) at zoom levels - Leaflet

texts:
I am adding a geoJSON to a map using a custom divIcon (from svg) for point markers with Leaflet.  The divIcon needs to resize on different zoom levels.  Currently I can see the icons and the popup works, but they do not resize (and actually, they are not anchored correctly either).  I have tried to reset the iconSize property on a zoom event, but that is not working; the zoom event is hit but the icon does not change.  I also tried to add multiple icons to each feature (points) in the geoJSON using the 'pointToLayer' function as described in the Leaflet docs, thinking that I could toggle the opacity at different zoom levels, but this does not work either (only one icon is added). How can I use divIcons that will resize?

some code snippets:

geoJSONLayer = L.geoJson(newFeatures, {
//onEachFeature: onEachFeature,
pointToLayer: function (feature, latLng) {

    return new L.Marker(latLng, {
        icon: aDivIcon, clickable:false, opacity:1

      }).bindPopup(popupContent, popupOptions),
    new L.Marker(latLng, {
        icon: aDivIcon2, clickable:false, opacity:0

      }).bindPopup(popupContent, popupOptions); /// this doesn't work

  }
}).addTo(map);

var aDivIcon = L.divIcon({"iconSize": [12,12], "iconAnchor": Anchorvar,  "popupAnchor": Popupvar, html: '<svg class="marker-anIcon"><use xlink:href="#marker-anIcon"/></svg>'});

(event example)

map.on('zoomend', function(e) {
   if (map.getZoom() === 15)  {
        iconSize = [12,12];
        Anchorvar = [12,12] ;
        Popupvar = [ 0, -12];
    }
...
}


texts: Is it possible to create a geocoding and reverse geocoding service based on the shapefiles data

texts:
We have some shapefile data and we have export to them to a MySQL table.

The table have these columns:

id name description fulladdress code lng lat ...

Now I wonder if it is possible to create a geocoding/reverse geocoding service based on these data?

If yes, what jobs will to be done?


texts: Displaying loaded project file name in QGIS title bar when opening via QgsProject.instance().read()

texts:
On Windows, if I open a project in QGIS then the project file name appears in the title bar.

If I do the same thing in PyQGIS using the following code the title bar does not change.

from PyQt4.QtCore import QFileInfo
QgsProject.instance().read(QFileInfo("C:/FilePath.qgs"))

But the project loads fine, I can see my project data and if I click the save button then the file name appears in the title bar.

I was wondering if there was a bridge or some hack needed to get the file name to appear. I need this because my client has similar looking projects and uses the file name to know that my code has loaded the right one for them. Of course I could show it to them another way like in the message bar but seeing as opened projects normally get their file name in the title bar I thought I should be keeping to this 'standard' way of doing it.

Things I've tried...

Saving the file via:

QgsProject.instance().write()

Does save the file but doesn't update the title bar.

Canvas bridge code

from qgis.gui import *
canvas = QgsMapCanvas()

bridge = QgsLayerTreeMapCanvasBridge( \
     QgsProject.instance().layerTreeRoot(), canvas)

Also didn't help.

How can I get the file name to appear in the title bar of QGIS when using .read() ?


texts: How to get layer list of layers that are published?

texts:
I'm stuck.

I've been using geoserver with python extension gsconfig library for a while now. Almost everything works as it should be. While there are little things like publish raster images are problematic, there is still workaround that isn't written into the catalog.py class.

But there is this where I have no idea how to get around.
How do you get a list of layer resource objects that unpublished? catalog.py provides get_layers and get_resources that return layers that are published but not the unpublished layers in a datastore.

I could list all the layers (publish and unpublished) in a datastore through geoserver interface but wms/wfs/wcs and even rest does not

The wms/wfs/wcs services can return layer information in xml but still no layers that are unpublished. I could understand why that is in this case since it is publish facing but I don't know why REST services have this problem as well.

I have searched a lot of forum with no luck.

Did anyone manage to find a solution for this?


texts: Qgis - Projection Problem with EPSG:2272

texts:
I work in both Qgis and ArcGIS.  Most of my data has been created using ArcGIS and all of my data should be in the EPSG:2272.  When I add my data to Qgis, it does not recognize my projection and gives it a user:100000 projection.  Below is a comparison from the first part of projection files from ArcGIS and Qgis.  As you can see the beginning section of each file is slightly different.  Is this causing Qgis to not recognize my projection correctly?  How can I fix it so my projections are recognized in both ArcGIS and Qgis?

from ArcGIS
PROJCS["NAD_1983_StatePlane_Pennsylvania_South_FIPS_3702_Feet",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137.0,298.257222101]],

from GIS
PROJCS["NAD83_Pennsylvania_South_ftUS",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137,298.257222101]],


texts: How to get land_area polygons with osm2pgsql?

texts:
I want to get the polygons of several administrative land areas. For example  https://wiki.openstreetmap.org/wiki/DE:Grenze#Innerstaatliche_Grenzen contains the relation numbers for the land areas of Germany and of five german states that differ from boundary=administrative and are tagged as land_area=administrative. I get entries in osm_polygon for that six boundary relations 51477, 62718, 62782, 28322, 62771 and 51529, but entries for the land_area relations 62781, 2833343, 451087, 62774, 454192 and 62775 are missing.

For osm2pgsql I copied empty.style as fds.style and added a line

way         land_area               text    polygon,nocolumn

(but I notice empty.style has no entry for "boundary" keys).

Data is for Germany and some surrounding areas, downloaded from geofabrik.de and merged with osmosis.

This is the command line I tried for osm2pgsql:

C:\Daten2\osm2pgsql-bin\osm2pgsql.exe --host ubuntu.forplan.local
--username postgres --password --slim --cache 30000 --proj 32632
--hstore --hstore-add-index --keep-coastlines --multi-geometry --prefix osm
--style C:\Daten2\osm2pgsql-bin\fds.style
C:\Daten2\osm-2019-11-06.pbf


texts: Using Batch Project Script to read from text file as Input to ArcPy?

texts:
I was wondering how I would modify the code of the Batch Project script shown below, to read from a text file of shapefile names as input. A sample of the textfile is shown under the code:

#Import required modules
import ConversionUtils, time

#Define message constants so they may be translated easily
msgWorkspace = ConversionUtils.gp.GetIDMessage(86109) # Message "Output workspace does not exist: "
msgCoordinateSystem = ConversionUtils.gp.GetIDMessage(86110)  #Message "Must Enter a Spatial Reference or Template Feature Class."
msgFail =  ConversionUtils.gp.GetIDMessage(86111) # Message "Failed to project "

#Set the input datasets
inputs = ConversionUtils.gp.GetParameterAsText(0)
# The list is split by semicolons ";"
inputs = ConversionUtils.SplitMultiInputs(inputs)

#Set the output workspace
output_workspace = ConversionUtils.gp.GetParameterAsText(1)

#Set the spatial reference
output_coordinate_system = ConversionUtils.gp.GetParameterAsText(2)

#Set the template dataset
template_dataset = ConversionUtils.gp.GetParameterAsText(3)

#Set the transformation
transformation = ConversionUtils.gp.GetParameterAsText(4)
#Message 86112   "Projecting multiple datasets ..."
ConversionUtils.gp.SetProgressor("step",ConversionUtils.gp.GetIDMessage(86112), 0, len(inputs))

if (output_coordinate_system == "" or output_coordinate_system == "#") and (template_dataset == "" or template_dataset == "#"):
    raise ConversionUtils.GPError(msgCoordinateSystem)
elif (output_coordinate_system != "") and (output_coordinate_system != "#"):
    sr = output_coordinate_system
elif (template_dataset != "") and (template_dataset != "#"):
    dsc = ConversionUtils.gp.Describe(template_dataset)
    sr = dsc.SpatialReference

for input in inputs:
    try:
        outdata = ConversionUtils.GenerateOutputName(input, output_workspace)
        #Message 86113 "Projecting "
        ConversionUtils.gp.SetProgressorLabel(ConversionUtils.gp.GetIDMessage(86113) + input)
        ConversionUtils.gp.Project_management(input, outdata, sr, transformation)
        #Message 86114 "Projected %s to %s successfully."
        ConversionUtils.gp.AddMessage(ConversionUtils.gp.GetIDMessage(86114) % (input, outdata))

    except Exception, ErrorDesc:
        msgWarning = msgFail + "%s" % input
        msgStr = ConversionUtils.gp.GetMessages(2)
        ConversionUtils.gp.AddWarning(ConversionUtils.ExceptionMessages(msgWarning, msgStr, ErrorDesc))

    ConversionUtils.gp.SetProgressorPosition()

time.sleep(0.5)

Text File Sample:

  C:\Users\revresources\scratch\shapefiles\Test_Output\Geochem
  Merged\Pyramid_PanCon.shp
  C:\Users\revresources\scratch\shapefiles\Test_Output\Geochem
  Merged\Pyramid_Rocks.shp
  C:\Users\revresources\scratch\shapefiles\Test_Output\Geochem
  Merged\Pyramid_Silts.shp
  C:\Users\revresources\scratch\shapefiles\Test_Output\Geochem
  Merged\Pyramid_Soils.shp
  C:\Users\revresources\scratch\shapefiles\Test_Output\test2\mod\interpretive
  alteration_point.shp
  C:\Users\revresources\scratch\shapefiles\Test_Output\test2\mod\interpretive
  alteration_polyline.shp
  C:\Users\revresources\scratch\shapefiles\Test_Output\test2\mod\Interpretive
  Pyramid geology_polyline.shp
  C:\Users\revresources\scratch\shapefiles\Test_Output\test2\mod\Interpretive
  Pyramid geology_region.shp


texts: How to use classified styles with joined layers?

texts:
I am trying to set a color range based on this tutorial

In LAYER PROPERTIES for "new", in JOINS, I have
joinlayer - data-join
joinfield - GEO_ID2
targetfield - GEOID10

Where "new" is my shapefile and "data-join" is my DBF file with the data I want to use for the graduated color. I want to join column "geo_id2" in "data-join" with "geoID10" in "new."

After I am set up, I click "Apply," then "OK."

When I toggle over to Style, I see the column I wish to use --- "total-hous" I have selected a color ramp "blues" and 5 classes at equal intervals.

But I'm not getting a range like in the Knight example, just five categories all with a range from 0 to 0.

Suggestions?


texts: I know it is possible to use natural Breaks legends with vector files in QGIS. Is this possible also with raster files? Actually I want to reclassify raster files by using Natural Breaks.


texts: How do you change the titles displayed inside infowindows in Carto?

texts:
The titles in my infowindows are copied straight from my column headings in the underlying data... So it says "STREET_ADDRESS" when I'd like it to display "Address" or "YEAR" when I'd like it to display "Year registered"
Usually, Carto is point and click for these kinds of options, but I'm hitting a wall.


texts:
Added a sidebar to my Leaflet map that is toggled by a click on a feature and contains all feature properties.

On local devices (Firefox/Chrome) everything works fine. My problems shows up only on mobile devices: If the sidebar is toggled and displayed on the map, I can't click on any other buttons on my map. They are blocked. As long as the the sidebar is shown the only thing I can do is panning through the map. When the sidebar gets closed the buttons are working again.

I tried stuff like:

sidebar.BringToBack();

and

examplebutton1.BringToFront();

but this didn't work.

Maybe it's a simple stylesheet problem because of mobile usage? For example an option that allows me to click on clickable items (marker, buttons) that are covered by the sidebar. I already set opacity to 0.5 so that I can see covered markers and features but they are not clickable.


texts: Distance from one city to others by a mutual ID

texts:
I don't understand spatial.data at all. I have been studying but I'm missing something.
What I have: data.frame enterprises with the columns: id, parent_subsidiary, city_cod.
What I need: the mean and the max distance from the parent's city to the subsidiary cities.
Ex:
    id         |     mean_dist     | max_dist
 1111          |         25km      |     50km
 232           |        110km      |    180km
 333           |          0km      |      0km

What I did :
library("tidyverse")
library("sf")
# library("brazilmaps")   not working anymore
library("geobr")

parent <- enterprises %>% filter(parent_subsidiary==1)
subsidiary <- enterprises %>% filter(parent_subsidiary==2)

# Cities - polygons
m_city_br <- read_municipality(code_muni="all", year=2019)

# or shp_city<- st_read("/BR_Municipios_2019.shp")

# data.frame with the column geom
map_parent  <- left_join(parent, m_city_br, by=c("city_cod"="code_muni"))
map_subsidiary <- left_join(subsidiary, m_city_br, by=c("city_cod"="code_muni"))

st_distance(map_parent$geom[1],map_subsidiary$geom[2]) %>% units::set_units(km)
# it took a long time and the result is different from google.maps
# is it ok?!

# To do by ID -- I also stucked here

distance_p_s <- data.frame(id=as.numeric(),subsidiar=as.numeric(),mean_dist=as.numeric(),max_dist=as.numeric())

id_v <- as.vector(parent$id)

for (i in 1:length(id_v)){

  test_p <- map_parent %>% filter(id==id_v[i])
  test_s <- map_subsidiary %>% filter(id==id_v[i])
  total <- 0
  value <- 0
  max <- 0
  l <- 0

  l <- nrow(test_s)

      for (j in 1:l){

         value <- as.numeric(round(st_distance(test_p$geom[1],test_s$geom[j]) %>% units::set_units(km),2))

         total <- total + value
         ifelse(value>max,max<-value,NA)
      }

  mean_dist <- total/l
  done <- data.frame(id=id[i],subsidiary=l,mean_dist=round(mean_dist,2),max_dist=max)
  distance_p_s <- rbind(distance_p_s,done)

  rm(done)

}
}

Is it right?
Can I calculate  the centroid of the cities and than calculate the distance?
Data example: structure(list(id = c("1111", "1111", "1111",  "1111", "232", "232", "232", "232", "3123",  "3123", "4455", "4455", "686", "333", "333",  "14112", "14112", "14112", "3633", "3633","77172","77172"), parent_subsidiary = c("1","2", "2", "2", "1", "2", "2", "2", "1", "2", "1", "2", "1", "2",  "1", "1", "2", "2", "1", "2","1","2"), city_cod = c(4305801L,4202404L, 4314803L, 4314902L, 4318705L, 1303403L, 4304507L, 4314100L,  2408102L, 3144409L, 5208707L, 4205407L, 5210000L, 3203908L, 3518800L,  3118601L, 4217303L, 3118601L, 5003702L, 5205109L,4111407L,4110102L)), row.names = c(NA,  22L), class = "data.frame")
PS: this is Brazilian cities
https://github.com/ipeaGIT/geobr/tree/master/r-package


texts: How to use stored procedure in Geoserver?

texts:
I have a MS SQL Server 2008 R2 and have a stored procedure in here. I want to create a layer in Geoserver using this. I try:

EXECUTE [dd_base].[dbo].[Request BaseSration]
  'param1',
  'param2',
  'param3',
  'param4'
GO

And get error near key word EXECUTE.

In Geoserver logs:

ERROR [geoserver.web] - Failed to create feature type
java.io.IOException: Error occurred building feature type

What am I doing wrong?


texts: Saving data of original layer while buffering in QGIS?

texts:
I'm using QGIS 2.18.27.

I have the following problem: I am trying to create some buffers in the city center, as you can see. And I need to calculate some statistical parameter on every buffer.

As you can see, the Buffer5 area should be used to calculate some statistical parameter, but after I save these Buffers into shapefile, I found that these buffer.shp do not have the original data.

However, I found some other data by using "Identify features" tool, as the picture shows.

My question is: Is there actually a way to keep these original data in buffer.shp, that I can use this data to calculate? Or can the original shapefile be split into these buffers, so that I can keep these data?


texts: When running gdalinfo on an SRTM file, I get the corner coordinates. However, when I attempt to get the elevation at those coordinates using gdallocationinfo -wgs84 file.hgt {corner coordinates} I get that they are off the file, why is this? It seems to me that if I got the extents I should be able to query those points for their corresponding elevation. How would I go about getting the inclusive corner coordinates?


texts: How to edit feature in SQLServer geodatabase

texts:
I am using ArcMap 10.2.1 and ArcObjects 10.2.1, I just want to connect to a geodatabase and edit, add, delete a feature from a feature class using ArcObjects, I tried to connect to database and it works fine but I can't start editing and adding a feature to it.

Here is my code:

public void open_Workspace()
    {
    Type factoryType = Type.GetTypeFromProgID("esriDataSourcesGDB.SqlWorkspaceFactory");
    IWorkspaceFactory workspaceFactory = (IWorkspaceFactory)Activator.CreateInstance(factoryType);
    IWorkspaceFactory2 workspaceFactory2 = workspaceFactory as IWorkspaceFactory2;
    // Build a connection string.
    String[] connectionProps =
{
    "dbclient=SQLServer", "serverinstance=CRAZYVIRUS-HP",
        "database=land_gis", "authentication_mode=OSA"
};
    String connString = String.Join(";", connectionProps);

    IWorkspace workspace = workspaceFactory2.OpenFromString(connString, 0);
    ESRI.ArcGIS.Geodatabase.IFeatureWorkspace featureWorkspace = (ESRI.ArcGIS.Geodatabase.IFeatureWorkspace)workspace; // Explict Cast
    ESRI.ArcGIS.Geodatabase.IFeatureClass featureClass = featureWorkspace.OpenFeatureClass("DBO.test");

    ESRI.ArcGIS.Carto.IFeatureLayer featureLayer = new ESRI.ArcGIS.Carto.FeatureLayerClass();

    featureLayer.FeatureClass = featureClass;
    featureLayer.Name = featureClass.AliasName;
    featureLayer.Visible = true;

    IWorkspaceEdit workspaceEdit = (IWorkspaceEdit)workspace;
    IWorkspaceEdit2 workspaceEdit2 = (IWorkspaceEdit2)workspaceEdit;
    IMultiuserWorkspaceEdit muWorkspaceEdit = (IMultiuserWorkspaceEdit)workspace;

    Make sure that non-versioned editing is supported. If not, throw an exception.
   if (!muWorkspaceEdit.SupportsMultiuserEditSessionMode
       (esriMultiuserEditSessionMode.esriMESMNonVersioned))
   {
       throw new ArgumentException(
           "The workspace does not support non-versioned editing.");
   }

   // Start a non-versioned edit session.
   muWorkspaceEdit.StartMultiuserEditing
   (esriMultiuserEditSessionMode.esriMESMNonVersioned);

    featureClass.CreateFeature();

    // Stop the edit session. The saveEdits parameter indicates the edit session
    // will be committed.
    workspaceEdit.StopEditing(true);

}

When I run above code I got this error Unable to cast COM object of type 'System.__ComObject' to interface type 'ESRI.ArcGIS.Geodatabase.IWorkspaceEdit2'. This operation failed because the QueryInterface call on the COM component for the interface with IID '{199A1F9D-435A-4118-9B3F-4E0B8F984AA0}' failed due to the following error: No such interface supported (Exception from HRESULT: 0x80004002 (E_NOINTERFACE)).

I tried to enable versioning or disable it but there is not any chance to make it work...  Can you help me find the solution into this error?
Is there any thing wrong with my connection?


texts: MapFish PDF: Geoserver WMS layers not generated as tree structure

texts:
I am using the mapfish printing module in my web map application.In that everything is working fine but layers not generated as tree structure in output PDF file and "ODC_Graphics_Data" wms layer not printed.  I have used below code to generate pdf.

System: Windows XP / Tomcat7 / Geoserver 2.4.0/printing-2.1-SNAPSHOT.jar

Kindly help me to fix this. Thanks in advance.

  http://localhost:8080/geoserver/pdf/print.pdf?spec={
  "layout":"A3 portrait",
  "srs":"EPSG:404000",
  "units":"m",
  "dpi":300,
  "outputFilename": "ggg",
  "layers":[
    {
        "baseURL":"http://xxx.xxx.xx.xx:8080/geoserver/cite/wms",
        "opacity":1,
        "singleTile":false,
        "type":"WMS",
        "layers":["Oxford_Map"],
        "format":"image/png",
        "styles":[],
        "overview":true
    },
    {
        "baseURL":"http://xxx.xxx.xx.xx:8080/geoserver/cite/wms",
        "opacity":1,
        "singleTile":true,
        "customParams":{"CQL_FILTER": "User_id = 234"},
        "type":"WMS",
        "layers":["ODC_Graphics_Data"],
        "format":"image/png",
        "styles":[],
        "overview":false
    },
    {
        "baseURL":"http://xxx.xxx.xx.xx:8080/geoserver/cite/wms",
        "opacity":1,
        "singleTile":true,
        "type":"WMS",
        "layers":["Plymouth_Buildings"],
        "format":"image/png",
        "styles":[],
        "overview":false
    },
    {
        "baseURL":"http://xxx.xxx.xx.xx:8080/geoserver/cite/wms",
        "opacity":1,
        "singleTile":true,
        "type":"WMS",
        "layers":["Plymouth_Hard_Surface_Block_Paving"],
        "format":"image/png",
        "styles":[],
        "overview":false
    },
    {
        "baseURL":"http://localhost:8080/geoserver/cite/wms",
        "opacity":1,
        "singleTile":true,
        "type":"WMS",
        "layers":["Plymouth_QUERIES"],
        "format":"image/png",
        "styles":[],
        "overview":false
    },
    {
        "baseURL":"http://xxx.xxx.xx.xx:8080/geoserver/cite/wms",
        "opacity":1,
        "singleTile":true,
        "type":"WMS",
        "layers":["Plymouth_Street_Lights"],
        "format":"image/png","styles":[],
        "overview":false
    }
],
 "pages":[
    {
        "center":[248744.93342196,55350.945],
        "mapTitle":"",
        "company":"testuser2_portrait",
        "comment":"",
        "scale":439,
        "rotation":0
    }
] }

Please find screenshot for output PDF:


texts: Editing Python toolbox template in ArcPy is not working?

texts:
Using ArcGIS 10.4 for Desktop, I am, for the first time, trying to work with a Python Toolbox (*.pyt).

According to Creating a new Python toolbox documentation I am trying to edit The Python toolbox template, but it seems not working.

I have edited the template like this:

import arcpy

class StackProfile(object):
    def __init__(self):
        """Define the toolbox (the name of the toolbox is the name of the
        .pyt file)."""
        self.label = "TOOLBOX LABEL PROPERTY"
        self.alias = "TOOLBOX ALIAS PROPERTY"

        # List of tool classes associated with this toolbox
        self.tools = [SP]

class SP(object):
    def __init__(self):
        """Define the tool (tool name is the name of the class)."""
        self.label = "Stack Profile for Multiple Elevations"
        self.description = ""
        self.canRunInBackground = False

I tried to rename the toolbox name, label and alias. But in toolbox Properties it seems only the name is changed not label, or alias or even not the tool name.

I have a screen shot of my properties window:

It also did not change the tool name.

What could be the possible reason for this problem?


texts: Translating one raster's contents closer to another’s with each at different scales using QGIS?

texts:
I'm working with multiple raster images, each in the Albers Equal Area Conic projection: one of the mainland US, one of Alaska, and one of Hawaii.

I'd like to translate the Hawaii layer to be closer to mainland US, reduce its size, and then output the resultant image as a .tiff file.

Is this is possible, is there a way to do it in QGIS?


texts: Table to point/raster in Google Earth Engine

texts:
I am trying to visualize the following table in google Earth Engine:

What I need is either:

Visualize it as a point feature, drawing the value 'tipEst_2012' with different colors. Radius value around 180 m.
Export it to a raster with raster value = tipEstiu_2012 and pixel size=180 m.

I've been looking for scripts to do so but I couldn't find any. You can find my data here:
https://code.earthengine.google.com/?asset=users/oriolbaena/tipestiu

My script so far is:

var tipestiupoints=ee.FeatureCollection("users/oriolbaena/tipestiu");
var tipestiuraster = tipestiupoints
  .filter(ee.Filter.notNull(['tip_Est2012']))
  .reduceToImage({
    properties: ['tip_Est2012'],
    reducer: ee.Reducer.first()
});

Map.addLayer(tipestiuraster,{},'tipestiu')


texts: gdalwarp to change block size also changes pixel values

texts:
I'm working with a GeoTIFF from http://www.fao.org/geonetwork/srv/en/main.home?uuid=ba4526fd-cdbf-4028-a1bd-5a559c4bff38

The original has data values from 0-11:

$ gdalinfo -hist ./data/FAO/glc_shv10_DOM.Tif
Driver: GTiff/GeoTIFF
Size is 43200, 21600
Metadata:
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  TIFFTAG_SOFTWARE=Arc/Info
  TIFFTAG_XRESOLUTION=119.999
  TIFFTAG_YRESOLUTION=119.999
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0,21600.0)
Upper Right (43200.0,    0.0)
Lower Right (43200.0,21600.0)
Center      (21600.0,10800.0)
Band 1 Block=43200x3 Type=Byte, ColorInterp=Gray
  256 buckets from -0.5 to 255.5:
  181470400 977689 29454097 33129685 71015243 18337497 3857498 122159 20544229
  28053726 11242330 534915447 0 0 0 ... (more zeros)

I tried to reorganize it with a 256x256 block size for subsequent steps in my processing pipeline, using:
gdalwarp -ot Byte -co COMPRESS=DEFLATE -co TILED=YES -s_srs "EPSG:4326" -t_srs "EPSG:4326" -to SRC_METHOD=NO_GEOTRANSFORM glc_shv10_DOM.Tif new_file.tif

I had used "-to SRC_METHOD=NO_GEOTRANSFORM" because the original file lacks information about its coordinate system but looks like EPSG:4326.

The new file has very different pixel values: 0, 19, 26, 30, 34, 37, 39, 41, 43, 45, 47, 48.

Driver: GTiff/GeoTIFF
Size is 43200, 21600
Coordinate System is:
GEOGCRS["WGS 84",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["geodetic latitude (Lat)",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic longitude (Lon)",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    USAGE[
        SCOPE["unknown"],
        AREA["World"],
        BBOX[-90,-180,90,180]],
    ID["EPSG",4326]]
Data axis to CRS axis mapping: 2,1
Origin = (0.000000000000000,21600.000000000000000)
Pixel Size = (1.000000000000000,-1.000000000000000)
Metadata:
  AREA_OR_POINT=Area
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  TIFFTAG_SOFTWARE=Arc/Info
  TIFFTAG_XRESOLUTION=119.999
  TIFFTAG_YRESOLUTION=119.999
Image Structure Metadata:
  COMPRESSION=DEFLATE
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (       0.000,   21600.000) (  0d 0' 0.01"E,Invalid angle)
Lower Left  (   0.0000000,   0.0000000) (  0d 0' 0.01"E,  0d 0' 0.01"N)
Upper Right (   43200.000,   21600.000) (Invalid angle,Invalid angle)
Lower Right (   43200.000,       0.000) (Invalid angle,  0d 0' 0.01"N)
Center      (   21600.000,   10800.000) (Invalid angle,Invalid angle)
Band 1 Block=256x256 Type=Byte, ColorInterp=Gray
  256 buckets from -0.5 to 255.5:
  181470400 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 977689 0 0 0 0 0 0
  29454097 0 0 0 33129685 0 0 0 71015243 0 0 18337497 0 3857498 0
  122159 0 20544229 0 28053726 0 11242330 534915447 0 0 (more zeros)

I don't understand why the pixel values changed. Other questions where data values changed involve a scale_factor in the metadata, but my file does not have a scale_factor.

I also need to figure out what I'm missing about adding a coordinate system, the Corner Coordinates are still pixels but now being interpreted as degrees.


texts: When editing feature-linked annotation with font color white, the annotation disappears

texts: When moving or editing the annotation, the annotation disappears. We are editing on a black background, so the text isn't just blending in to the background. The attribute table does not show any abnormality, but the "Zoom to Selected" option is grayed out. It can be selected in the table, but does not display in the display window. Anyone have this problem? Found any solutions?