texts: Preprocess a thematic raster to produce multi-band raster of proportions?

texts:
I've been trying to figure out a way to pre-process a single-band thematic raster dataset (e.g., land cover) to produce a multi-band raster of proportional cover.

For example, consider a land cover dataset with two land cover types: water and land.

What I would like is to calculate for each cell the proportional coverage by a given land cover type within some distance from the center of that cell.

The final product would be a 2-band raster where band 1 is proportional cover by water, and band 2 is proportional cover by land.

You can conceptualize this as a raster representing the output of a moving window exercise.

The window is a buffer of a given radius, and the centroid of the buffer is the center of the target cell.

The target cell is assigned the proportional cover by the land cover type, then the window moves to the next cell and repeats the process.

The full process would calculate a raster band for each land cover type in the dataset.

I've yet to figure out an efficient workflow for this; my solutions have all involved iteration (for loops) in R, which takes an absolute eternity for large raster datasets.

Is there a more efficient approach, preferably using ArcGIS or R?


texts: Merge connecting lines which have the same attributes?

texts:
In QGIS, there is a road map shape file, and I would like merge connecting lines which have the same attributums ("name" and "fclass").

I tried the "Dissolve" processing tool, but it doest't care about the geometry (it merges those lines which are not connecting, too). Can I split those which are not connencting?


texts: GeoServer GetCoverage with a specific resolution

texts:
I have a GeoServer serving GeoTIFFs as WCS. I would like to get a slice (bbox) of a coverage at a specific resolution.

MapServer seems to support this:

WIDTH=output_width: Width in pixels of map picture. One of WIDTH/HEIGHT or RESX/Y is required.

HEIGHT=output_height: Height in pixels of map picture. One of WIDTH/HEIGHT or RESX/Y is required.

(source: https://mapserver.org/de/ogc/wcs_server.html)

If I add this to my WCS request, it is ignored:

http://192.168.37.10:8080/geoserver/wcs?request=GetCoverage&coverageId=ml__Duesseldorf&height=512&width=512&service=WCS&version=2.0.1&subset=E(355017.6407251246,355761.0022115205)&subset=N(5679871.219933263,5680544.231752177)

This returns the TIFF in the same spatial resolution as in the source data.

adding scalefactor=0.25 works:

http://192.168.37.10:8080/geoserver/wcs?request=GetCoverage&coverageId=ml__Duesseldorf&scalefactor=0.25&service=WCS&version=2.0.1&subset=E(355017.6407251246,355761.0022115205)&subset=N(5679871.219933263,5680544.231752177)

This works but I can only scale the whole image, while I would like to set a specific resolution.
for example: height=51 width=512

Also I wonder if there is any documentation for GeoServer WCS that explains all the parameters that can be used?

(bbox for example also doesn't work as well)


texts: Is there a way to filter vector layer with JSON source in OL3?

texts:
Is it possible to filter the ol.layer.Vector with a JSON source and display only specific features?

My vector layer looks like this:

var vectorLayer = new ol.layer.Vector({
      source: new ol.source.Vector({
        format: new ol.format.GeoJSON(),
        url: 'http://vagrant-test.local.com/json/countries.js'
      }),
      style: function(feature) {
        return style;
      }
  });

I check the ol3 documentation and didn't find some filter functionality.


texts: How to call arcpy.RefreshActiveView() every 5 seconds using ArcGIS Python?

texts:
I have an Excel macro that generates random values after every 5 secs and I am reading it in ArcGIS for Desktop... to refresh the view I have to press f5 and the view gets updated. Now I want to automate it and run the loop infinitely. I found the command for refreshing the view is

arcpy.RefreshActiveView();

I want to call this function after a delay of 5 seconds.. I tried this code

import arcpy;
import time;

while True:
    arcpy.RefreshActiveView();
    time.sleep(5);

But the problem is when I run this the main thread is blocked and the program is not responding. So I tried to move it to a new thread.

import arcpy;
import time;
from threading import Thread;

def refView() :
    while True:
        arcpy.RefreshActiveView();
        time.sleep(5);

try:
    Thread(target=refView, args=()).start();
except Exception, errtxt:
    print errtxt;

This is also not working and the application becomes not responding, please let me know what I am doing wrong ...


texts: I recently completed a project that involved slope, aspect, and solar radiation to find optimal rooftop locations for solar panel installations. For my final map I combined the rasters that showed the most promising rooftops based on my ranking. I get the result I want but I am also left with little single cell or a cluster of cells that I don't want. Is there a method of removing cells that don't have a neighbor or cells that are in small clusters? my goal is to clean up the rooftops and be left with the larger raster areas.


texts: Proxy pass for GeoServer not working properly

texts:
I have GeoServer 2.18.2 hosted on tomcat 9.
http://xx.xx.xx.xx:8080/geoserver/web
With the IP it is working perfectly fine.
To use it within OpenLayers I was needed it behind SSL. So through proxypass I have hosted it in apache2 to give it a domain name.
Below are the configurations
<VirtualHost *:80>
    ServerName geoserver.example.com
    ProxyRequests Off
    ProxyPreserveHost On
    <Proxy /geoserver>
        Order deny,allow
        Allow from all
    </Proxy>
    ProxyPass /geoserver http://127.0.0.1:8080/geoserver
    ProxyPassReverse /geoserver http://127.0.0.1:8080/geoserver
RewriteEngine on
RewriteCond %{SERVER_NAME} =geoserver.example.com
RewriteRule ^ https://%{SERVER_NAME}%{REQUEST_URI} [END,NE,R=permanent]
</VirtualHost>

https://geoserver.example.com/geoserver/web/
Now I can access the layers in OpenLayers and login into admin successfully but when I click on navigation and use search and save other settings it gives me 400 Bad Request Error.
Apparently, it seems like there is an issue in proxy pass configurations because when I use it along with port 8080 it works prefectly. But I needed this way along with SSL https://geoserver.example.com/geoserver/web/


texts: OSM getting ways and nodes of a relation instead of an area

texts:
I am new with OSM queries, I understand that areas are not real entities but they work very well, let's suppose I call the following query:

[out:json][timeout:25];
area["name"="Tenerife"]->.a;
(
 way(area.a)[amenity=parking];
 node(area.a)[amenity=hospital];
 rel(area.a)[admin_level=6]["is_in:country"="Spain"];
);
out geom;

So I thought that if I have a relation number instead of using an area would be faster and more specific, but how can I do that? because I have look for it many hours and I cannot get it work, I understand it must be something like:

  [out:json][timeout:25];
  rel(2108882)->.a;
   (
   way(XXXXX)[amenity=parking];
   node(XXXXX)[amenity=hospital];
  );
  out geom;

Any clue on what can make it work? How I set the filtering to use that relation?


texts: How to clip/symmetrical difference a shapefile and ArcGIS Mapserver layer?

texts:
I have 2 layers. A portion of a shapefile and an ArcGIS Mapserver layer. I would like to subtract the satellite image area outside the shapefile (the rectangle) from the view. IOW, leave only the area inside the shapefile with the satellite view. screenshot is attached. My problem is when I try Vector --> Geoprocessing Tools --> the satellite layer is not seen. I have tried exporting the sat layer to a tiff. No luck.
I have tried everything I can think of and Googled this to death. Can anyone offer suggestions?


texts: Python code executes but nothing happens

texts:
I'm trying to write a script that looks into a database and lists all the datasets that begin with zone. The script then enters into a for loop and copies a map into each on the datasets listed. The code returned exit code 0 and it even prints the zone in the interactive window, but the map isn't copied into the respective datasets. This is the code I have so far:

import arcpy
import os
PATH = r"C:\Users\Daimon Nurse\Desktop\DFMPROJECT\DFM.gdb"
PATH2 = r"C:\Users\Daimon Nurse\Desktop\DFMPROJECT"
arcpy.env.workspace = PATH
arcpy.env.overwriteOutput=True

zones = arcpy.ListDatasets("Zone*")
mxd = arcpy.mapping.MapDocument(r"C:\Users\Daimon Nurse\Desktop\DFMPROJECT\Complete Final Capstone Map\Complete Z2 Campus Map_Landscape1.mxd")
for zone in zones:
    print zone
    mxd.saveACopy(PATH2 + zone +  ".mxd")


texts: How to link SQL file with relative path in Processing Scripts

texts:
I wrote a processing tool, which exports some tables from PostGIS to package, and as part of the export, I am calling an external SQL file to add some triggers, etc. The SQL file is located as my script (C:\Users\XXXX\AppData\Roaming\QGIS\QGIS3\profiles\default\processing\scripts)
Currently, I linked the SQL file as follows.
path_to_sql = str(Path.home()) + "\\AppData\\Roaming\\QGIS\\QGIS3\\profiles\\default\\processing\\scripts\\bwk.sql"
sql_file = open(path_to_sql)

But I was wondering, why it does not work when I wrote  path_to_sql = "\\bwk.sql" or path_to_sql ="bwk.sql" (No such file or directory: 'bwk.sql')


texts:
Is there a tool / library / plugin / API that recognizes geographic place names in their various forms (for example: "USA", "United States", "United States of America" would all map to the same object)?

I seem to run into this problem all the time. I have a file containing data about some geographic places (e.g. GDP by country) and I want to join that data to the corresponding map features (e.g. a shapefile of world countries). But the names of the places do not match up (one file says "United Kingdom", the other says "UK").

The only solution I've found so far is to go through and change the names manually, but it seems like there has to be a better solution. Any ideas?

I should also mention, I expect to be doing this for a large number of places, so a paid geocoding service like Google's would end up being very expensive. I'm hoping to find a free alternative.


texts: Writing band names into header (ENVI file type) using GDAL?

texts:
I am writing numpy arrays out as ENVI images using gdal in python. Until today, I always accessed the header the hard way by reading in all lines, changing some information and writing it out again. This is not only coding intensive, but also yields multiple error sources.

A better solution would be by using SetMetadata. I tried:

myArray = np.full(shape=(4,3), fill_value=10)

driver = gdal.GetDriverByName('ENVI')
destination = driver.Create("test.bsq", ncols, nrows, 1, gdal.GDT_Float32)
destination.GetRasterBand(1).WriteArray(myArray)
destination.SetMetadataItem('band names', 'band hello world', 'ENVI')

but in the output file, the band is still called "band 1". Is my information overwritten by the built-in function for header creation?

Or maybe I have to set the name for each band individually like

b = destination.GetRasterBand(1)
b.SetBandName("band hello world")
b.WriteArray(myArray)

But "SetBandName" does not exist and I cannot find how else it is supposed to be called.

How do I supply band names for ENVI-images then?


texts: Adding symbology template to geodatabase using Arcpy

texts:
I have a geodatabase with ~1000 rasters that I want to render with the same symbology.  I have a script that (1) imports the rasters into an mxd and (2) applies a template symbology to the raster layers within the mxd.

However, I want to be able to send the geodatabase to other people so that when they add the raster data from the geodatase into an mxd, they will have the same symbology. From what I have read, layer files are associated with mxd or outside geodatabases.

Is this possible, and if so, how can I do it?


texts:
Currently trying to use the following code to add and remove Legends as a result of an overlay being added to the map. Getting no error messages in browser but still doesn't seem to be working.

        map.on('overlayadd', function (eventLayer) {
    if (eventLayer.name === 'Bombing') {
        BombingL.addTo(map);
        map.removeControl(ArmedAssaultL);
    } else {
        map.removeControl(ArmedAssaultL);
    }
    });

in this case I'm trying to say that when a layer is added to the map called 'Bombing', also add the bombing layer legend (BombingL) and remove the Armed Assault Layer legend (Armed AssaultL)

If there is nothing wrong with that code i'll try and get a jsfiddle working to evidence it better.


texts: Mapping underground Utilities, like pipes, gas lines and telephone lines?

texts:
I just got hired on to a small local government and I am sill learning the ropes as a GIS Analyst.

What is an inexpensive way to map underground Utilities, like pipes, gas lines and telephone lines?

What equipment is used to do this?


texts: Using alias name instead of field name in GetFeatureInfo in GeoServer

texts:
I've created a Leaflet map and use GetFeatureInfo to send request to GeoServer. Everything is ok but I want to use Alias field name instead of field name. Indeed I want to replace real field with alias name when GeoServer is creating the return JSON file.

Is there any way to handle this issue?


texts:
I would like to use SQLite for storing some information in my plugin for QGIS.

I have a problem: if I put the db, let say test.db, among the resources as file, after compiled the resorce.qrc file, QGIS find the db, but if I insert some row in the db, QGIS does not find any more I have to compile again the file,  I tried to recompile the file every time I commit something in the db, it seems to work, but is a dirty solution and I do not think I can use it on every system different from the mine, is there a better solution?


texts: Forbidden error connecting to WMS with QGIS

texts:
FR :

  j'essaie de me connecter à un service WMS
  http://geod-prod-ws.cg94.loc/WMS_ASSAIN_DSEA/service.svc/get?SERVICE=WMS&REQUEST=GetCapabilities
  mais ça bloque

  Par contre ça fonctionne sous MapInfo ! Avez vous une idée d'ou cela
  peut provenir ?

EN :

I try to connect to a WMS service :

http://geod-prod-ws.cg94.loc/WMS_ASSAIN_DSEA/service.svc/get?SERVICE=WMS&REQUEST=GetCapabilities

but it blocks, but it works under Mapinfo!

Do you have any idea where this might come from?


texts: Hiding sublayer with toggle button using LayerActions?

texts:
I have a toggle map (like this example) with sublayers. I want all the sublayers to be visible when you open the map, but I would like to hide some sublayers when the toggle buttons are clicked. I've tried adding a hide sublayer code into the LayerActions section of my code, but it's not hiding the sublayers (the button does the same function no matter whether the 'layer.getSubLayer(1).hide()' is there or not). This is the code I'm using:

var LayerActions = {
 all: function(){
  layer.getSubLayer(1).hide();
  sublayers[0].setSQL("SELECT * FROM spreadsheet");
  return true;
 },
 2011: function(){
  layer.getSubLayer(1).hide();
  sublayers[0].setSQL("SELECT * FROM spreadsheet WHERE year IN ('2011') AND year IS NOT NULL");
  return true;
 },
 2012: function(){
  layer.getSubLayer(1).hide();
  sublayers[0].setSQL("SELECT * FROM spreadsheet WHERE year IN ('2012', '2012/2013') AND year IS NOT NULL");
  return true;
 }
}

 $('.button').click(function() {
 $('.button').removeClass('selected');
 $(this).addClass('selected');
 LayerActions[$(this).attr('id')]();
});

Could someone point me to where I'm messing up with my code? Or is it not possible to hide a sublayer using the LayerActions?

I'm now loading my layers differently, and this is my new code:

cartodb.createLayer(map, url1)
    .addTo(map, 0)
    .done(function (layer) {
        console.log('added url1 ',url1);
        cartodb.createLayer(map, url2)
        .addTo(map, 1)
        .done(function (layer) {
            console.log('added url2 ',url2);
        })
        .error(function (){
            console.log('problem adding 2nd layer');
        });
    })
    .error(function () {
        console.log('There is something wrong with requested data!');
    });

So if I wanted to hide the first layer that's loaded (url1) and select only certain points from the second layer (url2), how would I code that in my LayerActions section? I've tried the following, but nothing changes when I click the buttons:

var LayerActions = {
all: function(){
  layer.getLayer(0).hide;
  layer[1].setSQL("SELECT * FROM spreadsheet");
  return true;
},
2011: function(){
  layer.getLayer(0).hide;
  layer[1].setSQL("SELECT * FROM spreadsheet WHERE year IN ('2011') AND year IS NOT NULL");
  return true;
},
...


texts: How to simulate tsunamis and other natural disasters?

texts:
How do you go about simulating natural disasters like tsunamis? We're mainly looking for a place that might have more too.

Currently, we're trying out GrassGIS's Tsunami tool, but it seems to have been discontinued since grass 7, and we're unsure to it because of it

We're also looking into ClawPack as well because it was mentioned on ResearchGate. It looks promising as well, but complicated.

We're familiar with QGIS, PostGIS, and prefer to use those tools if possible, but are open to trying new things


texts: Printing NDWI and NDVI comparision graph

texts:
Trying to print the chart for both NDVI and NDWI

  function addNDVI(image) {
            var ndvi = image.normalizedDifference(['B4', 'B3']);
        return image.addBands(ndvi);
            }
    var filtered = L7.filterDate('2000-01-01', '2019-12-31')
     .filterBounds(mumbai);
    var with_ndvi = filtered.map(addNDVI);
    var rgb_vis = {min: 0, max: 0.1, bands: ['B3', 'B2', 'B1']};
    //Map.addLayer(filtered.median(), rgb_vis, 'RGB (median)');
    var greenest = with_ndvi.qualityMosaic('nd');
    //Map.addLayer(greenest, rgb_vis, 'RGB (greenest pixel)');
    //print(ui.Chart.image.series(with_ndvi.select('nd'), mumbai));

    /********** NDWI***************/

    var dataset = ee.ImageCollection('LANDSAT/LE07/C01/T1_8DAY_NDWI')
                          .filterDate('1999-01-01', '2019-12-31')
                          .filterBounds(mumbai);
        var colorized = dataset.select('NDWI');
        var colorizedVis = {
          min: 0.0,
          max: 1.0,
          palette: ['0000ff', '00ffff', 'ffff00', 'ff0000', 'ffffff'],
        };

        Map.setCenter(72.88, 19.154, 10);
       // Map.addLayer(colorized, colorizedVis, 'Colorized');
        print(dataset.select("NDWI"));
        print(with_ndvi.select("nd"));

       var y1 = ee.Array(with_ndvi.get('nd'));
       var y2 = ee.Array(colorized.get('NDWI'))
        //concat
        var yValues = ee.Array.cat([y1, y2], 1);

       // code to print chart of ndvi and ndwi
        var chart = ui.Chart.array.values(yValues, 0, mumbai)
        .setSeriesNames(['B5', 'B6'])
        .setOptions({
          title: 'LC8 TOA B4 vs. {B5,B6}',
          hAxis: {'title': 'B4'},
          vAxis: {'title': '{B5,B6}'},
          pointSize: 3,
    });

    // Print the chart.
    print(chart);

Error:

  Error generating chart: Array: Parameter 'values' is required.

  Array: Parameter 'values' is required.


texts: Real time changing status of a point, web mapping

texts:
I’m new to web mapping and developing and I try to make a web map with points. The points are bank ATMs. The points will change color depending on their condition. There will be three different colors and conditions. The condition and the color of a point will change in real time. I'll use the Leaflet For Web platform map software.

Do I need a database and if yes how to connect a point with a database. Is any plugin for changing the status in real time?

Edit:
Web platform = Web map software
The data will be in a database. The data that relate with the status of each point not the Lat Lon information.


texts: Can I simplify geometry on the fly in a ArcGIS JavaScript map viewer?

texts:
We've got a polygon feature class of parcels around a lake boundary that is incredibly accurate. The editor did a really thorough job of following all the inlets and meanders along the edge. From just a cartographic perspective, it's a great set of polygons.

Unfortunately, the result is an incredibly slow layer when we need to add it to a map service. Even at a fairly close zoom scale, it takes up to half a minute to load some of the more complex areas.

The obvious solution is to simplify the geometry. Unfortunately, the data owner doesn't really want to do that on his original layer (he's concerned about losing complexity when making static maps). There is department reluctance to make a new simplified copy because of recent pushes to reduce duplication of data.

Is there a way to simplify geometry on the fly when displaying highly complex polygons in an ArcGIS JavaScript API map viewer?


texts:
I am attempting to invoke HMS to run headlessly by prompting a batch file containing the following:
cd C:\Program Files\HEC\HEC-HMS\4.4
hec-hms -s C:\Users\abudhabikid\Desktop\TESTING\Project\BasinName_batch.script

The script itself contains the following:
OpenProject("Project", "C:\Users\abudhabikid\Desktop\TESTING\Project")
Compute("Run_Title")
Save
Exit(1)

These files are created in accordance to the section in the HMS manual regarding command line operation.
However, when my code invokes the script (or if I run the script myself through the command line, I receive the following error:

Note that this does not occur every single time that the script is run and I cannot tell that there is any pattern.
My requirements are that HMS run multiple times to satisfy a loop and I can't have users sitting at their desks clicking "OK".
I can see three directions here:

There's some additional code I should put in the .script that tells HMS to OpenProject and Compute and so-on
There's some code I can add to my VBA program to listen for this error message window and close it automatically
There's some edits to the .config file I need to make (possibly regarding the java.exe and javaw.exe paths); I have played around with this to no avail and I'd really rather not have to edit every user's .config file before they could use my program
(Bonus) Maybe I should limit users to a certain version of HMS?

I really hope somebody knows what do to here.


texts: How to use a polygon grid to calculate statistics on line and point features

texts: What seems a simple task but I can't find the right tool to do it. I have a polygon grid cells and say a point feature class. I want to be able to create a table that contains the grid cell ID and for example the sum of points within that cell or the average of the value for those points. For a line feature I would like to simply have a 0/1 value if the line falls into a grid cell. I am doing this in ArcMap 10. This is like Zonal statistics but on a vector layer instead of a raster.


texts: ArcGIS API for JavaScript - Featurelayer problem with CORS

texts:
I created simple simple web site:

        require([
      "esri/map",
      "esri/config",
      "esri/dijit/InfoWindowLite",
      "esri/InfoTemplate",
      "dojo/dom-construct",
      "dojo/domReady!"
    ], function (
        Map,
        esriConfig,
        InfoWindowLite,
        InfoTemplate,
        domConstruct
       ) {
        esriConfig.defaults.io.proxyUrl = "http://111.111.111.111:6969/agsproxy/proxy.ashx"
        esriConfig.defaults.io.alwaysUseProxy = true;
        var customExtentAndSR = new esri.geometry.Extent(-437000, -1244500, -409000, -1239500, new esri.SpatialReference({ wkid: 5514 }));
        map = new Map("mapDiv", { logo: false, extent: customExtentAndSR, nav: false, slider: false, fitExtent: true, fadeOnZoom: true, showAttribution: false });
        var infoWindow = new InfoWindowLite(null, domConstruct.create("div", null, null, map.root));

infoWindow.startup();
 map.setInfoWindow(infoWindow);
 var template = new InfoTemplate();
var featureLayer = new esri.layers.FeatureLayer("http://111.111.111.111:6969/ags/rest/services/VSLP2014/jprl_plochy/MapServer/0", {
 mode: esri.layers.FeatureLayer.MODE_SNAPSHOT,
 infoTemplate: template,
 outFields: ["*"]
 });
 map.addLayer(featureLayer);
    });

When I use "MODE_SNAPSHOT" for my featurelayer, layer is displayed and is possible to display infotemplate. But when I run my web site again, It seems to be problem with CORS (No 'Access-Control-Allow-Origin' header is present on the requested resource. Origin 'http://localhost:12345' is therefore not allowed access.). Once is everything OK, once not. When I run my web site 10 times, approximatelly 5 times is problem with CORS. When I use "MODE_ONDEMAND" it´s OK but when I zoom my layer, problem with CORS is displayed and zoomed layer is not displayed. When I used this sample code with the same layer but as dynamicservicelayer with labels, It seems to be OK. Could someone help me with this problem ? Is it problem with settings in proxy ?


texts: Drawing boundary through DEM of glacier valley using ArcGIS for Desktop?

texts:
I want to calculate/draw the boundary extent in Glacier valley with the help of DEM.

Can you please guide me which is the best tools in ArcGIS and steps to calculate.


texts: openlayers 3, osm and custom projection display issue

texts:
I have successfully combined mapserver and openlayers-3 in a custom projection, using proj4 library. As soon as I add osm layer, I can't see the mapserver layer anymore, even if it is the last to be called.
Is it a problem of projection ? I have read that osm data can be reprojected now since OL3 3.11.

<!doctype html>
<html lang="en">
 <head>
 <link rel="stylesheet" href="http://openlayers.org/en/v3.17.1/css/ol.css"     type="text/css">
<style>
   #map {
    height: 256px;
    width: 512px;
   }
</style>
<script src="/proj4js-2.3.14/dist/proj4.js" type="text/javascript"></script>
<script src="http://openlayers.org/en/v3.17.1/build/ol.js" type="text/javascript"></script>
<title>OpenLayers 3</title>

proj4.defs("EPSG:2154","+proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m   +no_defs");

  var map = new ol.Map({
    target: 'map',
    controls: ol.control.defaults({
        attributionOptions: ({
            collapsible: false //ceci ne fonctionne QUE si un fond OSM est appelé !!
                })
            }).extend([
        new ol.control.ScaleLine(),
        new ol.control.MousePosition({
                coordinateFormat: ol.coordinate.createStringXY(4),
                projection: 'EPSG:2154' //juste à l'affichage
            })
        ]),
    layers: [

      new ol.layer.Image({
        title: 'operations Inrap',
        source: new ol.source.ImageWMS({
          url: 'http://10.3.1.51:8080/cgi-bin/mapserv.exe?map=C:\\OSGeo4W\\apache\\htdocs\\fond_mapserver_ol3.map',
          params: {'LAYERS': 'Interregions',
                    //'CRS': 'EPSG:2154',
                    'FORMAT': 'image/png'},
          serverType: 'mapserver',
          extent: [93501.1, 6044342.00, 1248133.25, 7115784.50]
        })
      }),
      new ol.layer.Tile({
        source: new ol.source.OSM()
        })
    ],
    view: new ol.View({
      projection: new ol.proj.Projection({
      code:'EPSG:2154',
      units: 'm'
        }),
      //extent: [93501.1, 6044342.00, 1248133.25, 7115784.50],
      center: [681845, 6606361],
      zoom: 6
    })
  });
</script>


texts:
I want to make many processes for multi raster , and i want to make a loop so it will continue , i want a general case so i can add on it for other . like Clipping as example , this is the code i get .

gdalwarp -q -cutline C:/Users/user/Desktop/IMG_DATA/AOI.shp -tr 10.0 10.0 -of HFA C:/Users/user/Desktop/IMG_DATA/T36PVC_20170329T080601_B08.jp2 C:/Users/user/Desktop/IMG_DATA/Clipped/b8.tif


texts:
If I have this vector layer:

var vectors = new OpenLayers.Layer.Vector("vector", {isBaseLayer: false});
var feature = new OpenLayers.Feature.Vector(
                OpenLayers.Geometry.fromWKT(
                    "POLYGON ((-18 -3, -18 14, 15 14, 15 -3, -18 -3))"
                )
            );
vectors.addFeatures([feature]);

How I can zoom to a selected feature of the layer by a single click event?


texts: How to automatically generate useful scales for data driven pages in ArcGIS desktop

texts:
I have been tasked with creating several hundred maps of properties that intersect a long pipeline.  The maps will zoom to the part of the property which shows the intersection.  The lengths of pipeline that intersect the properties  is highly variable -from two metres to over 16,000 metres.  Therefore, each map will be at very different scales.

I am trying to work out a way to automate the updating of a SCALE field for use in data driven pages.  I have found that in general, for my maps, there is between a 1:4 to 1:5 ratio between the length of pipe segment to the scale.  That is, if the pipe segment is 1100 metres, then a useful scale would be about 1:5000.  I would prefer to round the scale up rather than down so that I can see the full feature on the map.

I want to avoid using scales that are not commonly used such as 1:7000 or 1:19,000, for example.  I have come up with the following solution for use in the field calculator:

import bisect
def scale(n):
  scaleMult=n*5
  scaleList=[500,1000,1250,1500,2000,2500,4000,5000,7500,10000,12500,15000,20000,25000,30000,40000,50000,75000,100000,150000,200000]
  roundedScale=scaleList[bisect.bisect(scaleList,scaleMult)]
  return roundedScale

The above code multiplies the length by some number (I chose 5) then finds out where in my predefined list of scales it should go and returns the next one higher.  The code works pretty well for about 90% of the properties.  I played around with the scaleMult variable and multiplied n by a few different values between 4 and 5 but nothing really improved the maps.

I was thinking that there might be a better solution.  If the pipe is not very straight and doubles back on itself slightly or has a major bend in it, then the feature didn't fit well on the map.  It was usually zoomed in too far.

One idea I had was to use the features bounding box.  Maybe there is a relationship between the feature's bounding box area and the scale.  But, if the bounding box is long and thin, the scale could be different to a bounding box that was a square of the same area.

Does anyone have any ideas on how to implement such a process?

Are there any better solutions that I'm not thinking of?

Looking for solutions in ArcGIS 10.1 and 10.2.


texts: I want to add a new table inside my OSM database with extra data for every restaurant, and then when retrieving the places join the table by the "osm_id", is this possible and if yes, any tips or instructions?


texts:
So I have a point shapefile, with each point representing a bus stop. The shapefile contains a separate stop point for both sides of the road, which I want to consolidate into one point.

So first, I used Integrate to snap the points together so that the points representing the same bus stop are in the same geographic location.
Integrate doesn't get rid of duplicate points though, so to do this I thought Dissolve would work - and it does - but I can't seem to figure out how to use Dissolve while retaining the shapefile's attributes.

From what I understand this should be achievable by adding the fields you want to keep in the "Statistics Fields" box. So I did this with the Stop ID, and for statistic I chose "FIRST" - from what I understand this should retain the first field value in each pair. I also tried "MIN" since the stop ID is a numeric value. Both times, the points dissolve just fine, but in the attribute table the Stop ID field contains the same value for every point.

Here is a screenshot of my inputs for the Dissolve tool:

Unfortunately there is no field in the data that marks whether points represent the same stop, so I can't Dissolve by a certain field to achieve the desired outcome.

So how can I make it so that it properly retains this field when I use Dissolve?


texts: GeoTIFF reprojection to polar stereographic: Can I speed it up?

texts:
For a GeoServer project I need to reproject my high resolution (~40000 pixels wide) GeoTIFF files to polar Stereographic (EPSG 3031 and EPSG 3413).
I manage to do this with gdalwarp, but this is very slow (>20min per file), even if I use intermediate VRT files. I tried to play with several cache configuration options of gdalwarp, but it does not help very much. The reprojection will be part of an operational system, so speed is important.
Are there other tools/ways to do these reprojections? With GDAL or otherwise?

For one of my files, here is the gdalinfo output:

Size is 36000, 12166

Coordinate System is:

GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]

Origin = (-180.000000000000000,82.840000000000003)

Pixel Size = (0.010000000000000,-0.010000000000000)

Metadata:
  AREA_OR_POINT=Area

Image Structure Metadata:
  COMPRESSION=DEFLATE
  INTERLEAVE=BAND

Corner Coordinates:
Upper Left  (-180.0000000,  82.8400000) (180d 0' 0.00"W, 82d50'24.00"N)
Lower Left  (-180.0000000, -38.8200000) (180d 0' 0.00"W, 38d49'12.00"S)
Upper Right ( 180.0000000,  82.8400000) (180d 0' 0.00"E, 82d50'24.00"N)
Lower Right ( 180.0000000, -38.8200000) (180d 0' 0.00"E, 38d49'12.00"S)
Center      (   0.0000000,  22.0100000) (  0d 0' 0.01"E, 22d 0'36.00"N)

Band 1 Block=256x256 Type=Float32, ColorInterp=Gray
  NoData Value=-9999
  Overviews: 18000x6083, 9000x3042, 4500x1521, 2250x761, 1125x381, 563x191, 282x96, 141x48


texts:
I have time-enabled data in ArcGIS Pro. The time field is strictly YYYY-MM-DD, but the time slider shows YYYY-MM-DD HH:MM:SS.
Is there a way to configure the time slider display to strictly show YYYY-MM-DD?
The HH:MM:SS values really aren't meaningful in my data.


texts: How to add a filter to “QFileDialog.getOpenFileName”

texts:
I would like to add a filter to the following code in order to select only the “CSV files”, but I did not succeed:

def filePath(title):
    qfd = QFileDialog()
    path = "D:\ennine\SIG HTB\BGN"
    f = QFileDialog.getOpenFileName(qfd, title, path)
    return f


texts: Creating an outline around groups of polygons

texts:
I am trying to find a tool or workflow for creating a type of minimum bounding polygon or envelope for complex land areas. (Not just a rectangle, but a simplified perimeter)

Can anyone recommend a way to create the simplified black boundary lines around the irregular colored checkerboard areas in the attached image without having to manually digitize.

I’m using ESRI ArcMap Advanced, but open to other tools.


texts: How do you create upstream buffers at a specific distance from another polygon?

texts:
ArcGIS Pro / ArcMap
Licenses: Spatial Analyst, 3D Analyst
I'm working on a project in the Daniel Boone National Forest to create critical habitat areas in order to keep roadside pesticide applicators away from sensitive areas.
As part of this analysis I need to create buffers around streams a specific distance upstream from another stream buffer polygon layer that overlaps with a stream line. There are 15,000+ different stream segments (blue lines) and about 150 existing stream buffer polygons (red lines) I'm working with here.
Essentially, I'm trying to use these existing 150 stream buffers and extend them upstream 1/4 mile. This has to include any lower stream order branches that split within that distance.

I have done some research that led me to the idea this workflow will involve the Flow Direction tool and potentially some stream network analysis.
However, I've never done an analysis like this before and I'm not sure how to stitch these different tools together to get the output I'm looking for. I'm hoping someone here might be able to help.


texts: Geoserver raster legend in reverse order

texts:
If you get a legend from a geoserver raster the minimum value will be on top, the maximum at the bottom.

<sld:ColorMapEntry color="#FFFFFF" quantity="0" opacity="0" label="0 mm"/>
<sld:ColorMapEntry color="#FFFF81" quantity="1" label="1 mm"/>
<sld:ColorMapEntry color="#B5F557" quantity="2" label="2 mm"/>
<sld:ColorMapEntry color="#6FEB36" quantity="3" label="3 mm"/>
<sld:ColorMapEntry color="#3AD62D" quantity="4" label="4 mm"/>
<sld:ColorMapEntry color="#3DB868" quantity="5" label="5 mm"/>
<sld:ColorMapEntry color="#299E98" quantity="6" label="6 mm"/>
<sld:ColorMapEntry color="#216D9E" quantity="7" label="7 mm"/>
<sld:ColorMapEntry color="#1F3F8C" quantity="8" label="8 mm"/>
<sld:ColorMapEntry color="#0D1079" quantity="9" label="9 mm"/>
<sld:ColorMapEntry color="#000351" quantity="10" label="10 mm"/>

But I would like to change the order and to display the maximum value on top of the legend. I tried to change the order of the sld:ColorMapEntrys but now all the values at the raster has a opacity of 0 and no color.

<sld:ColorMapEntry color="#FFFFFF" quantity="10" label="10 mm"/>
<sld:ColorMapEntry color="#FFFF81" quantity="9" label="9 mm"/>
<sld:ColorMapEntry color="#B5F557" quantity="8" label="8 mm"/>
<sld:ColorMapEntry color="#6FEB36" quantity="7" label="7 mm"/>
<sld:ColorMapEntry color="#3AD62D" quantity="6" label="6 mm"/>
<sld:ColorMapEntry color="#3DB868" quantity="5" label="5 mm"/>
<sld:ColorMapEntry color="#299E98" quantity="4" label="4 mm"/>
<sld:ColorMapEntry color="#216D9E" quantity="3" label="3 mm"/>
<sld:ColorMapEntry color="#1F3F8C" quantity="2" label="2 mm"/>
<sld:ColorMapEntry color="#0D1079" quantity="1" label="1 mm"/>
<sld:ColorMapEntry color="#000351" quantity="0" opacity="0" label="0 mm"/>

Is there a solution to show the standard legend in a reverse order?


texts: Converting Cartesian to Spherical/Polar coordinates in PostgreSQL

texts:
I am trying to convert Cartesian to spherical coordinates in PostgreSQL.

Illustrations:

Input: x, y, z
Output: elevation, azimuth

I tried making a function by a code given in the link https://stackoverflow.com/questions/24714604/how-to-convert-spherical-to-cartesian-coordinates-in-postgresql.

I have a similar problem that my function is created by executing the query gives me the error. I have a fundamental knowledge of PostgreSQL/PostGIS.

CREATE OR REPLACE FUNCTION ST_cart2spr(point_i float, point_j float,
point_k float)
  RETURNS geometry AS
$BODY$ SELECT ST_MakePoint(
SQRT(POWER(i,2) + POWER(j,2) + POWER(k,2)),
DEGREES(ACOS(k / (SQRT(POWER(i,2) + POWER(j,2) + POWER(k,2))))),
DEGREES(ATAN2(j, i)))
FROM (SELECT
 $1 AS i,
 $2 AS j,
 $3 AS k
) AS f$BODY$
LANGUAGE SQL IMMUTABLE;

and Query I am trying to run is:

SELECT ST_AsText(ST_cart2spr(SELECT ST_x(geom) As X, ST_y(geom) As Y, ST_z(geom) As Z FROM (
    SELECT (St_DumpPoints(geometry)).geom FROM surface_geometry JOIN building ON surface_geometry.id = building.id
                               WHERE geometry IS NOT NULL
                               AND cityobject_id = 95) as foo)

Any suggestions to understand where I am doing it wrong?
I guess the error part is in the Function call.

But if I run the subquery without putting it in function, it runs fine — the result of the query separately. (See the [part in picture highlighted in grey)


texts: Interactive Web Maps That Update as Users Change Attributes

texts:
I want to create an interactive map where users can go to the website, click on a polygon in the map, and have a chioce of three options. When they click one, it will change the color and update the attribute.

I want to have a split screen showing images from Google Earth (the image would be the same location but a different time period) and a outline of polygons overlaid on it. The user would ideally click the polygon and have a choice of red, yellow, or green pop up. They would then chose the color, and in the background the attribute for that polygon will change.

How would I go about doing this? I am new to web mapping.Thanks.


texts: How to make two attributes display simultaneously on one shapefile?

texts:
I am trying to make a map that displays both an incidence rate of a type of cancer and also display at the same time a percent of the population per state that has a certain vaccination. I've added both the data sets as CSV files and have joined them so that they all show up in the attribute table. But I can't figure out if there is a way to display both at once or if this isn't possible.

The bottom layer of incidence rate is already displayed as a graduated layer using color and ideally I could display the vaccine rate on top as a transparent layer using dots, etc.

If this is not possible, just let me know.


texts: Update license for enterprise geodatabase with GP/Python?

texts:
I am attaching dozens of database files to a SQL Server instance and then I have to run sdesetup command to update the geodatabase license by pointing to my ArcGIS Server license. I am trying to find a way to automate this process.

Two ways to update the license (as to my knowledge are):

1.Connect to geodatabase instance as administrator user in ArcGIS Desktop (ArcMap for instance) and then use the Update ArcGIS Server License window to browse to the .ecp file or keycodes file.

2.Use sdesetup –o update_key ArcSDE command line to update the license. These command line utilities are going to be deprecated soon, though.

Neither of these alternatives provide any way to automate the process. In the latest Esri blog post on this topic, there was no information on this apart from what I have outlined above. I have not found any related ideas on ArcGIS Ideas either.

Is there any other way I am not aware of? Do you know if Esri is going to provide tools for automating license update process?


texts: Visualising clustered features in web maps?

texts:
I'm looking for the best symbol to display clustered features in the ArcGIS Server JS API.

I don't like the "fly-out" animated cluster symbol, as used in the Silverlight API clustering or the JS clustering sample (they feel a bit gimmicky).

So far the best option I've found is the Google Maps marker clusterer symbol.

I won't be showing a number in the centre of the cluster symbol, but will be varying the size of the symbol to show concentrations of features, as in this mockup:

What other clustering symbol options are there? Can you point me to a map showing a nice implementation of marker clustering?


texts: How to move an existing line to a parallel offset of another existing polyline

texts: I'm currently trying to modify an existing line to be offset x distance from another existing polyline, is there any tool/addon that does that? I currently use Copy Parallel and create a dummy line that I move/snap my line to. I was just wondering if there was any way to directly specify that I want an existing line to be offset parallel from another line x units.


texts: Making map by individual with neighborhood shapefile

texts:
I am working with a table* and a Shapefile in ArcGIS. In the table I have individuals living in a neighborhood and a series of characteristics, at the individual level, say weight and height for the example. The table includes a unique code by individual and a code for each neighborhood. The shapefile refers to the neighborhoods, identified with the same code as the table. What I would like to do is to create a map where the element of the table are represented as points within their neighborhood (scattered randomly) and for me to be able to create symbology based on the personal information. I have been able to do the following in ArcGIS:

Join the table to the shapefile (one:many) but only the first person in the neighborhood is considered.
Join the shapefile to the table (many:one)  but this only brings the information of the table and does not allow me to do a map representation. I have also been able to add the X and Y coordinates of the centroid of the shape of the neighborhood. However, mapping based on this is basically a bunch of points stacked unto one another.
Get neighborhood-level averages for weight and height and do the symbology of the shapefile based on these averages.

As I mentioned I would like to do a map with the total number of individuals represented as one point each within the neighborhood boundaries, and with symbology levels that I can select based on their individual info. I use Stata and ArcGIS Desktop.

*Created from a data base originally in Stata.


texts:
I am using SNAP 7.0 to explore Sentinel 2 images for a certain area,,, when I started importing the image using the xml file inside the folder (L1C_T36SYA_A005695_20160725T082012) for example, there was no errors, but when I try to load a product using xml file inside the folder (S2A_OPER_MSI_L1C_TL_SGS__20160625T084235_20160625T134158_A005266_T37SBS_N02_04_01), there is an error.
How do I solve this issue?


texts: Error when importing vector layer into GRASS

texts:
I am trying to import a vector layer into GRASS, but I keep receiving the same error.  The values in the vector layer are important because I would like to use them to create a heatmap.  What could I be doing wrong?

Yesterday I managed to import the layer into GRASS, but it had discarded all the values.


texts: UK Polygon Property Boundaries and Addresses

texts: I am working on a project where i require the mapping of many UK property boundaries individually. Rather than going through each property and manually creating a polygon for its boundary, I was hoping to obtain a layer that connects the property address to the property boundary and automatically produces the polygon. Do you know if such a layer exists or is this achievable to create?


texts: Create Polyline from Coordinates Array in OL OpenLayers 5.3.2

texts:
I'm using OpenLayers 5.3.2, and would like to create my own marker animation using almost exactly the example at https://openlayers.org/en/latest/examples/feature-move-animation.html?q=marker

All I'd like to do is create my own polyline to replace the one used in the demo. I'll also change the center point, but I can do that on my own. There is no information on how to generate a polyline for OL 5.3.2. I've read all examples for past versions, and none of them work.

From what I've read, a polyline is an array of coordinate arrays, such as:

[
    [lon,lat],
    [lon,lat],
    [lon,lat],
    [lon,lat]
]

So I've generated my array to match exactly that. Example of the data I created, stored in a global variable:

var currPath =
    [
        [
            -98.907449,
            30.241299
        ],
        [
            -98.907477,
            30.241337
        ],
        [
            -98.907557,
            30.241449
        ],
        [
            -98.907636,
            30.24156
        ]
    ]

But from here, I have no idea how to proceed.


texts: pyqgis - Editing a feature - Is it possible to set a listener?

texts:
is it possible to set a listener that is fired when a feature is edited?
What I imagine is this: A geometry is edited (e.g. cut operation) then a listener is fired. This event gives me back the altered geometry. Purpose: A geometry is edited and the current area of the geometry is displayed in a window.

The class QgsVectorLayer provides several listeners (e.g. layerModified). Unfortunately I just get the message that something has been edited, but not the edited feature.

Is there any way to set a listener for a specific feature?

All ideas and suggestions are welcome. (:


texts: How to generate features to test different implementations of LinearRing?

texts:
I have some code that throws the following exception:

Points of LinearRing do not form a closed linestring

and I want to check if my code is calling this correctly by supplying some features that either violate or do not violate the constraints of
LinearRing.validateConstruction. Is there an existing Java utility that can be used to generate appropriate feature data for this?


texts: QGIS Build error : zzz-db_manager-5-depend error - qgis master 2.99

texts:
I am compiling  QGIS from master in VS 2015 on x64.

I am getting error as :

Error MSB6006 "cmd.exe" exited with code 1. zzz-db_manager-5-depend 7 C:\Program Files (x86)\MSBuild\Microsoft.Cpp\v4.0\V140\Microsoft.CppCommon.targets

and output followed by :

Build started: Project: version, Configuration: RelWithDebInfo x64
Build started: Project: ui, Configuration: RelWithDebInfo x64
Build started: Project: zzz-db_manager-5-depend, Configuration: RelWithDebInfo x64
Building Custom Rule C:/QGIS/python/plugins/db_manager/db_plugins/postgis/plugins/versioning/CMakeLists.txt
CMake does not need to re-run because C:\build3\python\plugins\db_manager\db_plugins\postgis\plugins\versioning\CMakeFiles\generate.stamp is up-to-date.
Generating ui_DlgVersioning.py
Traceback (most recent call last):
File "C:\QGIS\scripts\pyuic-wrapper.py", line 26, in  import qgis.PyQt.uic.pyuic
File "C:\build3\output\python\qgis__init__.py", line 71, in  from qgis.PyQt import QtCore
File "C:\build3\output\python\qgis\PyQt\QtCore.py", line 1, in  from PyQt5.QtCore import *
ImportError: DLL load failed: The specified procedure could not be found.
C:\Program Files (x86)\MSBuild\Microsoft.Cpp\v4.0\V140\Microsoft.CppCommon.targets(171,5): error MSB6006: "cmd.exe" exited with code 1.
========== Build: 2 succeeded, 1 failed, 10 up-to-date, 0 skipped ==========
Can any one help me to get out of this?


texts: How to remove LayerName.FieldName join format in model builder?

texts:
I am running a join in model builder based on a calculated field in a geodatabase table and a created layer.  After the join I am selecting features and copying selection to a new Shapefile.

My issue is that during the join the resulting joined layer has all field named using the LayerName.FieldName naming convention.  I have tried both checking and unchecking the environment setting of "Maintain Fully Qualified Field Names," however, neither are working.  Only way I can seem to fix it is to manually delete and rename the fields during the copy features process.

Wondering is there a specific modelbuilder setting/tool/process I can use to only keep actual field names.


texts: Using Openlayers Protocol.HTTP GeoJSON format

texts:
just took Openlayers GeoJSON Example and inserted the features to file json1.json instead using a variable, then tried to preform HTTP call:

 var epsg4326 = new OpenLayers.Projection('EPSG:4326'),
        epsg900913 = new OpenLayers.Projection('EPSG:4326');
        function init() {
            var map = new OpenLayers.Map({
                div: "map",
                projection: epsg900913,
                displayProjection: epsg4326, //Is used for displaying coordinates in appropriate CRS by MousePosition control
                layers: [
                    new OpenLayers.Layer.WMS(
                        "WMS", "http://vmap0.tiles.osgeo.org/wms/vmap0",
                        { layers: "basic" }
                    )
                    ,new OpenLayers.Layer.Vector("GeoJSON", {
                        strategies: [new OpenLayers.Strategy.Fixed()],
                        protocol: new OpenLayers.Protocol.HTTP({
                            //readWithPOST: true,
                            url: "json1.json",
                            format: new OpenLayers.Format.GeoJSON()
                        })
                    })

                ],
                center: new OpenLayers.LonLat(-112.169, 36.099)
            });
        }

The json1.json looks  identical to OL example:

 { "type": "FeatureCollection",
    "features": [
      { "type": "Feature",
        "geometry": {"type": "Point", "coordinates": [102.0, 0.5]},
        "properties": {"prop0": "value0"}
        },
      { "type": "Feature",
        "geometry": {
          "type": "LineString",
          "coordinates": [
            [102.0, 0.0], [103.0, 1.0], [104.0, 0.0], [105.0, 1.0]
            ]
          },
        "properties": {
          "prop0": "value0",
          "prop1": 0.0
          }
        },
      { "type": "Feature",
         "geometry": {
           "type": "Polygon",
           "coordinates": [
             [ [100.0, 0.0], [101.0, 0.0], [101.0, 1.0],
               [100.0, 1.0], [100.0, 0.0] ]
             ]
         },
         "properties": {
           "prop0": "value0",
           "prop1": {"this": "that"}
           }
         }
       ]
     }

Found that the server not finds it:

But its there!


texts: Is SRTM 3 second arc free of non-terrain features?

texts:
I made an elevation map of my study area in Macapa, Brazil. But the region has the average elevation lesser than 30 m as I can see from my map derived from SRTM 3 second arc.

According to my knowledge some parts of the region has trees taller than 30 meter and I am afraid that my elevation map shows tree canopy elevation instead of bare land elevation.

My questions are:

Is any correction implemented on SRTM 3 second arc images to remove trees and non-terrain features?
I found a tree height map from NASA. Can I do: SRTM 3 second arc minus elevation height map?


texts:
I have someone how changed the view that the Attributes are displayed when right clicking on a layer > "Open Attribute Table".

I don't see the individual cells but instead see a view that appears to be the structure of the table. How to I change it back?


texts:
Can someone help me to install an add-on for GRASS 6.4.3RC3, on Windows 7? I tried several approaches, but wasn't able to install the add-on. I found several topics here and on the web but the approach doesn't work.

Using g.extension, I install r.basin and it seems that everything is OK, but when I try to use r.basin it doesn't exist.


texts: How to get attribute data based on the GeoJSON you clicked Leaflet?

texts:
I'm trying to pull the attribute data from a GeoJSON layer.

var neighLyr = new L.GeoJSON.AJAX("data/neighbourhood.geojson",{
    }).addTo(mymap);

neighLyr.on('click', function (e) {
    console.log(e.layer.properties.id);
    });

Why doesn't this code work ? It keeps saying the id is undefined.
I'm using https://unpkg.com/leaflet@1.5.1/dist/leaflet.js


texts: PostGIS to Find Dupes and Segregate Geometric Dupes from Originals?

texts:
I've written some code in Python and PostGIS to find duplicate geometric features within the same feature layer.

That part is a success.  The problem I have is the dupes as well as the original feature IDs are part of the result of my cursor.execute command.

What's the best way for me to segregate the first occurrence of a feature geometry from subsequent occurrences of the same geometry?

#Select all features within current layer (i.e. HydrographySrf)
cmd = "SELECT * FROM " + layer
cursor.execute(cmd)
numLayer = cursor.rowcount

#find all duplicates in the current layer
cmd = "SELECT b.ogc_fid FROM " + layer + " AS a, " + layer + " AS b WHERE ST_Equals(a.geom, b.geom) AND a.ctid <> b.ctid"
cursor.execute(cmd)

if (cursor.rowcount > 0):

    #Store Suspect Feature IDs for Reporting at End of Logfile
    for row in sorted(cursor):
        if row[0] not in discoveredFeatureIntegerIDs:
            discoveredFeatureIntegerIDs.append(row[0])


texts: What is a OGRGeometryShadow and why can't I use it in Intersects()-method?

texts:
I'm trying to script something with Python and GDAL/OGR and Shapely.

What I have:

two LineString geometries: "geom1" and "geom2"
I made an Intersection: "intersection = geom1.Intersection(geom2)"

Now I wanted to check, if the intersection-Point is either the starting- or endingpoint (or none of both) of one of the input geometries.

What I did:

from shapely.geometry import Point
point_start_1 = Point(geom1.GetPoint_2D(0))
point_end_1 = Point(geom1.GetPoint_2D(geom1.GetPointCount()-1))
point_start_2 = Point(geom2.GetPoint_2D(0))
point_end_2 = Point(geom2.GetPoint_2D(geom2.GetPointCount()-1))

The following code won't work:

if intersection.Intersects(point_start_1):
            print "intersection and Starting Point of geom1 is the same Point"

The Error Message is:

"TypeError: in method 'Geometry_Intersects', argument 2 of type 'OGRGeometryShadow *'"

My Question is:

What is a OGRGeometryShadow and why can't I intersect it with another Geometry?
When I print the start_point_1 it says: "POINT (6.331836999999999 53.6154012)"
So it looks like a normal geometry. But as I said "Intersects"-method doesn't work.


texts: Is it better to work with non-contiguous or contiguous polygons?

texts: I would like to find width-height of a non-contiguous polygon in order to draw paths within it. My question: Is it better to use a bounding box (contiguous polygon) and then clip keeping my real area or to work with the non-contiguous polygon from the beginning?


texts: Creating mapping application with user input to dynamically weight/shade map layers?

texts:
Does anyone know of any examples of web-apps that allow user-input to change variable weights used to calculate symbology?

My organization is looking to build a dynamic mapping application and while I have the tech-savvy to envision and outline the project, I don't have the coding experience to implement it myself (I've been working with GIS for close to 10 years, but I am new to development). I have access to both GeoServer and ArcGIS Server, and I'd appreciate any input on how to go about this task (FLEX vs REST vs JS vs MYSQL vs POSTGIS vs GeoServer vs ArcGIS, etc)

We are creating a series of opportunity maps based on neighborhood indicators. Essentially, we compile a set of data at the census-tract level for an MSA and convert the variables to Z-scores. Then we aggregate Z-scores to create a composite index for a certain category. We then use the composite index score to break tracts into 6 quantiles and shade them accordingly. For example:

Environmental Quality Index (heading)

- Cancer Risk
- Air Quality
- # of Hospitals accessible in 30 mins
- average ambulatory response time

In this example, we'd calculate Z-scores for each variable (cancer risk, air quality, etc) then add them together to create an additive index for environmental quality. Then we simply classify by quantile to create a map showing areas ranging from "Very High" to "Very Low" Environmental Quality. This would create a composite index where each variable is weighted equally. However, We would like to create an app that renders a general map using equal weights, but allows users to define variable weights as they choose. For example:

Environmental Quality

 - 30% | Cancer Risk
 - 15% | Air Quality
 - 10% | Hospitals
 - 45% | Ambulatory Response Time

In this case, the percentages are multiplied across the Z scores to create the desired weights, the quantile ranges will be recalculated, and the map shading will be updated accordingly.

Using databases, I believe this should be a fairly straightforward coding task, but I have been unable to find any examples of similar projects.

This example using ArcGIS server 10.1 looks to be promising:

http://blogs.esri.com/esri/apl/2012/03/20/mapping-health-indicators-using-dynamic-layers-in-arcgis-server-10-1/

Though I would still welcome advice on the most appropriate API (FLEX, REST, JSON, etc) and how to begin coding the app.


texts: Crosstable two fields of pointfile using ArcMap?

texts:
I'd like to produce a simple crosstabulation from two fields of a point feature class. These fields contain land use class information which I'd like to compare.

One option is to export a dbf and do this with Excel or with R.
Another option that might work is to use ArcMap TabulateIntersection -tool.
(Although it returns a flat table with matching classes and observations rather than a real crosstable.)

This is how it looks in arcpy:

TabulateIntersection_analysis (in_zone_features, zone_fields, in_class_features, out_table, {class_fields}, {sum_fields}, {xy_tolerance}, {out_units})

# My version: Cross-table class1 values and class2

 arcpy.TabulateIntersection_analysis(Plots,"class1",Plots,path+OutTable,"class2","#","#","UNKNOWN")

For what I understand of this, this should compare the values of two fields of the same file - as the file is given here twice. And it sort of does, except the result is different from the first method I told (I tabulated this in R). Also there are same 'extra' observations. There are 4560 plots in my data and the latter method returns c. 4700 observations.  The results are close to each other, but there are a some discrepancies. Does anyone have an idea what I might do wrong here?

Pivot table code snippets for crosstabulation:

 #1.calculate frequencies:
 arcpy.Frequency_analysis(Plots,path+OutFreq,"class1;class2","#")
 #2.make a pivot table
 arcpy.PivotTable_management(OutFreq,"class1","class2","FREQUENCY",path+OutPivot)

This video might also be helpful: https://www.youtube.com/watch?v=9dGjuEQie7Y


texts: Disable customized context menu for some cases using python

texts: I create a python add-in tool to add new item to the context menu, but I want the new item enabled only when user right click on map when within a buffer of polyline feature. I haven't figure out how to make it, which means for the area qualified (within a buffer of polyline) can see the new item on context menu, others won't show the new item (or the new item is disabled in grey). I am using ArcGIS 10.2.2 and python.


texts: QGIS 3.8 model to Python script, not working in 3.4

texts:
I am using QGIS 3.4 but it does not have the option of saving a model to a Python script. So I installed QGIS 3.8, which has that option, loaded my model and saved it as a Python script.
I loaded the script into the Python editor, clicked "run" and expected the script to work in Version 3.8, but it didn't.
The Python console just gives back

exec(open('P:/.../QGIS_Import.py'.encode('utf-8')).read())

With loading it as a script it worked in QGIS 3.8:
The model and the python script do not generate the same results, however, that is not the real problem.
When I run the python script in QGIS 3.4 it gives this error:

Traceback (most recent call last):
File "C:\Users\heubacha\AppData\Roaming\QGIS\QGIS3\profiles\default\processing\scripts\__QGIS_3.8_Import_Original.py", line 73, in processAlgorithm
'OUTPUT': QgsProcessing.TEMPORARY_OUTPUT
AttributeError: type object 'QgsProcessing' has no attribute 'TEMPORARY_OUTPUT'

The code might not be downwards compatible...what is the correct coding in QGIS 3.4?
This is the code snippet where the error occurs:

alg_params = {
            'DISCARD_NONMATCHING': True,
            'FIELD': 'STP',
            'FIELDS_TO_COPY': None,
            'FIELD_2': 'OZ2',
            'INPUT': outputs['FbRefactorFields']['OUTPUT'],
            'INPUT_2': parameters['formstp'],
            'METHOD': 0,
            'PREFIX': '',
            'OUTPUT': QgsProcessing.TEMPORARY_OUTPUT
        }
        outputs['FbfsJoinAttributesByFieldValue'] = processing.run('native:joinattributestable', alg_params, context=context, feedback=feedback, is_child_algorithm=True)

        feedback.setCurrentStep(50)
        if feedback.isCanceled():
            return {}


texts:
I'm trying to work with the shapefile on python,
the file " donnees.shp" contains 3d points
I need to extract attributes from this file
but my code returns an error:

>>> import shapefile
>>> r = shapefile.Reader("donnees.shp")
>>> r.shapes()[0].shapeType
... 1
>>> r.shapes()[0].points
...[[841001.0, 6399999.0]]
>>> r.shapes()[0].z

... Traceback (most recent call last):
File "", line 1, in
r.shapes()[0].z
AttributeError: _Shape instance has no attribute 'z'

what is the nature of this error?
what am I missing?


texts:
I currently try to mask clouds in raster files with a mask that I generate with the help of another band. In my currenty attempt I use python with rasterio with the following code:
src = rio.open(io_raster, mode="r+")

with rio.open(clouds_raster) as clouds:
    np_clouds = clouds.read(indexes=1)

    np_bool = np.where(np_clouds >= 7, 0, 255)

    # TODO adapt resolutions
    # TODO increase number of masked pixels
    # with rio.open(io_raster) as src:
    msk = src.read_masks()

    cld_msk = (msk & np_bool)

    src.write_mask(cld_msk)
    src.close()

Then I tried to open the masked raster file in qgis. Is there somethong wrong with my code, or doesn't qgis support this feature?


texts:
I'm using Qgis 3.10.

I have two layers, let's say A and B. In layer A, there is informations about construction work (one row by construction). Each construction has an ID (ID-A). In layer B, there informations about the roads affected by each constructions (one row by road). Each road is linked to the construction by an ID (ID-B).

I'm trying to get the maximum value from fieldC in layer B in a virtual field FieldA in layer A, sorted by construction work. I'm using this expression :

aggregate('LayerB', 'max',
 expression:= "fieldC",
 filter:= attributes( @parent )['ID-A']="ID-B"
 )

The problem is that this expression is working only for some construction.
I've tried using different filters, different aggregate function, I didn't manage to get the correct result. Is there something I did wrong ?


texts: Calculating mean difference using ArcGIS focal statistics?

texts:
For each pixel of my raster, I would like to calculate the mean difference with the neighbour pixels on a given radius.

The tool "focal statistic" looks like exactly what I need, however, there is no way to calculate the mean difference with it, only the sum, the mean.

Does anyone have an idea about how to solve this by using an approach equivalent to the focal statistic but that allows me to calculate the mean difference? It could be in several steps, like first calculating the sum or the differences from which I could, then calculate the mean... or whatever that works to reach my goal.


texts: Dissolving lines with same directions using ArcGIS Desktop?

texts:
I have a feature class containing road polylines which has to be dissolved by some fields. But recently noticed that if some road has same atrtributes but different directions, Dissolve tool doesn't keep it separate.

Tried Unsplit lines parameter but it didn't help.
Any thougths how to keep lines like that untouched?


texts: Bring buildings from one CityGML-dataset to same planar level

texts: I need an idea of how to bring the buildings of a raw CityGML-dataset to one planar level. What I can is writing a javascript, which goes through the xml and substract the relevant heights but... maybe there are already solutions/tools/workflows?


texts: Length of the boundary intersection between two polygons in R

texts:
I am trying to find the length of the border shared by two polygons. I have run cleangeo on these shapefiles. There were no flags that came up. I have looked through several other posts on the topic such as https://stackoverflow.com/questions/33924895/intersect-spatial-lines-in-r and https://stackoverflow.com/questions/45338384/calculate-the-length-of-shared-boundaries-between-multiple-polygons. But I was unable to resolve the errors I was getting.

I have provided three separate pairs of polygons that represent some of the errors that I am encountering. They can be found at: https://www.dropbox.com/sh/kgqrbiojye2tu7o/AACe0AAqqpGmMyzVfKTiwVula?dl=0

When I run the following code, I get an error saying there is no slot of name "lines" for this object of class "SpatialCollections":

    library(sf)
    library(sp)
    library(raster)
    library(rgeos)
    library(rgdal)
    library(maptools)

    setwd("C:/Users/jasplund/Downloads/shapefiles/Georgia/")

    Shapefile <- readOGR(dsn=".",layer="test_polygons1")

    xx <- as(Shapefile, "SpatialLines")
    line <- rgeos::gIntersection(xx[1,], xx[2,],byid = TRUE)
    l_line <- sp::SpatialLinesLengths(line)
    l_line

By using the same code and replacing test_polygons1 with test_polygons2, you will see the second example. When I run

    plot(Shapefile)

I see the following image for test_polygons1.

When I view the intersection between these two polygons, I get the following figure.

This "line" has a marker on it, but some of these intersections have no markers and simply look like a line. Or, in the case of test_polygons2, the intersection forms disconnected lines.

I'm fairly confident I know what is wrong, but I am unable to resolve it. I think when the intersection of the two polygons is found, there is a very slight difference between these two polygons causing the intersection to be not exactly a line or missing some lines or simply points. I even tried using gBuffer as suggested in https://stackoverflow.com/questions/33924895/intersect-spatial-lines-in-r, but that gave me a different error, though I cannot remember what exactly this error said.

How exactly do I calculate the length of shared boundary between these polygons?


texts:
I am working on a project of mine in which I have to render geospatial data.  The data is in the form of geoTIFF images and I have to render that using OpenGL. I am good at C++ and OpenGL, but I don't have any knowledge about geospatial data.

I want to render the images in the form of texture but the GDAL library just makes the hillshade of the image and I don't know how to use the query in C++.


texts:
I have 7 separate heatmaps representing 7 years of data, all styled the same in QGIS. The shapes of the heatmaps change slightly between the different years and it would be nice to quickly toggle 2016 off and immediately view 2015 without having to take my eyes off the map canvas to click the little boxes in the layers panel.

I could just take a print screen of each layer and paste them into something like powerpoint and use my arrow keys to have the same effect. Is there a faster approach?


texts: Auto Fill Field with Value based on Location

texts:
My Question is the same as this: except I am using ArcGIS 10.5.1
Autofill field based on location in QGIS

I do not want to use field calculator, is it possible to have certain fields auto-populate when creating a new feature, with values based on its location?

For example, if a new point is within polygon A, then field 1 = True
if not within polygon A field = false

Default value is not ideal, it does fit my need. Any suggestions?


texts: QGIS Graphical Modeler not finding output files

texts:
I am going through the GST 102 Module 5 tutorial which uses the graphical modeler. I keep getting this message when trying to use the output of an Extract by attribute algorithm as a clipper for a point layer.

  Prepare algorithm 0: Extract by attribute
  Running Extract by attribute [1/10]
  Parameters: LAYERNAME =C:\Projects\qgis\data2\Counties.shp, ATTRIBUTE =COUNTY, COMPARISON =0, COMPARISONVALUE =Nueces
  Converting outputs
  Outputs: RESULT =C:/Projects/qgis/data2/Output/NuecesCounty.shp
  OK. Execution took 0.946 ms (1 outputs).
  Prepare algorithm 1: Clip
  Running Clip [2/10]
  Parameters: INPUT =C:\Projects\qgis\data2\Airports.shp, OVERLAY =C:/Projects/qgis/data2/Output/NuecesCounty.shp
  Failed
  Error executing algorithm 1 [Errno 9] Bad file descriptor See log for more details

I have checked and the C:/Projects/qgis/data2/Output/NuecesCounty.shp layer is created during the first process as it should. I am using Windows 7.


texts: Calculating colors of the squares in a vector grid in QGIS 2.18.5

texts:
I have a layer with many buildings as points on a map. One of the attributes of each building is "download_speed".

I want to put a square grid (100 x 100 meters) on top of the map. The squares should behave as follows:

The square should only be visible if there is at least one building in the square.
The square should be red if none of the buildings in the square have a "download_speed" above 10 (Mbit/s).
The square should be gray if some of the buildings in the square have a "download_speed" above 10 (Mbit/s)
The square should be black if all of the buildings in the square have a "download_speed" above 10 (Mbit/s)

I'm a complete newbie on QGIS (and GIS-software generally), but an experienced user of Python in data science.

The end result should look somewhat like the picture below:


texts:
I'm writing a plugin in PyQGIS 3, and I've to test for vector joins before run a code, so I've found here about  vectorJoin() and write this:

def joinTables(self):
    d, f = self.lyrPair()  # d = {target;layerToJoin}
    for k, v in d.items():
        target = QgsProject.instance().mapLayer(k)
        layerToJoin = QgsProject.instance().mapLayer(v)
        v_lyr = QgsProject.instance().mapLayers().values()

        for lyr in v_lyr:
            joinsInfo = lyr.vectorJoins()

            if not joinsInfo:
                QMessageBox.information(iface.mainWindow(), "Success",
                                "No previously joins!")
                # do some stuff here
            else:
                QMessageBox.critical(iface.mainWindow(), "Error",
                                "Previously Joins exists!")
                # ask for remove join

        fieldToJoin = QgsProject.instance()
        symb = QgsVectorLayerJoinInfo()
        symb.setJoinFieldName('id_feature')
        symb.setTargetFieldName('id')
        symb.setJoinLayerId(layerToJoin.id())
        symb.setUsingMemoryCache(True)
        symb.setEditable(True)
        symb.setDynamicFormEnabled(True)
        symb.setUpsertOnEdit(True)
        symb.setPrefix('')
        symb.setJoinLayer(layerToJoin)
        target.addJoin(symb)
        layerToJoin.startEditing()
        target.triggerRepaint()

but in both case, when has joined and when hasn't joined, the code above returns that there is. First, vectorJoins() returns a list of joined layers and I've to test every single layer, the second I've missed something in this code?


texts:
Getting an error when trying to load a shapefile via readOGR. readOGR works fine on my personal computer (mac) but I am running into issues on my work computer (windows).

When loading a shapefile using the typical syntax, I get the following error:

    Warning in ogrInfo(dsn = dsn, layer = layer, encoding = encoding, use_iconv = use_iconv,  :
      ogrInfo: all features NULL
    Error in readOGR(dsn = ".", "shp_example")  : no features found

For something reproducible, when I run the following code from the example from help(readOGR):

    library(rgdal)

    ogrDrivers()
    dsn <- system.file("vectors", package = "rgdal")[1]
    ogrListLayers(dsn)
    ogrInfo(dsn)

I get the following warning after "ogrInfo(dsn)", telling me that all features are null. Note that I only get this error on my work computer; this error does not occur on my personal computer:

    Warning in ogrInfo(dsn) : ogrInfo: all features NULL

Any ideas as to what is going on?


texts: ogr2ogr - Convert geojson to shapefile create offset in coordinates

texts:
I want to convert this geojson file below into a shapefile using ogr2ogr command ogr2ogr output.shp input.geojson.

However, an offset is appearing in the output coordinates that make the shapefile not correctly aligned with the input.

input.geojson:

{
   "type":"FeatureCollection",
   "features":[
      {
         "type":"Feature",
         "geometry":{
            "type":"LineString",
            "coordinates":[[-7551132, 6293816],[-7560916, 6260489]],
         },
         "properties":{
            "id":36
         }
      }
   ],
   "crs":{
      "type":"name",
      "properties":{
         "name":"urn:ogc:def:crs:EPSG:3857"
      }
   }
}

Is there a problem with my geojson? What could be the cause of the offset?


texts: Performing least cost path analysis on QGIS using permeability scores for 60 points in network?

texts:
I have a map that is made up of polygons that apparently cover the whole area without leaving gaps, each representing a field, a wood, a road, a house etc.  Each polygon has a score in the attributes column representing how difficult, relatively, it is for an animal to move through that habitat type.  Some of the polygons are Home habitat patches, with a permeability score of 0. I would like to calculate the least-cost-path to the nearest or most easily accessible neighbouring Home habitat patch for each of the Home habitat patches.  My attempts so far result in a single small squiggle going in no particular direction from one point to nowhere.  All the instructions I can find either assume an understanding of programming language and/or the overall order of steps needed; and they are for only two points - a start and a finish; or they are for calculating paths through mountain ranges or similar, which seems more complicated than I need as it goes by relative rather than absolute scores for grid squares (for mountains the path depends on how much higher your neighbouring patch is, not how easy it is to move through in itself).

I am currently using QGIS on an elderly laptop, but I have access to ArcGIS in a computer lab if necessary.

I think the steps go like this, but I would welcome any clarification:

1) Turn the vector layer into a raster layer based on the permeability score column.  I can do this easily in QGIS, and it looks correct; although I am worried that I might need to make sure that all the polygons actually join up with no gaps first, as I guess that there might be some little hidden "no data" areas in there otherwise.

2) Define starting and stopping points. I have seen this done as points on separate layers, and as grids squares on the same layer, and as shapes, in the examples - I'm not sure which is best in this case.  From the ArcGIS manual it looks as though you can get it to make least-cost-paths connecting all points in a network, rather than just A and B, but I don't know how to do this, and I'd imagine that it multiplies up the computer power needed at a massive rate, as the more points you have the more possible paths and possible nearest points to get to there are. Or that most algorithms for networks are designed to amalgamate the paths to make roads for people, rather than animals going any old way.  But perhaps not?

3) Make some sort of intermediate layers that show how hard it is to get to any given square from a home patch?

4) Add a direction indicator layer on top of that layer?

5) Run an algorithm that chooses between the possible paths over these layers. Which is the right one?

Am I on the right track?
There are 60 patches. Am I attempting the impossible?
If not, are there any fool-proof instructions for beginners?


texts: Cannot "Move Feature" in QGIS for PostGIS layers

texts:
I cannot seem to Move Feature for PostGIS layers in QGIS 3 (2.99 43a9bc1) after I enabled editing as shown below.

Even though I can add features with the Add button, the Move Feature button is greyed out.

Does anyone know what might be the cause or how to fix it?

--EDIT--

Regarding the table, I do have the privileges as I am the owner of the db, and I just verified that I can create features in it. I also created a field called ogc_field because of prior issue with updating the table in question. Unable to update features in tables without FID: How to use gid field as the FID column with ogr/gdal python? .

ogc_fid | integer              | not null default nextval('cfl.match_opl_prmatch_ex_pl62__p5_ogc_fid_seq'::regclass) | plain    |

I don't know if this is sufficient for QGIS.


texts:
I'm trying to create a plunker to find out a solution when adding projection list in Openlayers. But when I add a sample map, it's not showing the map and there is no error in the console either.

var olview = new ol.View({
  center: [0, 0],
  zoom: 2
});
var map = new ol.Map({
  layers: [
    new ol.layer.Tile({
      source: new ol.source.OSM()
    })
  ],
  target: 'map',
  view: olview
});

Plunker


texts: Programmatically modifying layer datasource on project save event?

texts:
I have a QGIS plugin which allows users to add WFS/WMS layers to the map.

I would like to be able to scramble the user ID part (token) of the uri when the project is saved/closed and then use my own instance of QgsProjectBadLayerHandler to fix the uri when a project is opened. The result being that only users with the plugin and correct access levels are able to to view the WMS/WFS layers, so if a qgis project file is opened by someone without the plugin they will be unable to loads the layers.

I am having no luck changing the uri when the project is saved. Here is what I have so far:

self.project.projectSaved.connect(self.project_saved)
self.project = QgsProject.instance()

...

def project_saved(self):
    reg = self.LoginAndStore.check_reg()
    token = str(reg[1])
    toc = self.iface.legendInterface()
    toc_layers = toc.layers()
    for layer in toc_layers:
        try:
            dp = layer.dataProvider()
            src = str(dp.dataSourceUri())
            if token in src:
                no_token = src.replace(token, "NOTOKEN")
                dp.setDataSourceUri(no_token)
                #dp.reloadData()
        except AttributeError:
            pass
    f = QFileInfo(r"C:\Users\me\Documents\testing.qgs")
    self.project.write(f)

How do I get this working or are there other methods to prevent unauthorized access?


texts:
I´m trying to upload geoTIFF files to GeoNode 2.6.3, installed in Ubuntu 16.04.
When I try I get this error:

Could not save the layer layer_name, there was an upload error: Error while storing uploaded file::Unable to create /usr/share/geoserver/data/data/geonode/layer_name

It used to work in GeoNode 2.4, and the geoTIFF files I´m uploading are the same ones I had in 2.4 version. I tried the upload from the web portal interface and from the QGis plugin "Geoserver-explorer". Had same error in both cases.
What could be the problem?

edit: add parts of geoserver.log

2017-10-19 15:26:19,059 ERROR [geoserver.rest] - No such datastore: geonode,landsatcr_3_re8_58
2017-10-19 15:26:19,060 ERROR [geoserver.rest] -
org.geoserver.rest.RestletException

2017-10-19 15:26:19,078 ERROR [geoserver.rest] - No such coverage store: geonode,landsatcr_3_re8_58
2017-10-19 15:26:19,078 ERROR [geoserver.rest] -
org.geoserver.rest.RestletException

2017-10-19 15:26:19,096 ERROR [geoserver.rest] - No such wms store: geonode,landsatcr_3_re8_58
2017-10-19 15:26:19,096 ERROR [geoserver.rest] -
org.geoserver.rest.RestletException


texts:
When exporting a chart to CSV is there a way to configure the date to day month year (dd-mm-yyyy) format as it is difficult to convert the date in Excel for graphing?

Sample of export:

when using the sample code example from Earth Engine:

// Define a FeatureCollection: regions of the American West. var regions = ee.FeatureCollection([   ee.Feature(    // San Francisco.
    ee.Geometry.Rectangle(-122.45, 37.74, -122.4, 37.8), {label: 'City'}),   ee.Feature(  // Tahoe National Forest.
    ee.Geometry.Rectangle(-121, 39.4, -120.8, 39.8), {label: 'Forest'}),   ee.Feature(  // Black Rock Desert.
    ee.Geometry.Rectangle(-119.15, 40.8, -119, 41), {label: 'Desert'}) ]);

// Load Landsat 8 brightness temperature data for 1 year. var temps2013 = ee.ImageCollection('LANDSAT/LC8_L1T_32DAY_TOA')
    .filterDate('2012-12-25', '2013-12-25')
    .select('B11');

// Create a time series chart. var tempTimeSeries = ui.Chart.image.seriesByRegion(
    temps2013, regions, ee.Reducer.mean(), 'B11', 200, 'system:time_start', 'label')
        .setChartType('ScatterChart')
        .setOptions({
          title: 'Temperature over time in regions of the American West',
          vAxis: {title: 'Temperature (Kelvin)'},
          lineWidth: 1,
          pointSize: 4,
          series: {
            0: {color: 'FF0000'}, // urban
            1: {color: '00FF00'}, // forest
            2: {color: '0000FF'}  // desert }});

// Display. print(tempTimeSeries);


texts: C# Custom Dockable window: grouping them together

texts:
I'm adding in some custom windows for our clients. I want to group my dock windows together on load. The end result would look like they are tabbed (from the bottom).

I have the following code to my button (x 3 ):

ESRI.ArcGIS.Framework.IDockableWindow aaView;
...
if (aaView == null)
            {
                //Get dockable window.
                UID dockWinID = new UIDClass();
                dockWinID.Value = "myView";
                aaView = ArcMap.DockableWindowManager.GetDockableWindow(dockWinID);

                if (aaView == null)
                {
                    ErrorBox.DoModal("Launchmyview Error", "", "", "", 1);
                    return;
                }
            }
            aaView.Show(!aaView.IsVisible());

I have the following in the esriaddinx(x3):

 <DockableWindow id="MyView" class="CustomEvents.ImageMetaData+AddinImpl" caption="Raster Metadata">
          <InitialPlacement height="300" width="300" state="unpinned" position="group" neighbor="{E1F29C6D4E6B11D2AE2C080009EC732A}" />
        </DockableWindow>

I know that I can move the windows over eachother manually once the UI is open, but I want to group them programmtically.


texts: Executing functions for multiple directories using PyQGIS?

texts:
I would like the user to select 2 directories, each directory having their own defined algorithms. For example: First_Directory will add a column to all shapefiles in the chosen directory; Second_Directory will delete a column to all shapefiles in the chosen directory.

When selecting one directory and leaving the other blank, i.e. only First_Directory and leaving Second_Directory blank, the algorithms defined for First_Directory work fine and vice versa.

However when selecting for both, it seems that only First_Directory is processed. What I would like is, ofcourse, for both directories to be processed sequentially.

The code below is merely trying to show the logic, it's the IF statement that I am concerned about (whether or not it's even required in the first place!)

##Test=name
##First_Directory=folder
##Second_Directory=folder

path_1 = First_Directory
path_2 = Second_Directory
path = path_1 or path_2

os.chdir(path)

def First(path_1):
#Do some processing

def Second(path_2):
#Do some other processing

if path_1:
    score_1(path_1)
elif path_2:
   score_2(path_2)
elif path_1 and path_2:
    score_1(path_1)
    score_2(path_2)
else:
    pass


texts: Where is my model? (saved in project, QGIS 3.10/3.16)

texts: I saved a model (Graphical Modeler) into my project and I don't know where it was saved. How can I get it back?


texts: .tif image (made in QGIS 3) imported into Google Earth Pro is blank/has no symbology options

texts:
I am trying to import .tif images into Google Earth Pro, which have been created in QGIS 3.

The image has two colours, black and white, with only two values, 1 and 0, which was represented by each colour.

When this is brought into Google Earth Pro the image aligns in the right location but is just black. I cannot find any way to change its appearance or any information about the data which is held in it. How can I do either of these things?


texts: Rescaling coordinates of rasters, shapefiles and spatial objects from meters to km in R?

texts:
All my shapefiles, rasters and Spatial objects were projected to UTM NAD83. However, I need to convert them from meters to km for calculation reason. For rasters I can easily convert the "values" by dividing them with 1000. However when I plot them the coordinates were still in meters. It didn't look right if I convert the extent though (by dividing the extent with 1000). (the resolution was 25 meters).

Are methods for converting meters to km different for shapefiles, spatial objects and rasters?

Eventually I would like to create a discrete distance-to-road surface in which the unit and the coordinate system were both in km.

##Read road shp.file
roads<-readOGR(dsn="E:/Doccuments/GIS", layer="Roads_2005_INDOT")

#coerce SpatialLinesDataFrame object to SpatialLines object
rd.spl<-as(roads, "SpatialLines")

#create an empty raster with the same extent and grid size=25
r<-raster(ext=extent(rd.spl), res=25)
rd.rs<-rasterize(rd.spl, r)          #rasterize SpatialLine
#assign 0 to all roads values
rd.rs<-reclassify(rd.rs, c(0,108,0))
rd.rs

class       : RasterLayer
dimensions  : 125, 104, 13000  (nrow, ncol, ncell)
resolution  : 25, 25  (x, y)
extent      : 590847.7, 593447.7, 4235915, 4239040  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=16 +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0
data source : in memory
names       : layer
values      : 0, 0  (min, max)
plot(rd.rs)

#create discrete "distance-to-road" surface
d<-distance(rd.rs)/1000   #distance in km
plot(d)
plot(locations, add=TRUE)
plot(movements, add=TRUE)


texts: I am interested in the determination of lithology by using ASTER band ratio. I try to a lot of band ratio (for example RGB 4/8, 4/2, 8/9) but I didn't get a good result. Do you think what is the best band ratio to distinguish sedimentary rocks and other lithologies in ASTER data ? Especially, the discrimination of sedimantary rocks is important for me.


texts: ArcGIS JavaScript manually render a layer with supplied JSON

texts:
Given the following JSON which is retrieved outside the ArcGIS/dojo framework, how would I go about rendering this layer?

    {
        "currentVersion": 10.01,
        "serviceDescription": "",
        "mapName": "Layers",
        "description": "",
        "copyrightText": "",
        "layers": [
            {
                "id": 0,
                "name": "Transportation Stations-unexpired",
                "parentLayerId": -1,
                "defaultVisibility": true,
                "subLayerIds": null,
                "minScale": 400000,
                "maxScale": 0
            }
        ],
        "tables": [],
        "spatialReference": {
            "wkid": 102100
        },
        "singleFusedMapCache": false,
        "initialExtent": {
            "xmin": -10398329.13791846,
            "ymin": 2687052.938617364,
            "xmax": -7541451.81403804,
            "ymax": 5379667.868452799,
            "spatialReference": {
                "wkid": 102100
            }
        },
        "fullExtent": {
            "xmin": -10152568.449939176,
            "ymin": 2795685.894150028,
            "xmax": -7787215.334097848,
            "ymax": 5828927.975739868,
            "spatialReference": {
                "wkid": 102100
            }
        },
        "units": "esriMeters",
        "supportedImageFormatTypes": "PNG24,PNG,JPG,DIB,TIFF,EMF,PS,PDF,GIF,SVG,SVGZ,AI,BMP",
        "documentInfo": {
            "Title": "test.mxd",
            "Author": "test",
            "Comments": "",
            "Subject": "",
            "Category": "",
            "Keywords": "",
            "Credits": ""
        },
        "capabilities": "Map,Query,Data"
    }


texts: I would like to know the path of a feature class that is in my TOC in ArcMap, how can this be done with ArcPy?


texts: Identifying coordinate system of San Francisco data?

texts:
I just found a dataset of trees in San Francisco, but I'm unable to determine the system used for the coordinates. Here are some samples:

5998894.045593,2119007.035323
6018697.2048701001,2097076.1109207

Does this look familiar to anyone?


texts:
We have some layer which is creating data every day. On every 3 months or if data is more than 100K geometry we need to archive the layer.

How can we archive this in GeoServer and any best practice how realworld is doing it?


texts: Is there a way to connect 2 tiles of downloaded landsat data in ArcMap?

texts: I have downloaded two Landsat 5 images because my study area is divided between the 2 tiles. I have created a composite of all the bands and then tried the mosaic tool in ArcMap but the result I have contains the black border that the raw images have and it is not correct. Is there a way how I could connect the 2 images in order to have only one image in the end to continue with my processing?


texts: How to begin with tool development for GRASS gis or QGIS GRASS

texts: I wanted to begin with customized tool development in Grass GIS or QGIS grass. How should i begin in to develop the tool. I search a lot in google search engine for begin with tool development in Grass GIS or QGIS grass. Can any one please help me in the tool development process in Grass.


texts: How to query distict row with the same label from different linestring with multiple segment?

texts:
I have a multilinestring street table in postgis with two columns: label and prefix.
This table has 2 different roads with the same name and each road contains multiple segments.
Roads were splitted on each intersection during authoring in QGIS so that the table has total 14 rows. How to query the two so that i have two unique street rows with prefix and label column in the result?


texts:
I have installed Grass 7.0.5 but I cannot open a mapset via QGIS 2.18. I receive this error message:

Cannot open the mapset. Mapset lock failed (exitStatus=1, exitCode=0, errorCode=1, error=Process crashed stdout=, stderr=dyld: Library not loaded: /usr/local/lib/libintl.8.dylib
  Referenced from: /Applications/GRASS-7.0.app/Contents/MacOS/lib/libgrass_gis.7.0.5.dylib
  Reason: image not found
)

If I create a new mapset instead of opening one, the error reads "New mapset successfully created, but cannot be opened:" and then the rest of the error as above.

I am running Mac OSX El Capitan, and I have already disabled SIP. I also encountered, and solved, Problems finding GRASS module after upgrading QGIS (Mac).

Can this latest problem be solved, or should I revert to Grass 6?


texts: QGIS - connect to ArcGIS FeatureServer hosted behind Azure proxy

texts:
I am using the QGIS (v3.16.0) ArcGisFeatureServer tool in the Browser window to try and connect to the following ArcGIS Feature Server: https://gojdippmaps.azurewebsites.net/proxy.ashx?https://maps.gov.je/arcgis/rest/services/Historic_Buildings/Historic_Buildings/FeatureServer
Note that this server requires an http referrer, such as: https://www.gov.je//citizen/Planning/Pages/HistoricEnvironmentDetail.aspx
Note also that the host is gojdippmaps.azurewebsites.net and the https://maps.gov.je.. url is treated as a parameter in the query string (preceded by ? after the path /proxy.ashx).
If I add the url as shown above I get "Connection failed: Network error" from the ArcGisFeatureServer tool. Hovering over this yields the following message:

Error transferring https://gojdippmaps.azurewebsites.net/proxy.ashx?https://maps.gov.je/arcgis/rest/services/Historic_Buildings/Historic_Buildings/FeatureServer&f=json - server replied: Bad Request

This appears to be because QGIS is trying to add f=json as a second parameter to the query string (preceded by &). If I instead make the request in a web browser (using the referrer header above) but replace the & with ? the request is successful.
Question: Is it possible (using the ArcGisFeatureServer tool or otherwise) to connect to this server and add the layer it contains to a QGIS project?
I have no expertise in Azure-hosted web applications and this is the first time I've seen a proxied url in a query string like this treated as a parameter. It appears QGIS cannot correctly deal with this scenario.


texts: Trying to update separate table with SearchCursor

texts:
Learning Python to use in ArcGis 10.2.

I want to get value from a certain field from one attribute table into another attribute table. IN PT1, If Meas of Object1 > Meas of Object2, then  in SPL1 End = Punkt_ID of Object1 and Start = Punkt_ID of Object2

Like this:

So I have written this bit (based on existing scripts) where individual segments are used to select points in order to update the Split_Lines attribute table

## Update field End with max value
OIDField = arcpy.Describe(SplitLines).OIDFieldName
Lines = arcpy.UpdateCursor(SplitLines)
for Line in Lines:
    LineId = OIDField
    Where_clause = '"' + OIDField + '"=' +str(LineId)
    arcpy.MakeFeatureLayer_management (SplitLines, "Current_Line", Where_clause)
    arcpy.SelectLayerByLocation_management ("Points_Lyr","WITHIN_A_DISTANCE","Current_Line", "0.02 METERS")
    pcur = arcpy.SearchCursor("Points_Lyr")
    prow = pcur.next()
    p1val = prow.getValue("MEAS")
    p1HH = prow.getValue("HH")
    p1ID = prow.getValue("Punkt_ID")
    prow = pcur.next()
    p2val = prow.getValue("MEAS")
    p2HH = prow.getValue("HH")
    p2ID = prow.getValue("Punkt_ID")

if p2val > p1val:
    Line.setValue("End", p2ID)
    Line.setValue("Start", p1ID)
    Line.setValue("HH_Stromab", p1HH)
    Line.setValue("HH_Stromauf", p2HH)

else:
    Line.setValue("End", p1ID)
    Line.setValue("Start", p2ID)
    Line.setValue("HH_Stromab", p2HH)
    Line.setValue("HH_Stromauf", p1HH)
Lines.updateRow(Line)
del prow
del pcur

I am almost there, only the Search Cursor on the Selected Points doesn't seem to update for each new selection and copy the Value of Punkt_ID =1 and Punkt_ID =2.

Like this:

Where is the mistake?


texts: Cannot open Esri Geodatabase files on QGIS 3.12 and 3.13 on MacOS Catalina

texts:
I've been struggling for the last two days to get QGIS 3.12 running on macOS Catalina to open ESRI Geodatabase (GDB) files. I've tried using OpenFileGDB with no luck - I just get the following error "Invalid Data Source: /Users/xxxx/Documents/GIS/1. Geology/Australia/National/1M Geology - SHP/geodatabase/geology_1M.gdb is not a valid or recognized data source." I've tried with multiple different GDB files with the same error message.
I've also tried Mark Egge's solution for Mojave Installing File Geodatabase (*.gdb) support in QGIS? with no success - Steps 1-3 work ok but when I run ogrinfo --formats FileGDB doesn't appear.
Any thoughts?


texts: How to get an attribute of the selected feature and pass it into a URL?

texts:
As a total Novice with Python, I thought 2 weeks would be enough time to learn how to build some tools for our GIS...  I was wrong and I'm running up on a deadline.

I am attempting to get a value from a selected feature in geodatabase, and pass that value into a URL link / query after the equals sign.

Opening the link to the image viewer works, if I don't try to pass in the field value.

After 9 days of reading training manuals my eyes are raw. How can I fix this...

import arcpy
import pythonaddins
import webbrowser
from threading import Thread

def GetParameterAsText (str.Unit_No):
"This prints a passed through string into this function"

import str
row.getValue(Unit_No)

def OpenBrowserURL ():
    url = 'http://fcweb/caviewer#unit='
    webbrowser.open(url,new=2)

class CAVButton(object):
"""Implementation for PythonButton_addin.CAVbutton (Button)"""
def __init__(self):
    self.enabled = True
    self.checked = False
def onClick(self):
        t = Thread(target=OpenBrowserURL)
        t.start()
        t.join()