texts:
I am trying to create a CSV file of all spatial data in a directory structure using arcpy walk as in the below code.
However, when the code is run on a large file system it keeps stopping after a limited number of records. The problem was occurring in different places in different directory structures. After running it numerous times on the same dir structure, it was always occurring in the same place. There is nothing wrong with the files themselves.

def fileListing():
    try:
        with open(os.path.join(outputLocation, fileName), "w") as csv:
            # write out the column headers
            csv.write("Layer Name" + "," + "Layer Type" + "," + "Location" + "," + "Spatial Reference" + "\n" )

            walk = arcpy.da.Walk(inputLocation, datatype="Any", type="All")

            for root, dirs, files in walk:
                for f in files:
                    desc = arcpy.Describe(root + "\\" + f)
                    descDataType = desc.dataType
                    print(descDataType)
                    try:
                        descSR = removeChars(desc.SpatialReference.name)
                    except AttributeError:
                            continue
                    descName = removeChars(desc.name)
                    descCatPath = removeChars(desc.catalogPath)
                    #write info to text file
                    csv.write(descName + "," + descDataType + "," + descCatPath + "," + descSR + "\n")
        csv.close()
    except Exception as e:
    print("type error: " + str(e))
    print(traceback.format_exc())

Having run the script on the same directory structure numerous times it keeps failing at record number 225. There doesn't seem to be anything wrong with that particular record. It is a shapefile which loads into ArcMap fine.

The following error keeps occurring in PyScripter.


texts: ArcGIS Geoprocessing - how to set a condition for data

texts:
I would like to write a script in Python which takes data and does some things with it.

So far, I have written the following code :

import arcpy
arcpy.env.workspace = "C:\wody.gdb"
arcpy.env.overwriteOutput = True
fcList = arcpy.ListFeatureClasses()
bufferList = []
for fc in fcList:
   if fc == "RZEKI":
       arcpy.Buffer_analysis(fc, fc + "Buffer", "200 meters")
       bufferList.append(fc + "Buffer")

Now, I would like to take RIVERS (RZEKI) only if river is longer than SOMETHING (they do not give data in meters but in units, so I just want to see an example seems like river_long > WHATEVER)

Thank you in advance :)


texts: Integrating standalone GeoWebcache with client?

texts:
I have installed GeoServer on machine 1 (windows) and GeoWebCache on machine 2 (windows). What I need to achieve is tile caching so that the client can be served maps from cache there by enhancing the performance of the client.

I am still not sure if My client should connect to GeoWebcache or Geoserver for fetching maps? In either of the case, what configurations are required and at what place to have my client leverage the cache feature.

Requesting some assistance/guidance on the above queries which is purely from integration side of it.


texts: Visualizing raster files stored in postgresql DB in Openlayers or Openstreetmap Layers

texts:
I am relatively new to this domain. I want to make a web app that will show raster data on openstreetmap layers or openlayers. So far I have successfuly managed to upload the files to postgresql.
My main question is how are these files supposed to be visualized on an interractive map. Looking forward to your replies.


texts: Display EPSG:3844 ArcGIS WMTS service in OpenLayers: origin axes in wrong order

texts:
I'm trying to display a EPSG:3844 WMTS service from an ArcGIS server in OpenLayers 3.20.1. In order for it to be displayed correctly, I need to manually swap the origin coordinates found in the capabilities file. This is my code:

proj4.defs('EPSG:3844', '+proj=sterea +lat_0=46 +lon_0=25 +k=0.99975 +x_0=500000 +y_0=500000 +ellps=krass +towgs84=2.3283,-147.0416,-92.0802,-0.30924979,0.32482188,0.49730012,5.68907711 +units=m +no_defs');

const response = await fetch(serviceUrl);
const parser = new ol.format.WMTSCapabilities();
const capabilities: any = parser.read(await response.text());
const options = ol.source.WMTS.optionsFromCapabilities(capabilities, {
    layer: capabilities.Contents.Layer[0].Identifier,
    matrixSet: capabilities.Contents.TileMatrixSet[0].Identifier,
});

const originalTileGrid = options.tileGrid;
options.tileGrid = new ol.tilegrid.WMTS({
    origins: originalTileGrid.getResolutions().map((_, z) => {
        const o = originalTileGrid.getOrigin(z);
        const transposedOrigin: ol.Coordinate = [o[1], o[0]];
        return transposedOrigin;
    }),
    resolutions: originalTileGrid.getResolutions(),
    matrixIds: originalTileGrid.getMatrixIds(),
    tileSizes: originalTileGrid.getResolutions().map((_, z) => {
        return originalTileGrid.getTileSize(z);
    }),
});

and this is the relevant portion from WMTSCapabilities.xml:

<TileMatrixSet>
    <ows:Title>TileMatrix using 0.28mm</ows:Title>
    <ows:Abstract>
        ...
    </ows:Abstract>
    <ows:Identifier>default028mm</ows:Identifier>
    <ows:SupportedCRS>urn:ogc:def:crs:EPSG::3844</ows:SupportedCRS>
    <TileMatrix>
        <ows:Identifier>0</ows:Identifier>
        <ScaleDenominator>3657607.3152141706</ScaleDenominator>
        <TopLeftCorner>754000.0 134000.0</TopLeftCorner>
        <TileWidth>256</TileWidth>
        <TileHeight>256</TileHeight>
        <MatrixWidth>4</MatrixWidth>
        <MatrixHeight>3</MatrixHeight>
    </TileMatrix>
    ...
</TileMatrixSet>

From the WMTS specification (pdf):

  CRS will be inherited from the supportedCRS parameter of the parent TileMatrixSet. The order of these axes, shall be as specified by the supported CRS. These are the precise coordinates of the top left corner of top left pixel of the 0,0 tile.

So my questions are:

Who's doing something wrong here?

Is it ESRI, should they be sending the axes the other way around? Is it OpenLayers, should they switch the axes?

Or is it me and I'm missing something?


texts: Generating polygons with streets as borders using QGIS?

texts:
I have searched online for almost a week and can not find an answer to my issue.

I need to generate polygons from street borders in a given city.

I have downloaded data from openstreetmaps and imported it into QGIS.

Is there a way to generate polygons and use streets as borders?


texts: Modifying File GDB Feature Class Features Using QGIS

texts:
I've installed the File GDB API using the OSGEO4W advanced installation when installing QGIS. I've had a bit of trouble with it initially and have posted a few questions with some helpful Stack members. I can now open File geodatabase feature classes in QGIS and start an edit session, but the only thing I can do is create new features using the 'Add Feature' tool. I cannot edit existing features because the 'Move Features', 'Node Tools', 'Delete Selected Features' and the 'Cut, Copy and Paste Features' Tools are all greyed out.

I've tested this using both ArcGIS 10.0 and 10.1 file geodatabases with the same result. The tools are greyed out except the 'Add Feature' tool even though the layer is selected and active.

I tested to see if I can get the tools to be active by using an ESRI Shapefile and this works which leads me to believe that the API isn't working as I think it should.....or maybe I'm just being a bonehead. Has anyone else experienced this behavior? Should I be able to modify existing features or should I only be able to create new ones?


texts: Merge images first, or georeference first?

texts:
I'm wanting to georeference and vectorise some features from scans of the OS 1896 Town Plans that I have purchased from the NLS (http://maps.nls.uk/os/london-1890s/info.html).  The scans provided are not georeferenced and when I try to crop the borders they are not quite square so will need rotating first.  The problem is the features lie at the intersection of 4 separate sheets and I'm uncertain of the best / most accurate approach to take to merge them:

1) Treat each image separately: Rotate and crop each image in Gimp, georeference in QGIS, then "nudge" the individual GCP points to try get them to line up properly, then merge with gdal_merge?

2) Merge first: Rotate, crop and merge the images together in Gimp, then georeference the whole image in QGIS

3) Some magic automated solution I'm not aware of :-)

Thanks!


texts: Manually install WFS-T plugin to implement with Leaflet

texts:
I'd like to use this Flexberry WFS-T Leaflet plugin, and am required to manually install it.
I've put it in my assets folder in this way:

Imported leaflet and three js files:
    <script src="https://unpkg.com/leaflet@1.7.1/dist/leaflet.js"></script>
    <script src="assets/src/WFS.js"></script>
    <script src="assets/src/WFST.js"></script>
    <script src="assets/src/GML.js"></script>

Initialized this way:
var rooturl = 'http://localhost:8080/geoserver/ged/ows';

var wfst = new L.WFST({
            url: rooturl,
            typeNS: 'ged',
            typeName: 'GED_IDV',
            crs: L.CRS.EPSG4326,
            geometryField: 'geom',
            style: {
              color: 'blue',
              weight: 2
            }
          });

But, something isn't working, and I suspect the installation, but can't find anything in the docs.
Uncaught TypeError: Cannot read property 'GML' of undefined
    at initialize (WFS.js:43)
    at initialize (WFST.js:7)
    at new i (Class.js:22)
    at HTMLDocument.<anonymous> ((index):590)
    at e (jquery-3.5.1.min.js:2)
    at t (jquery-3.5.1.min.js:2)


texts: Why do 2 instances of `CreateUniqueName` behave differently?

texts:
The looping process I'm working on outputs two files (table and image), and I'm trying to name them with matching numbers. I'm using arcpy.CreateUniqueName to generate new filenames into a GDB each time the loop runs. The problem is, while the first CreateUniqueName usage produces new filenames in the expected sequence each time I run it, the second usage just keeps overwriting the original file, with filenames never progressing past the starting filename.

Here's the code, followed by the resulting print statements from 3 sequential runs:

out_prof = r'F:\Joe_School\Thesis\scripts\Jerry\scratch.gdb\testprof'
patch = r'F:\Joe_School\Thesis\data\grid_extract\mask400_1'

import arcpy, os
from arcpy import env
import arcpy.ddd as DDD
import arcpy.management as DM

env.overwriteOutput = True
env.outputCoordinateSystem = arcpy.Describe(patch).spatialReference
env.workspace = os.path.dirname(out_prof)

arcpy.CheckOutExtension('3D')
tablename = arcpy.CreateUniqueName('prof_table')
graphname = arcpy.CreateUniqueName('prof_graph')
prof = DDD.StackProfile(out_prof, patch, tablename, graphname)
DM.SaveGraph(graphname, graphname + '.gif')
arcpy.CheckInExtension('3D')
print "Created files:" + "\n    " + tablename + "\n    " + graphname

>>>
Created files:
    F:\Joe_School\Thesis\scripts\Jerry\scratch.gdb\prof_table
    F:\Joe_School\Thesis\scripts\Jerry\scratch.gdb\prof_graph
>>>
Created files:
    F:\Joe_School\Thesis\scripts\Jerry\scratch.gdb\prof_table0
    F:\Joe_School\Thesis\scripts\Jerry\scratch.gdb\prof_graph
>>>
Created files:
    F:\Joe_School\Thesis\scripts\Jerry\scratch.gdb\prof_table1
    F:\Joe_School\Thesis\scripts\Jerry\scratch.gdb\prof_graph

Why are these seemingly identical statements behaving differently?


texts: how to prepare tiffs to create a netcdf file

texts:
I´d like to create a netcdf file containing a mean and stddev for some variables for 12 months.

At this moment each geotiff (270MB) is either a mean or stddev of a month.

Now my problem is that I have no idea how to proceed from here, I have read plenty of documents about how to handle netcdf (How to extract, adapt and convert to other formats).

Is there a way to either use gdal or tools like ncdump to create a netcdf base file and then add each geotiff with extra metadata (variable, month and stddev or mean)? Or should the geotiff already have this metadata before using gdal to create a netcdf?

Here something similar is being done but I am not sure where in the python code the tiffs are imported and how each tiff is set to a certain date (My python skills are not up to par).


texts: Generating ground control points with GRASS?

texts:
I want to find out matching ground control points in two images for auto georeferencing and grass manual says that i.points.auto will do the job.

I am using Grass 7.2 but couldn't find this plugin.


texts: Converting / interpolating delimited text csv points grid to surface

texts:
I have 12 delimited text csv files with regular pseudo-gridded points for temperature and salinity which look like this when zoomed in (inc distances) and this when zoomed out. I'm trying to convert them to a smooth surface layer - assumedly interpolated - but using the interpolation plugin on a vector layer saved from one of the csv's produced this kind of output, irrespective of how many columns I use.

How do I turn my grids of dots into a nice smooth surface that I can use colour interpolation on?


texts: Add multiple WFS layers to Leaflet map from GeoServer using JQuery Ajax calls

texts:
I am attempting to add multiple WFS layers from GeoServer into a Leaflet map using JQuery ajax calls.  It works if I set async: false, but this is deprecated and should not be used.  If I remove async: false, not all layers load.  I've enabled jsonp in GeoServer and on the webserver, and that part seems to work. How can I load multiple WFS layers?

var map = L.map('map', {
    center: [<lat>, <lon>],
    zoom: 12,
    minZoom: 11,
    maxZoom: 21
});

var tileLayer_OpenStreetMap_Mapnik = L.tileLayer('http://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png', {
    maxZoom: 22,
    attribution: '&copy; OpenStreetMap'
});

var baseMaps = {
    "Base Map (OSM Colour)": tileLayer_OpenStreetMap_Mapnik

};

tileLayer_OpenStreetMap_Mapnik.addTo(map);

var popup = L.popup();
function onEachFeature(feature, layer) {
    var popupcontent = [];
    for (var prop in feature.properties) {
        popupcontent.push('
' + prop + ': 	' + feature.properties[prop] + '
');
    }
    layer.bindPopup('
' + popupcontent.join("
") + '
');
};

var owsrootUrl = 'http://myurl/geoserver/ows';

function doGetJson(layerName) {
    var defaultParameters = {
        service : 'WFS',
        version : '2.0',
        request : 'GetFeature',
        typeName : layerName,
        outputFormat : 'text/javascript',
        format_options : 'callback:getJson',
        SrsName : 'EPSG:4326'
    };

    var parameters = L.Util.extend(defaultParameters);
    var URL = owsrootUrl + L.Util.getParamString(parameters);
    var WFSLayer = null;
    return $.ajax({
        url : URL,
        dataType : 'jsonp',
        async: false,
        jsonpCallback : 'getJson',
        index: i,
        success : handleJson
    });
};

function handleJson(data) {
    WFSLayer = L.geoJson(data, {
        style: function (feature) {
            //style the objects based on type of geometry
            switch (feature.geometry.type)
            {
                case 'Polygon' : return {
                        color: '#000000',
                        fillColor: '#000000',
                        fillOpacity: 1,
                        weight: 2
                };
                case 'LineString' : return {
                        color: '#000000',
                        opacity: '#000000',
                        weight: 2
                };
            }
        },
        crs: L.CRS.EPSG4326,
        onEachFeature: onEachFeature
    });
    map.addLayer(WFSLayer);
};

var layers = "<WORKSPACE NAME>:<LAYER NAME1>,<WORKSPACE NAME>:<LAYER NAME2>,<WORKSPACE NAME>:<LAYER NAME3>".split(",");
for (var i = 0; i<layers.length ; i++ )
{
    doGetJson(layers[i]);
}


texts:
i've been working on my project to create QGIS plugin for raster analysis purpose.

It will be easier if we use gdal library.

Are we able to use gdal library to develop QGIS plugin?

If so, how do we set up the gdal package inside the qgis2 plugin environment?


texts: Using Overpass Turbo to return all peaks over certain elevation

texts: I would like to use Overpass Turbo to return all peaks over a certain elevation. I tried using the wizard, but I cannot figure out how to use the greater than (>) operator.


texts: OpenGeo Suite installed on Ubuntu 14.04, unable to proceed/delete issues with dependencies

texts:
I have installed the OpenGeo Suite on my "Trusty Tahr". I am able to access my Tomcat server as well as Geoserver on http://localhost.

The issue I'm dealing with currently stems from dependencies related PostGIS and PostgreSQL versions that didn't install properly.

As a result I am unable to install any additional applications nor have all my dependencies updated as such.

I have tried several methods from trying to remove the applications at issue through the Software Center
,

I've attempted running apt-get install -f as well as other troubleshooting techniques that I found search through the web.

Any thoughts how to rectify this?


texts:
I have a Landsat-NDVI-time series from 2013 - 2017 with 23 observations per year (115 scenes in total). My goal is to get a smooth time series for a selected single pixel (maybe by using a Savitzky-Golay filter) but also get a trend map for all pixels in my research area.
I first calculated the missing dates, created a raster with NA values for each of these scenes, and stacked all of the files together.

Example of the output for one pixel:

    Date          NDVI
 2013-01-02       NA
 2013-01-18       NA
 2013-02-03       NA
  [...]
 2017-08-25  0.4551642537
 2017-09-10  0.3888225853
 2017-09-26  0.6013333797

Now I would like to calculate a smooth time series for this pixel. If I interpolate before applying the time series I get an irregular time series because the first year than only contains 19 observations. I tried following this tutorial, but could not get it to work.
How can I handle this problem and more importantly, how can I calculate a pixel-wise regression for the whole scene? I found this, which seems to cover my needs, but can I apply it to an irregular time series?


texts: Setting style for categorized vector in PyQGIS?

texts:
I'm writing a plugin in QGIS3  which produces a vector file from a thresholded raster (i.e. polygons all have a value of 1 or 0). So far so good. I'd now like to load it in QGIS and display the results so that the 1 values in the vector are displayed (as a filled color) and the 0 values aren't displayed at all. This should just be a display of a categorized vector but I can't figure out the best approach in QGIS3.

Can anyone point me in the right direction, please? I don't have any code to share on this point. I got as far as

processing.run("gdal:polygonize", {'INPUT': self.threshold_raster_obj, 'BAND':1, 'FIELD':'DN','EIGHT_CONNECTEDNESS':False,'OUTPUT': self.out_vectorFilename})
current_layer = self.iface.addVectorLayer(self.out_vectorFilename, "Polygonize_result", "ogr")

i.e. I can polygonize my raster, save it and load it into the main window but I've no idea where to start with the display options. I've looked at some of the examples but I've not been able to understand the right approach for polygons. Unfortunately some things are not intuitive for novices like myself.


texts: How to create an image of mean values from multiple NDVI images?

texts: I have downloaded MODIS NDVI images of 12 months. I would like to create just 1 image with the mean values. I tried using IDRISI's 'NDVI compositing utility', but I get an error. Is there a way to do this in QGIS?


texts: Keep / append attributes to buffered PolyLines in QGIS

texts:
I have a shapefile containing multiple line files with key attributes. I buffered each using cap_style=3 to get the correct buffer, but all of the attributes are lost. I'm either looking for added on code to retain the attributes from the original shapefile or a new script to append the old attributes to the new shapfile
lines= gpd.read_file(r"D:\brett\Documents\CAFE\Folder_Test3\Lines\Lines_comp\Lines_comp.shp")
Line_buffer = lines.geometry.apply(lambda g: g.buffer(1, cap_style=3))
line_buffer.to_file(r"D:\brett\Documents\CAFE\Folder_Test3\Lines\Lines_comp\Lines_buffer3.shp")


texts: Different numbers between count points in polygon and select by location

texts:
When I'm running the points in polygon in QGIS I get one result (125) and when I try to select the points to extract as a layer, select by location, I get another result (213).

If I extract the points and re-map them, they all seem to fall within the polygon boundary, any reason as to why they might be different, and which one is more of a credible source?

Running QGIS 3.4.9


texts: Unexpected DotSpatial reprojection results - OSGB36 to 4326

texts:
I am trying to reproject from OSGB1936 to 4326 using DotSPatial.
I'm doing this as follows ('shapefile' is my loaded shapefile) pulling the source projection from the shapfile which I believe it picks up from the .prj file:
sourceProjection = shapefile.Projection;
destinationProjection = ProjectionInfo.FromEpsgCode(4326);
.
.
double[] test = new double[2];
test[0] = 360216.44430417975;
test[1] = 176543.59223759052;
Reproject.ReprojectPoints(oi, z, sourceProjection, destination, 0, 1);

I get the following results:
Before (OSGB36): 360216.44430417975, 176543.59223759052
After (4326): -2.5730259372697137, 51.4859574135158411
When I use online converters I get slight variations but about this:
After: -2.5743517, 51.4864342
This puts my points out by a 100 yards or so.
I've tried ignoring the .prj file that comes as part of the shapefile and forcing the projection to OSGB4326, but I get the same result.
Does anyone have any ideas what might be causing this offset?
I tried defining my projection direct from ESRI strings:
double[] z = new double[2];
ProjectionInfo src = ProjectionInfo.FromEsriString("PROJCS[\"OSGB 1936 / British National Grid\",GEOGCS[\"OSGB 1936\",DATUM[\"D_OSGB_1936\",SPHEROID[\"Airy_1830\",6377563.396,299.3249646]],PRIMEM[\"Greenwich\",0],UNIT[\"Degree\",0.017453292519943295]],PROJECTION[\"Transverse_Mercator\"],PARAMETER[\"latitude_of_origin\",49],PARAMETER[\"central_meridian\",-2],PARAMETER[\"scale_factor\",0.9996012717],PARAMETER[\"false_easting\",400000],PARAMETER[\"false_northing\",-100000],UNIT[\"Meter\",1]]");
ProjectionInfo des = ProjectionInfo.FromEsriString("GEOGCS[\"GCS_WGS_1984\",DATUM[\"D_WGS_1984\",SPHEROID[\"WGS_1984\",6378137,298.257223563]],PRIMEM[\"Greenwich\",0],UNIT[\"Degree\",0.017453292519943295]]");
z[0] = 0;
z[1] = 0;
double[] points = new double[2];
points[0] = 360216.44430417975;
points[1] = 176543.59223759052;
Reproject.ReprojectPoints(points, z, src, des, 0, 1);

And I'm still getting the same offset as described above.
I'm thinking maybe there is some setting in DotSpatial that is throwing off the conversion here.


texts: How To Change Google earth engine Base Map Language

texts:
Google Earth Engine docs have an Article on how to Change the base map attributes
Yet not a word about changing the language of this base map.
Is there a way with the Maps API ?


texts:
I installed QGIS 3, and tried to execute tools like "sieve" "merge" and "polygonize". All of them gave me the error message:

GDAL command output: /bin/sh: (scriptname).py: command not found

I found the scripts in both /Username/miniconda3/bin and /Username/miniconda3/lib/python3.6/site-packages/GDAL-2.2.2-py3.6-macosx-10.9-x86_64.egg-info/scripts.

Checking env in Terminal shows that miniconda3/bin is in the path:

PATH=/usr/local/bin:/usr/bin:/bin:/usr/sbin:/sbin://Users/Username/miniconda3/bin

What could be causing the failure of executing those scripts?


texts: WPF: How to change the Map Point coordinate(lat/long) to screen point (x-y) on ArcGIS

texts: My project required me to perform a close traverse to plotting by placing a marker on that particular point. My problem now is that the traverse calculation are in x-y plane, but to plot a point using marker on ArcGIS require Lat/long. So my idea just know to get the current x-y point of current coordinate point then calculate, but how can i change from x-y point to map coordinate?


texts: Integrate geodatabase with other nonspatial databases in oracle 11g using ArcGis for Desktop 10.2

texts:
I have an oracle database namely xe which  I have spatially enabled. using the  Create Enterprise Geodatabase geoprocessing tool.

The application I am going to implement has an already created database/schema containing entitities which have gps data. For example it contains stations that I Have to project on the map.

I  want my the application to use, through the arcSDE, this non spatial information, as well as spatial information, layers, geometries, etc.

I have also created an sde and one schema-user (owner of a geodatabase).

How will I integrate these two schemas? should I create my non spatial enterprise schema inside of the schema-user schema and connect them with a sort of foreign keys?
Or use the CREATE DATABASE LINK, CREATE MATERIALIZED VIEW, and CREATE VIEW  privildges provided to the schema-user?

I am trying to understand the mentality behind the use of a geodatabase in accordance with a non spatial database.


texts:
In select statement how to avoid divide by zero in MapBasic:

select (col1/col2) "result"

If col2 is 0 it returns null.

How to use ISNULL or something in mapbasic?


texts: Load multipolygons from GeoJSON into GeoDjango model

texts:
I have a simple model:

from django.contrib.gis.db import models

class CountryShapes(models.Model):
    geonameid = models.PositiveIntegerField(primary_key=True)
    geom = models.MultiPolygonField()

and a GeoJSON file from here:
http://download.geonames.org/export/dump/shapes_simplified_low.json.zip

I am attempting to load the data into my model using this:

        import json
        from django.contrib.gis.geos import GEOSGeometry

        with open('shapes_simplified_low.json', 'r') as fd:
            data = json.load(fd)
            for feature in data['features']:

                geom = GEOSGeometry(str(feature['geometry']))
                geonameid = feature['properties']['geoNameId']

                country = CountryShapes( # throws error here
                    geonameid=geonameid,
                    geom=GEOSGeometry(geom))
                country.save()

and I get this error:

  TypeError: Cannot set CountryShapes SpatialProxy (MULTIPOLYGON) with
  value of type: class 'django.contrib.gis.geos.polygon.Polygon'

GEOSGeometry seems to create the geometry fine, its when I pass the geom argument to my model it throws the error.


texts: How to get directory from filepath in QGIS action

texts:
I have in QGIS layer that contains filepath (e.g. D:/file/28238/pic/D4D3-06_1.png) and I would like to get directory path only like D:/file/new/28238/pic.

Is it possible to do by expression in field calculator?


texts: How to calculate (summarize) atribute values of each point included in certain polygons?

texts:
in ArcMap 10.0 I have one point shapefile (some address layer with number of citizens) and one polygon shapefile (geomorphological segmentation). I need to add up citizens who live in certain area.

I can do select by location of course, but how can I automatically count all citizens from attribute in selected features (point) and add it to the polygon shapefile attribute to each polygon. Without manually summarizing fields from selected feature (there are plenty of polygons).

I would appreciate some hint on which tool should I use.


texts: matplotlib does not work in ArcGIS 10 and Eclipse (Python)

texts:
I have installed matplotlib (v. 1.0.1) which works with iPython quite good.
But if I insert the same script in ArcGIS 10 or Eclipse 3.7.0 (with PyDev) I get always errors. As I consider to use matplotlib to create graphs with ArcGIS Python Scripts developed in Eclipse thats both pretty unfortunate - so far eclipse works quite well for ArcGIS 10.

F.i. I tested this one:
http://matplotlib.sourceforge.net/examples/api/unicode_minus.html
In iPython it works perfect.

In Eclipse I get this error:
import matplotlib.pyplot as plt
ImportError: No module named pyplot

And in ArcGIS I get this one:
Running script Script1...
: bad screen distance "640.0"
Failed to execute (Script1).

Do I still have something to configure?

I have tested several other things and posted in other forums as well (like ESRI). As there seem to be no solution to the problem maybe someone does have a recommendation about which other tool might be suitable to create graphs while running python with ArcGIS 10. Besides regular pie and column charts I need in particular the possibility to create star plots/ radar charts / amoeba diagram - most likely something like:

http://matplotlib.sourceforge.net/examples/pylab_examples/polar_bar.html


texts:
I'm trying to build a web mapping app that retrieves data from our WFS server. I'd like to avoid JavaScript for the main search scripts.

From what's I've read the best way to send the XML POST request to the server is to use CURL.

I have tested the following very simple XML in GeoServer and it's fine. Hwever, when I send the request as below I get an error.

I have trawled Google but can't seem to find any examples of doing this. I can't believe I'm the only one trying to do this!

The PHP Code:

<?php
//Set variables.
$url = "http://squirrel.bis.local:8080/geoserver/wfs";
$post = '<GetCapabilities
    service="WFS"
    xmlns="http://www.opengis.net/wfs"
    xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
    xsi:schemaLocation="http://www.opengis.net/wfs http://schemas.opengis.net/wfs/1.1.0    /wfs.xsd"/>';

$post = str_replace("\r", "", $post); //Clean string.
$post = str_replace("\n", "", $post);
$post = str_replace("\t", "", $post);
$post = str_replace("> <", "><", $post); //Clean any leading spaces between elements.

//What are we sending?
echo "Post Request:
";
echo htmlentities($post);
echo "

";

//Send post request.
$ch = curl_init();
curl_setopt($ch, CURLOPT_HEADER, 0);
curl_setopt($ch, CURLOPT_RETURNTRANSFER,1);
curl_setopt($ch, CURLOPT_URL, $url);
curl_setopt($ch, CURLOPT_POST, 1);
curl_setopt($ch, CURLOPT_POSTFIELDS, $post);
$content=curl_exec($ch);

//Collect result.
$result = curl_exec($ch);
curl_close ($ch);

//What did we get?
echo "Post Response:
";
echo "
" . htmlentities($result) . "
";
echo "

";

?>

The error is:

<ows:ExceptionReport version="1.0.0"
xsi:schemaLocation="http://www.opengis.net/ows
http://squirrel.bis.local:8080/geoserver/schemas/ows/1.0.0                   /owsExceptionReport.xsd"
xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"     xmlns:ows="http://www.opengis.net/ows">
<ows:Exception exceptionCode="MissingParameterValue" locator="request">
<ows:ExceptionText>Could not determine geoserver request from http request         org.geoserver.platform.AdvancedDispatchFilter$AdvancedDispatchHttpRequest@19f0b0a6
</ows:ExceptionText>
</ows:Exception>
</ows:ExceptionReport>

Thanks
Steve


texts: Changing legend labels in ArcGIS Desktop?

texts:
I don't know how to change legend labels using attribute table. I want to use values from "litologia" field instead of "kod".

Is this even possible or do i have to switch value field in symbology properties and create all symbols again?


texts:
I need to copy many elements (and their attributes) from a shapefile and paste to another one. If I use copy/paste functions it pastes only the polygons but not the attributes.

How can I solve it?


texts: Integrate Python script with QGIS plugin

texts:
I am new to QGIS and trying different things. A Python script (link), which is running perfectly in QGIS Python console but I want it to run as a QGIS plugin. I am trying to do but getting confused about where to write my script logic in Plugin file. I copied the code below for reference which I want into Python plugin.

class TestEditTool(QgsMapToolEdit):

    def __init__(self, canvas):
        self.canvas = canvas
        QgsMapToolEdit.__init__(self, self.canvas)
        self.rb = QgsRubberBand(self.canvas, QgsWkbTypes.LineGeometry)
        self.rb.setStrokeColor(QColor('Red'))
        self.rb.setWidth(2.0)
        self.fixed_points = []

    def canvasReleaseEvent(self, event):
        click_point = event.snapPoint()
        if event.button() == Qt.LeftButton:
            if len(self.fixed_points) < 1:
                if not event.isSnapped():
                    iface.messageBar().pushMessage('Edit Error', 'Start point of line feature \
                    must be snapped to point layer', level=Qgis.Warning, duration=3)
                else:
                    self.fixed_points.append(QgsPoint(click_point))
            else:
                self.fixed_points.append(QgsPoint(click_point))
        elif event.button() == Qt.RightButton:
            if self.fixed_points:
                if not event.isSnapped():
                    iface.messageBar().pushMessage('Edit Error', 'End point of line feature \
                    must be snapped to point layer', level=Qgis.Warning, duration=3)
                else:
                    self.fixed_points.append(QgsPoint(click_point))
                    new_line = QgsLineString(self.fixed_points)
                    geom = QgsGeometry().fromPolyline(new_line)
                    ### Add Features here
                    layer = iface.activeLayer()
                    if layer.isEditable():
                        feat = QgsFeature()
                        feat.setGeometry(geom)
                        feat.setFields(layer.fields())
                        iface.openFeatureForm(layer, feat, False)
                        layer.addFeature(feat)
                        layer.triggerRepaint()
                    ###
                    self.clean_up()

    def canvasMoveEvent(self, event):
        if self.fixed_points:
            self.rb.reset()
            if len(self.fixed_points) == 1:
                pt1 = self.fixed_points[0]
                pt2 = QgsPoint(event.snapPoint())
                rb_line = QgsLineString(pt1, pt2)
                self.rb.setToGeometry(QgsGeometry().fromPolyline(rb_line), QgsProject().instance().crs())
            else:
                rb_line = QgsLineString(self.fixed_points)
                rb_line.addVertex(QgsPoint(event.snapPoint()))
                self.rb.setToGeometry(QgsGeometry().fromPolyline(rb_line), QgsProject().instance().crs())
            self.rb.show()

    def clean_up(self):
        self.rb.reset()
        self.fixed_points.clear()


texts: Applying a mask when running MESMA using VIPER tools plugin for ENVI

texts:
I am running MESMA on a landsat image that has large areas masked out, but I keep having all the background (value = 0 in all bands) pixels included in the analysis. They all end up unmodeled (no surprise there) so it's not a problem in regard to outputs, but computationally it's slowing the process down a LOT by trying to model all those background pixels. Is there any way I can somehow instruct VIPER tools to ignore the background pixels when applying MESMA? I have looked for an option to include a mask band but I can't find one.

I am running VIPER tools 2.0 plugin on ENVI 5.4.


texts:
I am looking to filter LiDAR data (.las format) by scan angle specifications (e.g. keep all points < 15). From a basic Google search, it appears that the libLAS library has the ability to filter by scan angle using las2las.

What other tools exist to filter .las files by scan angle, either open-source or commercial?  I am particularly interested in tools that are suited for bulk processing.


texts: ogr2ogr: export intersected gml using shapefile

texts:
I'm trying to do the same as the ogr2ogr example given by https://gis.stackexchange.com/a/119404/6458

The difference is that my inputs are a GML and a shapefile. Export.shp should only contain the data in input.gml intersecting with buffer.shp. I don't receive an error, but I also don't get an output:

ogr2ogr -f "ESRI Shapefile" -sql "CREATE VIRTUAL TABLE Buffer USING VirtualOGR('buffer.shp',
0, 'buffer');SELECT ST_Intersection(A.geometry, B.geometry) AS geometry, A.*, B.* FROM
Layername A, Buffer B WHERE ST_Intersects(A.geometry, B.geometry)" -dialect SQLITE export.shp
input.gml

'Layername' here in the SQL refers to the layer name in the input.gml

Does anyone have an idea of what is going wrong?


texts: An image on top of the marker with Leaflet

texts:
How can I achieve the same result (image below) with Leaflet?

L.imageOverlay() only can supply one image, and the image size is based on (lat, lon).

But what I need is the image size the same as the marker even when I zoom in/out.

Maybe someone can guide me through?


texts: Selecting minimum number of neighbouring polygons required to give a fixed total area

texts:
I have a polygon layer. I would like to select the minimum number of neighbouring polygons to give a total area greater than a fixed value. In reality the polygons are land parcels and I would like to know which parcels I need purchase in addition to create a continuous area of, for example 100 hectares.

Is this possible in QGIS using virtual-layers or by querying the spatialite database in the DB Manager? I need to repeat this process for potentially hundreds of sites, hence the need for an automated process. Is there some other way that could be used to achieve this within QGIS such as the processing modeler?

I had thought of creating a buffer of the appropriate area and then seeing which polygons intersect. However, the resulting area could be too small (not acceptable. Which would mean lots of trial and error for each site until I got it right. An area larger than my target value would be acceptable, however, the excess should minimised.


texts: Synchronize data edits between ArcMap and ArcGIS Online

texts:
What, if any, workflows are available to round trip data between ArcGIS Desktop and an ArcGIS Online for Organizations (AGOL) web map with editable feature layers?

For example, we have a web map with an editable feature service layer. Some staff might add features to the map using ArcGIS.com Map Viewer, others with ArcGIS Online Explorer, and still others with ArcGIS Desktop. Nothing needs to be done for the first two, as soon as the map is saved the other application will see it. However for ArcGIS Desktop, the so-called powerhouse in this equation, there doesn't seem to be a way to synchronize or upload changes back to ArcGIS Online. Is this really the case?

There aren't any SDE or databases in the mix, just ArcGIS Desktop with file geodatabase and AGOL.

Update: ArcGIS Server is not involved at all (disregarding that's what Esri uses power AGOL behind the scenes, as there is no meaningful interaction with it from a user perspective).


texts: How large is the QGIS community?  I'm interested in understanding the shift from traditional proprietary software to open sourced solutions.


texts: Why would Select By Location appear not to work properly after Multipart to Singlepart?

texts:
I am having issues using ArcGIS Desktop 9.3 and the multipart to singlepart feature and then selecting by location.

Basically, I am able to run the multipart to single part feature to split up some of my polygons. I then use Select By Location, setting the parameters to within a polygon in another layer.

The problem being that it does not pick up all the features within that polygon, missing most of the ones that were originally multiparts out.

Has anyone come across anything like this before?


texts: SQL queries in shapefiles using python/gdal return none result

texts:
I want try use Python and GDAL to create some SQL queries in my shapefiles.
I try something and it looks good, but doesn't work. I get a none result.

Here is the code :

from osgeo import ogr
ogr_ds = ogr.Open('line.shp')
TEST=3
sql = "SELECT id FROM {} where id='{}'".format(ogr_ds,TEST)
layer = ogr_ds.ExecuteSQL(sql)
print layer


texts: Utility GIS - Type of data cleanup needed & Effort

texts:
I am trying to get into Utility GIS and I have been asked by my boss to size up the effort of what it takes to clean up Pipeline, Utility asset GIS data and figure out if there is a way think of a factor like a per linear mile man days of effort.

The data is in the form of ESRI File
Geodatabase.

Typical clean up related tasks I have are :

Data consistency by using domains (Manual Effort) - Per Record Multiplier
Topology checks (ArcToolbox) Per Featureclass Multiplier
CAD to GIS conversion (Per Feature Class Multiplier)
Visual QCs

Can you suggest any other tasks typically for cleaning pipeline data and any assessment of what you would think of a man-day multiplier?


texts: Somebody knows a test to my application using WMS 1.1 (as a client) with OGC Standard?

texts: I need a test to certificate my application with the OGC Standard. Please if somebody knows a site or something in which i can do this please contact me or answer this topic.


texts: QGIS crash when using args=0 with qgsfunction

texts:
I am trying to define a custom function in QGIS 2.8.1 that takes no arguments.

from qgis.core import *
from qgis.gui import *

@qgsfunction(args=0, group='Custom')
def test(value1, feature, parent):
    return 'test'

This results in a QGIS crash. What am I doing wrong?


texts: Aggregating in FME and subsequent proportional split problem in MapInfo

texts:
I have recently aggregated 6 spatially identical MapInfo TABs in FME.  The data has Lotplan and PPN attributes that are then used to group by in the aggregation:

The polygons and figures in the resultant output appear completely fine with all the Floor Area (float)attribute data  being correct. Then these polygons are taken into MapInfo and split by another polygon layer with Floor Area figures split proportionally based on the new split area. The problem is that the resulting Floor Area data is divided by 6 (the number of input datasets) consistently across the dataset.

I have worked around this problem by multiplying the figures by 6.


texts:
I performed an object-based image analysis using ENVI EX, and each of the classified features were exported as a shapefile. I transferred those shapefiles into ArcMap 10, and tried to intersect them with other shapefiles of the same area that were acquired from GeoBase.

When I try to run the intersect tool, I get an error message stating that the topology of the shapefiles extracted from ENVI EX are invalid. I tried to turn on the topology tool to repair the topology of those shapefiles, but the buttons are greyed out, and I can't seem to click on them. The topology tool does not work for any of the shapefiles that I have.

Any solutions? How I can I fix the topology and run the intersect?


texts:
How do you increase Cache size? think QGIS is falling-over as haven't given it enough...

image to export 1080 x 800 mm!

Would like at least 120 dpi or 200 : )


texts:
I have two rasters, one with NDVI values, and one with greyscale values (=(red band + green band + blue band) / 3).
The rasters have the same name.

I need to reclassify them in arcpy with the following formular:

if ndvivaluefrompixel >=-0.000005*greyscalepixelvalue+0.314367:
     return 2

if ndvivaluefrompixel <=-0.000005*greyscalepixelvalue+0.314367:
     return 1

I thought about using map algebra in arcpy but I'm not sure if it is possible to use such a function for reclassification?

Is this possible using map algebra/reclassify in arcpy? Of course this has to run in a loop then, so I have to use python.

Any idea?


texts: Speeding up processing of geotiff background for QGIS?

texts:
Originally I had an image as a mrsid file.  That was very convenient as it was smaller and fast to pan etc.  The down side was that every time I upgraded QGIS I had to go through the ritual of reinstalling the mrsid stuff (I am on a Mac if that matters).

So I recently made a GeoTIFF of the image -- this works fine but doing anything that involves redrawing the image is slow.  I know that I can use gdal to produce a pyramid tiles and that this is the way web apps deal with the issue but I am unsure if that is the best approach for QGIS.

If it is the best approach then pointers to a tutorial would help greatly -- there are a lot of parameters to play with!


texts: Copy features from several layers into one layer

texts:
I am trying to merge many layers to one and the only working method I found is by copying features from one layer and paste in the desired layer. I tried to use MMQGIS - Merge tool but I get some geometry error (for e.g. polygon 2.5 ! = Unknown WKB 3003).

I tried with another plugin called multiplelayerselection but something strange happened: you can select all features from visible layers but it will only copy current layer.

What can I do to speed up the process? I have about 10000 layers, so doing one by one is very slow.


texts: Grid with data table attached to each polygon

texts:
I'm trying to wrap my head around a request a colleague gave for some data visualization.

She is performing field observations and for every observation, she is assigning a grid # that corresponds to a pre-existing grid I've created with Alpha/Numeric (row/column) ID's for each grid polygon.

Her question was if it is possible to take the data table she is creating, and attach it to the grid work, so that when she clicks on for instance grid polygon G5, it will pull up the table and show all observations that have the grid ID attribute of G5.

using ArcGIS for Desktop 10.2.2


texts: Extract by location not working in qgis modeler

texts: I am trying to extract lines to a new layer that don't fall 100% within a polygon.  I want to leave out anything that does. I can not get extract by location to work in the modeler, it does work when run by itself.  Is there another option or some information on how the plugin works.


texts: What is the most suitable method for filling voids (flaws) in raster DEM?

texts:
I am trying to delineate sub-basins in my area of study, and I have been faced with the following dilemma: currently what would be the best method for filling voids (flaws) in DEM raster using free software? Would it be using QGIS, Whitebox, GRASS or SAGA GIS?
My PC has 4GB of RAM and 2.3GHZ of processing.
I have an ALOS PALSAR DEM, but it is too heavy for my PC. I also have DEM SRTM 30 meters.


texts: Extracting "World Reference Base (2006) Soil Groups" data in R and searching for the legend of the extracted values

texts:
I downloaded the "World Reference Base (2006) Soil Groups" data per the SoilGrid website as tile1 and processed the data in R (xmin = -116, xmax = -114, ymin = 30, ymax = 32):
### Soil properties of the localities
soil = raster("C:/Users/Besitzer/Documents/MEP/data/soil/tile1.tif")
example_Test@data$soil = raster::extract(soil, example_Test)

I obtained the following values (> example_Test):
[1] NA NA NA NA NA NA NA NA NA 24 NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
[37] NA NA NA NA 24 NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA  4 NA NA NA NA NA 24  4 NA 24 24  4 NA NA NA NA

I expected values like: HS or AT or TC. Where do I find a key on how to interpret the extracted numbers? --> 4, 24, see above


texts: Extract polygons of PNG with GdalOgrInJava

texts:
I have PNG pictures with polygons of different colours on them. What I plan to do is to extract these Polygons using the polygonize function (see here)
Then I will try to encode the polygons in GeoJSON (maybe using this? Or is there also a java tool for that?)

I am quite new to GIS so I don't quite understand how to use the polygonize function mentioned above. I don't know which of the three polygonize functions I have to use.

I have very simple PNGs with some polygons on it. Would someone please show me a code snippet of how to use polygonize function with regard to my problem?


texts: Calculate Central Meridian Offset for Custom Projection of CAD File

texts:
I need to align an unprojected CAD file in ArcMap. I know how to adjust the False Easting and Northing, but this CAD file also has a rotation issue. I already adjusted the Central Meridian from -76.0 to -80.0, which helped, but it's still off. My adjustments of the Central Meridian parameter have all been trial and error.

Is there a way to measure the angular adjustment needed for the Central Meridian? See picture below:


texts:
I am trying to create OpenLayers vector features based on WKT strings that I retrieve from a MapGuide server.  I have successfully gotten the proper geometry to be rendered, but I do not know how to handle the SRS WKT string.  As a result, the vectors are being placed in the north end of the DR of Congo rather than eastern Pennsylvania.

I was able to retrieve the SRS WKT string:

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]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["False_Easting",1968500.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-77.75],PARAMETER["Standard_Parallel_1",39.93333333333333],PARAMETER["Standard_Parallel_2",40.96666666666667],PARAMETER["Latitude_Of_Origin",39.33333333333334],UNIT["Foot_US",0.3048006096012192]]

What would be the best way to handle this?  The base map is in EPSG: 900913.


texts: Finding roads in one particular province of The Netherlands using QGIS?

texts:
I have two shapefiles that I want to use to see all the roads in one particular province of The Netherlands.

So I have one layer with all the roads in The Netherlands and one layer with the boundaries of every province.

How can I 'select' or 'show' the roads in this one particular province? The attributes of the roads file doesn't contain data saying in wich province it is.

I am a beginner in using QGIS 2.0.1.


texts: Created isochrones are not displayed (QGIS QNEAT3 Iso-Area as polygons algorithm)

texts: I am trying to create isochrones using Iso-Area as polygons algorithm in QNEAT3 plugin. I choose the fastest past as the cost strategy. when I use the default speed to create isochrones the output would be displayed. But when I feed the algorithm with the speed field of the attribute table, the created isochrones are not shown and the attribute table of the output is also empty! what could be the cause of this?


texts: Clipping LAS data with multiple shapefiles

texts:
I am looking for a method to clip with multiple shapefiles (139 shapefiles polygon) a lidar.las file.

I'm using the software Fusion, as described in: How to clip LAS data using shapefile polygons and open source software?

However, because of the overlapping polygons, the process tend to misses around 18 tree plots. So, I split the polygons to become and individual polygon. There is about 139 separated shapefiles representing one polygon each.

The question is:

How am I going to write it in the batch file command to clip this 139 shapefiles data with one LAS file in one execution? I understand the step is to quote the location of the software, shapefiles and lasfile. But what if I have hundreds of shapefiles to clip with one las file.

I cannot dissolve the polygon because each of the polygons represent crown tree (overlapping crown in the dense tropical forest). I'm going to do analysis for individual tree modelling so that's why an alternative for individual polygon is necessary. Using lasclip from lastools only allows to clip one shapefile with one las file.

The pattern in the shapefile's names is: tree_0.shp, tree_1.shp,.....until tree_138.shp

Can you help me on this? In terms of how to write the command?


texts: Using temporary layer as input for other algorithm in Processing script

texts:
My script has a single point input layer that will be processed with some logic to extract just some points.

I'd like to pass these extracted points as input of another algorithm (simple buffer) within the same scripts.

I define the final output as QgsProcessingParameterFeatureSink:

self.addParameter(QgsProcessingParameterFeatureSink(
    self.OUTPUT,
    "Extracted points",
    QgsProcessing.TypeVectorPoint
    )
)

....

   (sink, dest_id) = self.parameterAsSink(
        parameters,
        self.OUTPUT,
        context,
        source.fields(),
        source.wkbType(),
        source.sourceCrs()
    )

running some logic to extract some points:

features = source.getFeatures()
for current, feat in enumerate(features):
    if feat.geometry().asPoint() == maximum:
        sink.addFeature(feat, QgsFeatureSink.FastInsert)

then I'm trying to call the buffer algorithm:

   buffer = processing.run(
        "native:buffer", {
            'INPUT': sink,
            'DISTANCE': distance,
            'SEGMENTS': 5,
            'END_CAP_STYLE': 0,
            'JOIN_STYLE': 0,
            'MITER_LIMIT': 2,
            'DISSOLVE': False,
            'OUTPUT': 'memory:'
        }
    )

but it does not accept the input:

Incorrect parameter value for INPUT


texts: Apply Join by Location to a Selection in QGIS Processing Modeller

texts:
In QGIS 3.10.2 I have created a process modeler file that takes inputs from 5 input layers and creates a new layer as an output.  The output layer is needed by an external program and the processing modeler is to ensure it can be updated frequently.

At the end of the file I want to "select by expression" based on two attributes, and then apply the "Join by Location" function to the selection and another point layer.  I tested it manually and it worked.  I copied the same expression and used this in the Processing Modeler.

However when I run the processing model, I get different result.  First I cant figure out how to restrict the Join to only the selected features - there is a tick box when doing this command manually. Secondly the Processing Model Join by Location doesnt seem to work correctly, giving different results than doing it manually.

My process model follows a linear pattern, where I have selected most of the layers.

The last two steps I want to "Select by Expression" firstly:

And then "Join by Location" on the Selection above

The results firstly dont return the "Selected" features, but also the Join doesnt seem to work as well as using the function manually.

I am not sure if I am doing something fundamentally wrong on this, but have tried many combinations and cant seem to get it function as if I were doing it manually.


texts: Importing World Map image changes color in ArcMap?

texts:
I am having a hard time finding a world relief map with "lighter" ocean colors that will actually display as such in Arc.

This is the image I download from the web

but this is how it displays in ArcGIS.


texts:
I have been tasked with finding out where an airplane was when it took the photo used in some aerial imagery that we have. I have been supplied with 'RAW' data and need to know how to interpret it towards determining the photo location for each image.

Here is what I have in file folders:

%image_name%_WM2.ads
%image_name%_WM2.hist
%image_name%_WM2.odf
%image_name%_WM2.odf.adj
%image_name%_WM2.prj
%image_name%_WM2.sup
%image_name%_WM2_x_x.tif

where x_x are a number of TIF files under the same base image name (e.g. 0_0 or 0_4 or 0_99).

There is a lot of data inside the SUP file that I feel is useful, but I don't know how to interpret it into my needs, if indeed it can be.

There is a lot of other information just outside this file structure within a file folder called CAM which contains a lot of files apparently to do with the sensors. I'm not sure if this data is useful in my current problem, but I have it if it helps.

With supposedly 'RAW' imagery data files, how can I extract out where each image was taken in space (X,Y,Z)?

Is there another way of determining photo location? What files should I request from the provider?


texts: How to generate landmark-based route description using ArcPy?

texts:
I am now using ArcGIS 10.1 and currently a beginner in ArcPy.

There is a project requiring to generate landmark-based route description. I am wondering how to use a buffer to find the landmarks in the neighborhood of the pre-defined route, in additional to provide a relationship between landmark and the route

are there any alternatives way to find out this "landmark-based route description" instead of using ArcGIS existing option?

The final product of the project requires a description of route in every turn/edge, such as:

"Turn right in McDonald to the Manhattan Road"

"Walk along the Manhattan Road, and there is a Park on your right"

It is not necessary to provide a distance in the description, as the target of the description is pedestrian and he/she doesn't really need this piece of information.


texts: Getting the compass point given longitude, latitude and direction in degrees

texts:
I have this data gathered using a GPS device

No.     Time       Longitude    Latitude Speed(Km/h)   Compass point      Direction (deg)
1   10/20/15 1:34   40.11724    -3.21032     11        Northwest by north  341
2   10/20/15 1:34   40.11714    -3.20999      6        Northwest by north  334
3   10/20/15 1:34   40.11713    -3.2098      10        North               352

My question is, how can I know the compass point given only the longitude, latitude and direction in degrees?


texts: Comparing field values to update new field using ArcGIS for Desktop?

texts:
I am currently doing an analysis in ArcMap 10 using python.

I have two shape layers, one is a point layer and the other one is a line layer. The point layer have attributes for point's unique ID, as well as X and Y coordinates (expressed as coord_X, coord_Y). The line layer include line's ID, and its has coordinates for starting point and ending point. Here they are expressed as FromID_X, FromID_Y, ToID_X, ToID_Y. There are also two empty column in line layer's attribute table as FromID and ToID for points.

All coordinates in both point and line layer are already know by field calculation, but the FromID and ToID of starting point and ending point in line layer is empty. I want to figure out the ID of point from point layer and fill them in the table of line layer, i.e. FromID and ToID.

I have an idea to compare the coordinates of FromID_X (line's layer) with Coord_X(point's layer), if they are the same, output the Point's ID to line layer as FromID and repeat this for ToID.

But now, I am not familiar with python in ArcMap 10, I am wondering if you can tell me some idea or grammar which I can look into so I could use python in ArcMap to do this work? Or you know some tools for such kind of work?


texts:
I am working in ArcMap 10.2 (and Modelbuilder) and I would like to solve the following problem:

I have Buffer A (single part) and Buffer B (multipart). I would like to extend Buffer A with Buffer B where there is an overlap between Buffer A and Buffer B or where both lines are within halve the distance of Buffer B. Though the extention used should be halve of the buffer distance used to create Buffer B (there are many different Buffer distance values used). So assume Buffer B1 has a buffer distance of 5m I would like Buffer A to be extended with 2.5m where there is an overlap between Buffer A and Buffer B1. If Buffer B2 has a buffer distance of 8m Buffer A should be extended with 4m where Buffer A and Buffer B2 overlaps. When Buffer B1 or any other part of Buffer B lies within halve the distance of the buffer distance used from Buffer A it should take that part with Buffer A as well.

The result (when dissolving the boundaries of the Buffer) would look something like the following:


texts: Inserting newlines in Rectangle Text elements via ArcPy causes overlap?

texts:
I ran across an issue the other day when I tried to use ArcPy's mapping module to edit rectangle text elements with newlines (\n) in an ArcMap document.  Here's what the output looked like:

Here's the code I used to generate that output.  The first column are rectangle text elements Text1, Text2, Text3 going down; the second column are "plain" text elements Text4, Text5, and Text6 going down.

import os
import arcpy

HomeDir = r"C:\Desktop"
arcpy.env.workspace = HomeDir

CurrentMXD = arcpy.mapping.MapDocument(r"C:\Desktop\TextTest.mxd")
OutputFilename = r"C:\Desktop\TextTest.pdf"
if os.path.exists(OutputFilename):
    os.remove(OutputFilename)

for TextElement in arcpy.mapping.ListLayoutElements(CurrentMXD, "TEXT_ELEMENT"):
    TextElementName = TextElement.name

    String1 = "The quick brown fox jumped over the lazy dog.\nShe sells sea shells by the sea shore."
    String2 = "The quick brown fox \njumped over the lazy dog.\nShe sells sea shells by the sea shore."
    String3 = "The quick brown fox jumped \nover the lazy dog.\nShe sells sea shells by the sea shore."

    if TextElementName == "Text1":
        TextElement.text = String1
    if TextElementName == "Text2":
        TextElement.text = String2
    if TextElementName == "Text3":
        TextElement.text = String3
    if TextElementName == "Text4":
        TextElement.text = String1
    if TextElementName == "Text5":
        TextElement.text = String2
    if TextElementName == "Text6":
        TextElement.text = String3

arcpy.mapping.ExportToPDF(CurrentMXD, OutputFilename)

So far, it looks like the presence of the messed up text depends on whether the line is longer enough to wrap, and whether the line before the newline is longer than the line after the newline.

Any ideas about what could be going wrong?  Is there a workaround?  I could use plain text elements and worry about wrapping lines using Python, but I'm hoping I can figure something out.


texts: How to add an item after an item in My Places in Google Earth Pro?

texts:
Here in My Places in Google Earth Pro we see we are about to add an item, e.g., folder to an item (Helicopter).
But how do we instead add an item after an item?


texts: Unknown coordinate system warning when programmatically adding layer

texts:
I have another interesting question for you.

I have data with unknown coordinate system.
When I add this data using a standard tool(File->Add Data) in ArcMap, a warning dialog appears.
It looks like:

But this dialog didn't appears if I add data programmatically using arcobjects.
Is there a way to show this dialog, if I add layers programmatically?

P.S. I found DatumCheckerClass, but that's not what I need.


texts: QGIS versions yielding different results on $area calculations?

texts:
When grading an exam, I discovered that QGIS 3.2 and 3.4 produce different values on computing $area over polygons.

I experimented using $area and area($geometry) to check the differences between planar and ellipsoid computations. It seems that QGIS 2.18 computes planar areas as default for $area, whereas QGIS 3.4 gives ellipsoidal areas for $area (as explained at the help section of the function).

In QGIS 3.2, $area produces different values (quite higher) than those for $area in QGIS 3.4, whereas area($geometry) produces the same values in the two versions.

As context, my data is projected on EPSG:32717 (WGS84 UTM z17 S).

Any insights on what's going on? Is this a bug? I'd like to explain the students these caveats and to be aware when computing area-derived variables.


texts: Zooming to ArcGIS JavaScript API graphic feature, but not to full extent?

texts:
I have a Python based geoprocessing service that gets called in my web app. The web app takes a few input parameters, and zooms to a returned buffered section of cable that meets the criteria. This all works really great.

I'm using the following code to zoom to the returned graphic feature:

map.setExtent(graphicsUtils.graphicsExtent(map.graphics.graphics));

But, I really don't want the map to zoom this close to the graphic:

I'd prefer something like this:


texts: Exporting each spatial bookmark to jpeg using arcpy gives IndexError?

texts:
I have some spatial bookmarks in my MXD file and i want to export each and every bookmark as one jpeg file. I have found below code but it is giving an error.

import arcpy

mxd = arcpy.mapping.MapDocument(r"C:\SOFTWARE\Python\test.mxd")

df = arcpy.mapping.ListDataFrames(mxd, "test")[0]

for bkmk in arcpy.mapping.ListBookmarks(mxd, data_frame=df):

    df.extent = bkmk.extent

    outFile = r"C:\SOFTWARE\Python\Out\\" + bkmk.name + ".jpg"

    arcpy.mapping.ExportToJPEG(mxd, outFile, df)

del mxd

After running above code, it is giving error as

df = arcpy.mapping.ListDataFrames(mxd, "test")[0]

  IndexError: list index out of range

How do I export each bookmark as jpeg with the above code?


texts: How to perform dynamic filtering of multiple check boxes using cql_filter and openlayers javascript?

texts: I have multiple checkboxes for example hospital checkbox where hospital datasets get displayed when checked, phase 1 check box where phase 1 hospitals get displayed when checked and grade A and B check boxes where when A is checked grade A hospital gets displayed. This is a kind of static filter.But I need dynamic filtering. For example, when phase 1 is checked and then grade A is checked I need phase 1 hospitals having grade A to get displayed. Likewise for other flters. Can anyone give suggestion on it using CQL_FILTERING and OpenLayers?


texts: QGIS Error- Layer ELRreal: OGR[3] error 1: REPACK of layer 'ELRreal' failed?

texts:
anytime I try to save my shapefille after editing it I get the next error messege-

  Layer ELRreal: OGR[3] error 1: REPACK of layer 'ELRreal' failed.

I don't know if it's relevant but the attribute table has data in hebrew.

Does anyone know what this error means and how to deal with it?


texts: How to stop ArcGIS Pro report from inserting cover page?

texts:
I have a simple report I've designed in ArcGIS Pro that outputs a handful of fields in a table for each selected feature in a point layer (see image below). Everything works as expected with the table itself and the formatting all looks correct.

However, Pro always inserts a cover page that contains only the Report Header no matter what I've tried. Turning off the Report Header has no effect on the layout; Pro just inserts a blank page instead and then starts the table on Page 2.

Is there a way to turn off this cover page by default in a report? It's obviously possible to tell Pro not to print this page in the PDF export options, but that's not ideal for reports that will be generated numerous times in the future.


texts:
My map works fine. Popup works fine. For some reason, the bottom of the popup window has some layout issues. I have all the resources included.

script:

dojo.require("esri.map");
dojo.require("esri.dijit.Popup");

And css:

<link rel="stylesheet" href="css/esri.css" type="text/css">
<link rel="stylesheet" href="css/Popup.css" type="text/css">

Do you have any suggestions on how I could debug this issue? Thanks.


texts:
I have an architecture plan with elevation measurements. When I display the labels (Bullet Leader style), the automatic scaling moves the label on top of the background layers, obstructing them.

I have tried defining the Feature Weight in Placement Properties as High in the architectural features, but that creates new problems with the zoom level at which labels are shown.

So my question is, is there a way to make labels behave in such a way so that their absolute size is constant (shrink and grow when I change the zoom) yet the placement properties are automatic?

I thought about using plain text labels and manually placing them in the correct position, but I have hundreds of elevation points, and I want to keep the plan neat (for publication) without different length.

Just to clear things, I have one points feature class of elevation points, and many polygon and line feature classes of the architecture.

Appreciate your help.


texts: I have seen a few DEM/hillshades in published journal articles that use a phong lighting scheme which gives surfaces a metallic luster. Does anyone know how/if this can be achieved in arcmap?


texts: Find lines that have specific distance between them in QGIS

texts: I have a rather large dataset that contains duplicated lines with an offset of 1m-5m and I want to find an efficient way to solve this problem using QGIS. Is there anyway in QGIS to find the lines that have a specific distance between them?


texts:
I am curious about making a thermal map to display where the hot spots and where the cool spots in a community are. Obviously the roads will be hot and the shaded areas with trees will be cool but this is more to brush up on my skills.

How would I go about collecting this data for an area where none currently exists?

It will be some form of raster data but I am lost from there on out.


texts: Center a selected point on map. Openlayers + Geoext1 + ExtJs3.4

texts:
I've got an ExtJs combobox next to the map with several points of the map that are provided from a Postgis database through Php/XmlHttprequest.

At the same time I've got a map with all points containing their respective attributes.

I need to select a point from the list and then this point must be searched on the map an centered on it.

Can somebody tell me which is the best way to proceed with this case? Relating combobox valuefield with map feature id maybe? Retrieving coordinates directly from database and centering the map on this coordinates?

Thank you


texts:
I'm trying to make a web application that is able to :

1) Upload raster data to GeoServer so when I use this command on cmd it works

curl -u admin:geoserver -v -XPUT -H "Content-type: image/tiff" --data-binary @"C:\Path\odm_orthophoto.tif"  http://localhost:8080/geoserver/rest/workspaces/raster_data/coveragestores/int_dec/file.geotiff

but on AJAX it doesn't work, giving me error 401 unauthorized here's my code

$.ajax({
            beforeSend: function (xhr) {
                xhr.setRequestHeader ("Authorization", "Basic YWRtaW46Z2Vvc2VydmVy");
            },
             url: "http://localhost:8080/geoserver/rest/workspaces/raster_data/coveragestores/int_decc/file.geotiff",
             contentType: "image/tiff",
             method : 'PUT',
             xhrFields:{
                  withCredentials: true
             },
             data : '--data-binary @"C:/Path/odm_orthophoto.tif"' ,
             success : function(data){
                  console.log(data);
             },
             error : function(httpReq,status,exception){
                 alert(status+" "+exception);
             }
        });

I also tried to use

 beforeSend: function (xhr) {
      xhr.setRequestHeader('Authorization', 'Basic' + 'admin:geoserver'
  }

but still getting the same error. I enabled CORS in web.xml in WEB-INF GeoServer folder, without the jetty.jar file like said here https://gis.stackexchange.com/questions/210109/enabling-cors-in-geoserver-jetty because adding it give me error 503 when opening geoserver.

2) Read layer's name from the coverage store also with REST API it works on Windows with this command

curl -u admin:geoserver -XGET http://localhost:8080/geoserver/rest/workspaces/raster_data/coveragestores/int_dec/coverages.json

My JavaScript code is :

$.ajax({
    url: "http://localhost:8080/geoserver/rest/workspaces/raster_data/coveragestores/int_dec/coverages.json",
    type: 'GET',
    xhrFields: {
        withCredentials: false
    },
    dataType: 'jsonp',
    success : function(data){
        console.log(data);
        },
    error : function(httpReq,status,exception){
        alert(status+" "+exception);
    }
});

but I get this error

  unexpected token: ':'

When I add

xhrFields: {
        withCredentials: false
}

to AJAX parameters I get the same error as 1) 401 unauthorized , which should eliminate CORS problem according to https://stackoverflow.com/questions/38953965/uncaught-syntaxerror-and-cors-error-using-ajax-with-jsonp/38956911#38956911 using json as datatype returns the same error

I'm using GeoServer 2.15 and JRE 8u221 32bits

UPDATE :

using withCredentials: true , worked for getting json list but not for the PUT method i get know HTTP/1.1 500 Server Error
here's the output of my curl command i think i need to change url to http:/localhost:8080 but i don't know what to do with the /geoserver/rest/workspaces/raster_data/coveragestores/int_decc/file.geotiff


texts: Using decimals in QGIS raster calculator

texts:
I am using the QGIS raster calculator to change the values of pixels of a specified value to a different value.  For instance, to change all pixels with a value of -32768 to a value of -999, I can use the following formula:

("RasterFileName@1"=-32768)*(-999)+("RasterFileName@1"!=-32768)*("RasterFileName@1")

This seems to work for all integer values, including positive, negative and zero.  However, it only works for some decimal values.  For example, the raster calculator can correctly identify and change a pixel with a value of -0.5, but not a pixel with a value of 0.1 or -0.6.

The raster calculator can change the value of a pixel to any of these decimal values, it just can't recognize these decimal values in an existing raster.

Any ideas on what is happening or is there a better way to change the values of pixels in QGIS?

The raster calculator seems to be a quirky tool.  It also doesn't recognize -3.40282347e+38 (sometimes used as a nodata value).


texts: I have some points in polygons green spaces and trees I want to intersect the data correctly, I try with Intersect tools but some data are missing and I do not why? the point layer is trees and polygon green spaces


texts: Install and configure ArcGIS Portal without Web Adaptor in Linux?

texts:
I have successfully installed and configured ArcGIS Server, ArcGIS DataStore, ArcGIS GeoEvent Server and I can access them without any problem.Now I'm trying to install and configure ArcGIS Portal in my server from past few days. I have successfully installed it but when I try to open the url, it redirects to internal ip/machine name.

Below are my server details.

OS             : Ubuntu 16.04 server
Web server : Nginx
RAM           : 16GB
Processor    : Quad core
FQDN          : Yes and have SSL installed
Open ports  : 80, 443, 6080, 6443, 6180, 6143, 7443, 2443

And all internal ports are open, so there's no problem with firewall.

ArcGIS Server, Arcgis GeoEvent Server, ArcGIS GIS Services, ArcGIS Portal Licenses authorized.

I know ArcGIS Web Adaptor requires Apache Tomcat as web server but I can proxy URLs using Nginx and I'm more comfortable with Nginx. I tried to use Apache Tomcat and again no instructions on how to deploy a war file (arcgis.war) and when I referred to other sources, the procedure seems very complex.

My Questions are

1.) Is a Web Adaptor necessary for ArcGIS Portal,  if I want to use only in single machine?
2.) Is it possible to use ArcGIS Portal from outside without Web Adaptor (Not from localhost)?
Because I accessed ArcGIS Portal without any problem from localhost. I had to install numerous GUI supported files and a browser in server to do that.

Please guide me with proper instructions on how to install and configure an ArcGIS Portal in Linux without Web Adaptor. I want to access it from public internet, not from localhost.

P.S: I have tried multiple times and I followed all the instructions from ArcGIS documentation. I think this redirection problem is a bug becuase ArcGIS Server and GeoEvent Server are not redirecting to internal IP. I can fix my problem if anyone can point me to the source code file where this redirection happens.


texts: Using arcpy.da.UpdateCursor to populate field with dictionary values

texts:
I want to populate a field, named FnotExist, with values from a dictionary. I want to check if the key of the dictionary matches the row value in field FI_PATH and if it does populate the corresponding row in FnotExist. The keys in the dictionary should match row values in FI_PATH because keys were created from this field. The code below is giving me the error:

  TypeError: value #1 - unsupported type: list

when it reaches the line cursor.updateRow(row)

import arcpy

arcpy.env.workspace = r"C:\TestFolder"

dbf = "ProcessControlTable"

dbf1 = "APExistOnlyEdit.dbf"

myfield = "FI_PATH"

filepathlist2 =[row[0] for row in arcpy.da.SearchCursor(dbf1,myfield)]
#print(filepathlist2)

#Dictionary one holds field names that are within shapefile
d1 = {}

for thisFile in filepathlist2:
    ContainedFieldNames = [f.name.upper() for f in arcpy.ListFields(thisFile)]
    d1[thisFile] = ContainedFieldNames

print(d1)

#Dictionary two holds field names that are thought to be in shapefile
d2 = {}

with arcpy.da.SearchCursor(dbf1,[myfield,'AP_FIPS','AP_BLDGCOM','AP_TYPE','AP_STATUS','AP_ID','AP_HOUSENU','AP_HALFADD','AP_PREDIR','AP_PRETYPE','AP_STNAME','AP_SUFTYPE','AP_SUFDIR','AP_UNITTYP','AP_UNIT','AP_BUILDIN','AP_FULLADD','AP_CITY','AP_STATE','AP_ZIP']) as rows:
    for row in rows:
        d2[row[0]] = [a.upper() for a in row[1:] if a.strip() != '']

print(d2)

#arcpy.AddField_management(dbf1,'APNOTEXIST','STRING')
##Gets the difference detween dictionary two and dictionary 1
##https://stackoverflow.com/questions/34298613/python-find-difference-between-two-dictionaries-containing-lists
diff = {}

for key in d1:
    diff[key] = list(set(d2[key])-set(d1.get(key,[])))
print(diff)

with arcpy.da.UpdateCursor(dbf1, ['FI_PATH', 'FnotExist']) as cursor:
    for row in cursor:
        keyValue = row[0]
        if keyValue in diff:
            row[1] = diff[keyValue]
            cursor.updateRow(row)

Below is an example of the dictionary, named diff, that I'm using to update the field FnotExist.

{'C:\\TestFolder\\Folder1\\A\\address_pt.shp': ['TYPEB', ' BUILDINGA'], 'C:\\TestFolder\\Folder1\\B\\ADDRESS_VALID_POINTS.shp': [' PRETYPEA'], 'C:\\TestFolder\\Folder1\\E\\Public_Safety_Addresses.shp': [' TYPEA']}

Below is a pic of the dbf table that I'm trying to update. The field FI_PATH was used to create the keys in the dictionary so keys should match values in row.

I'm not sure why I'm getting an error and after research I can't figure out how to structure the code within the UpdateCursor.


texts: QGIS - Python Error on Startup (Ubuntu 18.04)

texts:
After installing the 3.4.x Madeira LTR release, upon opening QGIS I am encountering the following error:

Couldn't load SIP module.
Python support will be disabled.

Traceback (most recent call last):
  File "", line 1, in
ModuleNotFoundError: No module named 'qgis'

Python version:
3.6.8 (default, Jan 14 2019, 11:02:34)
[GCC 8.0.1 20180414 (experimental) [trunk revision 259383]]

QGIS version:
3.4.10-Madeira 'Madeira', 3d97aaae88

Python path:
['/usr/share/qgis/python', '/home/wporter/.local/share/QGIS/QGIS3/profiles/default/python', '/home/wporter/.local/share/QGIS/QGIS3/profiles/default/python/plugins', '/usr/share/qgis/python/plugins', '/usr/local/lib/python36.zip', '/usr/local/lib/python3.6', '/usr/local/lib/python3.6/lib-dynload', '/home/wporter/.local/lib/python3.6/site-packages', '/usr/local/lib/python3.6/site-packages']

I suspect the issue lies with my Python path (or possibly anaconda?), but I did recently upgrade from Xenial to Bionic so I suspect that may also be a factor.

How can I resolve or begin to troubleshoot this error ModuleNotFoundError: No module named 'qgis'?

sudo vim /etc/apt/sources.list

...
deb https://qgis.org/ubuntu-ltr bionic main
deb-src https://qgis.org/ubuntu-ltr bionic main

Possibly related:

QGIS 2.18.3 path and sip module cannot be load in ubuntu 16.04

Couldn't load SIP Module

Update 7-23-2019

Digging further, I've noticed the qgis directory and PyQt5 directory seem to be missing from /usr/share/qgis/python

ls -l /usr/share/qgis/python returns...

drwxr-xr-x 2 root root 4096 Jul 22 16:16 console
drwxr-xr-x 5 root root 4096 Jul 22 16:16 plugins
drwxr-xr-x 2 root root 4096 Jul 22 16:16 pyplugin_installer
drwxr-xr-x 2 root root 4096 Jul 22 16:16 qsci_apis


texts:
In addition to https://gis.stackexchange.com/questions/173746/pgsql2shp-command-join-error?noredirect=1#comment260108_173746, the next command line works fine in order to join a table to a shapefile but the result is a dbf table and not a shapefile as I want. I see that exists a -g geom parameter but I try and it doesn't work. How do I get the geometry?

pgsql2shp -f "C:\....\teste2.shp" -h localhost -p 5432 -u postgres -P postgres database_name "SELECT table2.example, table2.example2 FROM table1 LEFT OUTER JOIN table2 ON table1.dicofre = table2.reside"


texts:
I'm new to GIS and tile servers. I've installed TinyOWS pretty successfully except for what looks like to be some dependency issues. I have installed PostGIS and also PostgreSQL. PostGIS is supposed to come with raster and topology support, so I assume it is already installed but not enabled. Anyway I can't find any good online resources for this so I need some help.

Here is my error:

NOTICE:  Function postgis_gdal_version() not found.  Is raster support enabled     and rtpostgis.sql installed?
NOTICE:  Function postgis_topology_scripts_installed() not found. Is topology support enabled and topology.sql installed?
NOTICE:  Function postgis_raster_scripts_installed() not found. Is raster support enabled and rtpostgis.sql installed?
NOTICE:  Function postgis_raster_lib_version() not found. Is raster support enabled and rtpostgis.sql installed?

Update:
I've tried creating an extension to these like as follows:

CREATE EXTENSION postgis;

Then I get the following:

ERROR: PostGIS is already installed in schema 'public', uninstall first


texts: How to generate a single kml file for each feature in a kml file?

texts:
I have a kml file containing features that I would like to parse out into separate kml files for each feature/header.

I'm hoping there is an easy way to do this using geojson.io or QGIS.

How can I generate separate kml files for each feature/header in a single kml file?


texts:
I'm learning how to use Python on Qgis. I'm using the 'Therminal Python' and my Qgis is the 2.16.3, 64 bits.

I loaded layers from Postgis. I'm trying to generate a PDF. I use the code on the chapter 9 on PyQGIS developer cookbook. Basically 'copy and paste'.

The class QPrinter() is not recognized. Then, the rest of the code has problem.

Note that I included 'import"'s as I read in many solutions on this forum. What I must to do?

from qgis.core import *
from qgis.utils import *
printer = QPrinter()

Traceback (most recent call last):
  File "<input>", line 1, in <module>
NameError: name 'QPrinter' is not defined"


texts: Error on GEE: No features contain non-null values of "system:time_start"

texts:
I'm calculating DVI over a region using MODIS. I am able to get the printout of the map properly but not the chart to retrieve the data. To calculate DVI, I did:

var dvi = img.expression(
    '(nir - red)',
    {
        red: img.select('sur_refl_b01'),    // 620-670nm, RED
        nir: img.select('sur_refl_b02'),    // 841-876nm, NIR
        blue: img.select('sur_refl_b03')    // 459-479nm, BLUE
      .rename('DVI')
    });

When I look at the inspector, I see the value for the DVI band, but the code for charting seems not to work. It says there's no system time.

var bandChart = ui.Chart.image.series(dvi, roi)

I also tried to add .filterDate ('2016-05-01', '2016-10-31') but it says the parameter is not a function. I.e.

var dvi = img.expression(
    '(nir - red)',
    {
        red: img.select('sur_refl_b01'),    // 620-670nm, RED
        nir: img.select('sur_refl_b02'),    // 841-876nm, NIR
        blue: img.select('sur_refl_b03')    // 459-479nm, BLUE
      .rename('DVI')
    });

I also tried to create a function (so as to specify timeframe) but that did not work either. How can I proceeed to get the chart?

var addDVI = function(image) {
  var dvi = img.expression(
    '(nir - red)',
    {
        red: img.select('sur_refl_b01'),    // 620-670nm, RED
        nir: img.select('sur_refl_b02'),    // 841-876nm, NIR
        blue: img.select('sur_refl_b03')    // 459-479nm, BLUE
      .rename('DVI')
      .filterDate('2016-05-01', '2016-10-31')

    });
    return image.addBands (dvi);
};


texts: How to input correct coordinates depending on various CRS with Python in QGIS

texts:
I used a script from the PyQGISDeveloperCookbook and another StackExchange post to render an image from a specific coordinate. The problem I have is that I am unsure how to input the correct coordinates for a custom CRS. The CRS I am using in my project is ETRS89 / GK25FIN, And the coordinate I want to render is 25500381,6677153. When I attempt to render, I only output garbage similar to the result Leolo got in the previously mentioned post.

I have tried using a smaller number but the odds of me brute forcing the correct location is practically zero. Being a noob I am unsure what the given CRS actually means, and what coordinates it is supposed to use. I have tried using this site to translate the coordinates with no luck.

Here is the my slightly modified script from the previous post:

from PyQt4.QtGui import QImage, QPainter
from PyQt4.QtCore import QSize
import os

# radius for image in meters
radius = 512
# center point
p = [25499771,6677153] # <- I don't know what to write here in this CRS
# bounding box
bbox = QgsRectangle(p[0]-radius, p[1]-radius, p[0]+radius, p[1]+radius)
# image
img = QImage(QSize(200,200), QImage.Format_ARGB32_Premultiplied)
## create painter
painter = QPainter()
painter.begin(img)
render = QgsMapRenderer()
## set layer set
layers = QgsMapLayerRegistry.instance().mapLayers()
lst = layers.keys()
render.setLayerSet(lst)

# Set destination CRS to match the CRS of the first layer
render.setDestinationCrs(layers.values()[0].crs())
# Enable OTF reprojection
render.setProjectionsEnabled(True)

## set extent
render.setExtent(bbox)
## set output size
render.setOutputSize(img.size(), img.logicalDpiX())
## do the rendering
render.render(painter)
painter.end()
## save image
img.save("K:\\QgisOutput\\TestHeight.png", "png")


texts: Displaying grid (TIFF) data on Plotly Dash mapbox with Holoviews and Datashader

texts:
Full Colab Notebook:
https://colab.research.google.com/drive/1_Rr71QbjlEjBmTW_r0hOoBLFVoxaMfBh?usp=sharing
Summary:
I'd like to display a TIFF image in a Mapbox layer in Plotly Dash, with tooltip showing the actual values.
I'm using Holoviews and I easily can get it work in the Bokeh backend:

As you can see, the image is overlaid on the map, and a tooltip shows up where the mouse pointer is, displaying the actual value at that location.
However if I use the plotly backend or I try to feed this exact result into Dash (which I assumes forces the plotly backend), I get an error:
Elements of type <class 'holoviews.element.raster.Image'> cannot be overlaid with Tiles elements using the plotly backend

It sounds weird, but ok.
Then I tried to feed the image into Datashader and display the result (as the live data will be much more detailed than the test data I'm going to use DS either way). Now I don't get the hover tooltip and weirdly enough also the image is vertically mirrored:

I've been trying to make this work for a week now, and I'm going nuts. Does anybody have an idea how to display a TIFF on a Mapbox properly while showing hover tooltips for the data? I know howto questions are usually frowned upon, but it seems like such a basic thing to do and yet with 20 years of developer experience and by reading through every documentation I could find and checking out every example I can't manage to do it.
Edit: After some research I concluded that as of now, Plotly itself is not capable to display raster data on mapbox (eg image that is not just a visual reprezentation but carries the underlaying data values too) at all, so of course it doesn't work through Holovies either. Thankfully I'm still at early stages of the project, so I will switch completely to the Panel-Holoviz-Bokeh tools and leaving Plotly-Dash behind.


texts: How do I create several concentric buffers in a loop in R

texts:
I am working with polygons and I need to create concentric buffers every ten kilometers until I reach a distance of 70 kilometers, in this sense I want only the 10 kilometers rings without the polygon that generated the buffer.

I want something like this:

I wanted to do a loop, because I have to do it for several shapefiles or a lapply.

I have eight minas sahpes and for each of them I will have to make the seven buffers (from ten to 70 kilometers)

this code works but would like to do a loop:

library(raster)
library(rgeos)
library(rgdal)

buff_min10km<- gBuffer(minas, width=10000, byid = TRUE)
buff_min20km<- gBuffer(minas, width=20000, byid = TRUE)
buff_min30km<- gBuffer(minas, width=30000, byid = TRUE)
buff_min40km<- gBuffer(minas, width=40000, byid = TRUE)
buff_min40km<- gBuffer(minas, width=40000, byid = TRUE)
buff_min50km<- gBuffer(minas, width=50000, byid = TRUE)
buff_min60km<- gBuffer(minas, width=60000, byid = TRUE)
buff_min70km<- gBuffer(minas, width=70000, byid = TRUE)

After taht I take out the center:

ring10km<- gDifference(buff_min10km, minas)
ring20km <- gDifference(buff_min20km, buff_min10km, byid = TRUE)
ring30km <- gDifference(buff_min30km, buff_min20km, byid = TRUE)
ring40km <- gDifference(buff_min40km, buff_min30km, byid = TRUE)
ring50km <- gDifference(buff_min50km, buff_min40km, byid = TRUE)
ring60km <- gDifference(buff_min60km, buff_min50km, byid = TRUE)
ring70km <- gDifference(buff_min70km, buff_min60km, byid = TRUE)

I tired this but didnt work :

bufferring<-seq(10000, 70000, by=10000)

outerBuffer<-function(x, dist){
 buff<-gBuffer(x, width = dist - 1, byid = T)
 e<-erase(buff,x)
 return(e)
 }

for (l in 1:length(bufferring) ){
ring<-outerBuffer(minas,bufferring)
}

Give me this error:

Error in gBuffer(x, width = dist - 1, byid = T) :
length(width) == n is not TRUE


texts: How to use a regression equation to create a map in QGIS?

texts:
I have been working on creating an exposure map based on measures taken at fixed points. I want to extrapolate the results to a grid to get an exposure map (air pollution data) for a region based on map data (road length, bus stops, etc.)

I have done the statistical work and have a regression equation from the collected data, but I am having a problem with the next step.

How do I use the regression equation to create the exposures for unmeasured areas?

I created 100mx100m grid over the area but am wondering if there is a tool I can use to input the regression equation for each grid section and have QGIS create a map with the calculated values?

Or do I need to calculate the values somewhere else and use that to create a map?


texts:
This is driving me mad. See the Leaflet tutorial for adding a layer control to your map:

L.Icon.Default.imagePath = "Scripts/images";

var mapUrl = 'http://{s}.tile.openstreetmap.se/hydda/full/{z}/{x}/{y}.png';
var secondMap = 'http://{s}.tile.osm.org/{z}/{x}/{y}.png';

var littleton = L.marker([39.61, -105.02]).bindPopup('This is Littleton, CO.'),
denver = L.marker([39.74, -104.99]).bindPopup('This is Denver, CO.'),
aurora = L.marker([39.73, -104.8]).bindPopup('This is Aurora, CO.'),
golden = L.marker([39.77, -105.23]).bindPopup('This is Golden, CO.');

var cities = L.layerGroup([littleton, denver, aurora, golden]);

var grayscale = L.tileLayer(mapUrl, { id: mapUrl }),
streets = L.tileLayer(secondMap, { id: secondMap });

var map = L.map('map', {
    center: [39.73, -104.99],
    zoom: 10,
    layers: [grayscale, cities]
});

var baseMaps = {
    "Grayscale": grayscale,
    "Streets": streets
};

var overlayMaps = {
    "Cities": cities
};

L.control.layers(baseMaps, overlayMaps).addTo(map);

Now, how do I use the removeLayer function to remove the "GrayScale" map from the control using a map.on('click') event? The code would be:

map.on('click', function(){
    //remove GrayScale
});


texts: Using the GDAL library, what would lng_var be and how to complete a polygon?

texts:
Just need some help understanding some code.

What I am trying to do here is create a Shapefile and not sure if I am doing that correctly or not using the library gdal. It should be a shapefile that contains a grid of square polygons. Where each polygon cell is going to be 10 degrees by
10 degrees and the grid should cover the entire world. With the shapefile I am also trying to include a variable that contains a decimal
column called ‘LAND’ that will be populated with 0.0 for all rows/cells. Not sure I am doing this correctly, I had some filler code and now I am trying to figure out the parts that are needed to make this work. In the end it should look like this (18x36 grid):

form osgeo import gdal
from osgeo import ogr
from osgeo import osr
import os
import random

# Specify that you want to work with a Shapefile
DriverName = "ESRI Shapefile"
# Not sure what this line does
driver = ogr.GetDriverByName(DriverName)

# Name of shapefile to create
FileName = 'gride.shp'

# simple check to see if that file already exists in the system, if so delete it
if os.path.exists(FileName):
  driver.DeleteDataSource(FileName)

# Now create an Empty shapefile Data source
shapeData = driver.CreateDataSource(FileName)

fieldName = 'LAND'

lng_min = -180
lng_max = 180
step = 10

lat_min = -90
lat_max = 90

# lng_var = ??
# lat_var = ??

for lng_var in range(lng_min, lng_max, step):
    ring = org.Geometry(ogr.wkbLinearRing)
    ring.AddPoint(lng_var, lat_var)
    ring.AddPoint(lng_var + step, lat_var)
    ring.AddPoint(lng_var + step, lat_var + step)
    ring.AddPoint()
    ring.AddPoint()
    poly = org.Geometry(ogr.wkbPolygon)
    poly.AddGeometry(ring)

    for lat_var in range(lat_min, lat_max, step):
        ring = org.Geometry(ogr.wkbLinearRing)
        ring.AddPoint(lng_var, lat_var)
        ring.AddPoint(lng_var + step, lat_var)
        ring.AddPoint(lng_var + step, lat_var + step)
        ring.AddPoint()
        ring.AddPoint()
        poly = org.Geometry(ogr.wkbPolygon)
        poly.AddGeometry(ring)

I haven't found a way to integrate my fieldName = 'LAND' yet. Where that should contain the decimal column called ‘LAND’, populated with 0.0 for all rows/cells.

I want to after this is done be able to open it on QGIS.

Not sure how to complete the polygon that is where there is two AddPoint()'s and also not sure what lng_var will be, or the lat_var when I add the inner loop.


texts: Creating irregular distance linear referencing hatches using ArcGIS for Desktop?

texts:
I am using ArcGIS 10.3 for Desktop.

I have stationing information which is placed at irregular distances from the starting point. These have a correlating mileage point. However, they are sometimes 0.1 mile apart and other times 0.5 miles apart. Therefore, although I have created routes from my layer, I cannot determine a way to bring in the stationing data table to create a hatching based on that information. How can I make my hatching irregular following a csv data file?


texts: Seeking C/C++ equivalent of gdal2tiles.py?

texts:
Is anyone aware of a permissively licenced open source C/C++ equivalent of gdal2tiles.py?

Want to generate OSM compatible tile sets from geo-referenced raster files, and have already confirmed that gdal2tiles.py will produce the output I need.

Am targeting a yocto generated embedded linux platform with relatively limited resources, so:

Want a native library or command line executable, rather than something that involves dragging in the 400lb gorilla of a python interpreter/JVM/etc.
Need something for which the source is available so I can cross-compile.

Am already aware of:

maptiler  Commercial tiling tool which I believe was created by by the original authors of the gdal2tiles.
mapnik Whole lot of stuff I don't need
simple-tiles.  Maybe, (haven't actually tried it yet) but note that it has a couple of extra dependencies above and beyond GDAL which I don't think I fundamentally need
gdalsplit A patch that has been sitting in the GDAL bug tracker for 10 years or more (effectively superseded by gdal2tiles I think).  Looking at source it appears to generate a simple tile grid with no depth, rather than the image pyramid  expected in an OSM/TMS tile set, but a useful reference.

Not too daunted by the prospect of mashing up my own implementation based on gdal2tiles & gdalsplit, or perhaps trying a compiled derivative of gdal2tiles using cython or whatever, but just wanted to confirm that there wasn't an existing solution somewhere out there that my google-fu had failed to identify


texts: Cannot create new spatial database connection in ArcCatalog?

texts:
I have searched a lot in the net, but i just found Add Spatial Database Connection fails? here.

However, I cannot see the confirmed answer.

I want to know if anyone has found the solution?
my ArcGIS desktop version and ArcSDE version is 10.0 and Sqlserver R2 2008.
The problem is that when i click test in spatial database connection wizard in ArcCatalog tree can see the message  "process succeed" but in fact after clicking ok nothing is created.


texts: In Python, reading a GDAL raster from memory instead of a file

texts:
I'd like to analyze some raster data in Python using NumPy.

This data doesn't necessarily exist as a file on the filesystem, it could be the result of a query against PostGIS, or something transferred over the network.

Is it possible to get GDAL to read data from a buffer in memory, instead of requiring a filename?

I realize one solution could be to open up a named temp file, write the data to that, and pass it as a filename to GDAL, but that seems inefficient, and more likely to incur I/O bottlenecks. Additionally, in the case of PostGIS, it would be nice to read the native format directly, as opposed to converting it to some other file format and reading that.

Edit: to clear up a bit of confusion, just using the PostGIS driver directly won't work either, since the query could be for a raster analysis, rather than just data stored in the DB.


texts: Has anyone here successfully exported from ArcGIS 10.5 to Adobe Illustrator?

texts:
If so, how did you do it?

I used the Export function within ArcMap 10.5 and while ArcMap could "see" the files it created, neither Windows nor my Mac could locate the files within the directory where the exported AI files were supposedly created. I am referring to both operating systems because I am running the emulator software Parallels on a Mac, which in turn enables me to run Windows and hence, ArcGIS 10.5.

Experience has taught me that the ArcMap 10.3 was the last and only version that exported smoothly to Illustrator and for that reason I continue to use it rather than opt for upgrading to 10.5. If not for the export problem, I would update my software to ensure that the version I use on my personal desktop computer is concurrent with the version being used at school.


texts:
I have a number of vectors (points, lines and polygons) loaded into PostGIS (not QGIS Server).  I would like to provide a colour schema based on attribute values which is saved as a default.  What I am hoping to be able to do is have my layers symbolized (ie: blue for lakes, black for roads and a star for capital cities) so that when users from other computers connect to the PostGIS database the layers are loaded with the symbology.

Is there a way to save all the symbols within PostGIS similar to a style file?

Thanks....