texts: After Feature Layer reference is added all dojo references are not recognized

texts:
I'm creating a simple javascript management system in which I have to query 6 feature layer services, edit, add or delete data from the services. Every time I want to add the Feature Layer reference all other functions stop working. It is like dojo and arcgis.js are not being loaded anymore.
Here is an example of the code I wrote:

require(["dojo/dom", "dojo/on", "dijit/form/ComboBox",
    "dojo/parser", "esri/tasks/query",
    "esri/tasks/QueryTask", "esri/tasks/GeometryService",
 "esri/layers/FeatureLayer", "dojo/domReady!"],

function (dom, on, ComboBox, parser, Query, QueryTask, GeometryService, FeatureLayer) {
    parser.parse();
    var link = '';
    var linkEN = "";
    var linkALB = "";
    esri.config.defaults.geometryService = new GeometryService("https://utility.arcgisonline.com/ArcGIS/rest/services/Geometry/GeometryServer");
    dijit.byId("id_category").on("change", function () {
        var category = dijit.byId("id_category").value;
        alert(category);
});
        function categorySelected(link) {
        var queryTask = new QueryTask(link);
        console.log("link " + link);
        var query = new Query();
        query.where = '1=1';
        query.returnGeometry = false;
        query.outFields = ["*"];
        queryTask.execute(query, populateList);
        on(dom.byId("id_category"), "change", execute);
    };
    function execute() {
        queryTask.execute(query, populateList);
    }
              .
             here goes the code for populating the list and displaying data
               ..

            dijit.byId("btn_save").on("click", function () {
        alert("start editing ");

        // create feature layer
        var pointsOfInterest = new FeatureLayer("link", {
           mode: FeatureLayer.MODE_ONDEMAND,
            outFields: ['*']
        });
        console.log(pointsOfInterest);
    });

All this code is working fine until the Feature Layer reference is added. After that none of the functionalities work, not even the alert from selected category.
I'm using the latest version of arcgis javascript api - 3.16


texts: Using Bing Map via NMF file in 2017 in ArcPad 10.2

texts:
Does anybody has recently experiences with using Bing Map Service on an Handheld Device with ArcPad 10.2?

I'm getting the Error message:

  Error

  Missing AuthenticationResultCode in response from Bing server

The Bing Map key I have set as described here.


texts: OpenLayers 3 incorrectly transforms coordinates into EPSG:900913

texts:
I've encountered a problem trying to transform coordinates from lat-lon (WGS84) format to EPSG:900913.

Consider the following code:

var olCoordinates = ol.proj.transform([55.75, 37.61], "WGS84", "EPSG:900913");

This is a coordinates of Moscow I try to transform and it yeilds the result [6206061.611725002, 4524477.448837536].
But on the EPSG site i get the other values altogether: [7508562.162509, 4188137.653801].

In other GIS system which we use I get the same values as EPSG site gives me, although it uses OL2, so I assume that there's something wrong either with my setup or with OL3 itself.

Any suggestions?


texts: Solar Potential roof-top analysis using GIS and Remote Sensing

texts:
3D Building models were generated by using Google Sketchup Pro and Google Earth.

Then shadow casted 2D orthographic raster models were exported as TIFF or JPEG format.

For further analysis with respect to shadow casted roof top areas, exported 2D raster models have to be classified or segmented. This process wanted to be executed for extraction of shadow polygons which casted over building roof tops and roof-tops open to sun light as well.

To achieving this objective unsupervised classification was done. But unfortunately shadows which casted over rooftops and ground area were being merged into one single class. Rooftop area shadow could not be extracted separately.
Supervised classification could not be done, because of unable to get training samples. (RGB values of training samples were same for shadow areas and open areas)
Therefore I figured out Matlab Application for image feature extraction.
I used some sort of extraction techniques such as Color threshold, RGB value extraction
But unfortunately they all were with lack of satisfaction.
Even though using Matlab with different approach may be result reliable output for this purpose. But up to now, I could not find out appropriate codes.

I have to extract features separately as output of image classification.


texts:
I have read a LAS file. If I import it in Global Mapper it says the projection, however when I read it with laspy (header.proj_id_1) there is no information associated.
I have edited the coordinates and updated the X,Y,Z values to new values (according to this) for the LAS file and now I would like to give a new projection after updating the coordinates.
What function do I use?


texts: Intersect/union two different polygon layers to produce a third polygon type

texts:
Can anyone instruct me on the correct way to union/intersect two different layers containing polygons so that a third layer is created containing the attributes of the non-overlapping regions from each original layer and a new 'class' or feature for any regions that do overlap?  I am using ArcMap 10.2.2

For instance, see in the image below:

I want the areas that overlap between the purple and green polygons to be considered a new dataset.  I know the "Intersect" tool in ArcMap would produce this "overlap" dataset, but only after stripping away everything else right?  On the other hand "Union" produces what we see below.  Is there anyway to do exactly what I am talking about?  Or do I have to do an intersect, generating an isolated "overlap" layer, and then intersect this "overlap" layer with the other two (one at a time) to finally end up with three different layers separated into the information I want - after which I would merge them back together.


texts: I am working on ArcGIS operation Dashboard. There is an widget on the top left corner “ Model Manager”. Here the Model manager will read the data from ForecastStage_MPT and ForecastTimePT database tables. These 2 database tables are in MS SQL server database. The data in these 2 tables are populated from NOAA. A script has been developed which will take the data from the NOAA website and populate the ForecastStage_MPT and ForecastTimePT database tables. My task is to read the data from the database tables into the ArcGIS operation dashboard model manager widget automatically. Do you have any idea how to read the MS SQL server database tables into the model manager widget? I am sending you 3 attachments and see if you can give me some directions.


texts: Classifying data as Natural Breaks using ArcPy?

texts:
Is there any chance I can create a simple python script to classify my data as Natural Breaks?

I have only created a script for quantile classification, but it is not enough for my work.

arrp = arcpy.da.FeatureClassToNumPyArray(in_features, field)
arr = np.array(arrp,np.float)

p1 = np.percentile(arr, 20)
p2 = np.percentile(arr, 40)
p3 = np.percentile(arr, 60)
p4 = np.percentile(arr, 80)
p5 = np.percentile(arr, 100)

if Method[0] == "5":
    with arcpy.da.UpdateCursor(in_features , [field,'Class']) as cursor:
        for row in cursor:
            if row[0] < p1:
                row[1] = 1  #rank 0
            elif p1 <= row[0] and row[0] < p2:
                 row[1] = 2
            elif p2 <= row[0] and row[0] < p3:
                 row[1] = 3
            elif p3 <= row[0] and row[0] < p4:
                 row[1] = 4
            else:
                 row[1] = 5
            cursor.updateRow(row)
    print Method

I need to process values in the field "Class" that are classified from the "field".


texts: Performing a spatial relationship in Oracle SQL with only lat/long.

texts: I currently have a query that pulls a latitude field and a longitude field (separate fields).  This is sufficient for one portion of customer product, but for a second table I need to add a second component that relates the coordinates to a polygon.  Is there a query to perform a spatial query or join using SQL (Oracle) without having the vector component, but instead the lat/long values in separate fields?


texts: Difference between Graphic() or a FeatureLayer() in ArcGIS

texts:
I've been adding paths and markers on my maps with the Graphic() service, but I know you can do it all with multiple FeatureLayer().

What is the main difference between doing it all with one or the other?


texts: Getting KML of external boundary of several areas from GeoServer

texts:
I have a GeoServer with PostGIS database.

PostGIS database has a table with small geographic areas (English LSOAs and output areas, to be precise), with codes and shapes.

I want to submit to GeoServer a list of area codes, and receive KML of their external boundary. The areas will be conterminous in most cases (not always).
GeoServer should merge the areas and dissolve the inner boundaries where possible.

I figured out how to do it in PostGIS:

SELECT ST_askml(ST_Union(ST_SnapToGrid(geom,0.0001)) )
   FROM geo_area
   where identifier = ANY( regexp_split_to_array('E01031338,E01031376,E01031377', ',') )
   group by level_identifier;

How do I pass this request to PostGIS via GeoServer?

I guess it can be made by a chain of WPS jobs

  [my_geoserver_url]/geoserver/web/?wicket:bookmarkablePage=:org.geoserver.wps.web.WPSRequestBuilder

http://docs.geoserver.org/stable/en/user/extensions/wps/install.html

but I am not sure yet.


texts: ArcMap Data Driven Pages Display Error - Grids Display Incorrectly

texts:
I am running ArcMap 10.4.1 with an Advanced Licence on Windows 7.
I have created a mapbook in ArcMap.

The mapbook consists of a mask, data-driven page grid (generated from a Fishnet).
The mask is a buffer of the Fishnet, and is symbolised so that only the information in the Data Driven Page is visible. The Data Frame has also been set up to only show the information within the Data Driven Page extent.

The issue I am having relates to how the Geographic Grids (coordinates bounding the page) are shown. The first page of the mapbook displays correctly, however as progress through the extents, the Geographic Grid appears to "leak" off the page. See examples below:

The position of the "leak" varies depending on the Data Driven Page position. The effect is carried through to the exporting of the mapbook to PDF.
I have tried repairing the geometry of my masks and Data Driven Page extents, with no success. I will be trying to reboot my system and trying to duplicate the issue on another machine. In the meantime, is there another potential solution, or something I should avoid doing when creating this type of mapbook?

Update

If I select "Measured Grid" the grids display perfectly, if I choose "Graticule" I get the same issues. Not sure if that is relevant.


texts: Populate more than one Silverlight Image from Bing Maps async request

texts:
I have the following code, which works fine if you just want to populate one Image with a response from Bing Maps. But if I try to do two then the variable _currentImage always ends up being "image1" because the calls are asynchronous. How can I pass the image variable along to the ImageryServiceGetMapUriCompleted method?

using System;
using System.Windows.Controls;
using System.Windows.Media.Imaging;
using BasicBingMapsImagerySvc.ImageryService;

namespace BasicBingMapsImagerySvc
{
    public partial class MainPage : UserControl
    {
        private const string BingMapsKey = "my key";

        private Image _currentImage;

        public MainPage()
        {
            InitializeComponent();

            GetMap(42.573377, -101.032251, image0, MapStyle.AerialWithLabels);

            GetMap(42.573377, -101.032251, image1, MapStyle.Road_v1);
        }

        private void GetMap(double lat, double lon, Image image, MapStyle mapStyle)
        {
            var mapUriRequest = new MapUriRequest();

            // Set credentials using a valid Bing Maps key
            mapUriRequest.Credentials = new Credentials();
            mapUriRequest.Credentials.ApplicationId = BingMapsKey;

            // Set the location of the requested image
            mapUriRequest.Center = new Location();
            mapUriRequest.Center.Latitude = lat;
            mapUriRequest.Center.Longitude = lon;

            // Set the map style and zoom level
            var mapUriOptions = new MapUriOptions();
            mapUriOptions.Style = mapStyle;
            mapUriOptions.ZoomLevel = 13;

            // Set the size of the requested image to match the size of the image control
            mapUriOptions.ImageSize = new SizeOfint();
            mapUriOptions.ImageSize.Height = 256;
            mapUriOptions.ImageSize.Width = 256;

            mapUriRequest.Options = mapUriOptions;

            var imageryService = new ImageryServiceClient("BasicHttpBinding_IImageryService");
            imageryService.GetMapUriCompleted += ImageryServiceGetMapUriCompleted;

            _currentImage = image;

            imageryService.GetMapUriAsync(mapUriRequest);

        }

        private void ImageryServiceGetMapUriCompleted(object sender, GetMapUriCompletedEventArgs e)
        {
            // The result is an MapUriResponse Object
            MapUriResponse mapUriResponse = e.Result;
            var bmpImg = new BitmapImage(new Uri(mapUriResponse.Uri));

            _currentImage.Source = bmpImg;
        }
    }
}


texts: Google maps not displaying until zoom in/out using openlayes2

texts:
I am using openlayers 2 for map displaying using following code,

    map = new OpenLayers.Map("mapDiv");

    map.addLayer(new OpenLayers.Layer.Google("Google Streets",{numZoomLevels: 20}));

    map.addLayer(new OpenLayers.Layer.Bing({
                                                name: "Bing",
                                                key: "AqTGBsziZHIJYYxgivLBf0hVdrAk9mWO5cQcb8Yux8sW5M8c8opEC2lZqKR1ZZXf",
                                                type: "Road",
                                                wrapDateLine: true
                                                }));
    map.addControl(new OpenLayers.Control.LayerSwitcher());

What happen is google maps isn't displaying until someone zoom in/out while bing maps get displayed. I don't have google key for map api. so currently using following JS.

<script src="https://maps.googleapis.com/maps/api/js?sensor=false" type="text/javascript" ></script>

I think its problem of not having google key. is it so?


texts: Why does the use of barrier lines significantly increase idw processing time?

texts: The IDW tool's documentation states that using barrier lines will significantly extend the processing time. It does not explain why this happen. Is anyone able to provide an insight into the inner workings of this tool?


texts:
I have a .CSV that is updated each time I make any changes to my input sheet in .XLSM.

When I add/display the .CSV by "Add Delimited Text Layer" I cannot then update the .CSV as QGIS is then using it.

Is it possible to have a link between the .CSV and QGIS that still allows me to edit the .CSV then changes automatically refresh in QGIS.

Do I need to do via a Join instead of delimited layer ?


texts: Use multiple colors in the same OpenLayers label

texts:
I made the following code that assigns a style to the field labels.
Now I would like each word or group of words in the label to be a different color previously selected by the user.
For example if we write T: 34 / P34 / S3400MP.
Then T: 34 should be reddish
P: 34 to be green
S3400MP. to be orange
Here is the code I tried with but not with much luck:
constructii.getSource().forEachFeature(function (feature) {

    console.log(feature.getProperties());
    var style = new ol.style.Style({
        fill: new ol.style.Fill({
            color: color.color_constructii
        }),
        stroke: new ol.style.Stroke({
            color: [0, 0, 230],
            width: 3
        }),
        text: new ol.style.Text({
            font: '15px Calibri,sans-serif',
            fill: new ol.style.Fill({
                color: color.color_eticheta
            }),
            stroke: new ol.style.Stroke({
                color: '#fff',
                width: 1
            })
        })
    });

    var data = "";

    var numar_const = "
" + feature.get('numar_const') + "
";
    var folosinta = "
" + feature.get('folosinta') + "
";
    var suprafata = "
" + feature.get('suprafata') + " MP.
";

    if ($('#etc_numar_const').is(":checked")) {
        data = data + "N:" + numar_const + "/";
    }
    if ($('#etc_folosinta_const').is(":checked")) {
        data = data + "F:" + folosinta + "/";
    }
    if ($('#etc_suprafata_const').is(":checked")) {
        data = data + "S:" + suprafata + "";
    }

    style.getText().setText(data);
    feature.setStyle(style);

})

The error is: The html :( tags also appear on the tag and the respective word tag is not colored in the given color :(
UPDATE:
my constructii color function is here :
features.forEach(function (feature) {

    feature.setId(undefined);
    feature.getGeometry().transform('EPSG:31700', 'EPSG:3857');
    //ol.extent.extend(extent, feature.getGeometry().getExtent());

    var const_stile = new ol.style.Style({
        fill: new ol.style.Fill({
            color: color.color_constructii
        }),
        stroke: new ol.style.Stroke({
            color: [0, 51, 204],
            width: 1
        }),
        text: new ol.style.Text({
            font: 'bold 11px "Open Sans", "Arial Unicode MS", "sans-serif"',
            placement: 'center',
            fill: new ol.style.Fill({
                color: 'white'
            })
        })
    });

    var data = "N:" + feature.get('numar_const') + "/F:" + feature.get('folosinta') + "/S:" + feature.get('suprafata') + " MP.";

    if (feature.get('layer') == "CONSTRUCTIE") {
        const_stile.getText().setText(data);
        feature.setStyle(const_stile);
    }

    vectorSource.constructii.clear();
    vectorSource.constructii.addFeatures(features);
    //console.log("features:");
    //console.log(features);
    map.getView().fit(extent);
})


texts: How to export selected records to a .txt file

texts:
I am working on a project where I need to create a tool/method to export the attributes of selected features to a .txt file.  I have point features containing addressing information, and the tool would need to export the attributes from one field to a text file.

I have been looking through the ArcPy documentation and found Copy Rows in the Tables Toolset - it can export selected features, but the only output option is exporting to a database table.  It seems like this should be a simple task, since exporting to a .txt file is one of the options when manually exporting from a table.

The ultimate goal here is to create a tool and publish it as a geoprocessing service to be used in a web application.  I would prefer to find a solution using Python, which I am just beginning to work with and have a little familiarity.

I am currently working with ArcGIS Desktop 10 SP1.


texts: How to synch Zoom control and MouseWheelZoom interaction in OpenLayers 6?

texts:
In OpenLayers v6.3.1 the map has a Zoom control and MouseWheelZoom is enabled. The map view has 10 resolutions.
In OpenLayers v5 prior to v6 using either the Zoom control or the MouseWheelZoom correctly incremented or decremented the map zoom level to one of the 10 possible zoom levels. Since upgrading the OpenLayers v6.3.1 the Zoom control works correctly for changing zoom level between the 10 levels. However, the MouseWheelZoom increments are greater, resulting in the map being zoomed to resolutions that are now part of the 10 resolutions defined for the view. In addition, the Zoom control does not move in synch w/ the MosueWheel.
How do you configure the MouseWheelZoom so it increments the same way as the Zoom control?


texts: Export .dbf maintaining column data type

texts:
I require a shapefile's .dbf to contain two types of text columns VARCHAR(n) and CHAR(n) to meet specs (attribute format) of gov receiving body. Have to send them the .shp, .dbf and .prj. Within arcmap or qgis add field only allows text and length for new fields. In postgres i can create the required format for example filed blah data type is character (10) and field blah2 data type is character varying(80), however if i export this as .dbf and then reimport it (to check data types are maintained) using the postgis plugin it changes the character(10) to character varying (10). Does anyone know a way to create a .dbf and export it maintaining set data types, i.e character (10) and character varying (80) for example.

Thanks


texts: Change in min/max values after Interpolation in QGIS

texts:
After rasterizing my shapefile with IWD Interpolation scheme, the min/max values change depending upon the distance coefficient and rows/columns values I choose, should I change the values in interpolated layer to the original min/max or would the new values generalize the original content and hence the output graph would be same?

What distance coefficient should be optimal if I need to point on exact locations in graph, how should one set the values for rows/columns (what should be their optimal range)?


texts: Simplest way to add built-in commands to an custom add-in

texts:
Currently in  Config.esriaddinx I have

<Button id="my_UniqueID" class="MyWrapperClass" caption="My Caption"
         image="Images\My_Image.png" ...> </Button>

and then I implement MyWrapperClass with:

protected override void OnClick()
{
    UID uid = new UIDClass();
    uid.Value = "esriArcMapUI.GraphicsToFeaturesCommand";

    ICommandItem cmd = ArcMap.Application.Document.CommandBars.Find(uid);
    if (cmd != null)
        cmd.Execute();
}

I'd like to just have the following in Config.esriaddinx

<Button id="esriArcMapUI.GraphicsToFeaturesCommand" />

and have the built-in command object do all the work for me.

Is this possible?


texts: GeoSever PointStacker Renders incorrect position

texts:
I've been trying the GeoServer Rendering Transformations - Point Stacker. According to this docs there's an optional cellSize parameter for setting the "Size of the cells in which to aggregate points, in pixels, default = 1".

I made the SLD based on the SLD example in the above OpenGeo docs.

I tried 2 cellSize (30 and 150). Here are the results:

Figure 1. GeoServer cluster (using Point Stacker) : there's displacement

Blue squares are points (exact, not clustered). Orange circles are clusters. Both layers (blue squares & orange circles) are WMS.

There are 41 points in the island.

The Graphic circles (defined in SLD as <PointSymbolizer>) seems to be positioned incorrectly. The above pictures show that displacement grow with cellSize parameter.

Why does this happen?

I saw PointStacker - Random locations at different zoom level asking about random position of Point Stacker, but no solid answer so far.

The below picture shows the same points (of another area). Again the blue squares are points (exact, not clustered) & orange circles are clusters.

Figure 2. GeoServer cluster (using Point Stacker) : there's displacement

Now because the zoom level is very high (more zoom-in than previous picture), the 4 points are not clustered. Ideally the orange points should be rendered in exact position as the blue point layer. But apparently GeoServer renders the orange points with some displacement. This is also asked in this question.

This picture also shows that displacement grow with cellSize parameter. Interesting enough that the displacements seem to be random.

Now let us see OpenLayers cluster. In contrast with GeoServer Point Stacker cluster (shown in Figure 1 above), OpenLayers cluster doesn't seem to have displacement. Here it is:

Figure 3. OpenLayers Cluster: no displacement

In the above figure, OpenLayers clusters GeoJSON vector points coming from GeoServer WFS. The right-side column is just the same with the left-side, only disabling the blue square for increasing visibility of the number label in OpenLayers clusters.

In the OpenLayers Cluster layer there's a circle with label "1" that is rendered in exact position as the blue square point, which I think is very make sense. To me, OpenLayers cluster is not perfect, but much-much better than GeoServer Point Stacker cluster (WMS).

Did I do anything wrong in my GeoServer Point Stacker cluster (WMS)?

The reason I am pursuing GeoServer Point Stacker cluster (WMS) over OpenLayers Cluster (WFS) is because I feel insecure with WFS. Serving my map as WFS means that I am distributing it in the raw format (text) that can easily be stored and reconstructed into vector map by (smart) user.

I'm using:

GeoServer 2.7.0 on Apache Tomcat 7.0.61 on Windows Server 2008 R2 Datacenter 64-bit on VMWare
JVM: Oracle Corporation: 1.7.0_80 (Java HotSpot(TM) 64-Bit Server VM)
OpenLayers 2.12
PostGIS - PostgreSQL 9.4.1 (on Windows Server 2008 R2 Datacenter 64-bit on VMWare)


texts:
I have a point shapefile and I create Thiessen (Voronoi) polygons programmatically using this scripting syntax:

CreateThiessenPolygons_analysis (in_features, out_feature_class, fields_to_copy)

However, each point is related with an area (i.e. the preferred size of each polygon) and I wish the Thiessen polygons to be weighted based on this field.

Is that possible and how?

Is there any relevant code in VBA?


texts: Earth Engine NodeJS API: Failed to load resource and CORS policy issue

texts:
I am developing a GEE nodeJS custom application. I have gone through the documentation and source(https://github.com/google/earthengine-api/tree/master/demos/server-auth-nodejs). I setup the service account and added the private key JSON file in the application. When I run the application it is showing:
i. Failed to load resource: the server responded with a status of 404 ()
and
ii. Access to XMLHttpRequest at 'https://earthengine.googleapis.com/map/projects/earthengine-......... from origin 'http://localhost:3000' has been blocked by CORS policy: No 'Access-Control-Allow-Origin' header is present on the requested resource.

Server.JS

const ee = require('@google/earthengine');
const express = require('express');
const handlebars  = require('express-handlebars');

const app = express()
  .engine('.hbs', handlebars({extname: '.hbs', cache: false}))
  .set('view engine', '.hbs')
  .use('/static', express.static('static'))
  .get('/', (request, response) => {
    var image = ee.Image('srtm90_v4');
    image.getMap({min: 0, max: 1000}, ({mapid, token}) => {
      console.log(mapid);
      console.log(token);
      response.render('index', {layout: false, mapid, token});
    });

  });

// Private key, in `.json` format, for an Earth Engine service account.
const PRIVATE_KEY = require('./privatekey.json');
const PORT = process.env.PORT || 3000;

ee.data.authenticateViaPrivateKey(PRIVATE_KEY, () => {
  ee.initialize(null, null, () => {
    app.listen(PORT);
    console.log(`Listening on port ${PORT}`);
  });
});

Script.JS

const EE_MAP_PATH = 'https://earthengine.googleapis.com/map';

/**
 * Initialize the Google Map and add our custom layer overlay.
 * @param {string} mapid
 * @param {string} token
 */
const initialize = (mapid, token) => {
  // Create an ImageOverlay using the MapID and token we got from Node.js.
  const tileSource = new ee.layers.EarthEngineTileSource({
    mapid,
    token,
    formatTileUrl: (x, y, z) =>
        `${EE_MAP_PATH}/${mapid}/${z}/${x}/${y}?token=${token}`
  });
  const layer = new ee.layers.ImageOverlay(tileSource);

  const myLatLng = new google.maps.LatLng(-34.397, 150.644);
  const mapOptions = {
    center: myLatLng,
    zoom: 8,
    maxZoom: 10,
    streetViewControl: false,
  };

  // Create the base Google Map.
  const map = new google.maps.Map(document.getElementById('map'), mapOptions);

  // Add the EE layer to the map.
  map.overlayMapTypes.push(layer);
};

index.hbs

<!DOCTYPE html>

    <title>SRTM</title>


texts: Using ArcPy dataSource Path with ArcSDE layers?

texts:
I have the following code:

When I run it:

conType = arcpy.mapping.ListLayers(imxd, "Construction_pt", df)[0]
print conType.dataSource
prints... "Database Connections\Edit - Editor.sde\Edit.SDE.Construction_pt"
and rest of program runs correctly

When another user runs it:

print conType.dataSource
prints... "\Edit.SDE.Construction_pt"
program crashes because cannot find path

My question is why is it only displaying the full path for me and not other users?

ArcGIS 10,
Windows 7


texts:
I'm trying to create a buffer but the layer I get is very large

Buffer settings:

Bad result:

What am I doing wrong? The idea is to have an area of influence of 10 meters


texts: Enabling multithreaded processing in QGIS?

texts:
I have found multiple links that appear to say that this now exists in QGIS 2.2, but it appears that it's only using one core to full capacity.

This is running on Ubuntu 14.04

The current process is creating a large amount of regular points. How can I make that process use more than the one core?


texts: Making stack from many images in Google Earth Engine?

texts:
Using the code below, I've computed 20 ndvi images for modis data.

How can I stack and save them in a single file?

var mod13 = ee.ImageCollection("MODIS/006/MOD13Q1")
.filterBounds(geometry)
.filterDate("2000-01-01","2001-01-01")
.select("NDVI");

var mod13ndvi = mod13.map( function(img){
  return img.multiply(0.0001)
  .copyProperties(img,['system:time_start','system:time_end']);
});

print(mod13ndvi);


texts: Export classified raster Class Values from DEM

texts:
I've classified an elevation raster into 4 classes. Is there a way to export these class values (0, 1, 2, 3) to a separate raster? I've tried exporting the raster to file, but this just gives the same elevation values.

In other words, I want to replace the pixel values in my DEM with their class values.


texts: How to animate .olControlLayerSwitcher with Jquery?

texts:
Not know much about css but some how got strucked in setting css property. I am using Openlayers .olControlLayerSwitcher class.  Here without using !important keyword, LayerSwitcher will not get move. Here is my code:

Its working.But very fast. I want to move it with some speed intervals like jquery animate.
Please suggest it.

$(".olControlLayerSwitcher").css("cssText", "top: 30px !important;");

$('.olControlLayerSwitcher').css("cssText", "top: 55px !important;");


texts: ArcObjects C++ - Has anyone successfully used AoToolbarAddCommand?

texts:
I am working with the C++ API for ArcObjects (Linux, QT), and I am running into problems with when creating a custom button on a toolbar. I have read up on the C++ specific implementation of creating a custom command (inheriting from AoCommandBase and adding the button with AoToolbarAddCommand).

After many hours of fiddling around, I have come to the conclusion that AoToolbarAddCommand simply isn't doing its job. I can add ESRI predefined commands with no problem, but my custom command shows up as a useless grey square (my icon doesn't even display). My debug output suggests that none of the functions defined in my command are being called during setup (or ever for that matter), and they definitely should be. If it helps, excerpts from my code are below:

From the add_toolbar_items function:

CComVariant varTool;

long itemindex;

varTool = L"esriControlCommands.ControlsOpenDocCommand";
ipToolbar->AddItem(varTool, 0, -1, VARIANT_FALSE, 0,
        esriCommandStyleIconOnly, &itemindex);

panRight = new PanRight(); /*panRight is global, PanRight is my command class*/

AoToolbarAddCommand(ipToolbar, panRight, esriCommandStyleIconOnly);

The PanRight class (my custom button) is stolen from a sample, but I changed the code for loading bitmaps because the provided code did not work. All the other functions remain unchanged, except for a cout statement letting me know if the function were ever to run.

PanRight::PanRight()
{
  m_ipHookHelper.CrreateInstance(CLSID_HookHelper);
  std::cout<<"panright constructor\n";

  //load bitmap
  IPictureElement3Ptr ipPict(CLSID_PictureElement);
  HRESULT hr = ipPict->ImportPictureFromFile(CComBSTR(L"<pathToPic>"));
  if(SUCCEEDED(hr))
  {
    OLE_HANDLE hBitmap;
    hr = ipPict->get_Handle(&hBitmap);
    if(SUCCEEDED(hr))
      {m_hBitmap = hBitmap;}
  }
}

Does anyone have any ideas as to why AoToolbarAddCommand might not be doing its job? Has anyone successfully added a custom command using the C++ API?


texts:
I recently shifted maps and waypoints from Delorme Xmap to QGIS on my Mac and wish to use iPhone/iGIS as the GPS (together with the Garmin GLO receiver). Everything is good in QGIS, but I need USGS topos as the base maps in iGIS, not the generic Google stuff used by iGIS. The app wants users to create tiled KML maps with MapTiler, which has a free version and pricier ones. I did create the tiles in MapTiler and they open just fine in Google Earth, but 1) the quality is very poor and 2) iGIS does not succeed in importing them. Seems like a very basic operation. Has anyone done this with any success?

[I should add that iGIS appears to be the only iPhone app out there that works with shapefiles and which reads / sends them to and from the Dropbox. This is important because QGIS only allows editing of waypoints in a shapefile layer (not CVS or GPX). If one creates the QGIS layer directly from the shapefile in the desktop Dropbox, the app communicates seamlessly with the desktop and waypoints are easily updated on both systems. If there is another iPhone app out there that does this, please let me know.]

The steps in the answer below worked well. Also, I was able to overcome the map quality issue in the 10KX10K file size limit in the free version of MapTiler by slicing my large geotiff file into 12 separate quadrants in QGIS and then processing each of those separately through MapTiler.


texts: Add / remove an ImageOverlay from Leaflet map

texts:
My goal is to add / hide an L.ImageOverlay image from my leaflet map on click. This code should work, but it only adds image and does not hide it on a following click.

var map = L.map('mainMap').setView([56.06, 92.500], 10);

function loadImage(e) {

    var imageUrl = 'data/kr1_cut.png',
        imageBounds = [[55.819224, 92.170872], [56.324808, 92.480904]],
        image = L.imageOverlay(imageUrl, imageBounds);
    // remove image layer if it already exists
    if (map.hasLayer(image)) {
        map.removLayer(image);
    }
    else {
        image.addTo(map);
    }
}

map.on('click', loadImage);

So what am I doing wrong?


texts: Effectively displaying a large amount of data on Google Maps APi

texts:
I have about 20 layers that I'd like to put on a Google maps API. I assume I'll need some form of a tile cache server and perhaps under mapserver or geoserver. These layers amount to about 400mb so the option of using static maps using kml/kmz is not an option. The data exists on two geodatabases, ESRI and PostGIS.

Can you provide insights on publishing large amounts of data to google maps api? especially in a way that allows users to easily do queries,load and unload layers at a relatively acceptable speed? Are there any good reads that improve the understanding of the issues and challenges?


texts: Transforming takes long time with st_transform

texts: I have two shape files contains considerably large data. I exported them to my geodatabase(using arccatolog), which takes 2 days. The problem is their srid must be changed. So i tried creating new table from the tables i exported with sde.st_transform method but the query did not end. At the end 5th day i stop the query. In previous versions of the arcgis server i use shp2sde function from linux command line and give the srid parameter directly. But in 10.2 version i could not find such a functionality. So i got only arcmap&arccatolog and SQL options. I need to set srid of this shape files before i exported them. So what can i do?


texts: Cannot change legend shape in print layout of QGIS

texts:
I am using QGIS 3.x.
I would like to customize the legend in my print layout. Unfortunately I am quite restricted, as I can't change my shape appearance.
According to the link below:
https://docs.qgis.org/2.14/en/docs/user_manual/print_composer/composer_items/composer_legend.html
it looks like I can only change the height or width of my shapes.

I have some symbols, which I would like to appear as circles instead of rectangles or thick lines.
There is some solution at Is it possible to change legend symbols in the QGIS print composer legend? but it looks like I need some Python skills for it.
Is there any simple way to customize the legend symbol one by one in the print layout area or should I go back to the geometry generator style and do something like at QGIS legend in layout manager?


texts: Proper way to assign CRS in rgdal 1.5-23

texts:
Perhaps my question has been already answered but if it did, I haven´t be able to find it. I´ve updated the rgal package to version 1.5-23 (GDAL 3.2.1 & PROJ 7.2.1) and when I run the following line
crs(res_data_PA) <- LAEA

I get the following warning message:
Warning message:
In showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj = prefer_proj) :
  Discarded datum Unknown based on GRS80 ellipsoid in Proj4 definition

I´m aware that I can get rid of it and since it´s not an error nothing should change in my data. I don´t like to suppress any type of message from R so, I would like to know how is it now the proper way to assign the CRS to my raster.


texts: How to separate a list of features from its shapefile_

texts:
I tried to copy all features in a list of dictionaries. The code looks as follows:

def iterate_features(self):
    feature_set = arcpy.SearchCursor(self.address)
    for feature in feature_set:
        yield feature

def create_new_dict(self):
    iterator = self.iterate_features()
    if u"Matched_to" not in self.fieldList[:]:
        #arcpy.AddField_management(input_data, u"Matched_to", "TEXT")
        for item in iterator:
            feature = dict()
            for f in self.fieldList:
                feature[f] = item.getValue(f)
            feature[u"Matched_to"]=[]
    else:
        for item in iterator:
            feature = dict()
            for f in self.fieldList:
                feature[f] = item.getValue(f)
        self.new_features_dict.append(feature)

Now that I want to fetch the Matched_to list, what it returns is a Unicode object instead of a list. Hence I cannot append anything to it. Moreover, the field Matched_to is added to the shapefile without being updated or edited. I think it has something to do with shallow and deep copy but I don't know how and where to implement it.


texts: Sphere rendering inside out using QGIS on Ubuntu?

texts:
I'm trying to use QGIS to show earthquakes around a 3d globe. I changed my CRS to 'Sphere_Van_der_Grinten_I' and it rendered like this:

Here's a report from Hardinfo on my system and display:

Processor:    12x Intel(R) Core(TM) i7-6800K CPU @ 3.40GHz
Memory:   65809MB (10399MB used)
Operating System: Ubuntu 16.04.2 LTS
Resolution    1920x1080 pixels
OpenGL Renderer:  Gallium 0.4 on llvmpipe (LLVM 3.8, 256 bits)
X11 Vendor:   The X.Org Foundation

Updating my OpenGL renderer doesn't seem to have done anything; this screenshot is unchanged before and after I did apt-get install mesa-utils

How do I make this sphere work?


texts:
I have an issue with merging several tiles. I am loading the original XYZ files into QGIS (works fine) and try to merge them afterwards yia Raster/miscellaneous/merge . If I only merge 2 tiles this works without any problem but as soon as I try to merge more than 2 tiles at once I get :
GDAL-Befehl:
gdal_merge.bat -ot Float32 -of GTiff -o "C:/Users/Jonas Wischnewski/AppData/Local/Temp/processing_ITBVYL/4c383c33f450475aa0b8886a7f7a8c2a/OUTPUT.tif" --optfile "C:/Users/Jonas Wischnewski/AppData/Local/Temp/processing_ITBVYL/4c6ee11824de43f9ab2d05a72c02f4b5/mergeInputFiles.txt"
GDAL-Befehlsausgaben:
**ERROR 1: Attempt to create 401x-1198 dataset is illegal,sizes must be larger than zero.**

The XYZ Files don't have anything suspicious IMO:
362000.0 5814000.0 39.315
362001.0 5814000.0 39.302
362002.0 5814000.0 39.304
...

Interestingly when I use other XYZ DEM Data from a different source with a similar input:
686000.00 5848000.00 33.59
686005.00 5848000.00 32.40
686010.00 5848000.00 33.69
...

this works perfectly fine.
Merging the rasters via Mosaic does work in this case, but it changes the z-values to an extend where it is not usable anymore for further processing.


texts: Query where one field includes the text string in another field

texts:
I am trying to query mailing address records that match physical address records, but the two fields sometimes have data in a different format.  For example, the property's mailing address may be shown as '123 NW Elm Street', while the physical address may be shown as '123 Elm St'.  Even if I can just query whether the address contains the same street name that would help me if I knew if 'Elm' was anywhere in the address.  There is a field with just street name, with no number or prefix or suffix.

I am trying to do a definition query to select nearly matching records with wildcards with something like this:

"StreetAddress" LIKE "%StreetName%"

I have tried numerous variations with no luck, including '%'+"STreetName"+'%' etc.

Is there a "FIELD_A" INCLUDES "FIELD_B" query or a "FIELD_A" CONTAINS "FIELD_B" command?


texts:
I want to reproject a raster from the rather common EPSG:4326 (WGS 84) to EPSG:5677 (DHDN / 3-degree Gauss-Kruger zone 3 (E-N)) :

gdalwarp -s_srs EPSG:4326 -t_srs EPSG:5677 -r near -of GTiff C:/.../raster_og.tif C:/Users/corne/OneDrive/gis/rheinauen/temp/raster_tg.tif

However, all I get is the following error:

ERROR 1: Too many points (8100 out of 8100) failed to transform, unable to compute output bounds.

Here is the complete gdal info of the source file :

GDAL command output:
Driver: GTiff/GeoTIFF

Files: C:/Users/corne/OneDrive/gis/rheinauen/wsel_hq1_IST.tif

Size is 4437, 8511

Coordinate System is:

GEOGCRS["WGS 84",

DATUM["World Geodetic System 1984",

ELLIPSOID["WGS 84",6378137,298.257223563,

LENGTHUNIT["metre",1]]],

PRIMEM["Greenwich",0,

ANGLEUNIT["degree",0.0174532925199433]],

CS[ellipsoidal,2],

AXIS["geodetic latitude (Lat)",north,

ORDER[1],

ANGLEUNIT["degree",0.0174532925199433]],

AXIS["geodetic longitude (Lon)",east,

ORDER[2],

ANGLEUNIT["degree",0.0174532925199433]],

ID["EPSG",4326]]

Data axis to CRS axis mapping: 2,1

Origin = (3434740.290000000037253,5421087.200000000186265)

Pixel Size = (1.000056344376831,-1.000064622253532)

Metadata:

AREA_OR_POINT=Area

Image Structure Metadata:

INTERLEAVE=BAND

Corner Coordinates:

Upper Left ( 3434740.290, 5421087.200) (Invalid angle,Invalid angle)

Lower Left ( 3434740.290, 5412575.650) (Invalid angle,Invalid angle)

Upper Right ( 3439177.540, 5421087.200) (Invalid angle,Invalid angle)

Lower Right ( 3439177.540, 5412575.650) (Invalid angle,Invalid angle)

Center ( 3436958.915, 5416831.425) (Invalid angle,Invalid angle)

Band 1 Block=4437x1 Type=Float64, ColorInterp=Gray

NoData Value=nan


texts: How to implement signals on add or modify feature?

texts:
I've created a test plug-in in QGIS successfully, but have had a hard time implementing a signal on a vector layer event (add or modify feature).  I've add this to my initGUI method:

QObject.connect(self.iface.activeLayer(),SIGNAL("featureAdded(QgsFeatureId)"),self.addedGeometry)

and this definition:

def addedGeometry(self, intValue):
        QMessageBox.information( self.iface.mainWindow(),"Info", "Hi!!!" )

I get the following error when I load my plug-in:

QObject.connect(self.iface.activeLayer(),SIGNAL("featureAdded(QgsFeatureId)"),self.addedGeometry)
TypeError: arguments did not match any overloaded call:
  QObject.connect(QObject, SIGNAL(), QObject, SLOT(), Qt.ConnectionType=Qt.AutoConnection): argument 1 has unexpected type 'NoneType'
  QObject.connect(QObject, SIGNAL(), callable, Qt.ConnectionType=Qt.AutoConnection): argument 1 has unexpected type 'NoneType'
  QObject.connect(QObject, SIGNAL(), SLOT(), Qt.ConnectionType=Qt.AutoConnection): first argument of unbound method must have type 'QObject'

If anyone could help I would be very grateful!  I can't find a lot of examples of implementing signals on QGIS objects.


texts: Hover on Airport Data using OpenLayers + Mapfile

texts:
I need to create a simple hover using openlayers for worldwide airport data. I'm thinking maybe the associated address, maybe a picture in regards to the airport, ...anything really. It's just an html file with some js, nothing too crazy but I've been out of the GIS/coding game for a while.

This is the current code for the html file. As you can see, there is an associated map file. Since there are so many airports, I'm guessing you would want to put the code in the mapfile (in regards to the hover(s))

Here is the map file. Apologies in advance for a little screwed up indentation.

So I guess my main question is: do I need to be messing with the map file or just the html? I found another example where they were referencing the hover in the map file but that example was not also a straight html file...if that makes sense. I looked through the openlayers documentation but it was pretty tough to grasp.

Any input is greatly appreciated.


texts: Display different symbols in QGIS using Point Displacement renderer

texts:
How to display different symbols (not different in size, different in shape and color) on the same point in QGIS 2.8 using Point Displacement renderer?

I want my layout with symbols around the centroid of areas, and symbols from rule-based renderer.

For some reason it seems QGIS forces all the symbols to be the same.

Worse case scenario I will do it manually.


texts:
I have installed Tomcat on a Windows server and have successfully deployed Geoserver. I am trying to deploy geoexplorer by building the war file from github and loading it into tomcat to deploy. It loads, but when I go to localhost:8080/geoexplorer, I get a 404 error. How do I load it into tomcat if not that way?

EDIT:
I did load Geoexplorer into the webapps folder of Tomcat and restarted Tomcat.
Here's what is in the log about Geoexplorer:

Jul 07, 2014 8:34:34 AM org.apache.catalina.startup.HostConfig deployWAR
INFO: Deploying web application archive C:\Program Files\Apache Software Foundation\Tomcat 7.0\webapps\geoexplorer.war
Jul 07, 2014 8:34:40 AM org.apache.catalina.util.SessionIdGenerator createSecureRandom
INFO: Creation of SecureRandom instance for session ID generation using [SHA1PRNG] took [240] milliseconds.
Jul 07, 2014 8:34:40 AM org.apache.catalina.startup.HostConfig deployWAR
INFO: Deployment of web application archive C:\Program Files\Apache Software Foundation\Tomcat 7.0\webapps\geoexplorer.war has finished in 5,844 ms


texts: findAndReplaceWorkspacePath ValueError: Layer: Unexpected error

texts:
Here is my code to replace the old SDE connection file with the new one. I am getting an error due to grouped layers, that I know. I have been trying to fix this but regardless of what I try I keep getting the ValueError: Layer: Unexpected error due to the grouped layer.

Can anyone let me know what I am doing wrong?

import arcpy
import os

for root, subFolders, files in os.walk(r"C:\Users\x\Desktop\lyrtest"):
    for filename in files:
        fullpath = os.path.join(root, filename)
        basename, extension = os.path.splitext(fullpath)
        if extension.lower() == ".mxd":
            print "------------------------------"
            print filename
            MXD = arcpy.mapping.MapDocument(fullpath)
            for lyr in arcpy.mapping.ListLayers(MXD):
                if not lyr.isGroupLayer:
                    if lyr.supports("DATASOURCE"):
                        source = lyr.workspacePath
                        print "%s -> %s" % (lyr, source)
                        basename, extension = os.path.splitext(source)
                        if extension.lower() == ".sde":
                            datapath = r"C:\Users\x\Desktop\lyrtest\xdbprod.sde"
                            lyr.findAndReplaceWorkspacePath(source, datapath)
                else:
                    print "exception"
            MXD.save()
            del MXD

Here is the error:

Traceback (most recent call last):
  File "C:\Users\x\Desktop\lyrtest\mxdchng5.py", line 20, in <module>
    lyr.findAndReplaceWorkspacePath(source, datapath)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.3\ArcPy\arcpy\utils.py", line 182, in fn_
    return fn(*args, **kw)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.3\ArcPy\arcpy\_mapping.py", line 696, in findAndReplaceWorkspacePath
    return convertArcObjectToPythonObject(self._arc_object.findAndReplaceWorkspacePath(*gp_fixargs((find_workspace_path, replace_workspace_path, validate), True)))
ValueError: Layer: Unexpected error


texts:
I've migrated a VB script that uses ArcObjects to a C++ Win32 DLL (using Visual Studio 2013).

Now I would like to test my DLL by launching it from my ArcMap 10.3 installation.

The problem is that when I try to add it through "Customize -> Customize Mode -> Add from file" it does not succeed and raises the following message: "Can't load type library from specified file", even when I run ArcMap as Administrator.
The same way, the "ESRIRegAsm.exe" I tried to use in command line also refuses to register my DLL.

I've made research but I could not find any answer to make it work.
But maybe I am not doing this the right way.

Do I need to implement a COM component to add my DLL to ArcMap?


texts: Calculating VRMS and HRMS from NMEA data?

texts:
I am trying to write code to parse the NMEA strings coming from survey-grade RTK/GNSS receivers. The requirements I have are to calculate the 3DRMS (distance RMS at 97.5% probability), VRMS (vertical RMS) and HRMS (horizontal RMS) values (amongst others).

Going by this article, I've been able to find a formula to derive the 3DRMS value from the GST NMEA message:

In pseudocode: 3DRMS = 3 * Math.sqrt(Math.pow(dx,2) + Math.pow(dy,2))

where dx is Latitude error and dy is Longitude error.

So, for example, if the GST message string is:

$GPGST,172814.0,0.006,0.023,0.020,273.6,0.023,0.020,0.031*6A

Then:

3DRMS = 3 * Math.sqrt(Math.pow(0.023,2) + Math.pow(0.020,2)) = 0.091

However, I'm struggling to find a definition or formula to calculate the VRMS and HRMS values, and therefore which NMEA strings I need to parse to calulate them.

I've found that the companion apps for several of the receivers display these, for example Eos Tools Pro displays them when used with the Arrow 200 receiver.

I've also found a PhD thesis which contains some equations, but I'm having trouble understanding how to abstract them into a generic formula. However, it looks like PDOP is involve, so I would assume I need to parse the GSA message in addition to the GST.

Where can I find the relevant formulas or how to calculate these values?


texts: Getting the Field Alias from a layer, from within a Map document

texts:
I've got a very simple script, that loops through all the layers in a mxd and loops through the fields, and prints out the Alias Name for the field, but I am having some problems, it seems to be printing the field alias of the feature class (in SDE) not the value I have changed it to  in the map document...

Here's the script

import os, arcpy

myMxdFolder = r"D:\Maps"
count = 0

for mxd in os.listdir(myMxdFolder):
    #print mxd
    if mxd.endswith(".mxd"):
        #if mxd.startswith("COMPASS") or mxd.startswith("COLGIS"):
        if mxd.startswith("COLGIS_City_Occupiers.mxd"):
            mxdoc = arcpy.mapping.MapDocument(myMxdFolder + os.sep + mxd)
            filePath =  mxdoc.filePath
            print filePath
            for lyr in arcpy.mapping.ListLayers(mxdoc):
                desc = arcpy.Describe(lyr)
                print lyr.name
                fields = desc.fields
                for field in fields:
                    print field.aliasName

The output is

OBJECTID
UPRN
ADDRESS
OCCUPIER
SHAPE

but the Aliases in the map document i have set it to are:

OBJECTID
HIGHLIGHT::UPRN
Address
Occupier
SHAPE

and that is what I was expecting, have a missed a trick somewhere in the code?


texts: Comparing ASTER/SRTM dataset in ArcGIS Desktop?

texts:
Currently I am working with two data set ASTER and SRTM. I want to check their accuracy with reference to toposheet.

How would I do that in ArcGIS Desktop 10.4?

I want to check accuracy in terms of area and length of my study area.


texts: Running processing toolbox outside QGIS in standalone app?

texts:
I am developing a standalone app in Python, using QGIS libraries including processing algorithms. I am developing in Linux environment and I put the libraries in a folder. Also I put the processing folder from QGIS in that same folder. I run the app but when I import the processing folder, it gives me nothing (not even error messages).

I import like that:

import processing.core.Processing

Can you help me?

CODE:

import sys
import qgis.core
import PyQt4.QtCore
import PyQt4.QtGui
app = qgis.core.QgsApplication([],True)
qgis.core.QgsApplication.setPrefixPath("/usr", True)
qgis.core.QgsApplication.initQgis()

(from here the program doesn't work)

from processing.core.Processing import Processing
inputLayer = 'c:....'
inputMask = 'c:....'

extent = ....
cellsize = ....
outPath = 'c:.....'

Processing.Initialize()
Processing.runalgorith(.....)

I downloaded the folder with the processing version 2.2.0-2 and I copied to my script folder. I run and it gives me this error:

Traceback (most recent call last):
  File "teste.py", line 12, in <module>
    from processing.core.Processing import Processing
  File "/home/p/Programas/processing/__init__.py", line 28, in <module>
    from processing.tools.dataobjects import *
  File "/home/p/Programas/processing/tools/dataobjects.py", line 33, in <module>
    from processing.core.ProcessingConfig import ProcessingConfig
  File "/home/p/Programas/processing/core/ProcessingConfig.py", line 30, in <module>
    from processing.tools.system import *
  File "/home/p/Programas/processing/core/processing.py", line 33, in <module>
    from processing import interface
ImportError: cannot import name interface

I got it!! What is missing in my case was to define the PYTHONPATH to the folder with the processing. Now I have the correct running of the code.


texts: GEE Export.image.toDrive output not matching bounding box?

texts:
I am new to GEE and javascript. Copy-pasting here and there, I managed to put together this code that extracts NDVI from cloud-free Sentinel-2 images on GEE. I created a fusion table from a KML with the polygon I need as bounding box. When I run the code, the result preview on GEE is exactly what I want, with boundaries matching the FeatureCollection geometry from the fusion table.

but when I export, the image does not match // Display NDVI results on map.

See screenshot: exported results cover an area that is not matching my geometry var

Export.image.toDrive({
  image: NDVI.clip(geometry),

Here I exported with maxPixels: 10E10,

Here with maxPixels: 10E11,

I have a 1TB GoogleDrive subscription, so space should not be the issue.

What am I doing wrong?

Here is the code:

//Import GEE Feature Collection (Somaliland kml)
var geometry = ee.FeatureCollection('ft:14f9y--LmmeDiBOkMFWA2eef5KvjrPRFDH68IRPfQ');

// Create image collection of S-2 imagery for the perdiod 2016-2018
var S2 = ee.ImageCollection('COPERNICUS/S2')
//filter start and end date
.filterDate('2018-08-01', '2018-08-31')
//filter according to drawn boundary
.filterBounds(geometry);
// Function to mask cloud from built-in quality band
// information on cloud
var maskcloud1 = function(image) {
var QA60 = image.select(['QA60']);
return image.updateMask(QA60.lt(1));
};
// Function to calculate and add an NDVI band
var addNDVI = function(image) {
return image.addBands(image.normalizedDifference(['B8', 'B4']));
};
// Add NDVI band to image collection
var S2 = S2.map(addNDVI);
// Extract NDVI band and create NDVI median composite image
var NDVI = S2.select(['nd']);
var NDVI = NDVI.median();

// Create palettes for display of NDVI
var ndvi_pal = ['#d73027', '#f46d43', '#fdae61', '#fee08b', '#d9ef8b',
'#a6d96a'];
// Display NDVI results on map
Map.addLayer(NDVI.clip(geometry), {min:-0.5, max:0.9, palette: ndvi_pal}, 'NDVI');

// Export image
Export.image.toDrive({
  image: NDVI.clip(geometry),
  description: 'NDVI_S2_Aug2018_2',
  scale: 10,
  maxPixels: 10E11,
  fileFormat: 'GeoTIFF',
});

Based on @Rodrigo E. Principe suggestion, I have changed the code to:

// Export image
Export.image.toDrive({
  image: NDVI,
  region: geometry.geometry().getInfo(),

The result is good, as it respect the minimum bounding geometry enclosing my input feature (as per fusiontable), but it is not clipped to the polygon boundaries. I guess I will have to clip the result with a GIS software outside GEE.


texts: Switch between ArcGIS and R changes my shapefile projection?

texts:
Not the best practice, but somethimes I need to change between ArcGIS 10.6 and R (version 3.6.1) to process my shapefiles. For example, I want to read my shapefile in R -> process it -> export it back to ArcGIS.

But, the projection of my data changes in between: I read the shapefile in one coordinate system, export it in same CRS, but ArcGIS reads it as different projection.

Let's see on watershed dataset:

I work on Finland watersheds. In ArcGIS, it has projection EUREF_FIN_TM35FIN, WKID 3067

After loading this file to R, the projection is a string:

library(rgdal)
library(sp)

setwd("U:/rawData/Finland/watershed/valuma")

watershed = readOGR(dsn = "U:/rawData/Finland/watershed/valuma",
                    layer = "Jako3")

proj4string(watershed)
# [1] "+proj=utm +zone=35 +ellps=GRS80 +units=m +no_defs"

This is how the this projection string should look like reffering to  spatialreference.org:

+proj=utm +zone=35 +ellps=GRS80 +units=m +no_defs

But, after exporting this file to new .shp:

# Add the stand geometry
writeOGR(obj    = watershed,
         dsn    = getwd(),
         layer  = "watershed_no",
         driver = "ESRI Shapefile",
         overwrite_layer=TRUE)

ArcGIS now reads this projection as Custom projection (even not listed within ArcGIS library!!!) as UTM_Zone_35_Northern_Hemisphere!!

I have tried to transform my data to correct CRS before exported in R:

watershed.t <- spTransform(watershed, CRS("+init=epsg:3067"))

But, ArcGIS still reads this as a Custom CRS. How to define my CRS in R, to be correctly read in ArcGIS? Or how can I avoid its change between softwares?

EDIT:

I have tried both suggestions:

Specify the watershed@proj4string@projargs by adding +init=epsg:3067

watershed@proj4string@projargs <- paste0('+init=epsg:3067 ', watershed@proj4string@projargs)

transform my data using spTranform and specify epsg:

watershed.t <- spTransform(watershed, CRSargs(CRS("+init=epsg:3067")))

but they both produce the same output: Custom projection UTM_Zone_35_Northern_Hemisphere
Authority: Custom in ArcGIS 10.6.

Do you have any other suggestions?


texts: Should resampling/downsampling a raster (using Rasterio) cause the coordinates to change?

texts:
I was able to successfully downsample a raster file, after getting Stack Exchange help. I loaded the new downsampled file into QGIS along with the original file, to check that things went well. The tricky thing was that the new image seems a bit shifted and smaller. So it seems like the coordinates in the image changed along with its size.

I was not sure what to expect after the operation, but I imagined that the two files should have the same coordinates--because the CRS was preserved in the downsampling. I figured that the two images would be the same size, but one would be blurrier than the other because its resolution was downsampled.

Is this image below the way that things should look after downsampling?

Here is the code I used to resize.

from contextlib import contextmanager
import rasterio
from rasterio import Affine
from rasterio.enums import Resampling

dat = 'original_image.tif'

@contextmanager
def resample_raster(raster, scale=2):
    t = raster.transform

    # rescale the metadata
    transform = Affine(t.a / scale, t.b, t.c, t.d, t.e / scale, t.f)
    height = raster.height / scale
    width = raster.width / scale

    profile = src.profile
    profile.update(transform=transform, driver='GTiff', height=height, width=width, crs=src.crs)

    data = raster.read( # Note changed order of indexes, arrays are band, row, col order not row, col, band
            out_shape=(int(raster.count), int(height), int(width)),
            resampling=Resampling.cubic)

    with rasterio.open('image-15cm.tif','w', **profile) as dst:
        dst.write(data)
        yield data

with rasterio.open(dat) as src:
    with resample_raster(src, 3.5) as resampled:
        print('Orig dims: {}, New dims: {}'.format(src.shape, resampled.shape))
        print(repr(resampled))


texts:
What software or template does Max Fisher use to create his maps for the Washington Post?

Do you recognise the software used to make this map?

As an example:

More examples here:
http://www.washingtonpost.com/blogs/worldviews/wp/2014/01/13/40-more-maps-that-explain-the-world/


texts: How do I remove the group name in the QGIS legend renderer

texts:
I am using QGIS 2.8.1 and the print composer to create a map for export. There is a lot of positive things to be said about the QGIS print composer. It gives the user good access to the rendering.

I want the legend to render without the group name. In this example I have two categories of protected areas I want to keep in the legend. The two categories are from the same spatial data source, in this case a shapefile with attributes making it possible to categorize styling. In short I do not want the text "Protected areas" in my legend.

In ArcGIS the below dialogue box alows for this functionality.

How can I do this in QGIS (2.8.1)? Is it possible?


texts: Finding Nearest Line Segments to Point using shapely?

texts:
Background

From a known Point, I require to establish the nearest surrounding "visible perimeter" against a table of MultiLineStrings, as shown on the diagram.

I've searched this site with a number of terms (e.g. minimum edge, minimum perimeter, nearest neighbour, clip, containing polygon, visibility, snap, cut nodes, ray-trace, flood fill, inner boundary, routing, concave hull) but cannot find any previous question which seems to match this scenario.

Diagram

The green circle is the known Point.
The black lines are the known MultiLineStrings.
The grey lines are an indication of a radial sweep from the known Point.
The red points are the nearest intersection of the radial sweep and the MultiLineStrings.

Parameters

The Point will never intersect the MultiLineStrings.
The Point will always be nominally centred within the MultiLineStrings.
The MultiLineStrings will never fully enclose the Point, therefore the perimeter will be a MultiLineString.
There will be a table containing approximately 1,000 MultiLineStrings (normally containing a single line of about 100 points).

Methodology Considered

Undertake a radial sweep by constructing a series of lines from the known Point (at, say, 1 degree increments).
Establish the nearest intersection point of each radial sweep line with the MultiLineStrings.
When one of the radial sweep lines doesn't intersect with any of the MultiLineStrings, this would indicate a gap in the perimeter which would be accommodated in the perimeter MultiLineString construction.

Summary

Whilst this technique will find the nearest intersections, it won't necessarily find all the closest perimeter node points, dependent on the resolution of the radial sweep. Can anyone recommend an alternative method to establish all the perimeter points or supplement the radial sweep technique with some form of buffering, sectoring or offsetting?

Software

My preference is to use SpatiaLite and/or Shapely for the solution but would welcome any suggestions which could be implemented using open source software.

Edit: Working Solution (based on answer by @gene)

from shapely.geometry import Point, LineString, mapping, shape
from shapely.ops import cascaded_union
from shapely import affinity
import fiona

sweep_res = 10  # sweep resolution (degrees)
focal_pt = Point(0, 0)  # radial sweep centre point
sweep_radius = 100.0  # sweep radius

# create the radial sweep lines
line = LineString([(focal_pt.x,focal_pt.y), \
                   (focal_pt.x, focal_pt.y + sweep_radius)])

sweep_lines = [affinity.rotate(line, i, (focal_pt.x, focal_pt.y)) \
               for i in range(0, 360, sweep_res)]

radial_sweep = cascaded_union(sweep_lines)

# load the input lines and combine them into one geometry
input_lines = fiona.open("input_lines.shp")
input_shapes = [shape(f['geometry']) for f in input_lines]
all_input_lines = cascaded_union(input_shapes)

perimeter = []
# traverse each radial sweep line and check for intersection with input lines
for radial_line in radial_sweep:
    inter = radial_line.intersection(all_input_lines)

    if inter.type == "MultiPoint":
       # radial line intersects at multiple points
       inter_dict = {}
       for inter_pt in inter:
           inter_dict[focal_pt.distance(inter_pt)] = inter_pt
       # save the nearest intersected point to the sweep centre point
       perimeter.append(inter_dict[min(inter_dict.keys())])

    if inter.type == "Point":
       # radial line intersects at one point only
       perimeter.append(inter)

    if inter.type == "GeometryCollection":
       # radial line doesn't intersect, so skip
       pass

# combine the nearest perimeter points into one geometry
solution = cascaded_union(perimeter)

# save the perimeter geometry
schema = {'geometry': 'MultiPoint', 'properties': {'test': 'int'}}
with fiona.open('perimeter.shp', 'w', 'ESRI Shapefile', schema) as e:
     e.write({'geometry':mapping(solution), 'properties':{'test':1}})


texts: Geoserver time dimension features not showing

texts:
I have a WMS layer from PostGIS data store. Setup is below on a figure. I have two attributes start and end of interval.

Temporal extent of layer is obviously calculated from start attribute. See two figures below:

Geometry viewer from PostGIS shows all the features. See figure below.

If I go to preview layer in Geoserver without time filter I have some features missing (with most recent end timestamp 2018-12-18T00:00:00).

If I include time filter which corresponds to temporal extent of the layer (i.e. &time=2018-10-05T00:00:00.000Z/2018-12-18T00:00:00.000Z/PT12H), than I get all the features:

My question is: Why are all the features not shown when I set that latest features should show. Why GS is ignoring features with latest timestamp? How can I get all the features without time filter?


texts:
I converted some GPS point data to LineStrings and I am sometimes getting these jagged jumps in coordinates. I am pretty sure my data is correct, this anomaly is just to perfect to be wrong...

Anyone seen this before and knows what is going on?

Detail:

Overview:


texts: Choosing graphics card for ArcGIS 3D Analyst?

texts:
in order to make the most of ESRI´s 3D analyst I would like to order an appropriate graphics card. I have been recommended the AMD FirePro V5900 2 GB by someone who doesn´t know much about GIS hence my question here.
Does anyone have any experience with 3D Graphics programs and decent hardware?
Would this card be suitable for 3D Analyst?
What is the most important factor to consider when buying graphics cards for 3D visualization?


texts: Automating Mastermap update in Oracle Spatial?

texts:
I've been set a challenge which is a bit out of my comfort zone. I'm an esri developer trying to break into the open source world.

I've been asked to automate the update of some oracle SDO datasets. I have a boundaries dataset (X) that needs to match topologically with OS Mastermap data. At the moment, the Mastermap update is a manual process (disk in the post!). If a boundary in the OS data is updated e.g a field is reduced as a result of a new road, the same field in dataset X needs to be updated accordingly. At first, I'm expected to flag up any changes that will need to be made to dataset X. Eventually, dataset X will be updated automatically. I also need to store the changes made for audit purposes. I've been challenged to achive the above with an open source stack.

My initial thoughts are that a workflow could look like – download update -> import dataset -> select all new and updated toids in OS dataset -> some kind of spatial query to identify effected features in dataset x-> export selected features to archive -> edit selected to snap to new OS boundaries.

Has anyone attempted something similar? What tooling would the community recommend?


texts:
I try to filter for catchwords by using the attribute table. I need to filter all attributes which starts with "industrial". The problem ist that the attributes are sometimes named as "industrialpark" or "industrialzone" or "industrialarea".

Is there a way to filter for just a special part of a catchword?

I could imagine that there must be a way using the "select by expression" tool and say something like "fieldname = 'indust...'" or "fieldname = 'industrial...'" , so that qgis prints out all catchwords which starts with the searched part of the word.


texts: Isolating utility pole points from small scale LiDAR scene?

texts:
I am working with some LiDAR data, and have a problem I just can’t seem to solve myself...

The image below shows a point cloud representation of a utility pole. In this point cloud, there is the pole itself, some portions of the power lines, and also a streetlight attached to the pole. I have labeled those accordingly. The area in the red box is the pole itself. It is not common for poles to have some amount of lean to them, as this pole demonstrates.

I am looking to isolate the points from just the pole (i.e. inside the red box). I have tried using both intensity and return index, but just don't get satisfactory results. To me, this seems like a nn (nearest neighbor) problem. If I could just filter out points that don't have another point X distance underneath them - problem solved.

I have a fairly simple approach in mind:

Divide LiDAR Extent into regularly spaced 3d grid. Each grid cell would span the full Z range.
For each cell - compute z distance between all points.  Remove all points that are > X distance from other point. Continue doing this until no points removed - at which point we move to next cell.

The problem I don't account for here is the pole lean. It is likely the lean will cause a pole to span multiple cells. Points near the boundaries of these cells will not appear to have points beneath them (in that same cell) - and thus I end up removing pole points...

I've been staring at this problem for 2 days now, and my ideas just keep getting more and more cockameny. I'm hoping some fresh eyes and different perspectives may have simpler approaches.

For what it is worth, these point clouds will be very small - never more than a couple hundred of points. So brute force is fine by me. I am familiar with PDAL and Laspy - but even rough pseudocode approach would be very welcome.


texts:
From How to extract elevation from .HGT file? I know how to get the elevation from SRTM1 and SRTM3 files. But how to get elevation from SRTM30 Plus ftp://topex.ucsd.edu/pub/srtm30_plus/srtm30/data/ ?
I tried to use the following code, but it returns the incorrect values:

def read_elevation_from_file_srtm30_plus(hgt_file, lon, lat):
with open('e060n40.Bathymetry.srtm', 'rb') as hgt_data:
    # HGT is 16bit signed integer(i2) - big endian(>)
    elevations = np.fromfile(hgt_data, np.dtype('>i2'), 6000*4800)\
                            .reshape((6000, 4800))

    lat_y = int(round((lat - int(lat)) * 6000, 0))
    lon_x = int(round((lon - int(lon)) * 4800, 0))

    return elevations[6000 - lat_y, lon_x].astype(int)

What am I doing wrong?


texts:
Formatted this CSV data and tried to add it to an existing QGIS project and for the life of me can't get the data to align. Project reference is Illinois East NAD 83. Have tried that for the CSV layer as well as EPSG:4326.

Data sets here: https://www.dropbox.com/sh/f6f53jvii6zgaky/AABLZlf0Fnclq3BQl9DkiGWna?dl=0


texts: Is there an open source library that offers geodesic buffering of large polygons?

texts: I need to buffer some large polygons and minimize distortion. Is there an open source library that offers negative and positive geodesic buffering of large polygons? I am aware ESRI has an implementation but it is closed source.


texts: Using EasyCostumLabeling tool with expression based label?

texts:
I'm using an expression based label (my label display 3 fields from my table on aqueduct line) but the resulting label create overlay. Is there a way to prohibition overlapping label from the same layer, but show ALL label?
Or is there a way to use EasyCustomLabeling tool with expression based label ?

My code look like :

"EG_DIAM_SECTION_TBD_SHORT_VALUE"  || 'ø - E.P.  '  ||    "EG_MATERIEL_SECTION_TBD_SHORT_VALUE"  || '\n' ||  "EG_SECTION_LONGUEUR"  ||  'm. @ '||  "EG_SECTION_POURCENTAGE_PENTE"  || '% - ' ||  "EG_SECTION_ANNEE_INSTALLATION"


texts: Batch export all layers in ArcMap document that have selected features?

texts:
I have several layers in an ArcMap Document. All of these layers are stored as shape files.

Specific features of the said layers are selected. I want to export the data as shapefile programatically. How can I do this using arcpy?


texts: Performance/Speed with different types of FeatureCollection filters?

texts:
Does anyone know if there is a speed/performance advantage to any of the following three approaches to filtering a large FeatureCollection?

OPTION 1:

myFilteredFeatureCollection = myLargeFeatureCollection
    .filterMetadata('START_N', "less_than", beginDateVar)
    .filterMetadata('END_N', "greater_than", endDateVar);

OPTION 2:

var predicates = ee.Filter.and(
  ee.Filter.lt('START_N', beginDateVar),
  ee.Filter.gt('END_N', endDateVar));
myFilteredFeatureCollection = myLargeFeatureCollection.filter(predicates);

OPTION 3:

var predicates = ee.Filter.and(
  ee.Filter.lessThan('START_N', beginDateVar),
  ee.Filter.greaterThan('END_N', endDateVar));
myFilteredFeatureCollection = myLargeFeatureCollection.filter(predicates);


texts: Downloading OpenStreetMap data at European scale?

texts:
I don't find a way to download specific data for Europe (E roads, i.e. roads from motorways/primary/trunk having the key "int_ref" not empty).

I tried several ways (using QGIS x64 on Windows 10), with the plugin QuickOSM:

with the script:

[out:json];
{{geocodeArea:Italy}}->.searchArea;
(
way['highway'='trunk'] ['int_ref'~'.'] (area.searchArea);
way['highway'='motorway'] ['int_ref'~'.'] (area.searchArea);
way['highway'='primary'] ['int_ref'~'.'] (area.searchArea);
);
out body;>;
out skel qt;

I can't do it for all the European countries (and moreover, often it blocks, saying "bad request", so I have to split into 3 scripts, one for each highway value).

opening the .osm.pbf file downloaded from Geofabrik, it work for a small file (tried with Luxembourg), but it take a lot of time (5 minutes), so I guess that for a 30 GB file, there is little hope...

Any hints on how to proceed?


texts: gdalwarp can't find libgrass_dgl or libgrass_vector files

texts:
Using gdalwarp I get the following errors:

  ERROR 1: libgrass_dgl.7.2.0.so: cannot open shared object file: No such file or directory
  ERROR 1: libgrass_vector.7.2.0.so: cannot open shared object file: No such file

The full command line I'm using is: gdalwarp f1.tif f2.tif -te -3.0 4.7 1.2 11.1

I'm using a Kubuntu 17.04 platform.
If I look in folder /usr/lib/grass72/lib these files are present. But gdalwarp can't find them. How do I correct this?


texts:
I am trying to use plot_usmap to map a data frame with two columns - fips and change2 - where the values for change2 range from -34.7 to 91.5. When I run this code:

plot_usmap(regions = "counties", values = "change2", data = map, color = "white") +
  scale_fill_gradient2(
    low = "red", mid = "grey", high = "dark green", name = "Change from SPLY", label = scales::percent) +
  theme(legend.position = "right")

The legend gives the values as -3000 to +9000.

Why doesn't the legend match the values of change2? Is there a fix to this that I'm missing?


texts: Interpolation error using QGIS 2.8.3 with nearest interpolation function

texts:
With several point imported from .csv file containing the information of longitude, latitude and VALUE.

To represent the 2-d distribution of value, I want to interpolate the sparse point into 2-d contour figure.

Here is the distribution of the points, the background is China which I want to use as the mask shapefile.

Using QGIS Kriging or nearest algorithm dealing with the point data

They all come with same error: The data source is invalid.

How to fix it. I uploaded my csv data and shapefile here:

(1) CSV containg lon,lat, value

(2) Shapefile of China


texts: Applying one color to many layers in QGIS 3?

texts:
I've been playing around with my QGIS Bonn (3.2.3) and I have a question regarding styles/symbology/colours.

Is it possible in QGIS to set one colour/style to multiple layers at once?


texts: Qgis Server Error while rendering postgis osm layer planet_osm_polygon

texts:
I have hosted QGIS server on a remote server. I want to publish OSM layer which is stored in PostgreSQL database. I have prepared one project file with .qgis extension and copied to the desired directory on hosted server. Now while accessing this layer as WMS service QGIS server is not able to generate the tiles for planet_osm_polygon layer and giving error in qgis-000.log as:

  Layer 'planet_osm_polygon20180803120119097' configFile:
  /usr/lib/cgi-bin/wms_prj/gpolygon.qgs not found in layer cache'

Note: I am able to generate the map tile if publising other OSM layer like planet_osm_point, planet_osm_line, planet_osm_roads.

Can anyone please guide me how to resolve it. I have made the same configration in my project file (.qgis) as I am doing it for other layers but only I am getting error for planet_osm_polygon layer.


texts: How to create polygon automatically with coordinates

texts: I have some bunch of points in an area and I want to draw connecting lines with the help of points to create a closed polygon. But all these needs manual efforts. Can any one explain to me how to create polygon automatically in QGIS with the help of outermost points.


texts: Display of zoomable image thumbnails instead of symbols on maps

texts:
sincere congratulations on your ingenous QGIS and the many extremely useful extensions. In my research on African archaeology I work a lot with geotagged pictures (drawings; photos; plates etc.) which are elements of self-defined GIS data records and which I display on my own maps, satellite imagery and the like using photo2shape. Now, what I am really dreaming of is to be able to look at preview thumbnails (that can be opened in larger size by a click) instead of symbols on my maps, such as it is possible, for instance, in ArcGIS or GoogleEarth.

Is this already possible in QGIS? Or would someone develop it perhaps? This feature would be extremely useful for all sorts of research and other applications.


texts: Inserting values with arcpy.da.InsertCursor()?

texts:
Building on Standalone Python script to split a polyline with a point layer, I don't think I am manipulating geometry objects correctly.  I am splitting an original line, preserving all the attribution, and inserting new lines with the original attribution but updated geometry.

import arcpy

mxd = arcpy.mapping.MapDocument("CURRENT")
layers = ["Primary OH Conductor", "Transformer Bank"]
line_lyr = layers[0]
point_lyr = layers[1]
points = [row[0] for row in arcpy.da.SearchCursor(point_lyr, "SHAPE@")]

line_fields = []
line_fields.append("SHAPE@")

try: arcpy.Delete_management("lines")
except: pass
arcpy.MakeFeatureLayer_management(line_lyr,"lines")
for rec in "lines":
    with arcpy.da.SearchCursor(in_table = "lines",
                               field_names = line_fields) as o_search_cur:
        feature_values = []
        for row in o_search_cur:
            feature_values.append(row)
            line = row[-1]
    for point in points:
        if line.distanceTo(point) < 0.001:
            snappoint = line.snapToLine(point).firstPoint
            if not (snappoint.equals(line.lastPoint) and snappoint.equals(line.firstPoint)):
                cutline1, cutline2 = line.cut(arcpy.Polyline(arcpy.Array([arcpy.Point(snappoint.X+10.0, snappoint.Y+10.0),
                                                                          arcpy.Point(snappoint.X-10.0, snappoint.Y-10.0)]),
                                                             arcpy.SpatialReference(2276)))
                a = feature_values[0][:-1] + tuple(cutline1)
                row_values = [a]
                with arcpy.da.InsertCursor(in_table = line_lyr,
                                           field_names = line_fields) as line_in_cur:
                    for row in row_values:
                        print("inserting into column " + str(line_fields))
                        print("values: " + str(row))
                        print("original feature values (for reference): " + str(feature_values[0]))
                        line_in_cur.insertRow(row)

And my results:

inserting into column ['SHAPE@']
values: (<Array [<Point (2336683.24, 7228439.51939, #, #)>, <Point (2336751.76429, 7228385.5766, #, #)>]>,)
original feature values (for reference): (<Polyline object at 0x5b409c90[0x8a000660L]>,)
Runtime error
Traceback (most recent call last):
  File "<string>", line 42, in <module>
TypeError: cannot read geometry sequence. expected list of floats
>>>


texts: How to find the nearest point by using PostGIS function?

texts:
table A:

lat | long | the_geom | code | sign
13.8433095 | 100.6360357 | 0101000020E61.... | ABC | start_point
13.7544738 | 100.5459646 | 0101000020E6..... | ABC | end_point
13.4124215 | 100.6232332 | 0101000020E61.... | DEF | start_point
13.2423438 | 100.2324426 | 0101000020E6..... | DEF | end_point

table B:

lat | long | the_geom | code
13.7546285 | 100.5458729 | 0101000020E.... | ABC
13.7546698 | 100.5458513 | 0101000020E.... | ABC
13.7547107 | 100.5458233 | 0101000020E.... | DEF
...

I would like to find the nearest point(s) of each points (start and end point of table A) compare with every points with the same 'code' in table B ?

What's the best PostGIS function/PostgreSQL query to solve this ?
What I have tried are ..

ST_Distance_Spheroid and ST_DWithin and ST_Distance

But I have no idea how to combine them within a single query.


texts: Filtering a Topojson on server-side. Which part of topology I need?

texts:
I have a big Topojson which I want to avoid loading it at page load. Instead I would like to parse on server-side at user request, sending back chunks of this original Topojson to be ploted.

So I want to filter by a list of attributes in topojson. Initially I thought that just sending its corresponding geometries with its arcs would be OK (sorry, i am used to GeoJSON) but then I see below the Topology section, which I am afraid is quite dependent on the previously selected arcs (in "Geometry").

Topology section is quite a very long part of this TopoJSON. The question is, as I cannot find a link between the arcs in "Geometry" and "Topology", do I have to load all Topology part, even if I only select a few of "Geometries"? it would not solve my problem at all… any alternative?

I want to filter by COD_POSTAL

  "geometries": [
    {
      "type": "Polygon",
      "properties": {
        "COD_POSTAL": "08630"
      },
      "arcs": [
        [
          0,
          1,
          2
   ]
]
},
    {
      "type": "Polygon",
      "properties": {
        "COD_POSTAL": "08256"
      },
      "arcs": [
        [
          8,
          9,
          10,
          11,
          12,
          13
        ]
      ]
    }


texts: Messed up frame annotation in map composer

texts:
I am trying to create a nice map figure with lat/lon frame annotation with the map composer. I tried following several tutorials on this and technically it works, however my problem is, that the frame annotation in both lat/lon is completely messed up: it is either centered on zero lat/lon or completely random mixed N+S for lat and E+W values for lon axis. In any way it does not match up with my actual map extent whatsoever and I have no idea why.
In the item properties tab the image extents are the right values of my map extent, I did not change the CRS in the "Grids" tab ( it is pseudo mercator) I chose Map unit for interval units and experimented with the "interval" and "offsets", to no success. Settings for grid frame are zebra, left/right divisions: latitude only, top/bottom divisions: longitude only. in the draw coordinates tab: format: degree, minute, second, with suffix and I only enabled left and bottom annotations.

I am running Qgis 2.12.1 on kubuntu 14.04.


texts: How do I process xy data to be compatible with census block data (so I can include them in my GWR model)?

texts:
I'm planning to do regression analysis (GWR) on several datasets, where all but one is based on census areas (UK OS lower level super output areas, lsoa) for which I have a shapefile that describes the polygons. Regression with the census-based variables is easy to set up and execute because both the dependent and explanatory variables are held within the same Input feature.

But I'm finding difficulty with my final explanatory variable - concentrations at x,y coordinates (UK OS co-ordinates). I believe I need to 'join' or 'relate' with my lsoa shapefile (which is what I've done with the other variables); I've tried all sorts of things such as kriging and extracting to vector, 'relates', 'joins', adding new fields, calculating new fields, but to no avail.

I'm sure there is someone out there who has done this before - I'd really appreciate your advice.


texts: Moving Multiple Points Into Single Polygon using ArcGIS Desktop?

texts:
I have several days worth of geotagged imagery to which I've attached a tree canopy mortality rating.
The goal is to eventually join all of these classifications to a shapefile that contains watershed catchments for a whole region.
The issue is that images were taken at an oblique angle, so the image location as reported by the camera may be different than the catchment with classified stands.
Determining which photos were in which catchment was done manually by using snapshot views in google earth. This was added to an excel file, from which XY coordinates were displayed. Each image has between 1-5 rows in the spreadsheet, each of which has following columns: Image Name, Catchment #, Easting, Northing, Elevation, and Mortality Classification. The same image may encompass up to 5 catchments, leading to the issue.
So far I've tried to select features by attributes to determine which photo points should be in which catchment.
Is there a way to move multiple points (different pictures) into a single polygon (the catchments) all at once? It's possible to manually move each photo individually into its respective catchment, but there are thousands of photos.

A catchment and 5 associated pictures. Both the catchment and the points have the same ID number, this ID number was used to determine which points were associated with which catchments. 3 Pictures fall outside of the catchment.

Where clause used to find a catchment of interest, as determined by that days spreadsheet.
Where clause used to determine which points should fall within the catchment of interest. The catchment ID for each attribute table (catchment polygons vs photo points) is different.


texts: Automatically create path between two points via roads in Google Earth (Pro)?

texts: Does Google Earth have the ability to automatically create a path between two points connected by roads?


texts: Converting SRTM to UTM makes pixels non-square

texts:
I have a SRTM originally set to EPSG 4326. I need to convert it to UTM (EPSG 31983), in order to generate viewsheds.

The original version (EPSG 4326) has square pixels (0.000833333 x 0.000833333) but, when I convert to UTM, the pixel size becomes 91.9372 x -91.8611.

This means that, upon conversion, pixels are becoming rectangular. This makes it impossible to use viewshed in the processing toolbox.

Has anyone encountered this issue? What may be causing it? Any suggestions on how to solve it?


texts:
QGIS version 3.2.1

Windows 10

Project File on GitHub

Shape files on Github

I am trying to create shapefile from coffee cup image.

The actual coffee should be the shapefile.

However I would like the background of shapefile to be original coffeecup.

Here is original coffee image

Here is how I'd wish shapefile to look like (the green polygon has shape properties). Remaining shapefile should have coffee cup and saucer as original image instead of additional polygon.

Both coffee image and created shapefile are the layers. When I unselect coffee image I get this

Perhaps there is a way to merge the image into the shapefile?


texts: using arcpy.UploadServiceDefinition_server to publish to portal

texts:
I'm using arcpy.UploadServiceDefinition_server() to publish SD files to a server that is federated with portal.  Using this ESRI Documentaion I am able to publish the service to a specific folder on my server with the arguments  in_folder_type="XXX", in_folder="xxx" Since my server is federated with portal those services are automatically shared with Portal which is what I want.  However I would like like the services to appear in a specific folder as I did with server.  I do not see a control to do this.  Is this possible?

I attached an image of an example structure a service that is published to the server needs to be told to reside in a certain folder in portal. Not the root directory of portal.


texts: How to make all the labels appear neatly in QGIS print composer?

texts: I have been trying to label a Virginia counties and cities map with their proper names, but when I try to do so, using all the different types of placements, I either get nice looking labels with not all of them showing, or really messy labels with all of them showing. Is there a solution to this problem?


texts: Bringing layout into another drawing using ArcMap?

texts:
I am using ArcMap 10.1.  I set up a layout view (a titleblock) and would like to bring it into another mxd file.

How do you do this?

If you are familiar with Autocad this would be the equivalent of importing a titleblock.


texts: ArcGIS 10 SP2 issue with AggregatePoints

texts:
Esri persuaded us to upgrade to SP2 and now we cannot use the code we have, as it doesn't work.

According to the download material, NIM063853 is when AggregatePoints doesn't create the derived table.

This worked in SP1, now doesn't in SP2.  Has anyone else experienced this?  Any workarounds?


texts: ESRI Shapefile --> Halliburton GeoGraphix .cdf via Data interoperability

texts:
I'm building a landgrid Public Land Survey System (PLSS) file for the state of Oklahoma using ArcMap and the data interoperability extension for use in LMKR's GeoGraphix. I've downloaded section grid for the state in the form of a shapefile, ran a quick export, generating the file format below.

I was able to import the file into GeoGraphix, but the relevant information needed to populate the GeoGraphix database is not imported with the file, and the features are drawn as quarries or mines and not sections.

After doing some research on the GeoGraphix file format, I came across this table which explains how the software reads the import file (attached below).

I then when into ArcMap and ran another quick export of a township, imported the .cdf into GeoGraphix and manually populated the database information for the imported sections. After getting everything populated and identified, I exported the file out of GeoGraphix as a .cdf resulting in the below format.

Where [B3001010] Tells GeoGraphix that it is a PLSS section feature [6] meridian [30.0]township [S]Township Direction [15.0] range [W] range direction [33] section number. followed by the coordinates for this feature.

Being that the shapefile contains all the information above as fields, its it possible to build a date interoperability model to convert the information stored in their respective fields into a string in the correct format so GeoGraphix can recognize the and import the file correctly.

I realize most people might not be familiar with GeoGraphix or the .cdf format, but I hope that people familiar with the data interoperability extension could lend enough advice to help me out.


texts:
I'm running QGIS 3.4.2 on Windows 7 64bit.
I need to import periodically a specific table from several different MS-Access 2010 mdb's.
The database and path changes quite often, so an odbc solution is  not the right solution for me.
I tried the solution "Opening Esri Personal Geodatabase (*.mdb) using QGIS?",
downloaded and installed the 64 bit ODBC driver with command line ( / passive),changed the qgis.bat file and added variable names in QGIS settings.  When I try to drag and drop an mdb-file onto the layer window, I am still getting the message "invalid data source".


texts: Syncing to a feature service fails after join and field calculator

texts:
I have access to a feature service, and I wish to make edits to it and then sync those edits.

Things I can do successfully:

Accessing the service:

Add the feature service to my local ArcMap instance
Create a local copy for editing

Editing:

I can edit a single value by hand
I can edit mutiple values by hand
I can join another layer and use field calculator to update values

Syncing:

I can sync a single edit done by hand
I can sync multiple edits done by hand
I CANNOT sync edits that were made using field calculator after a join

I get this error when trying to sync after edits were made using field calculator following a join:


texts: How to add Google Maps coordinates into ArcMap?

texts:
I have set the standard "World Imagery" as a basemap in ArcMap. Now I want to search for places from which I have the coordinates through Google Maps. Let's say I would like to display 51°38'27.1"N 6°33'53.9"E (51.639771, 6.566256) in ArcMap on my "World Imagery" basemap. I have been using the XY-Tool in ArcGIS and also tried to create a new point-shape class. However, whatever I do, it doesn't work. I checked the coordinate systems (GCS_WGS_1984), but it doesn't work.

How can I do this?


texts: Creating MBO files with MapBasic and Notepad++

texts:
[Note: I have been using Notepad++ (with some customisations as per this page), rather than MapBasic itself. This may make things a bit trickier, but so far the time saved in development has been worth the extra hassle during setup.]

I am trying to start using project files with MapBasic. As I understand it the process is as follows: compile your modules from .mb files into .mbo files. Then process the .mbp file which links the .mbo files together and, if there are no errors, creates a .mbx.

This makes sense to me except for one thing...how do i get my .mb files to compile to .mbo instead of .mbx? What am I missing?!


texts: How to properly draw and label overlapping polygons?

texts:
I would like to show overlapping polygons on a single layer.

The polygons are coloured based on a property and labelled according to another property. The problem is then twofold:

The polygons are rendered in random order in tiles where they overlap, resulting in some small ones to "disappear" under the larger ones (e.g. cities "under" regions), while I would like them to be rendered according to a given property (I was thinking about using Priority but it is available only for TextSymbolizer in geoserver).
The labels of larger polygons are shown "on top" of smaller polygons, while I would like to have them drawn ( repeated ) only at areas where the two (or more) polygons do not overlap. E.g. State names only to be shown where no other regional surface is being drawn (even if this might result in State names not being shown at all).

Is there a proper way of doing it?


texts:
I have a single band 32-bit integer Arc Binary Grid.  Using QGIS 2.12.3, when I try to add this Grid, it is listed in the table of contents, but nothing shows in the map window.
To add the Grid, I use the "Add Raster Layer" button and select the w001001.adf file from within the Grid's folder.  I've successively used this method in the past with previous QGIS versions on other Grids, and have confirmed the process at https://pvanb.wordpress.com/2010/03/15/esri-rasters-in-grass-or-qgis/
QGIS reports that the min/max values are both -1.14848e+06.  However, when I view the Grid in Arcmap 10.0, min and max are correctly reported as 1,000 and 1,000,000, respectively.
As a test, I converted the Grid to a Tiff using Arc.  QGIS reads this Tiff correctly.  Although I can use the Tiff for my project, I would like to understand why I cannot add the .adf into QGIS.


texts:
I have been using zonal_stats from the rasterstats module for quite some time in Python. And I have never come across this error, which I tried to Google to find some solution with no succes.
CPLE_AppDefinedError: PROJ: proj_create_from_database: Cannot find proj.db

It turns out when I tried to import the module it says:
In [1]: from rasterstats import zonal_stats
PROJ: proj_create_from_database: Cannot find proj.db

I am using Anaconda 3/Spyder with Python 3.6.10. My OS is Windows 10. I have tried to update the packages rasterstats and rasterio with no luck.
Can someone help me 'how to find proj.db'? Or have some idea of what I can do?
Thanks in advance!