texts: How to see what the equal interval breaks are in order to add them to the CartoDB Legend

texts: I am using median income for Oakland, CA, where each polygon represents a block group with a number associated with it for its median income. I've displayed the data as a choropleth and chosen equal interval for the data. How can I see what those intervals are in order to add them to the legend?


texts: How to zoom to selected Feature in ArcMap using ArcObjects?

texts:
By using the QueryFilter I am getting some features based on the query and I am placing that feature names in a listbox.
When the user selects a particular name in the listbox that feature has to be zoomed to in the map.

Below is my code:

IQueryFilter pQueryFilter = new QueryFilterClass();
            pQueryFilter.WhereClause = textBox1.Text + textBox2.Text;
            ICursor pCursor = pFeatureClass.Search(pQueryFilter, true) as ICursor;
            IRow pRow = pCursor.NextRow();
while (pRow != null)
            {
                    listBox2.Items.Add(pRow.get_Value(9).ToString());
                pRow = pCursor.NextRow();

            }

Here is what I tried:

IActiveView pActiveView = pMxDocument.ActiveView;
            IGeometry pGeometry = pFeature.Shape;
            IEnvelope pEnvelope = pGeometry.Envelope;

            pActiveView.Extent = pEnvelope.Expand(1.2, 1.2, true);

But I am unable to zoom.


texts: How to merge spatially separated polygons into one feature in QGIS?

texts:
As you can see below I have 11 polygons in one shapefile, so in the attribute table there are 11 features.

As they all have the same attribute, is it possible to merge them into one feature in QGIS?


texts: OpenLayers bbox strategy does not work when zoom in

texts:
I am using OpenLayers 6 vector source and bbox loading strategy.
this.source = new VectorSource({
  format: new EsriJSON(),
  strategy: bbox,
  loader: this.loader(layer)
});

So, when I open tha app first, it appears in 14 zoom level. The data is loading at first. But when I zoom in, the data does not refresh? My ArcGIS Server serves maximum 1000 feature in an extent. When I zoom in, the other unloaded data should be load and map data should be refresh. But does not work. Why?


texts: How to Lock Grid and Scalebar in QGIS Atlas

texts:
Have just started using Valmeira and I am looking at generating an Atlas.

I am using irregularly sized features as the index.

The challenge I have, is when I set up the template page and scalebar and grid, if I go to another page (with a different sized feature) the resulting scale change makes the scalebar push out, and the grid become too small too read.

There are two workarounds -

Use regularly sized index features - then Grid and Scalebar size are not an issue, as the scale doesn't change between pages.
Group the Features by their approximate scale to reduce the amount of changes in scale, and break the project into a number of atlases.

My preferred option would be if there were a 'Disconnect Object from Atlas' style of option. This would mean that I could set the Scalebar and Grid for that page, then disconnect them from the Atlas automaitc refresh once I have completed them.


texts: Running 'r.null' processing algorithm in Standalone Python script on Windows

texts:
I want to run this code outside QGIS as a standalone Python script.
import os
import processing

os.chdir('C:/Users/User/test/25')
alg_params = { '-c' : False,
'-f' : False,
'-i' : False,
'-n' : False,
'-r' : False,
'GRASS_RASTER_FORMAT_META' : '',
'GRASS_RASTER_FORMAT_OPT' : '',
'GRASS_REGION_CELLSIZE_PARAMETER' : 0,
'GRASS_REGION_PARAMETER' : None,
'map' : 'test.tif',
'null' : None,
'output' : 'Null/out.tif',
'setnull' : '0' }
processing.run('grass7:r.null', alg_params)

I searched for some solutions but got for Linux. I have QGIS 3.16. Please guide how to run it on Windows?


texts: How to change the WMS Background Color in OpenLayers?

texts:
I am using geoserver to serve up an image pyramid basemap.  I want to have all areas outside of the data extent be black (instead of white).

I have the map displayed inside of a div tag and changing the css on that tag or the container (another div) it resides in has no effect.  I have checked both the geoserver and openlayers documentation, but have not seen any way to change the background color.

Any suggestions on where to continue searching?  Do I need to generate a black tiles and serve them?


texts: Preventing map panning to last drawn marker using Leaflet?

texts:
I have a basic Leaflet map where I need to draw several markers and lines on.
How do I prevent the map from panning to the last drawn object?

I've read about setting autopan : False for popups from here.

Can the autopan : false option work for leaflet markers?


texts: Render Google Fusion Tables Layer Tiles in OpenLayers

texts:
I would like to know if there is a way now to include Google Fusion Table Layers Tiles in OpenLayers, and if yes, is it possible to use the dynamic styling just like in Google Maps?

I know you can use Vizualisation or FusionTableAPI to render Fusion Table Layers in geoJSON, but tiling is much more performant.

This example is it valid? I dont understand how it mixes Fusion Table Layers Tiles and OpenLayers.


texts: Which type and range has the component M in a LineStringXYZM for ESRI-Shapefile and the PostGIS context

texts:
I want to compress a bunch of sequentially recorded positions with timestamps into a linestring with XYZ an M dimension (LineStringXYZM). Each set of points is represented by a constant part (ship data from AIS like MMSI, length, width, type) and a dynamic part (x y z timestamp). I want to compress the repeated static info by creating a trajectory like LineString and the additional static attributes. The timestamps have a granularity of milliseconds (YYYY-MM-DD.HH:MM:SS.FFF).
I want to know of which datatype M is (int32, int64 or float32, float32) in the Esri-Shapefile and the PostGIS context to feed in that info into the geo-object.
Example operations posgresql:
INTEGER
 SELECT mmsi, lon, lat,
        extract(epoch FROM (dt||' '||tm)::timestamp(3))*1000 AS ts
 FROM   ais_pos
 WHERE  year = 2018 AND month = 1
 ORDER  BY mmsi, ts LIMIT 4;

   mmsi    |   lon    |    lat    |      ts
-----------+----------+-----------+---------------
 100000000 | 8.173233 | 53.573867 | 1514792874239
 100000000 | 8.173200 | 53.574033 | 1514792882455
 100000000 | 8.173183 | 53.574067 | 1514792884476
 100000000 | 8.173133 | 53.574283 | 1514792895864

FLOAT
 SELECT mmsi, lon, lat,
        extract(epoch FROM (dt||' '||tm)::timestamp(3)) AS ts
 FROM   ais_pos
 WHERE  year = 2018 AND month = 1
 ORDER  BY mmsi, ts LIMIT 4;

   mmsi    |   lon    |    lat    |       ts
-----------+----------+-----------+----------------
 100000000 | 8.173233 | 53.573867 | 1514792874.239
 100000000 | 8.173200 | 53.574033 | 1514792882.455
 100000000 | 8.173183 | 53.574067 | 1514792884.476
 100000000 | 8.173133 | 53.574283 | 1514792895.864

Which range has the measured value M (ordinal type/ float type) to prevent overflow errors (..see example).
 #!/usr/bin/perl
 # -- test overflow

 # sample from the sql query
 my $x = 1514792874239;

 # max unsigned 32 bit int
 my $y = 2**32-1;

 # make string
 my $s = $x>$y ? 'overflow' : 'in.range';

 # show result
 print "$x  $y $s \n"

 ! 1514792874239  4294967295 overflow


texts: Extracting Values to Points without ArcGIS Spatial Analyst?

texts: I have a raster DEM & multiple points.  I want to determine the DEM elevation at each point.  Other than manually examining the DEM at each point, (very time consuming with lots of points) is there any way to do this without access to the Spatial Analyst Extract Values to Points tool?  I do have access to 3d Analyst if that helps.


texts: How to add points (marker) in Mapbox real time geolocation with Cassandra database?

texts:
I am now using Angular in web and ionic in mobile with spring-boot for communicate between web and mobile and I use Cassandra for database.
I get longitude and latitude of mobile geolocation in real time, but I want put the points that I got in real time on the Mapbox in order to see the movements of the phone in my web page.
  var url = 'http://localhost:8081/api';
  map.on('load', function() {

    var iv=0;

    var request = new XMLHttpRequest();
    window.setInterval(function() {

    request.open('GET', url, true);
    request.onload = function() {
    if (this.status >= 200 && this.status < 400) {
    // retrieve the JSON from the response
    var json = JSON.parse(this.response);
    console.log(json);
    //// update the drone symbol's location on the map

  if(iv<json.length)
    {
      console.log(iv);

     var marker=new mapboxgl.Marker({color: 'red'});
     marker.setLngLat([json[iv].longitude_x, json[iv].laltitude_y])
     .setPopup(new mapboxgl.Popup().setHTML("<h1>"+json[iv].name+"</h1>")).addTo(map);
    }

  iv++;

 }

When i write this code the markers if seting this markers over some in map !
I want help, please !


texts: How to draw weighted voronoi polygons (thiessen) from points in QGIS?

texts:
I have a shapefile of polygons which have an amount of population, I need to divide each of these polygons into small ones with nearly equal population. So I thought attributing points on each area inside each polygon and attributing a weight depending on its density.

After that I have to draw voronoi polygons based on the points drawn before and their weights to get these small areas that should be equally populated.

I tried to draw them in qgis but all I get is polygons with equal weights.

How to draw weighted voronoi polygons in QGIS?


texts:
I have a bunch of Landsat 5 scenes that I want to mosaic together.

However Landsat 5 scenes have a zigzag pattern on the left and right edges as well as some bad data values along these edges.

Is there a program that can automatically trim the edges of the Landsat 5 scenes so that all the edges are nice and straight?


texts: Inputting formatted concatenate in field of attribute using ArcGIS Desktop?

texts:
I am new to working in ArcGIS Desktop.

I have feature point data with 3 fields in table attribute (objectID, shape and sequence).

I want to add another field lets say called "name_sequence)..

in the name_sequence field i want to add text which have format "x+yyy", ("x", "+", "yyy")

in loop every 100, start with 0+000..

so it become like this

0+000, 0+100, 0+200, .... 0+900

the point after +900 it will add "x+1" so in data no 11 it'll show 1+000 after that it will have same sequence as i mention above..
furthermore, after i searching the code i found we can use range function and concatenate function to do that..
so i come with the code like this in field calculator.. using VB

 for i = 1 to 10
     a = i
     b = i*100
     srt_join = concat(srt(a),"+", srt(b))
 end for

i know it might be totally wrong code or it is not complete yet.

so i want to know what should i add for this code?

or do you have another way to do this?

and i want know new understanding about coding in arcgis

for your understanding i attach the result of my data which i working in excel.


texts:
I have the geometry of the outline region as part of one shapefile (which is only one feature and there are other features in it apart from this). Other polygon regions' boundaries which are colored are present in another shapefile.
Is there a way to find out the geometry of non-coloured regions (altogether as one polygon or multiple polygons) of the image?
I tried with the code below but it was giving me the geometry of the entire region.

from shapely.geometry import shape, mapping, Polygon
import matplotlib.pyplot as plt
import geopandas as gpd
import fiona

schema = {'geometry': 'Polygon','properties': {'test': 'float'}}

outline_shape = fiona.open(shapefile1)
region_shape = fiona.open(shapefile2)

for feature in outline_shape:
    if feature['properties']['Name'] == 'Required field':
        schema = {'geometry': 'Polygon', 'properties': {'test': 'int'}}
        with fiona.open('diff.shp', 'w', 'ESRI Shapefile', schema) as e:
            for geom in [shape(feature['geometry']).difference(shape(j['geometry'])) for j in region_shape]:
                if not geom.is_empty:
                    e.write({'geometry': mapping(geom), 'properties': {'test': 1}})


texts:
There is a way to do a SUMIF, from EXCEL, in Model Builder? I tryed something like this in field calculator, but i can't iterate through an attribute from the attribute table.

def AreaTotal(A,B):
    for i in range(0,len(A)):
        aux=0
        for n in range(0,len(A)):
            if A[i]==A[n]:
                aux= aux+B[n]
            else:
                continue
        return(aux)


texts: Can i take the wms url from the user interface

texts:
If user wants to overlay his own wms layer overlay on the basemap of our portal, then how can user pass that url to openlayers method like below code enter code here

        overlayGroup.getLayers().push(new ol.layer.Tile({
        title: 'Africa Major Crop Lands',
        source: new ol.source.TileWMS({
        url: '',
        params: {'LAYERS': ''},
        serverType: 'geoserver'
        })
        }));
        overlayGroup.getLayers().push(new ol.layer.Tile({
        title: 'South Asia Crop Intensity',
        source: new ol.source.TileWMS({
        url: '',
        params: {'LAYERS': ''},
        serverType: 'geoserver'
        })
        }));

For this code I want to pass url and LAYERS from user side like below


texts:
I have multiple points and some of these points are sitting directly on a line which I do not want. You can see the image below. What I need is "to identify the points sitting on the line" I tried using the selection by location (used both intersect, touches boundary and contain) but no result was returned.
Can someone give me and idea of how to identify these points without having to manually search through over 200000 points and several lines?


texts: Breaking apart a polygon layer by unique values, creating many new layers

texts:
I have a large dataset of habitat information that I need to calculate distances from animal locations to the closest polygon of each habitat type. The way to do this, that I understand, would be to create Euclidean distance rasters for each habitat type. But in order to do this I need to break apart the original layer into many different parts. I can do this by "Select by attribute" and create a new layer for each one, but I was wondering if there was a tool that automates that and breaks apart a layer into multiple smaller layers by unique values in a column. If not, oh well, I'll get down to selectin' and creatin'.

Thanks in advance.

Ayden


texts: Q: Shapefile for NZ electorates excluding marine areas

texts:
I'm a complete newbie to this GIS stuff, but I've attempted some due diligence as follows -

It being election season I'm after shapefiles for the NZ electorates, however I'm just after a map of the land areas, not one that includes the marine margins (12 mile limit?), as the latter don't really look like the country, and how many fish vote anyway.

Related stuff I've found would include

Various administrative areas, but not general electorates, at http://www.diva-gis.org/gdata
General election shapes including marine areas, for example at http://m.stats.govt.nz/browse_for_stats/Maps_and_geography/Geographic-areas/digital-boundary-files.aspx
That same page (and elsewhere) has NZ coastline maps without electoral areas.

In theory I could clip (2) against (3), but after looking around a bit I couldn't figure out how to do that with the likes of ogr2ogr.

My output target is GeoJSON in d3.


texts: Usage of index in R raster setValues function

texts:
I have a RasterStack object and want to replace some values at a set of indices in all layers. That doesn't work as expected, however. setValues replaces all values even though I pass a set of indices. The following code illustrates my problem:
library (raster)

# Create sample object
r1 <- raster (nrows = 10, ncols = 10)
r2 <- raster (nrows = 10, ncols = 10)
r3 <- raster (nrows = 10, ncols = 10)
values (r1) <- runif (ncell (r1))
values (r2) <- runif (ncell (r2))
values (r3) <- runif (ncell (r3))
rst <- stack (r1, r2, r3)

# Pixel IDs to replace
indexvalues <- c (1, 2, 99)

# Values to substitute
newvalues <- c (10, -10, NA)

result <- setValues (rst, newvalues, index = indexvalues)

# All values are replaced, not just the ones at the indices
plot (result)

### This is what I would expect
values (r1) [indexvalues] <- newvalues
values (r2) [indexvalues] <- newvalues
values (r3) [indexvalues] <- newvalues
rst <- stack (r1, r2, r3)

plot (rst)

Any ideas how I can get setValues to produce my expected result or is there something else wrong? Thanks!


texts: Geometry vs geography vs spheroid for polygon-contains-point (PostgreSQL)

texts:
Many operations behave differently in geometry vs geography, e.g. the distance between two points.
But what if I only need to check whether a polygon contains a point?

Does that operation differ between geometry vs geography?
Does that operation differ between use_spheroid being false vs true?


texts: Flow accumulation map giving zero values

texts:
I have applied the flow accumulation tool in arcgis 10.5.

The results map has zero values in places that it should not, which makes the further processing of my images wrong.

Any idea why this happens?

If I set a threshold values before calculating, the streams will become smaller and I dont want that either.


texts:
I'm using OpenGeo Suite 3.1.I have stored large number of shapefiles(100) in a store()Directory of spatial files.

How to Publish all the files at a time in geoserver?


texts: Changing thickness of lines in GraduatedColorsSymbology using arcpy

texts:
I have been using the following code (works) to update the symbology of a layer in ArcGIS. Now, in addition to this, i need to change the thickness of the lines previously symbolized by using GraduatedColorsSymbology.

mxd = arcpy.mapping.MapDocument(mxdfile)
df = arcpy.mapping.ListDataFrames(mxd)[0]
lyr = arcpy.mapping.ListLayers(mxd, "Plot", df)[0]
lyr.visible = True
lyr.transparency = 0
sym = lyr.symbology
sym.valueField = "GrossTons"
sym.classBreakValues = [0, 1000000, 5000000, 10000000, 20000000, 40000000, 60000000, 600000000]
sym.classBreakLabels = ["<1M", "1M to 5M", "5M to 10M", "10M to 20M", "20M to 40M", "40M to 60M",">60M"]
arcpy.mapping.ExportToJPEG(mxd, "Export.jpg", "PAGE_LAYOUT", resolution=800)
del mxd

But when I try to add the color and thickness part, it gives an error saying 'GraduatedColorsSymbology' object has no attribute renderer.

Code added:

sym = lyr.symbology
sym.renderer.breakCount = 7

breakVal = 100000000000
cv = 0
lw = 1
for brk in sym.renderer.classBreaks:
  brk.upperBound = breakVal
  brk.label = "\u2264" + str(locale.format("%d", breakVal, grouping=True))
  brk.description = "Description " + str(cv)
  brk.symbol.color = {'HSV' : [cv, 100, 100, 100]}
  brk.symbol.outlineColor = {'HSV' : [240-cv, 100, 100, 100]}
  brk.symbol.size = lw

  breakVal +=100000000000
  cv += 40
  lw += 0.5
lyr.symbology = sym

I didnt find any reference to this anywhere else. Is there a way to change the thickness and colors of lines in GraduatedColorsSymbology? Please avoid using

arcpy.ApplySymbologyFromLayer_management(lyr, symbologyLayer)

because this will not allow me to be flexible with changing class breakvalues.


texts: ArcGIS JS API SimpleMarkerSymbol is not a subclass of Layer

texts:
I have polar ArcGIS map and want to display my own data on it. I will get data using AJAX REST queries to my server. As I understand, it isn't commonly used way to get data in ArcGIS, but presenting data in appropriate format for ArcGIS and loading it on ArcGIS server will be more difficult task and it may be not free.

I have the following code for my map:

require([
"esri/Map",
"esri/graphic",
"esri/symbols/SimpleMarkerSymbol",
"esri/views/MapView",
"esri/layers/TileLayer",
"esri/layers/FeatureLayer",
"esri/geometry/Point",
"dojo/domReady!"
],
function (Map, MapView, TileLayer, FeatureLayer, Graphic, Point, SimpleMarkerSymbol) {
    var map = new Map({
        layers: [
            new TileLayer({
                url: "https://gis.ngdc.noaa.gov/arcgis/rest/services/arctic_ps/arctic_basemap/MapServer"
            })
        ]
    });

    var point = {
        type: "point",
        longitude: 77.962148,
        latitude: 24.451656
    };
    var markerSymbol = {
        type: "simple-marker",
        color: [226, 119, 40]
    };

    map.on("load", function () {
        map.graphics.add(new Graphic({
            geometry: point,
            symbol: markerSymbol
        }));
    });

    var viewSpatialReference = 3995;
    var centerPoint = {
        x: 1000000,
        y: -1000000,
        spatialReference: viewSpatialReference
    };

    var viewOptions = {
        container: "map",
        map: map,
        viewingMode: "local",
        spatialReference: viewSpatialReference,
        scale: 12000000,
        center: centerPoint
    };

    var view = new MapView(viewOptions);
});

When I load map, I get error in console:

  [esri.core.Accessor] Accessor#set Assigning an instance of 'esri.symbols.SimpleMarkerSymbol' which is not a subclass of 'esri.layers.Layer'

The strangest thing is Google doesn't find anything relative with this error, that make sense, only link to SimpleMarkerSymbol documentation and so on.

Hope someone is familiar with this problem. If there is another way to display markers with certain lat long, I can change my approach.


texts: A couple of days ago, I went to this site heron-mc.org/examples.html and found some useful application templates. I viewed them, went to view source code and copied the source code to my editor to make it my code. So far the code is displaying the map and a couple of features but some aren't working as they are supposed to like on the web template. To be specific, the toolbar is working fine but I can't see the button's description(icon and text) and the sliding zoom is not complete. I'm asking for an advice on how to go about fixing my problem. I have also attached a picture of the output.


texts: How to avoid the pop-up callout in the toolbar in ArcGIS 10.3.1 for Desktop

texts:
I try to find a way to avoid from seen the pop up callout any time i stand on each of the toolbar functions.I mean this such of callout for example:

I try to change the default in the file "AdvancedArcMapSettings", but i didn't found any solution.


texts: How to create Shapefile on ArcGIS Server side?

texts:
I am developing a script which gets an JSON array with point coordinates. This array is encoded to a python dictionary. From this I filter all points (latitude and longitude) and an ID. This coordinates are converted into arcpy-points by using arcpy.Point(). Now in a next step I want to create a shapefile of these points by using
Code:

arcpy.CopyFeatures_management(ptList, r"pointTemp.shp")

My question is how to store this shapefile into a folder which is located on ArcGIS Server side. Do I have to give the skript a special path or location?

One folder that is on ArcGIS Server side in my case is "C:\temp", I have full read and wright rights, but the script doesn´t create the shapefile in this folder. If I run the script local on my pc, all works fine and the shape is created and stored locally on "C:\temp".

Any suggestions or tips?


texts: Why does ArcGIS map document (MXD) appear empty of data?

texts:
A colleague created a map for me using ArcGIS for Desktop and emailed it to me in another state. I was able to open the document using ArcGIS for Desktop but no map appears, why?

She also sent it to her boss using data from a central fileserver and again the map didn't appear why?


texts:
I want use the Python colorama module to add colour to the output text of my scripts in ArcMap so that they are easier to read. The problem I have is that while colorama works at command line, it does not work within the Python Window of ArcMap.
I don't get an error message but ArcMap simply ignore the ANSI sequences. Here is the output within ArcMap:

And here is the output at the command line:

Are there any alternatives to using colorama?


texts: What qualifications should auditor have to auditing GIS?

texts:
Background summary: Various systems are being audited and GIS is the next to be audited.  The initial kick-off for the audit was a meeting between me and a single auditor.  Initially, I thought it was going to be the typical security audit that is done for all RDBMS, which is fine.  That's good.  But then she said they would be executing a series of tests against our data (deletes, updates, inserts).

The geodatabase is versioned, so naturally I'm concerned.  When I received the meeting minutes from the auditor, I had to re-write it because the vocabulary was all wrong (lots of GIS words used in the wrong way/wrong order).

I'm not concerned about failing the audit--I'm worried about the auditors really mucking things up because they have no GIS experience.  Does anyone have experience with GIS audits and can provide some guidelines as to what qualifications the auditor (or audit team) should have? Do internal audit departments typically have GIS professionals on staff?

I added the esri-geodatabase and arcgis tags to this question because of the versioning--I might consider dropping those once a read answers.


texts: How to manage JSON data in a Leaflet to control layer

texts:
I want to add this json data connection to the map control.

This is my JSON data:

[
  {
    "id": "1234",
    "id_jenis": "15",
    "nama": "",
    "cabang": "01",
    "no_sp": "",
    "pelapor": "HERNI BUDIYATI",
    "alamat": "JL.SAMARINDA",
    "status_proses": null,
    "gps_lat": "-8.15901548898326",
    "gps_log": "114.33938219636661",
    "status": "0",
    "petugas_1": null,
    "tindakan": null,
    "tgl_tindakan": null,
    "keterangan": "PASANG BARU"
  },
  {
    "id": "1253",
    "id_jenis": "15",
    "nama": "",
    "cabang": "01",
    "no_sp": "",
    "pelapor": "AHMAD SULAIMAN",
    "alamat": "JL.BATUR",
    "status_proses": null,
    "gps_lat": "-8.20307616134188",
    "gps_log": "114.37318817745438",
    "status": "0",
    "petugas_1": null,
    "tindakan": null,
    "tgl_tindakan": null,
    "keterangan": "PASANG BARU"
  },
  {
    "id": "1264",
    "id_jenis": "15",
    "nama": "",
    "cabang": "01",
    "no_sp": "",
    "pelapor": "MASTIKAH",
    "alamat": "JL. KUTILANG",
    "status_proses": null,
    "gps_lat": "-8.244054483129496",
    "gps_log": "114.35282327992832",
    "status": "0",
    "petugas_1": null,
    "tindakan": null,
    "tgl_tindakan": null,
    "keterangan": "PASANG BARU"
  }
]

and this is my HTML script:

init();
function init() {
url2 = "urlfromjson";
var promise = $.getJSON(url2);
    promise.then(function(data) {
        var ob = L.layerGroup(data);
        pb = ob._layers;
        for (var i in pb){

            var lat= pb[i].gps_lat;
            var lng= pb[i].gps_log;
            var ps = L.marker([lat, lng]).addTo(map);;
        }
    });
}


texts: How to convert ECW to PNG with gdal_translate while preserving colors

texts:
I am trying to convert some ECW files to PNG using the command gdal_translate test.ecw -of png test.png, but the PNG file comes out grayscale.
How do I produce a colored PNG from an ECW file using this tool?
Output from gdalinfo test.ecw:
Driver: ECW/ERDAS Compressed Wavelets (SDK 3.x)
Files: test.ecw
Size is 16000, 24000
Coordinate System is:
PROJCRS["ETRS89 / UTM zone 32N",
    BASEGEOGCRS["ETRS89",
        DATUM["European Terrestrial Reference System 1989",
            ELLIPSOID["GRS 1980",6378137,298.257222101,
                LENGTHUNIT["metre",1]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4258]],
    CONVERSION["UTM zone 32N",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",9,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",0.9996,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",500000,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["easting",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["northing",north,
            ORDER[2],
            LENGTHUNIT["metre",1]],
    ID["EPSG",25832]]
Data axis to CRS axis mapping: 1,2
Origin = (448000.000000000000000,6283000.000000000000000)
Pixel Size = (0.125000000000000,-0.125000000000000)
Metadata:
  COLORSPACE=MULTIBAND
  COMPRESSION_RATE_TARGET=4
  VERSION=2
Corner Coordinates:
Upper Left  (  448000.000, 6283000.000) (  8d 9' 3.88"E, 56d41'17.50"N)
Lower Left  (  448000.000, 6280000.000) (  8d 9' 6.06"E, 56d39'40.48"N)
Upper Right (  450000.000, 6283000.000) (  8d11' 1.41"E, 56d41'18.29"N)
Lower Right (  450000.000, 6280000.000) (  8d11' 3.51"E, 56d39'41.27"N)
Center      (  449000.000, 6281500.000) (  8d10' 3.71"E, 56d40'29.39"N)
Band 1 Block=256x256 Type=Byte, ColorInterp=Undefined
  Description = Band #0
  Overviews: 8000x12000, 4000x6000, 2000x3000, 1000x1500, 500x750, 250x375
Band 2 Block=256x256 Type=Byte, ColorInterp=Undefined
  Description = Band #1
  Overviews: 8000x12000, 4000x6000, 2000x3000, 1000x1500, 500x750, 250x375
Band 3 Block=256x256 Type=Byte, ColorInterp=Undefined
  Description = Band #2
  Overviews: 8000x12000, 4000x6000, 2000x3000, 1000x1500, 500x750, 250x375
Band 4 Block=256x256 Type=Byte, ColorInterp=Undefined
  Description = Band #3
  Overviews: 8000x12000, 4000x6000, 2000x3000, 1000x1500, 500x750, 250x375
Band 5 Block=256x256 Type=Byte, ColorInterp=Undefined
  Description = Band #4
  Overviews: 8000x12000, 4000x6000, 2000x3000, 1000x1500, 500x750, 250x375


texts: If I have polygon data in shapefile, can I discretize the area into equidistant points?

texts:
If I have polygon data in shapefile, can I discretize the area into points?

e.g. a 3x3 polygon into 9 points. I tried MMQGIS - create grid layer and extract nodes, they produces messy things.

I use Python and QGIS


texts: Doing change detection of land cover between 2 landsat images in google earth engine?

texts:
I am doing land use/ land cover classification of 2000 and 2018 year . I had completed supervised classification now i wanted to do change detection analysis. Please guide me how to do change detection for 2 years .
Link to my code :https://code.earthengine.google.com/7627f681e4d50cb9cd8027d3ff34cf54

var image = ee.Image(ee.ImageCollection('LANDSAT/LC08/C01/T1_SR')    .filterBounds(geometry)    .filterDate('2018-01-01', '2019-1-31')    .sort('CLOUD_COVER')    .first());

print( image);

var Trivandrum = ee.Image ('LANDSAT/LC08/C01/T1_SR/LC08_144054_20180206');
var Trivandrum1 = ee.Image ('LANDSAT/LC08/C01/T1_SR/LC08_143054_20180911');

var mosaic = ee.ImageCollection.fromImages([Trivandrum,Trivandrum1]).mosaic();
Map.addLayer(mosaic, {bands: ['B4', 'B3', 'B2'], max: 0.3}, 'mosaic');
Map.addLayer(table);

var newfc = waterbodies.merge(forest).merge(settlements).merge(sandyarea).merge(agriculture_cropland).merge(agriculture_plantation).merge(fallowland);
print(newfc, 'newfc');

// Select the bands to be used in training
var bands = ['B2', 'B3', 'B4', 'B5','B6','B7'];

// Sample the input imagery to get a FeatureCollection of training data.
var training = mosaic.select(bands).sampleRegions({
  collection: newfc,
  properties: ['landcover'],
  scale: 30  // should reflect the scale of your imagery, the scale is 30 meters
});
var trained = ee.Classifier.cart().train(training, 'landcover', bands);

// Classify the image with the same bands used for training.
var classified = mosaic.select(bands).classify(trained);

var palette = [
  '0000FF', // wb
  '0F8C3A', // f
  'FF0000', // settlements
  'A52A2A',   // sandy area
  'FFFF00', // vegetation
 ' 0B6607',// forest

 '665A07' // wastland
];

var composite = {
  bands: ['B2','B3','B4'],
  gamma: 1,
  max: 0.80,
  min: 0.04
};

var clipped = classified.clip(table);

Map.setCenter(76.98244, 8.4686, 11);

Map.addLayer (mosaic, composite);
Map.addLayer (classified, {min: 0, max: 5, palette: palette}, 'Land Use Classification');

Map.addLayer(clipped, {min: 0, max: 5, palette: palette}, 'Land Use ');

var geometry = ee.Geometry.Polygon(
          [[[76.2593, 8.8932],
          [76.4076, 8.0402],
          [77.314, 8.915],
          [77.4897, 8.2686]]]);

Export.image.toDrive({
  image: clipped,
  description: 'Trivandrum_mosaic_17',
  scale: 30,
  region: table,
  crs: 'EPSG:4326',
  folder: 'GEE',
  maxPixels: 10000000000000
});
var options = {
  lineWidth: 1,
  pointSize: 2,
  hAxis: {title: 'Classes'},
  vAxis: {title: 'Area m^2'},
  title: 'Area by class',
  series: {
    0: { color: 'blue'},
    1: { color: 'green'},
    2: { color: 'red'},
    3: {color : 'brown'},
    4 : {color: 'yellow'},
    5 : {color : 'orange'},

     6 : {color : 'grey'},
    7: {color : 'peach'}
  }
};

var areaChart = ui.Chart.image.byClass({
  image: ee.Image.pixelArea().addBands(clipped),
  classBand: 'classification',
  region: table,
  scale: 60,
  reducer: ee.Reducer.sum()
}).setOptions(options)
  .setSeriesNames(['wb', 'forest', 'settlements','sandyarea','crop','plantation','fallow']);
print(areaChart);

// Reduce the region. The region parameter is the Feature geometry.
var meanDictionary = clipped.reduceRegion({
  reducer: ee.Reducer.mean(),
  geometry: table,
  scale: 30,
  maxPixels: 1e9
});

// The result is a Dictionary.  Print it.
print(meanDictionary);


texts: Merging TIFF files creates offset / deformation

texts:
I have a PNG library (524288 tiles of 256px * 256px) which I georeferenced via GDAL.

This first part works correctly.

ullr:
ulx = (360 * (column / 1024)) - 180
uly = 90 - (180 * (line / 512))
lrx = ulx + 0.3515625
lry = uly - 0.3515625
ULLR = ulx + '' + uly + '' + lrx + '' + lry

gdal_translate -of VRT -a_srs EPSG:4326 -co NBITS=1 Path/Of/Input.png Path/Of/Output.vrt -a_ullr + ULLR
gdalwarp -s_srs EPSG:4326 -t_srs EPSG:3857 -co NBITS=1 Path/Of/Input.vrt Path/Of/Output.tiff

The second part works correctly too

I merge the tiles of a row (1024 TIFF tiles) using a VRT then from it I create the TIFF file
#!/bin/bash

for i in `seq 0 1 511`;
do
        gdalbuildvrt Path/Of/lines_vrt/map_$i.vrt Path/Of/tiff/map_11_*_$i.tiff
        gdal_translate -co NBITS=1 Path/Of/lines_vrt/map_$i.vrt Path/Of/lines_tiff/map_$i.tiff
done

Here if I use gdal_merge.py, either I get an invalid file (1.xMo while the lines are approx 60Mo) or I observe a offset / deformation ... (In QGIS)

The last part doesn't work properly (and I don't know why ...)

Here I am merging all the lines (TIFF files) into a global TIFF file using the same method as before. (with a list.txt file)
#!/bin/bash

gdalbuildvrt -input_file_list Path/Of/lines_tiff/list.txt Path/Of/global_vrt/map_global.vrt
gdal_translate -co NBITS=1 Path/Of/global_vrt/map_global.vrt Path/Of/map_global.tiff

and I observe a offset / deformation. (In QGIS again)
I wonder what am I doing that is wrong
My final goal is to use the global file (TIFF) with gdal2tiles.py
(Tiles overlay on GMaps / OSM)


texts:
I'm looking for some suggestions about how to make my python code more efficient.  Normally efficiency doesn't matter for me but I am now working with a text file of US locations with over 1.5 million points. With the given setup it is taking about 5 seconds to run operations on one point; I need to get this figure way down.

I'm using three different python GIS packages to do a few different operations on the points and output a new delimited text file.

I use OGR to read a county boundary shapefile and get access to the boundary geometry.
Shapely checks to see if a point is within any of these counties.
If it is within one, I use the Python Shapefile Library to pull attribute information from the boundary .dbf.
I then write some information from both sources to a text file.

I suspect that the inefficiency lies in having a 2-3 tiered loop... not quite sure what to do about it.  I'm particularly looking for help with someone experienced in using any of these 3 packages, as it is my first time to use any of them.

import os, csv
from shapely.geometry import Point
from shapely.geometry import Polygon
from shapely.wkb import loads
from osgeo import ogr
import shapefile

pointFile = "C:\\NSF_Stuff\\NLTK_Scripts\\Gazetteer_New\\NationalFile_20110404.txt"
shapeFolder = "C:\NSF_Stuff\NLTK_Scripts\Gazetteer_New"
#historicBounds = "C:\\NSF_Stuff\\NLTK_Scripts\\Gazetteer_New\\US_Counties_1860s_NAD"
historicBounds = "US_Counties_1860s_NAD"
writeFile = "C:\\NSF_Stuff\\NLTK_Scripts\\Gazetteer_New\\NewNational_Gazet.txt"

#opens the point file, reads it as a delimited file, skips the first line
openPoints = open(pointFile, "r")
reader = csv.reader(openPoints, delimiter="|")
reader.next()

#opens the write file
openWriteFile = open(writeFile, "w")

#uses Python Shapefile Library to read attributes from .dbf
sf = shapefile.Reader("C:\\NSF_Stuff\\NLTK_Scripts\\Gazetteer_New\\US_Counties_1860s_NAD.dbf")
records = sf.records()
print "Starting loop..."

#This will loop through the points in pointFile
for row in reader:
    print row
    shpIndex = 0
    pointX = row[10]
    pointY = row[9]
    thePoint = Point(float(pointX), float(pointY))
    #This section uses OGR to read the geometry of the shapefile
    openShape = ogr.Open((str(historicBounds) + ".shp"))
    layers = openShape.GetLayerByName(historicBounds)
    #This section loops through the geometries, determines if the point is in a polygon
    for element in layers:
        geom = loads(element.GetGeometryRef().ExportToWkb())
        if geom.geom_type == "Polygon":
            if thePoint.within(geom) == True:
                print "!!!!!!!!!!!!! Found a Point Within Historic !!!!!!!!!!!!"
                print str(row[1]) + ", " + str(row[2]) + ", " + str(row[5]) + " County, " + str(row[3])
                print records[shpIndex]
                openWriteFile.write((str(row[0]) + "|" + str(row[1]) + "|" + str(row[2]) + "|" + str(row[5]) + "|" + str(row[3]) + "|" + str(row[9]) + "|" + str(row[10]) + "|" + str(records[shpIndex][3]) + "|" + str(records[shpIndex][9]) + "|\n"))
        if geom.geom_type == "MultiPolygon":
            for pol in geom:
                if thePoint.within(pol) == True:
                    print "!!!!!!!!!!!!!!!!! Found a Point Within MultiPolygon !!!!!!!!!!!!!!"
                    print str(row[1]) + ", " + str(row[2]) + ", " + str(row[5]) + " County, " + str(row[3])
                    print records[shpIndex]
                    openWriteFile.write((str(row[0]) + "|" + str(row[1]) + "|" + str(row[2]) + "|" + str(row[5]) + "|" + str(row[3]) + "|" + str(row[9]) + "|" + str(row[10]) + "|" + str(records[shpIndex][3]) + "|" + str(records[shpIndex][9]) + "|\n"))
        shpIndex = shpIndex + 1
    print "finished checking point"
    openShape = None
    layers = None

pointFile.close()
writeFile.close()
print "Done"


texts: Event not firing for map.graphic.add objects added

texts:
Using ESRI ArcGIS javascript api v 3.16

I programmatically add a graphic to my map, which works fine, for example :

var pnt = map.graphics.add(new Graphic(evt.mapPoint, stopSymbol));

I need to fire an event handler when this pnt is drawn on the map.  I have tried map.on("graphic-add"), but it does not fire for above.

I've also tried to fire a reposition event with map.centerAt(pnt.geometry), but also nothing fires .

Any advice ?


texts: Using terrain correction in SNAP sentinel 1 toolbox removes all the values

texts:
I am new to SAR image processing, I have a SAR sentinel 1 dataset downloaded from ESA site: S1A_IW_GRDH_1SDV_20190314T030905_20190314T030930_026328_02F196_2EC3

I am trying to do terrain correction via SNAP sentinel 1 toolbox using Range Doppler Terrain Correction.

My Region of Interest (is a coastal area)

I try to perform Range Dopler Terrain Correction with the default parameters.

The result generated is a blank map where there are no values:

Am I missing something? I understand the SRTM data is use correct the terrain. However, since my data is near the coast and considering the accuracy limitations of the SRTM data is it removing the entire remote sensed values considering they all are below mean sea level?

Below are the parameters for Terrain Correction using Range Dopler Terrain Correction:


texts:
I'm trying to find PostGIS Raster functions supports the Inverse Distance Weighted interpolation method (IDW) and Kriging interpolation method.
I'm especially trying to generate a temperature map from a geometry table with point locations in PostGIS.

Can someone please point me in the right direction?


texts: How to use arcpy.setValue for Title casing?

texts:
I am trying to figure out how to use fieldName.title() method while updating rows using updateCursor for multiple fields. I was trying to improve the code of this post to fit my purpose:

arcpy.UpdateCursor for multiple fields

def ChangetoTitleCase(sourceTable, fieldNames):
      try:
        rows=arcpy.UpdateCursor(sourceTable)
        for row in rows:
          for fieldName in fieldNames:
            uField='"%s"' % fieldName
            row.setValue(uField, "uField.title()")
            rows.updateRow(row)
        del row, rows

      except Exception, e:
        import traceback, sys
        tb = sys.exc_info()[2]
        print "Line %i" % tb.tb_lineno
        print e.message

I am working with Python 2.6.5, ArcPy 10 and Win 7 64. Please let me know where I am wrong?


texts: Unable to connect ArcGIS Server - Error: Proxy server got bad address from remote server

texts:
This is my problem. anybody know how to solve this?

  We were unable to connect to ...

  Error: Proxy server got bad address from remote server (verify ther server is running).


texts: Sea Routes shipment Dataset Caribbean South America, Free

texts: I would like to obtain gis data of cargo, container vessels as well as small cargo inter-islands ships and for an educational non commercial project/publication on food trade across the Caribbean Sea. Any information available.


texts: Convert DEM elevation points from centimeters to meters

texts:
I would like to convert a large DEM stored in PostGIS which contains elevation values in centimeters to elevation in meters.

I found an easy way by updating each point like this:

UPDATE mydem a
SET rast = ST_SetValue(rast, x, y, val / 100)
FROM (SELECT rid, x, y, val FROM (SELECT rid, (ST_PixelAsPoints(rast)).* FROM mydem) t) r
WHERE a.rid = r.rid;

It works but very slowly, much too slow to process large tables. Is there is any faster way to achieve this ?


texts: Percentages of web using mapping softwares

texts:
How accurate are the online mapping assessments on trends.builtwith.com?
http://trends.builtwith.com/mapping

Can ArcGIS really be only .33% ?
If so, what does this say about ArcGIS.

Also why aren't QGIS Server, GeoServer, or MapServer listed ?


texts: Taking only the largest polygon from multi-polygon layer in QGIS

texts: I have multi-polygon geometries in my file, and I want to take only the largest polygon within the multi-polygon geometry to further work with it. I tried using Geometry Tool Multipart to Single part, which allowed me to have single polygons. Now I suppose I could filter for non-unique IDs and have QGIS calculate the largest out of its geometries. I am not sure how to do this, and any alternative methods could be beneficial as well!


texts: Any way to have variable min and max in one visualization parameter?

texts:
What I mean is: I'm using the "3 standard deviations" value to visualize my data, but this changes seasonally so I need to update what the values actually are for each image. Can I create a "3 standard deviations" palette that will automatically calculate these new values for each image?

In pseudocode terms, I'm imagining something along these lines:

var Oa04_palette = {
bands: ['Oa04_radiance'],
min: sminus3,
max: splus3,
sminus3: standard deviation of 'Oa04_radiance' * -3,
splus3: standard deviation of 'Oa04_radiance' * 3
};

So when I add a bunch of layers all using this palette, it will always use the appropriate min-max range for that layer.


texts: Change Map Title Block Information through ArcPy

texts:
I am new to python, but I'm currently working on a script tool to loop through all MXDs within a folder and change their associated title block information(Map Title, Report Title, Date, Charge Code, etc...) then export the newly MXDs to PDF.

I have no problem changing some of the static information that is the same for all MXDs such as the Date/Client/Report through arcpy.GetParameterAsText, but my real question is ...

Is there a way to allow for multiple titles?

For instance I know you can allow for multi-values and the output string looks something like "Exhibit 1 ; Exhibit 2; Exhibit 3 " Is there a way to take that string and split it into individual titles for each based off the order of the mxds in the folder?

Here is my current script :

    env.workspace = r"E:/Python/ENV_TEST"
clientTitle = arcpy.GetParameterAsText(0)
reportTitle = arcpy.GetParameterAsText(1)
dateTitle = arcpy.GetParameterAsText(2)
PNTitle = arcpy.GetParameterAsText(3)
titleList = arcpy.GetParameterAsText(4)

for mxdname in arcpy.ListFiles("*.mxd"):
    mxd = arcpy.mapping.MapDocument(r"E:\Python\ENV_TEST\\" + mxdname)
    pdfname = mxdname + ".pdf"
    for titleItem in arcpy.mapping.ListLayoutElements(mxd,"TEXT_ELEMENT", "Title"):
        titleItem.text = titleList
    for ReportItem in arcpy.mapping.ListLayoutElements(mxd,"TEXT_ELEMENT", "Report"):
        ReportItem.text = reportTitle
    for ClientItem in arcpy.mapping.ListLayoutElements(mxd,"TEXT_ELEMENT", "Client"):
        ClientItem.text = clientTitle
    for DateItem in arcpy.mapping.ListLayoutElements(mxd,"TEXT_ELEMENT", "Date"):
        DateItem.text = dateTitle
    for PNItem in arcpy.mapping.ListLayoutElements(mxd,"TEXT_ELEMENT", "PN"):
        PNItem.text = PNTitle
    arcpy.mapping.ExportToPDF(mxd,pdfname)
    mxd.save()


texts: Satellite image is cut, when the composite tool is run on it- ArcGIS

texts:
Hi I am using satellite imagery in my project.  The images come in squares, after using composite to join the bands, one square shrinks or is clipped somehow. Does anyone know the solution?  Thanks.
See the before and after pics.
Before:
After: The small green rectangle is what remains of the upper left square


texts: What is the equivalent of arcpy "Con" in QGIS and/or R raster-package?

texts:
I'm looking for an equivalent of the conditional expression "Con" from arcpy in an open source spatial analysis package. I'd prefer to use QGIS or R.

To provide some additional detail, I need to evaluate some moderately complicated conditional statements in map algebra, e.g., if the value for one raster is above a certain threshold, use one value, but if the value is below that threshold, use a value based on whether another raster value is above or below another threshold, e.g.:

newraster = Con(raster1 < 4, 0, Con(raster2 > 3, 0.5, 1))

Caveat: I have not been able to confirm that the above code will execute, since the machine I'm currently using does not have any ESRI products installed.


texts: unable to calculate geohash from a large geometry

texts:
I am trying to insert a large buffer in the city table.

    INSERT INTO city (geom) (
    SELECT
    st_buffer(ST_multi(ST_union(geom)),
    10.0) FROM santa_clara_county);

In the city table I have a trigger defined which calculates geohash. During the insert I got the following error.

    [XX000] ERROR: Geohash requires inputs in decimal degrees, got
    (6.06788e+06 1.78535e+06, 6.35495e+06 2.00181e+06).
    Where: PL/pgSQL function centroid_perimeter_area_trigger() line 6
    at assignment

The same error remains when I remove the st_buffer and try to insert st_multi
Does anyone know what cause the problem?


texts: How to raster calc on rasterio with two image not on the same band

texts:
I want to compute for an NDVI using rasterio. I must use two raster images that are adjacent to each other. They do not have the same sizes. They have a great overlap area which I can calculate NDVI on. One image is slightly shifted from the other. How to calculate raster that are not of the same sizes and does not overlap entirely?
Data: https://drive.google.com/drive/folders/1J3gmJVObv0LVoiUM_JhziD8ip19IB3l3?usp=sharing
EDIT 1: I have this to start with:
red_filename = "data/red.tif"
nir_filename = "data/nir.tif"

red_dataset = rasterio.open(red_filename)
nir_dataset = rasterio.open(nir_filename)

red_band = red_dataset.read(1)
nir_band = nir_dataset.read(1)

fig, (ax1, ax2) = pyplot.subplots(1,2, figsize=(10,12))
show((red_dataset, 1), ax=ax1, cmap='Reds', title='red_band')
show((nir_dataset, 1), ax=ax2, cmap='Blues', title='nir_band')


texts:
I have a DEM, geo tif, Float32 whith some gaps.
Pixel value relative to the gaps is indicated in the metadata as -9999, but if I click on these pixel with the identify value tool -9999.96 is shown.
I tried to fill these gaps whith the tool raster > analysis > fill no data, but it doesn't produce any result.


texts: Geoserver WMS and filtering tables with cyrillic column names

texts:
I try to serve PostGIS tables with Geoserver WMS (both 2.1 and 2.2). Some column names are in cyrillic.

Problem is - service fails, when I try to filter by those cyrillic columns (both with CQL and regular, XML-based filters). For example:

FILTER=<Filter><PropertyIsGreaterThan><PropertyName>Магазин</PropertyName><Literal>1</Literal></PropertyIsGreaterThan></Filter>

or

CQL_FILTER=Магазин=10

As result, I get

<ServiceException code="internalError">
    Rendering process failed
    Property &apos;??????&apos; could not be found in magasiin
</ServiceException>

or

<ServiceException>
    Could not parse CQL filter list. Lexical error at line 1, column 1.
    Encountered: &quot;\u0430&quot; (1072), after : &quot;&quot; Parsing : ??????=10.
</ServiceException>

I can't find clear answer from Google. Closest posting is http://jira.codehaus.org/browse/GEOS-4177 that claims, it's bug and it's fixed, but apparently it's not working in my case.

When I get confirmation, it's problem of Geoserver, then onliest solution seems to renaming of those columns.

Also, WFS service seems to fail always (despite of filtering), when at least on column name is in cyrillic.


texts: Adding Horizon and Sky to a Pokemon Go-like Map?

texts:
I am trying to create a map similar to the Pokemon GO map. I have used mapbox studio to create a sytle which I am pretty happy with.

However I am wondering how to get a horizon and a sky on this? Not sure where to start with this.


texts: Calculating geodesic distance along a path (lat/lon points) at once?


texts: Visualizing LiDAR data in a detailed footprint of rooftops with ArcGIS

texts:
For the last couple of weeks, I have been doing a research for making a detailed footprint with LiDAR data. The detailed footprint is needed for visualizing a 3D model. With a detailed footprint I mean a footprint with information of the rooftops within blocks of buildings. This information are different levels/heights, walls, holes etc. Some questions arose during this research. The process consists of two steps: 1. Classifying the LiDAR data and 2. Creating a detailed footprint. The first step is done. The last step, creating the detailed footprint, cause me a lot of problems. For this step I tried a lot of different tools and workflows but none of them gives me the result I want.

The Classified LiDAR data for the footprint consists of blocks with flat roofs and different heights/levels:

The most tools for creating a building footprint are only creating the outline of a block/house. A good workflow for the solution of my example is LAS Point Statistics as Raster, Raster to Polygon followed by Regularize Building Footprint. But, there is always a problem with the first step. When I run this tool, I’ve got the next result:

In the Image above you can see a raster of the rooftops from a block of buildings. I want that raster to be visualized as a footprint (polygon). The footprint contains information of the rooftops with different heights/levels and components(walls, for example) in a block. In the image above is the output raster faced North. This is wrong, because the footprint needs to be straight. This is not the case in the image above. The edges and information within the blocks are ribbed.

In this case, the raster output needs to be rotated by 37 degrees. When it is rotated like this, the blocks and information within the blocks are straight. The raster output should look like this:

So, the end result is a detailed footprint/polygon of the rooftops based on elevation. With detailed I mean all the levels of the blocks are different polygons within the polygon of the block, like this:

This can be done with a raster. How can I rotate the output of the raster to get a straight raster that matches with the blocks

Or, is there another way to create polygons based on heights. For example, in a block are five buildings with different levels/heights. How can I visualize these different levels/heights in different polygons?

Or, has somebody another suggestion to visualize a footprint the way I want?


texts: How to identify township range and section from longitude and latitude

texts: I have several hundred coordinates in an Excel file. I would like to identify the township range and section they are located--what is the easiest way to do this? I have you QGIS and basic GIS skills (i.e. digitizing). I googled my question and came across the following website: http://www.esg.montana.edu/gl/trs-data.html. This allows you to go from TRS to longitude and latitude. I'm thinking GIS would provide the best solution but I'm not sure where to start. I would like to identify actual legal descriptions.


texts: Reading composite types of PostgreSQL in GeoServer or MapServer?

texts:
I'm creating a data model with user defined types (Composite Types) https://www.postgresql.org/docs/10/rowtypes.html to implement INSPIRE EU Data Model. And I want to publish my data with Web Feature Services (WFS). I tried that with GeoServer and MapServer.

GeoServer: doesn't read the column. When I tried, I got the follow warning:

09 mar 15:17:23 WARN [geotools.jdbc] - Could not find mapping for 'inspireid', ignoring the column and setting the feature type read only

MapServer: With MapServer I can have a WFS with all my data but not structured. For example, the column INSPIREID has 3 sub columns and it shows like this:

<base:inspireid>(vaedificio_100,https://catalog.inspire.geoportal.pt/BU/Building,"(2020-02-24,)")</base:inspireid>

insted of this

<base:inspireid>
    <base:localid>vaedificio_100</base:localid>
    <base:namespace>https://catalog.inspire.geoportal.pt/BU/Building</base:namespace>
    <base:versionid>2020-02-24</base:versionid>
</base:inspireid>

Possible workarounds in both cases is to give the query directly but I lose the structured with that to. MapServer have a interesting feature, "gml_groups" but it doesn't work for complex hierarchies. And with MapServer I have little control over the nampespaces, witch is another problem.

Is there alternative or know to can I resolve this problem?


texts: ArcObjects VB.Net - Looping through Layers from TOC and ListBox

texts:
I'm creating a tool in ArcMap 10 using a Windows form... its been a struggle (maybe some of you recognize the same themes in my questions) because this is the first time I've really messed with ArcObjects and VB.Net... either way I'm almost there!! (Largely thanks to a very helpful GIS community.)

The main chunk of the code works (grabs an attribute from one layer and assigns it to the "target" layers), and now I'm in the process of linking up the user parameters.

As it stands, the "NHD Reach Layer" combo box (cmbNHD) requires the source layer be in the table of contents (TOC).

The target layers, on the other hand, can either be in the TOC (through the "Target Layers" combo box (cmbTarget), or accessible through the GxDialog. Either way, the selections are placed in the list box (boxTarget).

    ...

    Dim pLayer As IFeatureLayer
    Dim pGxDataset As IGxDataset = pGxObjects.Next
    Do Until pGxDataset Is Nothing
        pLayer = New FeatureLayer
        pLayer.FeatureClass = pGxDataset.Dataset
        pLayer.Name = pLayer.FeatureClass.AliasName
        If Not boxTarget.Items.Contains(pLayer.Name) Then boxTarget.Items.Add(pLayer.Name)
        pGxDataset = pGxObjects.Next
    Loop
End Sub

OK, so what I can't figure out, is how to form some "loopable" list from the items in the list box. The original code worked only with the items in the table of contents, and so it could loop through them all using this type of loop...

    For i = 0 To pMap.LayerCount - 1
        Dim pLayer = pMap.Layer(i)

...and just skipping the index number of the source layer ("NHD Reach Layer").

AND SO, I'm pretty sure I just need to create an enumeration of the layers in the list box,  but I've only been able to dig up the IEnumLayers interface, and it seems as if that only works with the layers in the current map.

Any thoughts/ideas?

Also, I'm not sure if this is the place (though, I can't think of a better one), but when this is all finished, I was wondering if I could post the final code (the reason I ask is because it wouldn't really be a question). I've spent a good amount of energy trolling through forums and developer help files to find relatively simple answers in the form of examples, but usually they are snippets, which, when viewed out of context, can be confusing.

It's my belief that this code would be really helpful to anybody just getting started, because it offers few frills, and combines a lot of basic methods with some intermediates thrown in here and there.


texts: Using spatial clustering or statistical tools for attribute data

texts: I have a data-set of well test results with attributes like iron and manganese concentrations. What are the best tools to try and determine if there are any spatial trends for each attribute? I would like to know if there is a statistical trend and be able to display it if any trends exist. I've tried natural neighbors but that only tests the well location, I've also looked at Grouping Analysis but am unsure on how to select the proper number of groups.


texts: How does QGIS calculate distances when the data is not projected (WGS84)?

texts:
I am playing a bit around with Qgis 2.14.6 to understand how it deals with distance measurement.

I took a line-shapefile that I saved in 3 different coordinate systems (4326, 32642 and a national one that doesn't exist in Qgis yet).

Using 32642 and the national one, the length is quite similar (6442m vs. 6428m)

But using 4326 (unprojected WGS84), I get 8057m

In all cases, I just used the Field Calculator and introduced $length.

I was already positively surprised that the result using 4326 was in meters and not in degrees.

To calculate the distance on the dataset in 4326, I assume that Qgis does some kind of re-projection in the background. Does anyone know which one?

Shapefiles here if anyone wants to play around: http://bit.ly/2cwe4w9

EDIT - adding more info after more tests

As you can see, now both my shapefile and the project are in EPSG:4326, so unprojected:

Here the rough measurement of my shapefile, with the Qgis measure tool:

When I use the field calculator, the preview gives me a value in degrees, while the result seems to be in meters:

This would mean that the function $length, like the measurement tool, always gives a result in meters, but using a very coarse approximation of the degrees => meters transformation:

Is this a bug? or should it really behave like this? I find it at best counter-intuitive (but maybe it's just me) and, at worse, completely misleading (it makes very little sense to use, everywhere on earth, the length of 1 degree at the Equator). Any input or opinion here? Should I open a bug report (at least for the inconsistency between the preview and the actual result)?


texts: When symbolising by category can I make the transparency different for each category?

texts: Normally when I'm making something semi-transparent I edit the transparency in layer properties, but now im trying to make something in the layer (a symbol level) semi transparent. I want to edit the transparency of each fill symbol when I've symbolized by category. is this possible in ArcGIS 10.0?


texts: Atmospheric correction for Sentinel 2 level 1-c in the Google Earth Engine Java API

texts:
So there is a GitHub code by Sam Murphy for doing this: https://github.com/samsammurphy/gee-atmcorr-S2

But I have written all my codes in Java and I don't know python. Is there a way to implement the same code in java? Otherwise, is there a way to run my Java code on the atmospherically corrected images.


texts: ArcMap polyline renders incorrectly unless project reloaded

texts:
I am trying to get a graphic polyline to render over top of an imported shapefile.

I'm using C# with ArcMap 10.2.2.

The issue I'm having is that the first time the polyline is created, it merely draws a red dot in the top left corner of the map like so:

I believed this to be a projection issue, but after doing a File > Save As and saving it the .mxd file and re-opening it, I am greeted with an appropriately drawn polyline:

It would appear that something isn't getting initialized properly when I am initially creating my polyline.

Creating and adding polyline as element to graphics container:

            IMxDocument mxDocument = (IMxDocument)this.Document;
            IActiveView activeView = mxDocument.ActiveView;
            IScreenDisplay screenDisplay = activeView.ScreenDisplay;
            IGraphicsContainer graphicsContainer = activeView.GraphicsContainer;

            IGeographicCoordinateSystem geoCoordinateSystem =
                screenDisplay.DisplayTransformation.SpatialReference as IGeographicCoordinateSystem;

            // segments are just two geo coordinates
            // each a lat/long pair represented as
            // doubles, see the MakePolyline method
            // for more details
            var segments = dataProvider.GetSegments();

            IPolyline line = this.MakePolyline(segments, geoCoordinateSystem);

            var lineSymbol = new SimpleLineSymbolClass();
            lineSymbol.Width = 3;
            lineSymbol.Style = esriSimpleLineStyle.esriSLSSolid;
            lineSymbol.Color = new RgbColorClass { Red = 255, Green = 0, Blue = 0 };

            ILineElement lineElement = new LineElementClass();
            lineElement.Symbol = lineSymbol;

            IElement element = (IElement)lineElement;
            element.Geometry = line;

            graphicsContainer.AddElement(element, 0);
            activeView.PartialRefresh(esriViewDrawPhase.esriViewGraphics, null, null);

MakePolyline method:

    private IPolyline MakePolyline(IEnumerable<IGeoPolylineSegment> segments, ISpatialReference spatialReference)
    {
        IPointCollection4 pointCollection = new PolylineClass();
        object missing = Type.Missing;

        foreach (IGeoPolylineSegment segment in segments)
        {
            pointCollection.AddPoint(
                this.MakePoint(segment.Start.Longitude, segment.Start.Latitude, spatialReference),
                ref missing,
                ref missing);
            pointCollection.AddPoint(
                this.MakePoint(segment.End.Longitude, segment.End.Latitude, spatialReference),
                ref missing,
                ref missing);
        }

        IPolyline polyLine = (IPolyline)pointCollection;
        polyLine.SpatialReference = spatialReference;
        return polyLine;
    }

MakePoint method

private IPoint MakePoint(double x, double y, ISpatialReference spatialReference)
    {
        IPoint p = new PointClass();
        p.PutCoords(x, y);
        p.SpatialReference = spatialReference;
        return p;
    }


texts: Consolidating a one-to-many relationship table to one-to-one based on conditions in ArcPy

texts:
Below is an example subset of the attribute table for my data:

As shown, I have multiple Point_IDs, and I would like to group them based on their unique value (i.e. Point-ID-1 and Point-ID-2). However, I would like to add a condition, saying that for each group, I want the record with the highest score. If there are records with the same highest score, then I want the record with the lowest distance. The highlighted records show the records that I want to keep after consolidating the table.
I am on ArcGIS 10.6 (Python 2.7)
This is my ArcPy code that I currently have, but I'm only able to select the records with highest score, and I'm stuck on how to incorporate the lowest distance if scores are equal condition:
    fields = [FID_Point, FID_Line, Score]
    search_dict= {}
    with arcpy.da.SearchCursor(source_fc, fields) as scur:
        for row in scur:
            try:
                if row[2] > search_dict[row[0]][1]:
                    search_dict[row[0]] = [row[1], row[2]]
            except:
                search_dict[row[0]] = [row[1], row[2]]

    arcpy.AddField_management(out_fc, "JOIN_FID", "LONG")
    with arcpy.da.UpdateCursor(out_fc, ["ORIG_FID", "JOIN_FID"]) as ucur:
        for row in ucur:
            row[1] = search_dict[row[0]][0]
            ucur.updateRow(row)

    joinfields = [x.name for x in arcpy.ListFields(source_fc) if not x.required]
    arcpy.JoinField_management(out_fc, "JOIN_FID", source_fc, arcpy.Describe(source_fc).OIDFieldName, joinfields)


texts: Create Query either using Select by Attributes or Python

texts:
I am working on something for my job, where I have a list of water/sewer service customers. Some customers have one service or the other, while others have both which is stored as two separate entries. I need to figure out a way to create a table that separates the sewer customers, water customers, and combines those that have both. Does anyone have any ideas or scripts that would give me this outcome?

This is what I have:

This is what I need:


texts:
I have an arcmap project using an arcims service from NOAA ENC (http://ocs-spatial.ncd.noaa.gov) but maybe the service is unavailable right now and I can't access to imagery...

Is there an image service from GOES available for ArcGIS?


texts: Creating polygon feature class from UTM coordinates using ArcPy?

texts:
I am working to create a polygon from UTM coordinates.  The script I am using fails with a "must be a sequence, not float" error:

Traceback (most recent call last):
  File "F:/sample/clipTiles.py", line 104, in <module>
    arcpy.Array([arcpy.Point(*coords) for coords in feature])))
TypeError: type object argument after * must be a sequence, not float

I have been working off the ArcGIS Polygon help section.  My UTM coordinates are generated earlier in a script and take the following form:

[[734855.5142000001, 4744379.424500001], [740607.6685000001, 4744379.424500001], [740607.6685000001, 4736862.1753], [734855.5142000001, 4736862.1753]]

These are the min/max point pairs of a custom bounding box I need to convert to polygon.

                # feature_info eventually takes the following form: [[734855.5142000001, 4744379.424500001], [740607.6685000001, 4744379.424500001], [740607.6685000001, 4736862.1753], [734855.5142000001, 4736862.1753]]
                feature_info = [[xmin, ymax], [xmax,ymax], [xmax, ymin], [xmin, ymin]]

                # A list that will hold each of the Polygon objects
                features = []

                for feature in feature_info:
                    # Create a Polygon object based on the array of points
                    # Append to the list of Polygon objects
                    features.append(
                        arcpy.Polygon(
                            arcpy.Array([arcpy.Point(*coords) for coords in feature])))

                # Persist a copy of the Polyline objects using CopyFeatures
                outshp = os.path.join(outws, row[2] + "_shp")
                arcpy.CopyFeatures_management(features, outshp)

How can I write these UTM coordinates to a polygon featureclass?


texts: Efficient way to find polygons with no neighbours, eg, detached houses

texts:
I have a map of buildings, I'm looking to identify only those buildings which are detached (no neighbours). I can run select by location touches in QGIS and invert this but it does not seem the most efficient way. Is there a PostGIS or QGIS method that will be quicker?

Edit; I'm testing this approach. https://geospatial.commons.gc.cuny.edu/2013/11/21/finding-islands-or-the-converse/


texts: How to add Routing Table Query-Result on pgrouting?

texts:
I've followed the pgrouting latest workshop http://workshop.pgrouting.org. I wanted to know how to add the routing table query-result into my webpage directly after the route showed up?

regards


texts: HTML5 Geolocation with ArcGIS for JavaScript API?

texts:
My organization would prefer to use its own basemaps rather than those provided by ArcGIS Online due to the accuracy of information that is represented in them. The area we service is rapidly developing, and as such, there is no way for Esri to keep up with the data as well as we can.

As such, the basemaps are all created with a spatial reference based on a State Plane projection and hosted on our server this way. This causes problems whenever I attempt to use the HTML5 Geolocation feature in my code, as it returns a latitude/longitude value.

I have toyed with several methods for using these values in my code, primarily: using a Geometry Service to project them into my desired spatial reference or overlaying Esri's basemaps with our data in order to have our data converted to the spatial reference of their maps. However, the former method requires thorough error handling and the latter creates a significant loss of accuracy as well as some graphical issues due to labels overlapping.

Is there a better way to handle this issue?


texts: How to download monthly values for a time series? NCEP_RE/surface_temp in GEE

texts:
I want to download values of Air temperature from NCEP/NCAR Reanalysis Data, Surface Temperature in Google Earth Engine. I want the average temperature of every month from 1948 to present, but with my current script I get the monthly average, this is only one value for all the Januaries and so on..

What's wrong in the script?

// Importa una colección de imágenes
var TEMP = ee.ImageCollection('NCEP_RE/surface_temp')
//.select('precipitation')
//.filterDate('2015-01-01','2015-12-31');

print(TEMP)

//Funcion para sumar la PPT mensual (devuelve una lista de imagenes)

var months = ee.List.sequence(1, 12);

var TEMP_mensual = months.map(function(m) {
  // Filter to 1 month.
  return TEMP.filter(ee.Filter.calendarRange(m, m, 'month')).mean();
});

var TEMPmIC = ee.ImageCollection(TEMP_mensual)

print(TEMPmIC)

// Importa puntos para hacer la extracción desde una Fusion Table
var puntos =     ee.FeatureCollection('ft:1_53GKyY5xjt2EJwNLI9hw2LFh5D0sOmCgCdSA-Zj');

// Agrega los puntos al mapa
Map.addLayer(puntos);
Map.setCenter(-105, 19.6, 9);

// Función para extraer los puntos en un bucle
//Tiene el sacle fijado en 30 por que reduceregions no tabaja con WGS84
// que es la proj nativa de TRMM, pero en este caso el 30 no tiene mas
// sentido que darle un scale a la función

var datos = TEMPmIC.map(function(i){
  return i.reduceRegions(puntos, 'first', 30)
})

//var dat = collection.reduceRegions(puntos, 'first')

var extract = ee.FeatureCollection(datos.flatten())
print(extract)

//EXPORTAR

Export.table.toDrive({
  collection: extract,
  folder: 'ee_extract',
  description: 'airtemp_pruebamean',
  fileFormat: 'CSV',
});


texts: How to study the association between categorical rasterised environmental variables in R?

texts: I have a set of environmental variables in raster format where I want to study the dependencies between them. To do so I was thinking of using association tests such as chi-squared test. I want to know if there are any libraries in R dedicated to this calculus.


texts: Does a GPX file contain infomation of the track color?

texts:
I am trying to open a gpx file in Qgis, with works just fine.
the problem however is that it contains multply tracks.
When I walk my tracks in the app "Galileo", I give the tracks specific colors.
Can I somehow open it up Qgis with the color information?
Right now Qgis just give the tracks random colors,


texts: Selecting features of one layer based on another layer

texts: I have a polygon layer A, then a new layer is created by applying minimum bounded boxes algorithm and then I selected those features in which height and width of bounded boxes is less than 150 and saved it in another layer named B. Now, I want the corresponding original polygons from layer A whose bounding boxes are saved in B into an another layer. Number of features is around 41k in layer A and 40k in layer B. Can this be achieved? I am using qgis 3.8


texts: How to authenticate WFS requests being made to geoserver through openlayers application?

texts:
I posted this question How to access geoserver web services on a backend server through openlayers web application hosted on another server via IIS? and @IanT helped me to expose our web services running on a backend server so they load when logging into our openlayers web app. When you first load the app, you must supply your domain credentials and then the site and services load. However, there is a tool in the app that makes wfs requests to the geoserver, and these requests are getting blocked and I am getting the 401 unauthorized message in the console.

How can I authenticate these requests being made from the application?

Our web app is running on IIS and the credentials being supplied are windows domain credentials.


texts: Finding adjacent grid squares in MGRS/USNG

texts:
Given an MGRS grid reference, for example 4QFJ123678 (100m precision), is there a good way to obtain the adjacent squares (north, south, east and west) from that origin square?

Is it better to compute the neighbors within the MGRS system, or convert to lat/long, find the coordinates 100m away in each direction, and convert back to MGRS?


texts:
I created a Centerline object from a shapely polygon, which is basically a multiline string which contains the starting and endpoint of the centerlines. How can I transform this multiline string into a [start_x start_y end_x end_y] 4-column n-row numpy array?
x_coords = cloud[:,0]
y_coords = cloud[:,1]
array_2d = np.column_stack([x_coords, y_coords])
c_hull = alphashape.alphashape(array_2d, alpha)
kozep_vonal = Centerline(c_hull)


texts: OpenLayers Interaction - events are returning wrong features

texts:
Trying to extract the drawn features on openlayers. I have written the below code and it is working fine when we draw polygon on four objects. But, when we try to draw polygon on three objects it is returning four features instead of three, not sure why the function(evt.feature.getGeometry().getExtent()) is returning four features when polygon drawn on three objects on a map.

draw = new ol.interaction.Draw({
        source: sourceMeasure,
        type: /** @type {ol.geom.GeometryType} */ ('Polygon')
    });
    openlayerMap.addInteraction(draw);

    draw.on('drawstart',
            function(evt) {
        // set sketch
    }, this);

    draw.on('drawend',
            function(evt) {

        var extent = evt.feature.getGeometry().getExtent();

        vectorSource.forEachFeatureIntersectingExtent(extent, function(feature) {

            var layer = feature.get('LAYER');
            /* here my logic is written and it working fine with four coordinates */

        });

    }, this);


texts:
I found out that ArcGIS 10.2 was released yesterday, which is great since I am looking forward to it. Although so far I haven't been able to find any information on how the update will be installed.

I was hoping that it would be more of a patch update, but alas nothing to be found here.

So it's probably a full uninstall and re-install similar to 10.0 to 10.1. Does anyone have any concrete information (with s source) on this?


texts: Moving a Staging ArcGIS Server to the DMZ?

texts:
I am trying to move my EDN Arcgis Server with all the mapservices and caches to the DMZ ..

The move is from our intranet to the DMZ to provide access to mapservices to outside world using the internet port 80 or 8080..
Everything is fine when we try accessing applications and mapservices from within our network to this box in the DMZ!!

When it comes to external access it cannot access any ??

I have tried a similar setup using windows xp box in the DMZ! Both internal and external machines can access and works well.. This should account for Network configuration correctness, hopefully..

I am wondering if it has to do with ArcGIS Server on Windows 7 ??
Am I missing something wrt  settings with ArcGIS Server or IIS ??


texts:
I'm working on Ubuntu 12.04 server; I noticed this when I installed additional support for ECW format into the GDAL 1.9.1 which worked OK. Then I tried to add ECW files to QGIS but it said they still weren't supported. I checked Help->About and to my surprise, the GDAL version that QGIS was using was 1.7.3; no wonder they didn't work.
Through trial and error, I removed QGIS then reinstalled it via the synaptic installer, and also forcefully disabled the libgdal 1.7 as a dependency. Then QGIS loaded the version of GDAL I intended for QGIS from the beginning. But this seems like an overkill to do and I don't feel like I'm having much control here
Is there a better way to for example:

check which dependencies exist and don't need installing, or
is there a way, a config file, or something to tell QGIS which library version you want it to use

The packages for 12.04 seem to be quite old and have old dependencies installed with them. This same issue arises with qgis-mapserver for which I have no way of telling which version of GDAL it is using, but my guess is the old one because ECW imagery doesn't work while the other does.


texts: Getting Error generating chart: No features contain non-null values of system:time_start

texts:
I am getting the error "Error generating chart: No features contain non-null values of "system:time_start"

var CSmap = ee.FeatureCollection('USDOS/LSIB_SIMPLE/2017').filterMetadata('country_co', 'equals', 'CS');
var area = ee.FeatureCollection(CSmap);
var s2 = ee.ImageCollection("COPERNICUS/S5P/NRTI/L3_NO2").select('NO2_column_number_density')
.filterDate('2020-03-20', '2020-03-25').filterBounds(area);
var NO2 = s2.mosaic().clip(area);
var band_viz = {min: 0.0000489, max: 0.0001659,
  palette: ['black', 'blue', 'purple', 'cyan', 'green', 'yellow', 'red']};
Map.addLayer(NO2, band_viz, 'S5P N02');
Map.setCenter(67.5,30.4,5);
print(ui.Chart.image.series(NO2, area, ee.Reducer.mean() ,1000)
    .setSeriesNames(['Nitrogen Oxies Levels'])
    .setOptions({
      title: 'Nitrogen Oxies Levels',
      lineWidth: 1,
      pointSize: 3}));


texts: arcobjects determine when attribute dialog window closed

texts:
In arcmap, I choose to "Display the attributes dialog before storing new features" in editor options for several layers.

I intend to do further development upon user closing the attribute dialog window (for new created feature).

So I want to check how to determine when the attribute dialog window (see the picture below) is closed.

I am using arcmap and arcobjects 10.2.2, vb.net arcmap add-in.


texts: Error: can't display dash lines (broken lines) when zooming PageLayout

texts:
I am working with PageLayoutControl, its mapframe include elements: lines, texts. And I am getting an error when zooming: If I zoom the pagelayout to 100% ( smaller than 110%), it can display "Dash Lines". But, if zoom to 120% (or larger), all "dash lines" become "solid lines". Please see attached screenshots:

And another problem: when I print the pagelayout, all dash lines become solid lines.

Are there anyone know what's wrong ?

Thanks and regards,

Tai


texts:
I would like to map with palette and definition used in the example of "old" WorldPop data as in the image bellow (https://developers.google.com/earth-engine/datasets/catalog/WorldPop_POP)
the "new" world wide coverage data on WorldPop (WorldPop/GP/100m/pop at: https://developers.google.com/earth-engine/datasets/catalog/WorldPop_GP_100m_pop ) in Google Earth Engine

my current code does not plot as in the image attached ... but only black and white with all saturated ....

var worldPop = ee.ImageCollection("WorldPop/GP/100m/pop")
                          .filterDate('2017');
Map.addLayer(worldPop)


texts: Will Anaconda with Python 3.7 work with Python 3.6.5 environment of ArcGIS Pro?

texts:
I just uninstalled my Anaconda with Python 3.6 and need to reinstall, but the current version is now with Python 3.7.

Will that work for ArcGIS Pro (Python 3.6.5)?

I'll have to copy the ArcGIS Pro env and place it in my new Anaconda env in order to access my arcpy and arcgis mods. If it is not compatible, I tried looking in the Anaconda archives to find the old 3.6 install for Windows 10/64-bit, but I didn't see the correct one.


texts: QGIS OpenGeo Explorer plugin fails when Publishing with "'QgsRasterLayer' object has no attribute 'geometryType' "

texts:
I am trying the QGIS OpenGeo Explorer plugin tutorial at http://suite.opengeo.org/opengeo-docs/qgis/explorer/tutorial/.

When I try "Publishing a QGIS project", with the quickstart project supplied for the tutorial, I get the error message:

  'QgsRasterLayer' object has no attribute 'geometryType' (details below)

If I try publishing my own QGIS project I get exactly the same error, so it seems not be a problem with the project itself.

I do wonder if it is a version issue, as when I earlier created a new catalog I got a message something like

  the specified catalog seems to be running an older or unidentified
  version of Geoserver

But having searched for the issue, it seems that it can be a spurious message so I felt fairly safe continuing. It is a recent installation: QGIS 2.8.3-Wien, and OpenGeoSuite 4.8, on Windows 10.

Any ideas?

Details of the error:

Traceback (most recent call last):

  File "C:/PROGRA~2/BOUNDL~1/QGIS/qgis/apps/qgis/./python/plugins\opengeo\gui\explorer.py", line 111, in run
    command(*params)
  File "C:/PROGRA~2/BOUNDL~1/QGIS/qgis/apps/qgis/./python/plugins\opengeo\qgis\catalog.py", line 475, in publishLayer
    sld = self.publishStyle(layer, overwrite, name)
  File "C:/PROGRA~2/BOUNDL~1/QGIS/qgis/apps/qgis/./python/plugins\opengeo\qgis\catalog.py", line 148, in publishStyle
    sld, icons = getGsCompatibleSld(layer)
  File "C:/PROGRA~2/BOUNDL~1/QGIS/qgis/apps/qgis/./python/plugins\opengeo\qgis\sldadapter.py", line 199, in getGsCompatibleSld
    return adaptQgsToGs(sld, layer)
  File "C:/PROGRA~2/BOUNDL~1/QGIS/qgis/apps/qgis/./python/plugins\opengeo\qgis\sldadapter.py", line 61, in adaptQgsToGs
    if layer.geometryType() == QGis.Point:
AttributeError: 'QgsRasterLayer' object has no attribute 'geometryType'


texts: How to measure the area of polygons within a fishnet

texts:
I am currently trying to measure the area of polygons within each square of a fishnet.  The aim of what I am doing is to work out the proportion of different agri-environment schemes by 1km squares.  I have tried following the instructions in the pdf below but when it comes to the dissolve I get a 'warning empty output'.  Any suggestions?

Thanks!

http://www.cabnr.unr.edu/weisberg/downloads/polygons_within_polygons.pdf


texts:
I have a dataset with contours. Depression contours are separate of the normal contours and need to be digitized counter-clockwise. Still, there are usually some faulty ones with wrong direction and manually going through thousands of them is painful. I learned that it is possible to convert the contour lines to polygons and back to lines again since arcgis polygons outer edge is always digitized clock-wise - then I just reversed all the lines.
Unfortunately, if I have several contours inside each other it creates polygons with holes and hole edge is always digitized counter-clockwise. Now if I convert the polygons to lines and reverse them all the holes are with wrong direction.

Is there a way to create polygons of each closed contourline without it creating a hole if there is a smaller contour inside?
Or is there a possibility to somehow differentiate between lines that are the outer edge of the polygon and which have been the hole edge so I can only switch the direction of the outer ones?

Or maybe there is an even simpler solution in finding depression contours that are digitized in clock-wise movement and only flip them around?

If it can be done in some other software that keeps the elevation data, you can suggest that also.


texts: Finding lasview tool from LAStools within QGIS?

texts:
I am trying to make a DEM from a .las file.

I have installed the LAStools directory and have the plugin:

I am wanting to use the lasview algorithm, but I can't see it in the standard toolbox I installed. Is there a way to add more algorithms?


texts: Closest facility analysis using polygons

texts:
I have a polygon shapefile containing polygons representing green spaces and a point shapefile of postcode centroids and am hoping to carry out closest facility analysis in ArcMap to calculate the shortest distance from each postcode to a green space. However it doesn't seem that I can load locations as incidents or facilities which are polygons, as whenever I try a message comes up reading "the input features must have the correct geometry type or have location fields".

Because of this I was alternatively thinking of creating a point layer of access points by manually digitizing them for each green space, so that I could carry out the analysis with these instead, since points seem to load fine as locations. However, as I have a few hundred green spaces, it will be very time consuming to work out where all the access points are and then digitise them. So I was wondering whether there is a way I can just load the polygons as locations to use in the network analysis, or whether there's perhaps a way in which I can intersect the network dataset with the green space polygons to create a point shapefile of the points where they both intersect?

I suppose this would essentially generate a layer of access points I could use, rather than having to manually digitize them all, but I'm not sure whether it's possible to do this or not?


texts: Tabulate Area Using Multiple Shapefiles and Multiple Rasters

texts:
I currently have ~150 county boundary shapefiles and have already clipped the statewide landcover data set by each county. My end product needs to be a table that contains either the total number of cells of each landcover class (e.g., evergreen, hardwoods, etc.) or a percentage for each landcover class by county. For example,

I tried using ModelBuilder to iterate through multiple rasters (Iterate Raster tool) but wasn't sure if zonal statistics to table was the approriate tool. Furthermore, I was unable to add an additional iterate tool so I could loop through all raster images and all shapefiles plus add this into a table similar to the one above.


texts: Automatically apply symbology appropriate to each layer?

texts:
I was wondering whether there was a way to save and apply a symbology configuration.

I am not trying to use the same symbology from one layer on multiple other layers, rather apply the symbology configuration I am looking for.

For example:
I am trying to apply the "Quantile" classification method with specific numerical formatting options. It would be great if I could save this setting as a profile and reuse it somehow.

Applying the symbology from the same layer wouldn't work as the value classification in each layer is different, please correct me if there are workarounds to this limitation.


texts: Creating map series in tif format with ArcMap?

texts:
I would like to tile a rendered ArcMap 10.6 project in tif format(2km x 2km). The map is simply a colour topographic map rendered in black and white. Important is that the resultant images use the colormap configured in arcmap and not the raw data.

Which tools can do this or would this need to be scripted?


texts: Calculate Coordinates on the Beeline between two Coordinates on a Spheroid/Ellipsoid (Earth)

texts:
I have two coordinates of locations on the earth and want to calculate points in between these two coordinates on the Beeline, for example the coordinate in the middle or at 10%.

I am able to do this on a plane.  However, I am not sure how to do this on a spheroid or ellipsoid. It doesn't need to be 100% accurate, but the planar version is definitely not good enough.

As asked by Ian Turton, I tried a planar formula, in Java this could look like this (taken from jts):

  public Coordinate pointAlong(double segmentLengthFraction)
  {
    Coordinate coord = new Coordinate();
    coord.x = p0.x + segmentLengthFraction * (p1.x - p0.x);
    coord.y = p0.y + segmentLengthFraction * (p1.y - p0.y);
    return coord;
  }

I had a look at JTS as well, but I could only find a planar (LineSegment#pointAlong) version for this, not sure if they offer a different version?


texts:
I am about to import a raster (.tif with 50mb, size 3301x21676, pixel size 10x10) into postgres with raster2pgsql.

With small extracts of the image it works but ...

raster2pgsql [myfilename].tif -F [somename] > [somename].sql

My problem is that I end up with a huge sql file (about 1 GB) for the 50mb TIFF.
with psql \i  somename.sql
this ends up with an rollback and can't change 0 to 100000 byte error

I already tried to make tiles with raster2pgsql, but it doesnt change anything for the sql file.

Any advice?


texts: Can't export features from PostGIS table into ESRI FGDB when text fields contain more than 255 characters

texts:
I am loading the PostGIS table in ArcMap 10.4.1 using a database connection.  Then I export the loaded layer (Right click, Data, Export Data...) into a newly created File Geodatabase. The procedure is working, there are only a few of the features that are not exported.

Then I tried with QGIS 2.18.11 (Right click, Save, File GDB). The same happens. The same features are not exported to the FGDB.

I looked into the PostGIS table and the failing features all have something in common: they all have one or more text fields wider than 255 characters.

Those same fields in the destination FGDB have a limit of 255 characters.

Is there a way to make the text fields adjust automatically (and larger than 255 chars) to the size of the field in the table being imported?

It is possible to create text fields as wide as 2 147 483 647 characters according to the ESRI docs.


texts: How to access current job id from geoprocessing script running on esri server?

texts:
I want to programmatically access the current job id of a geoprocessing script running on ESRI 10.2.  The script is what needs to know the job id, not the caller of the script.

While searching, I have seen plenty of examples of how the submit job request returns with the job id.  This isn't what I want.

In my script that is running on the server, I want to include the job id in the logs, however I haven't found how to obtain the job in in the script.  I initially checked arcpy.env, but I didn't see anything.  Where else should I look?


texts: Is there hydrologically correct water-course data for Africa

texts:
I was wondering if there are any hydrologically correct river-channel data available for download for Africa?

By hydrologically correct, I mean

the topology should be correct. Connected river reaches should be connected in the GIS file.
there should be (ideally ) flow directions on the polylines reflecting real flow directions

I looked at VMAP0 data, the rivers are broken and disconnected like this:


texts: OpenLayers 3 getGetFeatureInfoUrl url is always undefined

texts:
I want to show the features of a layer from a GeoServer WMS at the point where I click with OpenLayers 3. I used this
tutorial but when I use my own server the URL that getGetFeatureInfoUrl returns is always unidentified.

Below is my code:

<!DOCTYPE html>

<title>WMS GetFeatureInfo (Tile Layer)</title>
<link rel="stylesheet" href="https://openlayers.org/en/v4.2.0/css/ol.css" type="text/css">

Info Bereich: &nbsp;

Maybe it is just something with the settings of my server but I don't know what I could have forgotten. I hope you can help my at this problem.


texts: Convert floating point raster array to integer raster using GDAL

texts:
I was searching a lot for everything related to this topic. Could someone give me an idea, a suggestion how could I convert the cell values of a float raster to integer in Python environment using GDAL?

I have tried this:

int(round(myarray))


texts: Make table views in ModelBuilder with Iterator

texts:
I'm trying to use the model as shown in the image to get 300 tables for the 300 points in the input feature class. I can get 300 viewshed rasters but can only get the table of the last point.

It seems the model did create 300 tables, however, each table was named the same (vshed_0_View) even though I had put the "_%n%" in the output name, so the created tables were always overwritten by the next one.

I had uncheck "Intermediate" for vshed_%n% so the obtained viewshed rasters were not deleted.

I don't know what caused the problem. Does "Make table view" not work with Iterator? Are there any other ways to avoid being overwritten? Many Thanks!


texts: qgis2web can't create map w/ large dataset, and potentially utilizing GeoServer

texts:
I have a very large dataset (about 120k records) and about a dozen layers, all focused on California for the moment. That is, all the shapes like zips, census tracts, counties etc. are selected using the query builder to narrow to just this state.

My issue is when I try to export a qgis2web map on my Win7, what ends up happening is my computer disk space gets maxed out to whatever I have available (currently 65GB) as there is a .tif file being created (ostensibly by the qgis2web plugin but, I'm not certain it's not a setting or another plugin conflict). Is that the norm and am I simply supposed to provide more disk space for such a gigantic image file?

I'm also using a hosted GeoServer to host my shapefiles and have my data in a postgres/postGIS-extended database. My belief is it's somewhat foolish to parse all that data in QGIS to create a map when in effect the data and shapes for the map I end up making (in OpenLayers or Leaflet) will be called from GeoServer.

Is there a way to minimize the load on QGIS w/ qgis2web so that it can do it's job more efficiently and I load the GeoJSON file with the 112k data records manually?

QGIS 2.16.1 | qgis2web v1.13.1 | GeoServer 2.9.0 | Postgres 9.5.4 | PostGIS 2.2.2r14797


texts: Splitting / Divide Polygon in same squares

texts: I have a polygon in a table with 100 meters to 100 meters ( 1 row) and I would like to splitt / devide it to same square parts for 10 meters. The result should be then 100 rows with 10x10 meters. How can I do this? Is there any function or Select in PostGIS?


texts: Does the esri javascript api Identify Task Return a FeatureSet?

texts:
I am trying to get the results of an identify task using ArcGIS JS API 2.7 into a featureSet so I can utilize this sample http://arcscripts.esri.com/details.asp?dbid=16528 that takes the results of a datagrid store and exports them to an excel file. I have it working where one feature is being exported, I just need to have all selected features export. Does an Identify Task produce a FeatureSet?

Thanks!

 function executeIdentifyTask(geom) {
  //clear the graphics layer
  map.graphics.clear();

  identifyParams.geometry = geom;
  identifyParams.mapExtent = map.extent;
  identifyTask.execute(identifyParams,function(response){

 var polygonSymbol = new esri.symbol.SimpleFillSymbol(esri.symbol.SimpleFillSymbol.STYLE_SOLID, new esri.symbol.SimpleLineSymbol(esri.symbol.SimpleLineSymbol.STYLE_DOT, new dojo.Color([151, 249,0,.80]), 3), new dojo.Color([151, 249, 0, 0.45]));
     var markerSymbol = new esri.symbol.SimpleMarkerSymbol(esri.symbol.SimpleMarkerSymbol.STYLE_CIRCLE, 20, new esri.symbol.SimpleLineSymbol(esri.symbol.SimpleLineSymbol.STYLE_SOLID, new dojo.Color([0, 0, 0]), 1), new dojo.Color([25,50,225,0.3]));

     var controlItems = [];
     var surveyItems = [];

     dojo.forEach(response,function(result){

       var feature = result.feature;

       if (result.layerName =="surveys"){
    showSurveysNameGrid();
        feature.setSymbol(polygonSymbol);
        searchType="selSurvey";
    surveyItems.push(attributes);
              var attributes = feature.attributes;
    var csvTextSurv = "DOCUMENT_N;TOWNSHIP_RANGE;SECTION\n";

    csvTextSurv += attributes["DOCUMENT_N"] + ";" +
                           attributes["TOWNSHIP_RANGE"] + ";" +
                           attributes["SECTION"] + "\n";

        document.getElementById("csvSurv").value = csvTextSurv;
        map.graphics.add(feature);
   }else{
     showPointNameGrid();
    feature.setSymbol(markerSymbol);
    var attributes = feature.attributes;
    controlItems.push(feature.attributes);
    searchType="selControl";
    var csvTextMons = "POINT_NAME;SECTION\n";

    csvTextMons += attributes["POINT_NAME"] + ";" +
                           attributes["TOWNSHIP_RANGE"] + ";" +
                           attributes["SECTION"] + "\n";
   }

     document.getElementById("csvMons").value = csvTextMons;
   //add selected feature to graphics layer
       map.graphics.add(feature);
 });

     if(surveyItems.length >0){
     showSurveysNameGrid();
     var surveysStore = new dojo.data.ItemFileReadStore({data:   {identifier:'DOCUMENT_N',items:surveyItems}});
        var grid = dijit.byId('grid5');
        grid.setStore(surveysStore);
     }
 if(controlItems.length >0){
    showPointNameGrid();
    var controlStore = new dojo.data.ItemFileReadStore({data:{identifier:'POINT_NAME',items:controlItems}});
    var grid = dijit.byId('grid4');
        grid.setStore(controlStore);

 }

 });

}