texts: Adding multiple points to the GraphicLAyer in aRcgis API for JavaScript 3.18

texts:
I am using ArcGis API for JavaScript 3.18 and I wanted to make a simple map with multiple points that display a short description when I click on them but I have issue with defining multiple points on one GraphicLayer, my code for one PictureMarkerSymbol:

  map.on("load", function() {
         var gl = new GraphicsLayer();
         var p = new Point(19.466172, 51.774539);
         var s = new PictureMarkerSymbol('hat.png', 32, 32);
         var i = { Wydział: "Ekonomiczno-Socjologiczny", Adres: "Polskiej Organizacji Wojskowej 3/5, 90-255 Łódź", Telefon: "42 635 53 56"};
         var it = new InfoTemplate("EK-SOC");
         var g = new Graphic(p, s, i).setInfoTemplate(it);
        gl.add(g);
        map.addLayer(gl);
      });

How I can add more points to the (gl) Layer? I have tried with addMany function but it seems that's only for the newer API (4.1)


texts: Optimizing number of printing pages in QGIS

texts:
When I try to print under QGIS, due to the framing, it prints me 4 pages. Looking more closely, and one by changing the angle of the 2 pages, I am at 2.
How do I automate this?
Current situation:

Hoped-for situation:


texts:
I was reading this answer by MikeT to Creating circle in PostGIS? which mentions that the following

SELECT ST_Area(geom)
FROM ( VALUES
  ('CURVEPOLYGON(CIRCULARSTRING(2 1, 1 2, 0 1, 1 0, 2 1))'::geometry)
) AS t(geom);

returns 3.14033115695475. I tried it myself and it's true. He says it's only a rough approximation of pi? Why is that?


texts: Creating paths within boundary using ArcGIS Desktop?

texts:
I am looking to create movement paths of fish based on observed locations.

Since I am looking at movements of fish within rivers and lakes, simply connecting the points to form a polyline will not work as many of the paths would be over land.
I need some way to restrict the movement paths to within the water boundaries.

I am not a programmer and rely on the tools in the Arc toolbox to perform the operations. A new column in the FAT with distance between sequential points would also be useful.

I am using ArcMap 10.

Does anyone have an suggestions as to how to generate the movement paths?

A bit more information on the data; the majority of the locations are collected from remote receivers which tend to generate a lot of observations (one of my projects has over 3 million detections and growing, each detection record includea ID, datetime, lat and long).  If a fish stays within the detection range of a receiver it is detected roughly once every two minutes which results in a lot of observations (points), so my first step is to average the lat and long over a time period in this case 1 day.  Averaging the detections in the lake habitat typically works well, but doing so in the river portions can result in detections outside the river/lake boundary.  So the first thing I need to do is "snap" the averaged locations into the River/Lake Boundary and then I would like to create a movement path that is restricted to within the river/lake boundary, ideally this path would include distance for each line segment.  My rivers are polylines and the lakes are polygons, but I can convert the rivers to polygons if necessary.  Alternatively, I could use the raw unaveraged data and all my detections would be in the boundary, although creating movement paths restricted to within the boundary is still problematic, but this would result in a lot of data.

The attached image shows the daily average positions for two different fish.


texts:
I am trying  to  visualize the full route of a vessel from start point to an  end point including  stops (vessel stops for some time within a port area). I have three separate shapefiles:
1. Route of Vessel including coordinates and timestamp for every point (lat,long, timestamp)
2. Check-In/ Check-Out  in a port (lat, long, type:Check-In/Out)
3. Port  basins

My  goal is to  full visualize the route of a vessel into an animation including  the stops. When the vessel enters a port area I would like to  have a label pump up or point  to change color which  will indicate a check-In/ out  accordingly.

I am quite new in GIS and now I am palying  around with time slicers and mapbox, trying to  figure out  how to  do this and which  is the best  way. Can s/one please guide me what is the best way and also  easiest  way?


texts: Merging many different feature classes into one but only taking certain fields from each?

texts:
I'm trying to merge ~50 feature classes into 1.  The problem is, these feature classes have many different fields, and some that i don't want in the final FC.  I can delete the fields one by one manually when the Merge geoprocessing window comes up in the Field Map area, but this really isn't ideal - i just deleted ~400 fields that took over an hour, and then the Merge ended up failing.

There must be an easier way to do this?


texts: How can I remove the thick blue boundary lines from the qgis2web output?

texts:
I created a map of Africa that shows all of Africa's administrative divisions. When I try to export it to a web map, it adds a thick blue line that represents the boundaries of these divisions. How can I remove the thick blue line?
Another problem this map has is that when I hover over any of the divisions, a pop up appears that contains info about that division. The popup causes the map to move lower in order to fit. How can I prevent it?


texts: Calculate geographic position based on azimuth, distance and height

texts:
How to calculate the geographic position of an object in sky (Point B)(like an aircraft) from Home Point (Point A) based on only four variables:

Lattitude and Longitude at Home Point, Azimuth at point B, Distance from home point to point B and height of point B?


texts: Re-docking/pining panel in layout window in QGIS?

texts:
Using QGIS v3.4.4 and accidentally closed the "Item Properties" and "Items" panels in the layout window.

Cannot now get them to re-dock/pin back into place. They will only float. Have googled everywhere and tried every suggestion I would find but nothing works.


texts: Is there "label wrap on character" option in GeoServer GeoCSS?

texts:
In GeoServer (using GeoCSS) I am importing layers which use the symbol _ to wrap characters in QGIS.
As I searched, technically it is possible to have such option, TileMill uses this code text-wrap-character: '_';.
Question is if this option is available in GeoServer? Can I do wrap on character?


texts:
I used QGIS 3.2 to load OSM data obtained from Geofabrik as

  Layer/ Add Layer/ Add Vector Layer

The problem is that QGIS does not show roads completely compared to the OpenStreetMap website. I converted the .osm format to .xml and searched for the ID of missing line segments from the OSM website. I found out the  tag with the specified ID exists in the XML file so the Geofabrik data is the same as OSM webpage. I also installed QuickOSM plugin in QGIS with no success of loading the whole data. Is it some bugs in QGIS or something else?


texts: Select border region from segmented image

texts:
I have several thousand aerial photographs which have been scanned in and orthorectified, but the image borders and fiducial markers remain on the image. I am looking to identify these parts of the images and later remove them.

To this point, I have attempted to perform k-means clustering on the original raster, but have found that the pixel values for the border are often so close to pixel values in the image that I lose quite a lot of the image.

I have also attempted to use https://github.com/gina-alaska/dans-gdal-scripts/wiki/Gdal_trace_outline which takes a global threshold value and applies this to the image. This is better but each image requires a different threshold value which cannot be calculated using the histogram (as the peaks often correspond to regions within the image) or a percentage-area approach as, due to scanning errors, I may only have a partial image, changing the percentage of the tif that is 'actual image'.

The best approach that I have found so far is to use mean-shift segmentation to break up the image into regions, and assign zonal statistics to each polygon based on the underlying raster. While this feels like a step forward, I still need to somehow select the regions which relate to the border. In the image below, I am talking about the darker red border areas not the blue outer region.

The images are grayscale, but I have colourised the segmented image based on mean pixel value to make it easier to differentiate between regions. To a human, it seems obvious which polygons relate to the border, but there are several issues in programmatically selecting the relevant areas:

The border is not always complete

The border polygons may be 'contained' within the outer polygon boundary

The mean pixel value for border polygons may be similar to the image value

This means I can't select polygons on the edge, or based on a threshold value, to identify the borders.

The test image is here

Is there another approach which may help in separating out the image and border regions?


texts: Creating regularly spaced points inside polygon based on attribute value using QGIS?

texts:
I have a polygon layer with an attribute value attached representing the count of individuals inside a polygon. I'm using QGIS 3.0.

I need to visualize the counts as a grid of points inside a polygon (i.e. 1 point = 100 frogs, colored by specie). The points do not necessarily need to occupy entire polygon, and preferably be drawn around a centroid with pre-defined spacing (i.e. 300 m), so that they can represent a neat grid at a specific resolution.

Current workaround I found is to use "Regular points", with counts parameter, to create grids inside extent drawn on a map, then delete the extra points created by algorithm (algorithm rounds to a grid and you might get 20 points instead of 17 in input). This produces desired result, but requires drawing the extent of area for regular points for each polygon, as well as manual attribute input and clean up of points created. Furthermore, since I'm running "Regular Points" with count specified, instead of spacing, each polygon gets differently spaced points.

In a nutshell:
Regular points (defined # of points, drawn extent) for each polygon (batch) > Delete extra points > assign attribute values to the # points

Alternative (avoids the irregular spacing, but requires even more manipulations):

Regular points (defined spacing i.e. 300m) in the extent of the entire layer >
Clip to the polygons extent > Delete extra points in each polygon until you get desired # of points > assign attribute values to points.

The main problem with my workarounds is the polygon-by-polygon processing, which makes it hard to update the data and work with larger number of polygons.

Analogous question was asked as Creating regularly spaced, defined # of points within polygon in QGIS, but specifics are different.


texts: How to extract coordinates from a raster in QGIS

texts:
I have the Harmonised Global Soil Database in raster format (tif). I want to extract the coordinates (cell centres in lon/lat) along with the field value in order to link them with the soil data containing in the mdb file.

Is there any way to extract coordinates from a raster file?


texts: Adding output parameter to Python script tool?

texts: How do you add an output parameter in a python script in a toolbox?


texts: Avoiding gaps of data on edges after clipping a raster?

texts: After using interpolation and clipping the grid I get the result you can see on the picture below. Next to the boundaries there are a lot of gaps caused by the raster cells. These gaps are problematic for me. Is it possible to buffer the result over the boundaries?


texts:
I have been trying to read in a stack of images for a timeseries and created an output raster that contains the max pixel values of the timeseries.

The arrays were created using numpy

img_stack = np.empty((image_width, image_length, len(image_files)), np.dtype('f'))
final_values = np.empty((image_length, image_width), np.dtype('f'))

I imported the image using the following loop:

for i, fname in enumerate(image_files):
    img = gdal.Open(os.path.join(fn + ("\{0}".format(fname)))).ReadAsArray()
    img_stack[:, :, i] = img

then I read the values into a 2D list using the following nested loop:

for x in range(image_length):
    for y in range(image_width):
        final_values[x][y] = 0
        for z in range(len(image_files)):
            if (final_values[x][y] < img_stack[y, 0, z]):
                final_values[x][y] = (img_stack[y, 0, z])
            else:
                next

then use

output_raster = dr.Create(final, image_width, image_length, 1, gdal.GDT_Float32)
output_raster.GetRasterBand(1).WriteArray(final_values)

the output contains the correct pixel values, but is rotated 90 degrees to the left.

I believe the issue is due to GDAL reading the image as an array first by width and then by length, but I am not sure.

Is there an issues with my code or is there a function I could use to rotate it to the proper orientation?

Thanks


texts: Hillshade function calculation: ValueError: Received a NULL pointer

texts:
I'm trying to get a raster from a DSM running the following code from python libraries osgeo-gdal in python qgis console:

import pandas as pd
from osgeo import gdal

df = pd.read_csv(azimut_altitude.csv')
src_ds = gdal.Open(my_raster.tif')

#array generation function
def shade_array (index):
    out = r'.\my_path\raster'
    azimut = df.azimut[index]
    altitude = df.altitude[index]
    data = gdal.DEMProcessingOptions(options = [],
    colorFilename=None,
    format="GTiff",
    creationOptions=None,
    computeEdges=False,
    alg='Horn',
    band=1,
    zFactor=None,
    scale=None,
    azimuth=azimut,
    altitude=altitude,
    combined=False,
    multiDirectional = False,
    slopeFormat=None,
    trigonometric=False,
    zeroForFlat=False,
    addAlpha=None,
    callback=None,
    callback_data=None)

    shade = gdal.DEMProcessing(out + str(index), src_ds, 'hillshade')

    return shade

shade_array(0)

The problem is that when I run the script it yields this message:

ValueError: Received a NULL pointer.

I've tried to solve introducing the input file as an array but then the error message is this:

TypeError: in method 'DEMProcessingInternal', argument 2 of type 'GDALDatasetShadow *'


texts: QGIS plugin creation using ArcPy Model in ArcGIS

texts:
I have created a ModelBuilder in ArcGIS 10.3 and exported the Python file.

Now I want to use that Python file in QGIS plugin. Can that be done?


texts: Convert ST_POINT subtype to ST_GEOMETRY supertype?

texts:
Oracle 18c; ArcGIS 10.7.1
I have a query that produces an ST_POINT column:
--Mimics the SDO_GEOMETRY get_vertices function
--https://docs.oracle.com/database/121/SPATL/sdo_util-getvertices.htm

create or replace view road_vertices as (
select
    road_id,
    sde.st_pointn(shape,n.num) vertex_st_point
from
    roads s
cross join
    numbers n
where
    n.num <= sde.st_numpoints(shape)
)

That query works fine in when I run it in an SQL client.

Problem:
I run into various issues when I try to use the query in ArcGIS Desktop. For example, when I try to use the query as a layer in ArcMap, I get an error from the ST_POINT shape column:

Warning

Could not add the specified data object to the map.
Wrong column type[Invalid Geometry Storage Type]

Note: ESRI has documented the problem here:
Unable to define a query layer in ArcGIS where the data source uses an st_geometry subtype in Oracle
But that article doesn't tell us what the actual fix is.

Ugly Fix:
If I convert the ST_POINT subtype to the ST_GEOMETRY supertype, then I can add the query to ArcMap without issue.
sde.st_geometry(sde.st_astext(   sde.st_pointn(shape,n.num)   ,26917) as vertex_st_geom

The only way I've found to do that is to convert the the ST_POINT to WKT and then to ST_GEOMETRY.
However, in my opinion, it seems pretty backwards to use sde.st_astext and sde.st_geometry to perform the conversion (and adds a significant cost).

Question:
Is there a cleaner way to convert the ST_POINT subtype to the ST_GEOMETRY supertype?


texts: How can I set the default layout export resolution in pyQgis?

texts:
I want to change the default export resolution (default value = 300 dpi) programatically when I open the layout designer
l = QgsPrintLayout(project)
l.initializeDefaults()
l.setUnits(QgsUnitTypes.LayoutMillimeters)
l.setName('my_layout')
self.iface.openLayoutDesigner(l)


texts: OpenLayers3: catch invalid CQL filter exception

texts:
My code is following:

var filterParams = {
    'FILTER': null,
    'CQL_FILTER': null,
    'FEATUREID': null
};
if (filter.length > 0) {
    filterParams["CQL_FILTER"] = filter;
    lyr.getSource().updateParams(filterParams);
}
else {
    alert ('Filter is empty')}

When the filter string is invalid, for example, PARAMETER VALUE (without =), lyr.getSource().updateParams(filterParams) silently fails and users don't know that their input is invalid. I'd like to alert the user if the filter is invalid

Open Layers 2 had a means to verify a CQL expression was valid as shown below:

var format = new OpenLayers.Format.CQL();
var filter;
try {
    filter = format.read(cql.value);
} catch (err) {
    output.value = err.message;

Does ol3, ol4, or ol5 have similar functionality?

Even if the expression is valid, lyr.getSource().updateParams(filterParams) silently fails if the fields referenced are incorrect This post describes how to return a list of fields in the layer with WFS or WMS.  From that, it would take some work to parse the CQL filter to validate the fields.

Is it asking to much to expect lyr.getSource().updateParams(filterParams) to return an error message so developers don't have to code substantial error checking?


texts: National Bridge Inventory (NBI) data coordinate system

texts: I am trying to import the National Bridge Inventory Lat/Long data into ArcGIS. But when I do it with WGS 1984 Geographic coordinate system, they don't match with the base maps. I don't where I am going wrong. To convert the Lat/Long values I am using the ddMMSSHH concept.


texts:
I have the following data:

The data have 5 level (in tones blue). I want to collect only the nearest polygons with the same level and sum the number (represent client number). When I collect only for level I lose the number or client.


texts: Shapefile for linked micromap is not working

texts:
I am trying to replicate the example provided in the micromap R package of a linked map with the poverty rate for the 50 US states. I am trying to achieve the same results with a Shapefile from Germany. Here is my code:

library(micromap)
library(ggplot2)

setwd ('C:/Users/Jesus/Dropbox/linked maps/Chapter4')
ger<-readShapePoly("germany3.shp")
edPov<-read.csv('gerpoverty.csv')

statePolys <- create_map_table(ger, 'VARNAME_1')
head(statePolys)
lmplot(stat.data=edPov,
       map.data=statePolys,
       panel.types=c('labels', 'dot','map'),
       panel.data=list('Id1','poverty',NA),
       ord.by='poverty',
       grouping=5, median.row=F,
       map.link=c('Id1','VARNAME_1'))

When I tried to parse it into a flat table for use with ggplot2 the procedure does not work (instead of having 16 observations it has 8000 plus observations). Also when I try to plot the linked map, the following error appears on screen:

Error in `[.data.frame`(DF, , ord.by) : undefined columns selected

Here is the link to the files in a zip file:

https://www.dropbox.com/s/c43k755aadvu2z6/germany.rar


texts: QGIS 2.14.3 export PostGIS layer to dxf issue

texts:
I am trying to export a PostGIS layer to dxf in QGIS (2.14.3) via Project>DXF-Export using the following settings (German GUI):

When I open the resulting dxf in DWG TrueView 2016 I see nothing, but the layer name is correct and the coordinates seem to fit too:

Re-opening in QGIS gives no error, but nothing is displayed as well:

Has this something to do with PostGIS or is there something wrong with my setting?


texts: Lag OpenLayers map when I add two layers

texts:
I'm writing a directive that controls OpenLayers 3 map. The problem arises when I want to add two layers. For example KML layer and GeoJSON layer. Map then does not work smoothly. You can only zoom but do not move it. Here my code, can you see something you can refactoring?

This is view:

<openlayer geojson-data="['file.geojson']" kml-data="['file.kml']"></openlayer>

and this is piece of code:

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

and after initial my map I check if scope exist and do some code, ex.:

if(scope.geojsonData) {
  setMapGeo();
}
if(scope.kmlData) {
  setMapKML();
}

and these are the two methods:

function setMapKML() {
  var KMLLayer = new ol.layer.Vector({
    source: new ol.source.Vector({
      url: scope.kmlData,
      format: new ol.format.KML()
    })
  });
  map.addLayer(KMLLayer);
}

function setMapGeo() {
  var geoJsonLayer = new ol.layer.Vector({
    source: new ol.source.Vector({
      url: scope.geojsonData,
      format: new ol.format.GeoJSON()
    })
  });
  map.addLayer(geoJsonLayer);
}

Maybe in the code I'm doing something wrong? The problem is that it adds two layers to one map?

EDIT
I do not think it depends on the size of the data. For example, if I add a layer marker is the situation also takes place.

var vectorMarker = new ol.layer.Vector({
          source: vectorSource,
          style: new ol.style.Style({
            image: new ol.style.Icon(({
              anchor: [0.5, 46],
              anchorXUnits: 'fraction',
              anchorYUnits: 'pixels',
              opacity: 1,
              scale: 0.07,
              src: 'marker.png'
            }))
          })
        });

map.addLayer(vectorMarker);

var point = new ol.geom.Point(ol.proj.transform([33.33, 19.22], 'EPSG:4326', 'EPSG:3857'));
var iconFeature = new ol.Feature({
                  geometry: point
                  });
vectorSource.addFeature(iconFeature);


texts: Saving Equations used in QGIS Field Calculator?

texts:
Im trying to calculate certain values of a newly calculated point. So I use the field calculator to write an equation to determine the values.

Is there a way to save that equation, so that if and when I start a new job, I don't have to retype it in?

Something similar to how you save styles and symbology for a array of points.


texts: Using replace function to replace value to one of value map using QGIS field calculator?

texts:
I kind of made a mistake, with a certain setting one column to a value map. Originnaly the "Regio Column' had text which consited of four specific value's (five if I count Null with it). To make thing simplier for field work, I decided to add a dropdown list, and for specific "numbering reason" of the object I gave each text a number, but in actuality they see only the text they need to choose. But something went wrong and it is kind of my own fault since I needed to replaced text. I cannot change this now.

Now because I have over 3000 records changing manualy those (centrum) and other values to the actual value of the value map choice by choosing down the dropdown list is kind of tedious. I tried to use the replace function but I cannot get it right. To illustrate see image below:

What I am doing wrong here?


texts: Drawing a road network: offsetting overlapping lines (directed graph)

texts:
Say I have a road network defined by a table, indicating, for each arc, its ID, and its head and tail node. In other words, my network is defined as a directed graph.

For example, we have two points, A and B, linked by a two-way road ; each way is thus represented by a specific arc. It gives us the following table :

Arc_ID | Head_node | To_node
Arc_01 |     A     |   B
Arc_02 |     B     |   A

Each node (A and B) has (X,Y) coordinates, so I drew the network using "LINESTRING" in WKT.

The problem is, Arc_01 and Arc_02 lines are overlapping. I am trying to find a way to draw them so that they do not overlap, and that they are drawn in respect with the way traffic flows on them (right-hand driving).

For that, I thought about calculating the normal vector of each arc, and forcing an offset of a certain number of pixels based on that (so that changing the scale of the map doesn't cause the lines to overlap again), but I have no idea on how to implement it, since I am fairly new to QGIS (I have version 2.0.1 installed) and programming in general.


texts:
I'm looking to implement a 3D/2D tactical display in a browser. This should be capable of plotting maybe 100 tracks and their histories with 0.5 second updates.

I've researched the available technologies and considering using the following stack.

GeoServer with PostGIS backend
CesiumJS with ol3-cesium
OpenLayers 3 - supports WFS

Most of these technologies seem to be based around serving static data rather than live dynamic data. I'd ideally like a high performance solution using WebSockets rather than HTTP polling, however I don't want to re-invent the wheel.

Can anyone provide any recommendations/experiences?


texts: Locate Features along routes expands to double the size of the original shapefile

texts:
I am attempting to run the Locate Features Along Route tool (Linear Referencing) and when I run it, the output is returning almost double the amount of data points as are in the input while leaving out some of the needed. I am attempting to locate line segments on an existing road network. Both the line segments and the road network have unique ID's.
Route Identifier Field: CDS_Num
XY Tolerance: 5m
Keep only the closest route location: FIRST
Include distance field on output table: DISTANCE
keep zero length line events: ZERO
Include all fields from input: FIELDS
Use M Direction Offsetting: M_DIRECTION


texts: geotools shapefile does not support utf-8 for dbf header

texts:
I create a ShapefileDataStore to create a new shapefile and call dataStore.createSchema(...) . But it does not work for utf-8 header fields.
I checked the geotools source code and realized that it does not support utf-8.

org.geotools.data.shapefile.dbf.DbaseFileHeader.java

public void writeHeader(...) throws IOException{

  ...
  line  844:
  buffer.put((byte)fields[i].fieldName.charAt(j));
  ...
}

How can I create a new shapefile with desired utf-8 header fields?

note: The shapefile plugin can read utf-8 header fields but it cannot write them in utf-8 !


texts: What ArcGIS interpolation tool would be most appropriate for price/square meter?

texts: I have the price per square meter of around 13000 properties ; for the other 8000 of my shapefile the price is missing. The easiest option is to transfer the average price per square meter to the missing properties. However, I believe it would be more realistic if a sort of neighbour effect would be considered in the calculation. I thought therefore to use interpolation tool. I am not very experienced in it, I have read the tutorials and explanations, but it is still to me not very clear what tool would be more appropriate for my case, and with what settings. Would you be able to advice on this?


texts: Intersection test for all records in a postgis table for each other

texts: I have a table "spheres". I want to check if there are any intersections between the records of these table in 3D. I want to use ST_3D Intersection. But I don't know how I can make using this function for intersection test between every geometry records in a single table. In examples intersections are made two between different tables and geometries. I want to make intersection tests between different geometries but they are in the same table. How can this be done??


texts: arcpy.SelectLayerByAttribute_management not selecting anything

texts:
I need to run a process on selected records only in a feature class stored in Enterprise geodatabase (SQL Server) that meet a certain condition:
WC = "Status = 'U' AND IN_City = 'Big City' "

But when I run the arcpy command:
arcpy.SelectLayerByAttribute_management(sde_featureclass, "NEW_SELECTION", WC)

Nothing gets selected
If I use the SQL builder in the select by attributes feature in ArcGIS Pro I get several hundred selected records.
What do I need to do in my script that will make a selection on the featureclass using this criteria where I can run a process that will run only on the selected records?


texts:
I'm trying to select POIs from around given location within some distance. I tried so far two approaches I found online, but unfortunately none seems to be working - I'm getting no rows despite the fact I can find records from that area by hand.

My current code is:

(1)

SELECT * FROM public.planet_osm_point
WHERE ST_DWITHIN(public.planet_osm_point.way,
   ST_TRANSFORM(ST_SETSRID(ST_MAKEPOINT(40.6978,-73.9792),4326),3857), 50);

(2)

SELECT *
FROM public.planet_osm_point
WHERE ST_DWithin(
  ST_SetSRID(ST_Point(40.6978,-73.9792), 3857),
  ST_Transform(public.planet_osm_point.way, 3857),
  500
)
LIMIT 200;


texts: Serving content from MXDs to ArcGIS for Server by publishing from ArcGIS Pro?

texts:
I have 10 ArcMap .mxd's serving content using ArcGIS Server. I have made the switch to using ArcGIS Pro for my all of my analysis operations and publishing of ArcGIS online content.

How do I now replace my .mxd's which are sharing ArcSDE data with .aprx files?
Is it possible to publish using arcgis server on ArcGIS Pro?

Originally I thought I would just share the data stored in ArcSDE using ArcGIS Pro in a webmap, that would mean the data is no longer upto date once we make changes to the data in ArcSDE.


texts: I have the data for several cities in a layer. I'd like to run a batch export, that goes through a list and zooms in on one city (coordinate + zoomlevel), then the next, etc and exports each results based into a pdf. Is this possible?


texts:
I am having trouble getting GDAL to read a file geodatabase with data projected in NAD83(2011) State Plane Arizona Central FIPS 0202.  I have developed the following code:

# standard imports
import sys

# import OGR
from osgeo import ogr

# use OGR specific exceptions
ogr.UseExceptions()

# get the driver
driver = ogr.GetDriverByName("OpenFileGDB")
print ("I made it to line 12")
# opening the FileGDB
try:
    gdb = driver.Open("CM_2009_DEM.gdb", 0)

except Exception, e:
    print e
    sys.exit()

The data can be downloaded here.  I am new to python and GDAL and I am simply trying to list out the items in the specified geodatabase.  When i run the code I get the following error:

EPSG PCS/GCS code 102988 not found in EPSG support files.  IS this a valid EPSG coordinate system?

The projection was made in house and there hasn't been an EPSG of ESRI code assigned to it yet.  I know its valid, but is there a way to make GDAL think its valid?

The projection file reads as follows:

PROJCS["NAD_1983_2011_StatePlane_Arizona_Central_FIPS_0202",GEOGCS["GCS_NAD_1983_2011",DATUM["D_NAD_1983_2011",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",213360.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-111.9166666666667],PARAMETER["Scale_Factor",0.9999],PARAMETER["Latitude_Of_Origin",31.0],UNIT["Meter",1.0]]


texts: Select features interactively in ArcGIS WebAppBuilder Using a Widget

texts:
I am looking to build a seemingly simple tool into my WebAppBuilder Map.

Some basics:

Content is created in ArcGIS 10.3.1
Content is published to ArcGIS 10.3.1 Server as a Cached Service
Map is developed on ArcGIS Online
Map is pushed to WebAppBuilder
Out-of-box widgets are used currently (have developer version of WebAppBuilder so I can work with code to custom create if required)

The Task:To replicate the function of Select Features Interactively in WebAppBuilder.

To date I have tried using variations of the Query Widget as it allows users to interactively select features by using the "Use Spatial Filter to Limit Features" and then the "Only Features Touching a User Defined Area" (see image below)

However, this has not proven successful as it is a variation of what I am seeking. I simply want the user to be able to select geometries (polygons, lines, points) and have them highlight on the map.

As a bonus - I am looking to have the user be able to control the color that highlights the selection - similar to the Selection Options. (see image below)


texts: Formula to calculate the next coordinate given a coordinate, distance,direction and probably speed

texts:
I have this coordinates 52.3545362,4.7638767 which I got from a Google Maps search of Amsterdam.

URL

https://www.google.com/maps/place/Amsterdam,+Netherlands/@52.3545362,4.7638767,11z/data=!3m1!4b1!4m5!3m4!1s0x47c63fb5949a7755:0x6600fd4cb7c0af8d!8m2!3d52.3666969!4d4.8945398

Is there a formula I can use to get the coordinates if I were to move one meter from my current coordinates in any direction for instance,

What would be the new coordinates if I moved 1 meter westwards from the coordinate 52.3545362,4.7638767

If there is a way I can do this along known routes, it would be an added bonus.


texts:
For example: you run the DesktopExtensionCS2010 sample that comes with the SDK, having it open arcCat at startup.

Now whenever I open ArcCat the extensions are there. I can uncheck them but I want to get rid of them (without hacking the registry).


texts: Selecting polygon features that are not fully contained against neighbors in QGIS?

texts:
In QGIS what would be the best method to select features that are not fully contained against it's neighbours?

For example I want to select local authorities that have a coastline so in the example below how would you select the ones highlighted in yellow.

Obviously you could manually select them but this can be time consuming on larger datasets and is prone to error.


texts: Calculate lat lon bounds for individual tile generated from gdal2tiles

texts:
I have multiple tile sources generated using gdal2tiles I'd like to present on the same map. So when presenting a tile I need to determine which source to serve it up from, checking against the bounds.

Anyone know how to calculate the lat lon bounds for a single tile based purely on zoom, x, and y (from the file structure generated by gdal2tiles)? BTW: I'm using Google Maps API v3 in case there's a need to call some functionality from the API to help with calculations.

The reason I'm restricted to purely zoom, x, and y is because the tiles aren't just getting called by an overlay on the map but also some custom print functionality that allows printing outside of the map.


texts: Interpolating geographically weighted regression (GWR) - fit to raster?

texts:
Using the spgwr package I followed this tutorial and found a way to interpolate the result to a raster layer (gwr on  pixel/raster  data using the meuse data set, further down the script). But the second step is taking hour by hour and I am not sure if R is hanging or if the code is wrong.

Here is what I did until now:

# calculating the bandwidth
bw<- gwr.sel(pb_train ~ gesteine + schnee_tag_1 + rs_hospso, data=env.spdf, adapt=FALSE)

# running first gwr model
fGo<- gwr(pb_train ~ gesteine + schnee_tag_1 + rs_hospso, env.spdf, bandwidth=bw, hatmatrix=TRUE)

spplot(fGo$SDF, "pred")

This is where R keeps hanging...

gwr.priclus<- gwr(pb_train ~ gesteine + schnee_tag_1 + rs_hospso, env.spdf, bandwidth=bw, fit.points=env.grid, predict=TRUE, se.fit=TRUE, fittedGWRobject = fGo)

env.spdf is a SpatialPointsDataFrame with my predictors at each location. env.grid is a SpatialPixelsDataFrame that contains the coordinates and the predictor values for the entire study area.

What am I doing wrong here? I would also like to know what exactly happens within the scond step.


texts:
Following the QGIS installation guide I have tried both UbuntuGIS and QGIS stable, but both fail with the error message:

The following packages have unmet dependencies:

python-qgis : Depends: python-qgis-common (= 2.6.1.1+trusty1) but it is not going to be installed

           Depends: libqgispython2.6.1 but it is not going to be installed
           Depends: libqgis-analysis2.6.1 but it is not going to be installed
           Depends: libqgis-core2.6.1 but it is not going to be installed
           Depends: libqgis-gui2.6.1 but it is not going to be installed
           Depends: libqgis-networkanalysis2.6.1 but it is not going to be installed

qgis : Depends: libgdal1h (>= 1.8.0) but it is not going to be installed

    Depends: libqgis-analysis2.6.1 but it is not going to be installed
    Depends: libqgis-core2.6.1 but it is not going to be installed
    Depends: libqgis-gui2.6.1 but it is not going to be installed
    Depends: libqgis-networkanalysis2.6.1 but it is not going to be installed
    Depends: qgis-providers (= 2.6.1.1+trusty1) but it is not going to be installed
    Recommends: qgis-plugin-grass but it is not going to be installed
    Recommends: qgis-plugin-globe but it is not going to be installed

E: Unable to correct problems, you have held broken packages.

I had QGISv.2.4 running on my machine with Ubuntu 12.04. After upgrading to Ubuntu 14.04 I was getting error message that pyqgis could not be loaded.

So I uninstalled QGIS using the Ubuntu Software Center and also run

sudo apt-get purge qgis

After that I tried to install the latest version of QGIS Debian Stable and UbuntuGIS unstable but getting the error messages described above.

The weird thing is that neither Ubuntu Software Center nor apt-get list have QGIS listed as installed, but if I run a search for qgis on my file system qgis and all its dependencies are still there. That's why I getting told to have broken packages.


texts:
I am new to QGIS and would like to convert some GeoTiff images to KML/KMZ for distribution to non-GIS users to look at in Google Earth. I have recently left the Windows OS and I am running QGIS in Linux-Ubuntu. I have not yet successfully installed GE Pro on this OS. In the past I would have used Arc, Global Mapper, or GE Pro to do this.

Can someone point me to a work around?


texts:
Following on from How to use `layer.getFeatures()` in a standalone script?, I would like my script to do the following (includes illustration):

Search for the first shapefile (Paris (O).shp) in Country_France > Original files and clip it onto Grid.shp.
Search for the first shapefile (Paris (P).shp) in Country_France > Processed files and take the first value in the Rank field.
Take Paris (O).shp and add a Rank field with the value from Paris (P).shp.
Output results in Country_France > Processed files with new filename, Paris (P)x.shp.
Search for second shapefile in Country_France > Original files (Toulon (O).shp) and repeat Step 1.
When all shapefiles are found in Country_France, move on to Country_Germany and repeat Step 1.

Problem:

ALL output shapefiles are based on the last shapefile to be clipped from Step 1. In other words, Paris (P)x.shp) has the exact same geometry with Toulon (P)x.shp. The only difference is in the Attributes Table where the Rank value is correctly applied:

I am expecting Paris (P)x to look like the image below (ie. Paris (O).shp clipped onto *Grid.shp, not Toulon (O).shp):

I'm quite certain the problem lies with the loops, the problem is I am unsure how to create the loop(s) which resembles the Steps above. I obtained the results using the following script:

root_dir = home + "\Desktop\Test\\"
original_dir = "Original files"
processed_dir = "Processed files"
country_dir = "Country_*/"
grid_path = root_dir + "Grid.shp"

for f in glob.glob(root_dir + country_dir + original_dir + "\*.shp"):

    outputs_1=general.runalg("qgis:clip", grid_path, f, None)

for g in glob.glob(root_dir + country_dir + processed_dir + "\*(P).shp"):

    pieces = list(os.path.splitext(g))
    pieces[0] = pieces[0][:-3] + "(P)x"
    newFile = "".join(pieces)

    vLayer = QgsVectorLayer(g,"any name","ogr")
    idx = vLayer.fieldNameIndex('Rank')
    Rank_value = vLayer.uniqueValues(idx)[0]

    general.runalg("qgis:advancedpythonfieldcalculator", outputs_1['SAVENAME'], 'Rank', 1, 10, 0, False, "value=%f" % Rank_value, newFile)


texts: Maxent can't get past "javalang.ClassCastExeption: [D cannot be cast to [F" error message

texts:
I have turned all stones I can think of how to solve this, and it feels like it should be something simple I am missing, but I just can't get Maxent to work.
When Maxent (v.3.4.3) reaches the step "Extracting random background and sample data" it stops and gives me the error:
Error reading file: C:\data.csv: javalang.ClassCastExeption: [D cannot be cast to [F

I formatted the csv files as specified in instructions to have only 3 columns: species, longitude, latitude. I have made sure that all rasters are in ASCII format, with the same projection, of the same extent, and with the same pixel size. I have even tried it with projected datasets and with unprotected geographic coordinate systems, both for the rasters and csv files. And I still get the same error.
Has anyone encountered a similar error or know why this happens?


texts:
Is there any reason why QGIS swap axis order if CRS (EPSG:3908) is defined inside GML, for example:

<au:geometry>
    <gml:MultiSurface gml:id="_42d05737-123d-4a52-a655-aef968cd9511" srsName="urn:ogc:def:crs:EPSG:3908" srsDimension="2">
        <gml:surfaceMember>
            <gml:Polygon gml:id="_ae5f0b7e-2315-4567-95af-dc5ba90183dc">
                <gml:exterior>
                    <gml:LinearRing>
                        <gml:posList>6614600.56 4825580.27 6614639.00 482571.91 6614785.47 4825616.26 6614600.56 4825580.27</gml:posList>
                    </gml:LinearRing>
                </gml:exterior>
            </gml:Polygon>
        </gml:surfaceMember>
    </gml:MultiSurface>
</au:geometry>

The representation of data, in this case, is not correct, because X and Y are swapped.

If I delete the srsName definition and manually choose CRS during import, the representation is correct, the same representation of data is returned by FME. GML file is produced using hale studio from data stored in PostGIS database. I am using QGIS 3.0.2.


texts:
Following Getting polygons from PNG images for the QGIS Shapefiles
I intended to digitize the following PNG image in QGIS 3.8 following this instruction:
https://www.qgistutorials.com/en/docs/digitizing_basics.html
which is valid for version 2.2 (2+)
In the result some steps don't match to the 3.8 version, therefore I have a few doubts/problems, that I already encountered:

I cannot manage with snapping options, as error occur:

The edit toolbar doesn't match with the tutorial for 2.2

I should have an option to drawing the polyline and polygon, which doesn't appear here. In turn I can only draw the straight line from the topleft corner of the map.

Settings don't provide the "Snapping options"

I have no option for drawing (redrawing) the existing lines on the PDF image

Is anyone able to give me some hint?

Even if I activate the Snapping toolbar I am still unable to draw the lines. I can only mark the polygons as per in the image below:


texts: Cutting down tiff file layer to smaller polygon using ArcGIS Desktop

texts:
I have a dynamic land map featuring open forest, closed forest etc of an area which is a tiff layer file. I want to be able to cut the layer file to my area of interest, the polygon I created.

How do I retrieve data from my dynamic land map just in that polygon region?


texts: What is faster for geoprocessing (intersect) big shapefiles?

texts:
I have 200 Mb shapefiles of polylines downloaded from OSM. I've been working with QGIS 2.0 in Ubuntu 64 Bit 8 Gb RAM to make some geoprocessing tasks and it's taking 2 days for processing.
Once I've got my final shapefile, i have to make further calculations/scripting for making a model, and i don't have that time.
So, in order to accelerate the process I've been thinking in:

Exporting my shapefiles to R
Converting them to SQL and process them by spatialite
Maybe gdal in shell?

But i don't know what is faster/better option, because I'm a beginner with big data and scripting, so python is not yet an option for me.  I have a few experience in R, but I do not know if it's the better option.

I have a database for Brazil. I'm interested in two shapefiles: "landuse" and "roads". I have another with the state of São Paulo "sp". So, i just need to intersect the landuse with sp = landuse_sp, roads with so = roads_sp, and later roads_sp with landuse_sp.
That's in order to have every the roads of the state of São Paulo with the landuse. Then I'll intercept with the municipalities and with another dataset with vehicle count I'll generate a model for vehicle count.
With this final shapefile, i need to perform case calculations adding fields. With field calculator, creating the field "count" as an example, would be something like this:
CASE WHEN roads IS 'primary' AND landuse IS 'residential' THEN exp(8 + 0.0033*2)
 ELSE 0
END
This is just an example, but it's quite the idea


texts: Multiple least cost paths - automated labeling of the calculated paths with formula starting-point_destination-point?

texts: I am using r.cost and r.drain to calculate a least cost path, using a single friction layer (slope), in a set of points, from each single point to the other. I would like to use the same approach on a significantly bigger set of points (up to 200). Is there a way to introduce automatic labeling of the calculated least cost paths (lines) with a formula starting-point_to_stopping-point? I would want to avoid to label each single calculated line manually, it would be time - consuming and ultimatelly would deter me from the whole process.


texts: Exporting feature class into multiple feature classes based on field values using ArcGIS Desktop?

texts:
I have a feature class with over 2,000 features, and I need to make them all individual feature classes based on a field.

Is there a way to do this?


texts: How do I expand the .SDE file specified in my current workspace?

texts: I'm storing my database connection in a .SDE file. When I output my feature class names in that database connection, I get the full path to the .SDE file. I want the details, like instance, for the database stored in the .SDE file. How do I expand (maybe not the right term) the .SDE file to expose the actual database details?


texts:
I have a few hundred invalid x-coordinates on points in my database, as a result of a GPS logger malfunction. Database is Spatialite, GPX-files were imported with gpx2spatialite.
As they are messing up my tracks, I'd like to replace these invalid x-coordinates with the average of the x-coordinates of the preceding and succeeding point.
The table I'm viewing the coordinates in is a virtual table:
CREATE VIRTUAL TABLE "idx_trackpoints_geom" USING rtree(pkid, xmin, xmax, ymin, ymax)

It seems this virtual table uses as its sources these tables:
idx_trackpoints_geom_node
idx_trackpoints_geom_parent
idx_trackpoints_geom_rowid
The idx_trackpoints_geom.pkid for my faulty row matches idx_trackpoints_geom_rowid.rowid, which gives me idx_trackpoints_geom_rowid.nodeno, which then lead me to idx_trackpoints_geom_node.data. However, it's datatype BLOB...
The faulty point:
xmin: 0.850792
xmax: 0.850792
ymin: 60.231220
ymax: 60.231228
Here's the output (1228 hex bytes) of HEX(idx_trackpoints_geom_node.data); this does not appear to adhere to the BLOB format described in the documentation.
0000002900000000000560E141C8189441C818964270EA374270EA39000000000006143F41C817AB41C817AD4270EA484270EA4A00000000000560E241C8199541C819974270EA484270EA4A00000000000560E341C81A9341C81A954270EA5A4270EA5C000000000006143E41C8178241C817844270EA604270EA6200000000000560E441C81B8741C81B894270EA6D4270EA6F000000000006143D41C8175541C817574270EA754270EA7700000000000560E541C81C7441C81C764270EA844270EA86000000000006143C41C8172441C817264270EA8B4270EA8D00000000000560E641C81D5A41C81D5C4270EA984270EA9A000000000006143B41C816F141C816F34270EAA64270EAA800000000000560E741C81E4441C81E464270EAB14270EAB3000000000006143A41C816B941C816BB4270EABC4270EABE00000000000560E841C81F2B41C81F2D4270EACB4270EACD000000000006143941C8167E41C816804270EAD34270EAD500000000000560E941C8201841C8201A4270EAE74270EAE9000000000006143841C8164441C816464270EAE94270EAEB000000000006143741C8160641C816084270EAFE4270EB0000000000000560EA41C8210B41C8210D4270EB074270EB0900000000000D54473DCCCCCB3DCCCCCD4270E9EF4270E9F100000000000D2B2B3D4CCCCB3D4CCCCD4270EB024270EB0400000000000E1F9441C3580641C358084270E9F54270E9F700000000000E1F9541C3578D41C3578F4270EA084270EA0A00000000000E1F9641C3571741C357194270EA1A4270EA1C00000000000E1F9741C3569E41C356A04270EA2B4270EA2D00000000000E1F9841C3562941C3562B4270EA3C4270EA3E00000000000E1F9941C355B641C355B84270EA514270EA5300000000000E1F9A41C3554A41C3554C4270EA624270EA6400000000000E1F9B41C354E441C354E64270EA774270EA7900000000000E1F9C41C3548441C354864270EA8D4270EA8F00000000000E1F9D41C3542B41C3542D4270EAA74270EAA900000000000E1F9E41C353E041C353E24270EAC24270EAC400000000000E1F9F41C3539941C3539B4270EADC4270EADE00000000000E1FA041C3536241C353644270EAF84270EAFA00000000001009B541C153A741C153A94270EA0F4270EA1100000000001009B641C153D141C153D34270EA3A4270EA3C00000000001009B741C153F741C153F94270EA664270EA6800000000001009B841C1541741C154194270EA934270EA9500000000001009B941C1542941C1542B4270EABE4270EAC000000000001009BA41C1543341C154354270EAEB4270EAED000000000019EC783F59CD813F59CD834270ECC54270ECC7000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
Is there any reasonably sane way to manually edit a point in Spatialite?


texts:
I am clustering points in a table based on their Mercator distance (~10km, it doesn't have to be that precise, so I'm fine with projection distortion). What I want to achieve is to only query the clusters which contain multiple points.

I have come up with this query:

select * from
(SELECT unnest(ST_ClusterWithin(ST_Transform(geom, 3857), 10000)) gc
FROM mils) f
where ST_NumGeometries(gc) > 1;

Which does exactly what I want. However, this will do the clustering first on the whole table, and only then filter on the results.

My question is: is this effective enough / the way it should be done, or is there a more efficient way of filtering clusters on the fly?

Please note that the table is not that big, its size is usually between 500 and 10000, tops. Also, I have added a spatial index on the column geom.


texts: Principal Component Analysis (PCA) issue on Spectral indices and transformations of Earth Engine

texts:
Following the Principal Component Analysis (PCA) section on Lab 4: Spectral indices and transformations from https://developers.google.com/earth-engine/edu I reached a point which its not clear, in fact its black, as seen in the picture

And I'm not sure on how to ask it. May I just ask for help on executing this Principal Component Analysis (PCA) section on Lab 4: Spectral indices and transformations?

What's the meaning of the black body? Or is there a bug around here?

May I provide further info on the issue?


texts: semi variogram analysis for buidings in raster dem with penta spherical model

texts:
I am carrying out semi-variogram analysis in arcgis using geostatistical tool
on a DEM. the plot gives the semi variance and lag distance. I want to know at what elevation of DEM the sill value is flatenning or is constant.


texts: POSTGIS intersection query is a extremely slow

texts:
I have two spatial postgis layers tablea(polygon) and tableb(multipolygon) and I want to create a postgis query to take a new intersection polygon layer.

here the postgis query :

CREATE TABLE myschema.new_table_int AS(
  SELECT myschema.tablaa.field1 AS field1,
         myschema.tableb.fieldb AS fieldb,
         (st_dump(
           ST_CollectionExtract(
           st_intersection(
              myschema.tablaa.geom,
              myschema.tableb.geom)::GEOMETRY(GEOMETRY,4326),3))).geom AS shape
  FROM  myschema.tablaa,
        myschema.tableb
  WHERE st_intersects(myschema.tablaa.geom,myschema.tableb.geom))

that code works fine but is extremely slow.

Is there a faster way to populate the table?

update

table a have over 250000 features

table b have over 10000 features

the calculate time cost is over 2 hours

here example :


texts: Disabling XML schema when using GML3 decodeFeatureCollection in GeoTools

texts:
I wrote some code based on http://docs.geotools.org/stable/userguide/library/xml/geometry.html#id6

GML gml = new GML(Version.GML3);

try {
    InputStream in = new FileInputStream("file.gml");

    SimpleFeatureIterator iter = gml.decodeFeatureIterator(in);

    while (iter.hasNext())
    {
        SimpleFeature feature = iter.next();
        System.out.println("Properties: "+feature.getProperties());
    }
}

This works well for most of my files but some other dataset comes with a XML schema that points to an offline host and a non-existent .xsd file. I trust the files so I would just like to ignore the schema.

How can I make GeoTools ignore the XML schema?

Searching online I found suggestions to do eg parser.setStrict(false), for example in https://github.com/boundlessgeo/geotools-2.7.x/raw/master/demo/gmlparsing/src/main/java/org/geotools/demo/gmlparsing/GMLParsing.java ,  but that uses a different approach in reading the file. If I created a Parser object, how do I tell gml.decodeFeatureIterator() do abide it?

I am very new to GeoTools so this question might be very confused.


texts:
I can't seem to access the plugin repository.

Is this an issue with my system or is there something going on at the server end?!

I get the following error as wells sometimes a window saying: "Error reading repository: QGIS Official Repository
HTTP request failed"

It's snowing like crazy here and that often messes with the internet but every thing else seems to be working fine.


texts:
I want to use PostGIS data to display it in Leaflet.
I've created a PostGIS view (v_baeume_GeoJSON) that delivers already valid GeoJSON. When I retrieve the data from this view using PDO I have to do a few steps to get again valid data to display it in Leaflet?

CREATE OR REPLACE VIEW bongert.v_baeume_geojson
AS SELECT row_to_json(fc.*) AS geojson
   FROM ( SELECT 'FeatureCollection'::text AS type,
            array_to_json(array_agg(f.*)) AS features
           FROM ( SELECT 'Feature'::text AS type,
                    st_asgeojson(st_transform(lg.the_geom, 4326), 15, 0)::json AS geometry,
                    ( SELECT row_to_json(t.*) AS row_to_json
                           FROM ( SELECT lg.gid,
                                    lg.id_sort) t) AS properties
                   FROM bongert.baeume lg
                   LIMIT 2) f) fc;

The data is valid GeoJSON and looks like this:

{"type":"FeatureCollection","features":[{"type":"Feature","geometry":{"type":"Point","coordinates":[6.47434138884473,49.8181346585368]},"properties":{"gid":21647,"id_sort":"ap_"}},{"type":"Feature","geometry":{"type":"Point","coordinates":[6.44536728156301,49.6813625945007]},"properties":{"gid":34501,"id_sort":"ap_starearl"}}]}

When I use this view as datasource in a PHP script, that is addressed by an AJAX call, I have to strip off parts of the resulting string, since I have problems to get the data 'as it is' out of the script. How can I get the blank data out of the PDO query without the need of JSON encode and the whole preg_replace to get back to the plain result of my view?

require("db_connection.php");

$query = "SELECT *
    FROM bongert.v_baeume_geojson";

$result = $conn->query($query);
$geojson = json_encode($result->fetchAll(PDO::FETCH_ASSOC));
$geojson = preg_replace('/\\\"/',"\"", $geojson);
$prefix = '[{"geojson":"';
$suffix = '"}]';
$geojson = preg_replace('/^' . preg_quote($prefix, '/') . '/', '', $geojson);
$geojson = preg_replace('/'.   preg_quote($suffix, '/') . '$/', '', $geojson);
echo($geojson);
?>


texts: Creating a GRASS location from a bash script

texts:
I would like to have a GRASS script that creates a new location of its own, instead of using an existing location. This facilitates the usage of the script by other users.

Going through the guide on Working with GRASS without starting it explicitly and the notes on creating locations form the command line, I arrived at the following script:

#!/bin/sh
# path to GRASS binaries and libraries:
export GISBASE=/usr/lib/grass74
export PATH=$PATH:$GISBASE/bin:$GISBASE/scripts
export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:$GISBASE/lib

# GRASS session variables
GRASSDB="/tmp/grassdata"
GRASSLOC="tmp_location"

mkdir -p "$GRASSDB/$GRASSLOC/PERMANENT"

# Create a new GRASS location
#grass74 -c $PATH_IN"ORCDRC_M_sl1_250m_clip.tif" "$GRASSDB$GRASSLOC"

grass74 -c /path/to/some/file.tif "$GRASSDB$GRASSLOC"

echo "Creation done"

# path to GRASS settings file
export GISRC=/tmp/grass7-${USER}-$GIS_LOCK/gisrc
# remove any leftover files/folder
rm -fr /tmp/grass7-${USER}-$GIS_LOCK
mkdir /tmp/grass7-${USER}-$GIS_LOCK
# set GISDBASE, LOCATION_NAME, and/or MAPSET
echo "GISDBASE: $GRASSDB" >>$GISRC
echo "LOCATION_NAME: $GRASSLOC" >>$GISRC
echo "MAPSET: PERMANENT" >>$GISRC
# start in text mode
echo "GRASS_GUI: text" >>$GISRC

# Test it
g.version

This kind of works, however, the grass74 -c command not only creates a new location, it also starts a new GRASS session and launches the GUI.

Is there another way of creating a location without launching the GUI (and creating a session)?


texts:
I'm using Leaflet JS with a modified GeoJSON file of world contries to display a choropleth, like this:

Is there any way to display the far-eastern part of Russia as a contiguous entity? Rather than it being split to due to its coordinates.


texts: How to link map to a timeline using Odyssey.js

texts:
Background info of our project:
we are trying to do a vertical timeline that is linked to a map and story for web. Visuals are attached.

We are using this plugin, http://cartodb.github.io/odyssey.js/ as of now. Is it possible to incorporate timeline functionality on top of this plugin? if so, how? (we have seen torque tool, we are looking for the map to move with the timeline) Also is it possible to remove the map while retaining marker placements?


texts: How do you construct a definition query to return only those features in the current month without hard coding the date?

texts:
I have a feature class containing a date field within a 9.3.1 File Geodatabase. There are "duplicate" features contained within this layer, one feature for each month over the past few months.

However, the final dataset may contain dates for the last few years.

| Feature 1 | '2012-05-01' |
| Feature 1 | '2012-04-13' |
| Feature 1 | '2012-03-23' |
| Feature 1 | '2011-03-23' |
| Feature 2 | '2012-05-12' |
| Feature 2 | '2012-04-14' |
| Feature 2 | '2012-03-30' |
| Feature 2 | '2011-03-30' |
... and so on

I want to apply a definition query such that the map will only display features for the current month and year and I do not want to update the definition query each month.

I have attempted to use a combination of date, CURRENT_DATE, EXTRACT, and CONCAT but without much success. Ultimately I think I need a definition query that looks something like this but with the hard-coded dates removed.

DateField >= date '2012-05-01' AND DateField <= date '2012-05-31'

CURRENT_DATE returns the current date and can replace the second component to give:

DateField >= date '2012-05-01' AND DateField <= CURRENT_DATE

However I cannot seem to format the first component such that it always starts at the beginning of the current month. I attempted to format the date string via the CONCAT but I am unable to form valid SQL and I am not sure what the error is. Ignoring the CURRENT_DATE and EXTRACT methods for now gives me this invalid SQL:

DateField >= date CONCAT(CONCAT(CONCAT('2012', '-'), CONCAT('05', '-')), '01')

Ultimately the year and month would be replace by the current year and month using EXTRACT and CURRENT_DATE. However the "simplified" SQL above fails with a generic "An invalid SQL statement was used." error. I know I am doing something wrong but I am not sure what it is.

So, how do you construct a definition query to return only those features in the current month without hard coding the date in ArcMap 9.3.1 in a 9.3.1 File Geodatabase?


texts:
I'm having a bit of an issue loading a WMTS into QGIS 2.8. The application can connect to the WMTS, correctly identify the layers that are available and load the layer into the layer pane. However, nothing draws. It seems to try to load something, but nothing appears on the mapping screen.

It is replicable using the following WMTS:

http://opencache.statkart.no/gatekeeper/gk/gk.open_wmts?Version=1.0.0&service=wmts&request=getcapabilities

I'll log a bug if need be, but wanted to see if anyone had a workaround or had found this as an issue already.


texts: GetFeatureInfo returns "no layers queryable"

texts:
I'm not able to get information with "GetFeatureInfo" request. I'm getting following reply:

"Either no layer was queryable, or no layers were specified using QUERY_LAYERS"

I made sure that my WMS (https://www.geoservice.apps.be.ch/geoservice2/services/a42geo/a42geo_hoehenwms_d_fk/MapServer/WMSServer?service=WMS&version=1.1.1&request=GetCapabilities&) is set to be "Queryable".

Using QGIS I manage to get attributes without a problem.

Are there any settings in my service that I may forget to enable in order to use the request "GetFeatureInfo"?


texts: Remove unwanted white bands when raster layer are overlapping

texts:
In my project when raster layers are overlapping a white band appears the edge of the top layer.

How can I remove it and obtain a seamless overlapping of the raster ?

(I have try to built a vrt with those layers but the with bands remained)

I am working with QGIS 2.4

Here is a example of the phenomena (each layer has 10% transparency)

Here is the gdalinfo output of one of the raster

Driver: GTiff/GeoTIFF
Files: /vdata/scot_009h_georef.tif
Size is 5985, 2349
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
Origin = (-3.873394354081865,58.119237085214102)
Pixel Size = (0.000601740178426,-0.000601740178426)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  -3.8733944,  58.1192371) (  3d52'24.22"W, 58d 7' 9.25"N)
Lower Left  (  -3.8733944,  56.7057494) (  3d52'24.22"W, 56d42'20.70"N)
Upper Right (  -0.2719794,  58.1192371) (  0d16'19.13"W, 58d 7' 9.25"N)
Lower Right (  -0.2719794,  56.7057494) (  0d16'19.13"W, 56d42'20.70"N)
Center      (  -2.0726869,  57.4124932) (  2d 4'21.67"W, 57d24'44.98"N)
Band 1 Block=5985x1 Type=Byte, ColorInterp=Red
  NoData Value=0
Band 2 Block=5985x1 Type=Byte, ColorInterp=Green
  NoData Value=0
Band 3 Block=5985x1 Type=Byte, ColorInterp=Blue
  NoData Value=0


texts: Use color map from .tif.aux.xml from ArcGIS in QGIS

texts:
I have a raster file (.tif) of different ecosystem types (from natureserve). Each value corresponds to a different type. I'm trying to get a table of those attributes and colors. The file includes (what looks like) an attribute table with a color map in a .tif.aux.xml file.

Opening the xml file in a text editor shows entries such as: GDALRasterAttributeTable tableType="thematic". and then rows and rows of the different raster values and their attributes.

How can I open this in QGIS or convert it to a format I can use?


texts: MapServer Runtime Substitution | Omitting WHERE (in DATA) if no variable present

texts:
I am using MapServer's RunTime Substitution.  They are a great feature.

In my case, I have layer setup pointing at a PostGIS table.  This layer is exposed via WMS.  I allow an optional query string variable in GetMap Requests (handled via Runtime Substitution) that ultimately is used in the Data SQL Query

VALIDATION
    "health"  "[a-zA-Z0-9,_]*$"
    "default_health" "'Good','Fair','Poor','Unknown'"
END

DATA "geom FROM (select * from foo where health IN (%health%)) as q using unique id using srid=4326"

In the above example, if no health variable is passed in via Runtime Substitution - I have a default that just builds an IN clause that encompasses all possible health values.  And this works.  My problem is I have other layers setup like this.  But these layers have many, many possible values.  So I end up having these monstrous IN clauses.

I am hoping someone can suggest a way that if a variable is not passed in via Runtime Substitution - I totally omit the where clause in my DATA line.  Is there some way to do a conditional IF in the MAPFILE?


texts: Storing GPS tracks with timestamps in PostGIS

texts:
I need to store a series of GPS points with timestamps in the database (traces of various vehicles).

I am wondering which of two approaches to choose for storing them in the PostGIS table:
1.) Each point in own row with separate column for timestamp and one for track ID

2.) Use Linestring or MultiLineString with fourth M coordinate which would contain UNIX timestamp of the time for the point on the track. This would mean single row for each separate track.


texts: Create a raster with georeferenced information in R

texts:
Referring to my previous question, how do I "geocode"/convert - using R - a generated raster image to a format that can be shown as a layer on a map? In short, I'm aiming to recreate the results of the commercial SpatialKey (screenshot) software.

I'm using a GeoServer instance to host the data for the maps, but how can I incorporate the required projection and coordinate info into the image so it will align correctly with its true geographical position. The density heatmap was produced using spatial data.


texts: Grouping neighbor points based on field sum in QGIS

texts:
I have about 14,000 points. These points have a specific field "time", an integer value, mainly between 1, 2 or 3.

I would like to split these 14,000 points into groups based on the "time" field totalling up to 1,000 hours between the neighbours.

So some groups may have 700 points and others may have 500 or whatever it takes to add up the "time" field value to 1,000.


texts: Distinguishing between table and feature class using ArcPy, when exporting from ArcSDE

texts:
I'm working on a project where I'm pulling some data out of ArcSDE and putting it into a file geodatabase for a group of people to look at with a custom map. I'm writing a script using the featureclasstofeatureclass function to pull out the features. Pretty basic stuff, however when trying to run the script I'm getting the " does not exist or is not supported" error.

After checking for all the basic errors the only thing I can figure is that the script is failing because in our main SDE connection, every feature has a corresponding table with the exact same name. I tried running the script on my edit version for the same feature, where there aren't any same-named tables, and it ran fine. All I changed was the SDE path and it worked.

Is the script failing because of these duplicate names or is there some other possible explanation?

If this is the issue can I use the describe function to get the script to select just the FC and not the table?

What code or function do I need to use to check for tables or select only features?

Right now in the code I'm just copying features because I'm not sure what I need to do next. I know its probably something to do with Describe or listfeatureclasses I'm just not familiar enough with these functions to know how to implement them.

Code snippet:

#############################

#### Update Map Data Script

#############################

#### Imports:

import arcpy

#############################

#### Local variables:

warehouse = "C:\\Users\\aaronmanuel\\AppData\\Roaming\\ESRI\\Desktop10.0\\ArcCatalog\\Warehouse.sde"
meters = "\\coagiswarehouse.coagis.coa_water_meter"
MapData = "G:\\WATERMNT\\AaronManuel\\MapData\\MapData.gdb"
MapData_old = "G:\\WATERMNT\\AaronManuel\\MapData\\MapData_old.gdb"

############################

### GDB Management: Create Backup GDB and create new one

try:

    if arcpy.Exists(MapData_old):
        arcpy.Delete_management(MapData_old)
        print( "Old GDB deleted..." '\n')

    if arcpy.Exists(MapData):
        arcpy.Rename_management(MapData,"MapData_old.gdb")
        print( "Renamed MapData GDB..." '\n')

    arcpy.CreateFileGDB_management("G:\\WATERMNT\\AaronManuel\\MapData", "MapData.gdb")
    print( "Created New GDB..." '\n' '\n')

except:
    print arcpy.GetMessages(2)

#### Processes: Export features from SDE

# Meters Export
arcpy.FeatureClassToFeatureClass_conversion(warehouse + meters, MapData, "meters_test")


texts: Import CSV: Is there a cap for number of row

texts:
I am trying to import several CSV files which contain more than 1,000 rows respectively.
However, but only 129 records were only imported each file.

Can't see any options that deal with this?

I am using Mac(10.7.5), QGIS 1.8.0-Lisboa, GDAL/OGR 1.9.1.
Anyone can advise me?

Many thanks.


texts:
Being aware that there were already similar questions concerning the output of hillshade/slope layer with an unwanted grid texture, but still I m confused, so I hope someone can explain it.
Working with ArcMap 10.7. / QGIS (Madeira), I made a SLOPE layer out of a DEM ( cell size 1x1 m) and got a slope raster with a grid texture , which I don't want because I need to do further analyses with MAXENt.
What could be the reason for such texture and how can I eliminate it?


texts: Calculating Volume differences between two DTM from different years with QGIS?

texts:
I need to calculate the volume of the differences between two DTM from different years with QGIS.  I am trying to figure out how to do it, but I do not see anything clear about that.

Could anyone give me a tip to be able to do it?


texts: Using QGIS calculate distance of point along a line

texts:
I have lots of points, and a road (lines) and I would like to calculate from point A to point B (the distance must route through the roads (lines).

I found this QGIS calculate distance of point along a line , but when I'm trying to do a v.distance, I don't see it saving the distance into the DB.

I need the DB to show:

Source      Destination       Distance
Point A     Point B          10.0011
Point A     Point C           101.12
point A     Point D            1212
Point B     point c           121.0
point B     Point D         31.0
...         ...                ...

Any idea?

Thanks!


texts:
Is there any way to use a map journal offline? That is assume all resources are local, like maps on an internal server and images hosted on the local pc.

I've downloaded the developer code of the Map Journals and been playing around for a while, but really, there is no way to alter the code of the actual Map Journal content outside of on ArcGIS Online.

Any idea?


texts:
I want to use a folder with cached bing map tiles for using it with openlayers 2.
After a google search, i haven't find any example on how to use the tiles locally. The tiles are saved as ".png" in my folder. For example the following image is saved as 1221002.png

The only example that i have found is this link from openlayers site
https://wiki.openstreetmap.org/wiki/OpenLayers_Local_Tiles_Example

with the following piece of code

// This is the layer that uses the locally stored tiles
 var newLayer = new OpenLayers.Layer.OSM(
             "Local Tiles",
             "tiles/${z}/${x}/${y}.png",
             {numZoomLevels: 19,
              alpha: true,
              isBaseLayer: false}
         );

What i have tried so far is to set my code as follows

var newLayer = new OpenLayers.Layer.OSM(
         "Local Tiles",
         "tiles/${z}${x}${y}.png",           <--- Slashes removed
         {numZoomLevels: 19,
          alpha: true,
          isBaseLayer: false}
     );

but with no luck.
Is this the correct approach?


texts: GeoServer WFS 2.0.0 get number of features

texts:
I am writing a Python script to fetch data from GeoServer WFS 2.0.0.  To avoid downloading too large an amount of data I am using pagination with count and startIndex. The problem I am facing is how to know the amount of data in feature layer. In the JSON sent by the server a field called totalFeatures is included but this figure change in the last page.

Is there any way to get the number of features in advance?

resultType=hits& does not work in version 2.0.0


texts: Exporting all features from file geodatabase to PostgreSQL database to specified schema

texts:
How can I export features from Esri file geodatabase (.gdb) to PostgreSQL database, to specified schema, with specified geometry column name?

Now I am doing this that way:

os module to get all gdb from specified dir,
ogr2ogr to get all the features from gdb and export them to PostgreSQL,
psycopg2 to update current schema ('public') to disired one and change geometry column name to 'geom' like all the other tables in PostgreSQL database.

Can I do the same using only ogr2ogr or it's another way?


texts: Problem with performance when loading PostGIS layer using pyQGIS

texts:
I'm having a question on loading a PostGIS layer into the project from my QGIS plugin.
The code is like this:

uri1 = qgis.core.QgsDataSourceURI()
uri1.setConnection(str(host), str(port), str(database), str(user), str(password) )
uri1.setDataSource(str(schema), "layer", "geom")
qgis.utils.iface.addVectorLayer(uri1.uri(),"layer", "postgres")

It actually works fine but the performance is causing me lots of troubles. It takes up to 10 minutes, since there are about a 100,000,000 entries in the table.

When I load the layer manually using the PostGIS button on QGIS, there is a checkbox in the settings which says: "Geschätzte Tabellenmetadaten nutzen". In english this would be: "Use estimated table metadata".
When checked, the layer loads in realtime.

My question is: Is there a way to load a PostGIS layer into QGIS using the estimated table metadata, and thus increase performance significantly?


texts: Joined MSSQL Tables are not editable in QGIS?

texts:
I am trying to join two tables (MSSQL tables) in QGIS and then match some fields and insert/update them.

However, when I join two tables from DB Manager and show them as a new layer, I cannot edit any data.

Does anyone know why I can't edit, and how to make it editable?


texts: Clipping several raster datasets with mask?

texts:
I would like to clip many images (raster datasets, input) with a mask (shapefile) for my study area. Is it possible to do with ArcGIS 10.3 for Desktop?

I tried to merge all the images but it is too much for the computer (152 GB); that is why I was thinking in clipping the rasters directly.


texts: ArcGIS: Locate features along routes, line vs polygon

texts:
I want to use Locate Features along Routes to locate the polylines below:

The routed feature is in red and I want to locate the green lines. I have buffered these but found it only locates the item when the buffer intersects the route. Would running this just on the liner items with a locate radius of ~9m still have the effect I'm after? I have >500 features I want to locate and I am asking because I need to be very accurate and do not trust manual validation alone to be sure.


texts: Using the function “directionalDistanceTransform” from Google Earth Engine on local data

texts:
I am trying to implement the function directionalDistanceTransform from Google Earth Engine API to run it on my local GeoTIFF, so I was looking at Google Earth Engine's GitHub project to find that function and copy it on my python code so I can run it locally but I can't find it.
Is there somewhere I can find the source code or a similar function from other library instead of Google Earth Engine API?
I found a similar function from OpenCV: distanceTransform but it doesn't receive an angle as parameter so it isn't the same.
I know I can upload the GeoTIFF to Google Drive and use Google Earth Engine API from there but I want to launch it in my python environment.


texts: Converting raster to polygon: vertical lines

texts:
I'm converting a reclassified raster (with a resolution of 25.18m) into a polygon. The reclassified raster has values of only 0s and 1s. However, the conversion results in a polygon with vertical lines of missing values. Would anyone happen to know why this is happening, and how to fix it?

A few things:

I've tried using geographic and projected coordinated systems
The original file is a land use raster
Tried this on ArcGIS and QGIS.


texts:
When importing data from OSM using osm2pgsql in slim mode, the resulting nodes table has two columns, named 'lat' and 'lon', which contain coordinates for each node. The values are large (~9 digit) integers, not lat-lon, values and the projection is not given anywhere. e.g.

lat,lon
395980942,-860674711
394871880,-860148511
394879027,-860153601
382614407,-857420524
382613507,-857418018

These coordinates do not appear to be effected by setting a projection in osm2pgsql, though that does effect geometry fields in other tables. The nodes table doesn't have a geometry column, and I am simply trying to create one.

I have tried assuming that these coordinates are in osm2pgsql's default projection, spherical mercator, EPSG:3857.

ALTER TABLE nodes ADD COLUMN geom geometry(POINT,4326);
UPDATE nodes SET geom = ST_Transform(
    ST_SetSRID( ST_MakePoint(lon,lat), 3857),
    4326
);

This doesn't throw any errors, but gives lat, lon values that are clearly wrong, and which fail to render in QGIS. I read somewhere that the coordinate values needed to be divided by 100 (I can't find the source now...):

UPDATE nodes SET geom = ST_Transform(
    ST_SetSRID( ST_MakePoint(lon/100.0,lat/100.0), 3857),
    4326
);

And this seems to get closer. QGIS renders the geometry now, and I can see the right shapes in the data, but it's still off by a few hundred miles from where it should be.

To the northwest is the points table projected correctly. The other dataset is the nodes table, as projected in the preceding code block.

What is the projection of these coordinates?


texts: Looking for TSP with Target different from Source

texts:
I know the TSP (traveling sales person) in osm2po is beta but it works quite well.

Is there a way to designate a final target location that is different from the initial source location. In other words, I would like to designate:

Start location
Several intermediate locations
End location (different from start location)
Determine the best way to begin at the Start location, travel to all the intermediate locations while hitting the End location last.

Currently, I see the TSP only does a round trip from an initial start location.


texts: Generate XYZ coordinates of DEM mesh/grid

texts:
Is there a way to generate a wire grid of plain text XYZ coordinates starting and ending points for each wire segment? I want to generate them from a digital elevation model. Preferably with gdal/qgis and python. I have experience with gdal and python, I just don't know how to get started in this case.

I'd like to import the ground into an antenna modeling software.


texts: QGIS Creating 4km point grid within polygons

texts:
I'm working in QGIS right now to create a grid of points within ~36 irregularly shaped polygons. The spacing on the of the points within in the polygons needs to be 3.8/4km.

I attempted using Vector -> Research Tools -> Vector Grid to create an evenly spaced point grid within the polygons but even on the larger polygons that are 30-40 km^2 the tool only creates one point in the very top left corner of the bounding box. Once I create these grid points I'll be exporting them out to KML or XML for use by a 3rd party program that doesn't support shp. The grid will create a "tiling" effect for spot images over an area for our third party software.

I've looked at several of the stack exchange discussions here and none seem to get it quite where I need it.

I was considering creating a 4km fishnet grid over the entire extent (semi-global) and then clipping it to the polygons, then creating centroids for those remaining polygons to produce my grid.

I looked at:

Creating grid constrained within polygon in QGIS?

Creating evenly distributed points within an irregular boundary

Creating regularly spaced, defined # of points within polygon in QGIS

And those have somewhat led me to using the fishnet/centroid solution. I'm wondering if there is a more practical, scriptable way to do this. I'm new to use python in general but am making myself learn it for obvious reasons!


texts: Can *.ssf files from Trimble GPS be imported as features to display in ArcMap?

texts:
I have not come across *.ssf files before but I believe they are a Trimble GPS binary format.  Someone has asked me if these files which he thinks have line/poly features can be imported for display in ArcMap 10.0.

I am keen to know whether this is possible because my research so far seems to indicate that they are for correcting feature accuracy rather than containing features.

Can *.ssf files be converted for display in ArcMap?


texts: What is the state of the art in html5 geospatial applications?

texts: I am very interested in the HTML5 canvas element for maps; work like Cartagen from an MIT Media Lab researcher looks very promising, for example.  There is some interesting SVG based work at carto.net.  WebSocket seems like a fantastic API for live geospatial data.  I have done a few HTML5 experiments, for example here, with Flickr and Picasa data.  What are people getting up to with these new technologies, or what have you tried?


texts: Finding matching rows across columns and populating a separate column

texts:
I am trying to write a script that recategorizes 70 attributes in to 1 of 5 attributes.

I have two tables, each with two columns that I am focusing on. The first table has one column with 20k rows that all contain 1 of the 70 strings.

The second table is the key used as a key. The first column has 70 rows, each listing one attribute type. The second column tells what that attribute converts in to. There are 5 options, so each of the 70 will fall in to one of those five.

What I need to script to do is read through the key one row at a time, then find all matching rows in the other table's column. For all of those matches, a second column in the original table needs to be populated with the corresponding attribute (from the 5 options).

I have trouble explaining this, but it is really just recategorizing 70 attributes into 5 attributes. These are all strings. The original table is a layer in ArcMap, and the key table is an excel file. I am not sure how to write this script.

Here is what I have so far:

import xlrd
from xlrd import open_workbook

xl_workbook = xlrd.open_workbook(conversion_table)
xl_sheet = xl_workbook.sheet_by_index(0)
num_cols = xl_sheet.ncols   # Number of columns
for row_idx in range(0, xl_sheet.nrows):    # Iterate through rows
    print ('-'*40)
    print ('Row: %s' % row_idx)   # Print row number
for col_idx in range(0, num_cols):  # Iterate through columns
    cell_obj = xl_sheet.cell(row_idx, col_idx)  # Get cell object by row, col
    print ('Column: [%s] cell_obj: [%s]' % (col_idx, cell_obj))

# print water land use category column
conversion_column = xl_sheet.col(0)
list1 = list()
list1.append(conversion_column)
print list1

# print land use field from roadway
fc = original_table
field = "LAND_USE"
cursor = arcpy.SearchCursor(fc)
row = cursor.next()
list2 = list()
while row:
    results = (row.getValue(field))
    list2.append(results)
    row = cursor.next()
print list2

I'm using version 10.3.1.


texts: Creating heatmap (concentration of points) in QGIS?

texts:
There was already How to build effective heat-maps? from 2010 where the answer was basically "use QGIS Heatmap plugin".

I've also read on QGIS webpage that it is now "implemented as a native C++ plugin"

What is the best way to generate a heatmap from a list of observations, in QGIS?


texts: Get data for ArcPad validation error on unfulfilled relationships

texts:
I have a file GDB with a handful of feature classes that have related tables that I am trying to check out to an AXF. I can check out an empty database just fine, no errors reported upon validation. However, when I add some sample data so that I can begin field testing, for every unfulfilled relationship (so any 1:0 relationships) I get a validation error. It's as if ArcMap is treating the relationship as a mandatory relationship.

Some additional information: I'm running ArcMap and ArcPad 10.2. All relationships are  1:M type and I've allowed Null values for each field involved in the relationship.

Anyone have any insight as to why this is happening? Are relationships mandatory?

UPDATE: I've narrowed my error issue down to this type of relationship (see attached) I'm trying to create relationship rules depending on a subtype in a table. So similar to the riverside data's poles and lights examples except for when the subtype is filament or tube I want it relate to poles using Polie_ID_PK and Pole_ID_FK and then when the subtype is "energy saver" I want the relationship to be based on Pole_ID2_PK and Pole_ID_FK.  This is the type of relationship that keeps throwing the validation error. I don’t seem to run into the same trouble if it is the origin class that has the subtypes. Does anyone know why or how this type of relationship class is not valid?


texts: Automating powerpoint slides creation with maps for bulk jobs

texts:
I have about 50 states in US which are subdivided into 200 territories. Often I have to generate powerpoint presentations consisting of map views for all or some territories by going back to each territory and taking a snapshot for each of them. Many a time I need to repeat this task for different layers. This is very time consuming.

Is there any way I could just specify the states and the territories QGIS should generate PowerPoint slides using given layers?


texts: How can I validate a table name on SDE with Oracle backend?

texts:
Second title: "Is arcpy.ValidateTableName broken on Oracle?"

I have an Oracle SDE database 9.3.1 (service connection). My arcGIS client tools are version 10.0 SP5). I have some names of feature classes(fc) that I would like create.

I want to make sure that the name of the feature class is valid before I create it so I use the arcpy.ValidaTableName function.

However when I try to validate a long fc I do NOT get a valid table name for Oracle (Oracle has a table identifier character limit on 30).

import arcpy
oracleWorkspace = "F:\shr016\all_is_sde_databases\esri_mogpros_5167_gisowner.sde"

print arcpy.ValidateTableName("very_long_table_name_with_too_many_characters_for_an_oracle_database", oracleWorkspace)

# prints: very_long_table_name_with_too_many_characters_for_an_oracle_database

Is there way I can validate table names on SDE with a Oracle backend?

Update

It seems like this feature is indeed not working correctly.

Esri support writes in an email:

I have tested with arcpy.ValidateTableName with an Oracle database and it does behave as you have mentioned in your previous email. Currently, I am not aware of any alternate way of checking if the feature class name is valid from the Oracle side. I will be investigating this issue further and will update you with my findings.

Update 2

Esri support has followed up:

This tool was designed to check for naming conventions, i.e. presence of any invalid characters or use of any reserverd keywords specific to the RDBMS being used. For checking the character length, I am contemplating logging an enhancement request. In the meantime you could include something like the following in your script to make sure the table name does not exceed 30 characters:

import string, os

test = "tableName"
oracleWorkspace = "F:\shr016\all_is_sde_databases\esri_mogpros_5167_gisowner.sde"

if len(test) <=30:
  print arcpy.ValidateTableName(test, oracleWorkspace)
else:
   print "Invalid table name"

Conclusion:

If you need a function that returns a valid table name for an Oracle SDE table name you have to write it yourself.

Update 3

A final update from esri support which confirms the conclusion:

I have logged an enhancement request for this incident: NIM090928-arcpy.ValidateTableName should check that table names do not exceed the workspace's maximum number of characters.


texts: ASP .NET and ArcGIS API for Javascript : FeatureTable won't show up

texts:
I"m a newbie with the ASP .NET, a have an issue. I want to show up a feature table in a viewer but I don"t get any results.

I think I"m missing something?

Here is example that I took from the ArcGIS API example and include it in my project :

    <link rel="stylesheet" href="https://js.arcgis.com/3.20/dijit/themes/claro/claro.css">
<link rel="stylesheet" href="https://js.arcgis.com/3.20/esri/css/esri.css">

The FeatureTable.


texts: Converting and clipping nc4 files into raster file using QGIS?

texts:
Currently I am working on organizing the NASA's Gravity Recovery and Climate Experiment (GRACE) data with QGIS 3.10 and wish to clip it with the state boundary and rasterize its data so that I can do raster calculation.

It cannot be opened as a raster so I have to use "add Mesh Layer" to open it. I tried Crayfish to edit it but it failed to rasterize (with oversized pixel).

How I can generate the raster with the state boundary of my choice?


texts: Crop a satellite image to 32 by 32 pixel

texts: I have a satellite image of resolution size 10 m. I want to get a cropped image of 32x32 pixels of that image using a polygon. I randomly generate a polygon no matter the size but to reduce the polygon to 32x32 pixel I should add 320m to the x-axis and subtract 320 m to the y-axis. Am I correct or there is another way to reduce the polygon to the size I want. I tried this but my width is always high but my height is really close to 320 m.


texts: Custom Buttons for Google Drawing Manager

texts:
can we design my own custom button for the Drawing manager?

According to API

https://developers.google.com/maps/documentation/javascript/reference#DrawingManagerOptions

we can customize only the symbols which is drawn out of the drawing manager.

for example i need to make this drawing manager look like Google Map maker

http://www.google.co.in/mapmaker


texts: Reshape polygon to boundary of other polygon in QGIS

texts:
I have two polygons.

The first one has a white background and brown boundary.
To other has the pink background.

In QGIS, how is it possible to reshape to first polygon to the boundary of the other polygon?


texts: Raster calculation with string in modeler

texts: In my model I have one raster file and also a string input (factor for further calculation). I tried to use the grass r.mapcalculator where I can find the raster and also the string input. I tried to calculate the following: string + A without success. I also tried 'string'+A or "string"+A. Is it possible to solve this issue?


texts: Measuring distance between two points on curve line in QGIS

texts:
I am using QGIS 3.10 version.
Image those two blue lines are points. How to measure the white line between those points?


texts: How to add multiple rasters to a layer file?

texts:
I'm trying to make a .lyr file to bring into ArcMap with arcpy.mapping. I can bring a single raster into a .lyr with the script below. I'm unsure how to proceed to get multiple rasters from the same folder into a .lyr. Any help appreciated.

Thanks

# Import system modules
import arcpy, os
from arcpy import env

# Workspace Directory
env.workspace = r'F:\Ortho_prelim\TEST'
workspace = r'F:\Ortho_prelim\OUT'

# Set local variables
inLyr = 'ImageLyr'
outLyr0 = 'ImageLyr'
outLyr = workspace + os.sep + 'Image.lyr'

# Make feature layer variables
fList = os.listdir(r'F:\Ortho_prelim\TEST')
for f in fList:
    if fList.endswith('.ecw'):
        print fList
        arcpy.MakeRasterLayer_management(f, outLyr0, "", "")
        print 'Made ' + outLyr0
        print outLyr0
        # Save Feature Layer
        arcpy.SaveToLayerFile_management(inLyr, outLyr)
        print 'Made layer ' + outLyr

So the pseudo code should be something like:

Look in folder
Look up all files with extension .ecw
Add all the .ecw files in the folder to a .lyr file...


texts: Cannot Uncheck Preserve Aspect Ratio of text box in ArcMap

texts:
I am using ArcMap 10.5.

I'm have dynamic text box I'm trying to add to the buttom of my legend. The text is about 350 characters so need to resize the text box to fit into my legend. But, the "Preserve Aspect Ratio checkbox is grayed-out.

Tried to look in ESRI help with no luck and even found a thread on a forum with the same question, but the answer did not help.

How do I make my text box size editable?


texts: What strategies exist to speed up point-in-multiple-polygon searches?

texts:
Consider a scenario where you have polygons representing the border of every country in Europe. You also have the Lat/Long of a user's mobile device. You now want to display which country that user is currently located in.
In essence, you have multiple polygons and a single point, and you'd like to find the polygon that contains that point.
Point-in-polygon algorithms are well known, but what are some strategies to reduce the search time?


texts:
I am working with spatial data UK.

The map is plotted here. I want to resize the map and make it bigger and clear!

here is my code

LA <- readOGR(dsn =".", "infuse_dist_lyr_2011_clipped")

plot(LA)


texts: Merging DEM Files into a Single DEM and Exporting Cropped DEM

texts:
I am looking to combine .DEM files into a single .DEM file and export a crop of the center.

My site lies across 4 DEM files, about a quarter in each. It only occupies a small portion of each. I am using an outside program to convert DEM files into 3D data for Rhino to use as a NURBS site map. The problem is the program will combine 3 of the files, but one of the DEMs will not combine with the others reading "Out of Memory." I am thinking that if I crop just the site it will hopefully upload- as it is a much smaller area than any one of the DEMs- but again it exists across 4 separate DEM files.

Essentially, combine all 4, export out a small portion of the middle as a .DEM.

The program only works with .DEM extensions and not GeoTIFFs.


texts: Avoiding ImportError: No module named numpy?

texts:
When I try to import numpy, I get the following message.

>>> import numpy
Runtime error
Traceback (most recent call last):
  File "<string>", line 1, in <module>
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\__init__.py", line 24, in <module>
    from arcpy.toolbox import *
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\toolbox.py", line 342, in <module>
    from management import Graph, GraphTemplate
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\management.py", line 22, in <module>
    import _management
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\_management.py", line 14, in <module>
    import _graph
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\_graph.py", line 27, in <module>
    import numpy
ImportError: No module named numpy


texts: How to extract subdataset from HDF raster?

texts:
I am trying to extract a subdataset from HDF raster. Later on, I intend to make operations such as mosaicking and extracting raster averages from the extracted subdataset. The following code is supposed to extract the subdataset, but it seems to me that all it does is extracting the subdataset NAME and not the object itself:

import gdal
from gdalconst import *
filename= "Image.hdf"
dataset = gdal.Open(filename,GA_ReadOnly)
subdataset = dataset.GetSubDatasets()[1]
print(type(dataset))
print(type(subdataset))

Now, if you run this, the type of the dataset and subdataset objects will be printed out:

<class 'osgeo.gdal.Dataset'>
<type 'tuple'>

So, the extracted subdataset is just a tuple. How can I use it to do mosaicking? Am I missing something here?
P.S. HDF rasters can be easily downloaded from this FTP site: ftp://ladsweb.nascom.nasa.gov/allData/5/MOD11A1/2012/193/


texts: I want to normalize dem images that min= -49 and max=229 to -1 and +1 in ENVI. what is that formula in band math?


texts: Is the geometry type of PostGIS represented as WKT or just the Point type?

texts:
Despite reading a lot of materials on the subject, including the PostGIS manual I still haven't been able to understand, I'm a beginner in this world of geographic objects. I have two tables in my database that have a column called the_geom. Table ways and ways_vertices_pgr.

ways: the_geom -> geometry (Point, 4326)

Example of the_geom data with a select in a column:
POINT (-48.2486036 -18.8790378)

ways_vertices_pgr: the_geom -> geometry (LineString, 4326)

Example of the_geom data with a select in a column:
LINESTRING (-48.2662325 -18.9026985, -48.2657083 -18.902022)

The geometry type is a PostGIS spatial data type used to represent a feature in coordinate systems. OpenGIS Consortium has defined two ways to represent spatial objects in a bank, WKT and WKB.
In my example, I am using WKT and not EWKT, correct? Because I'm using a Cartesian coordinate system, 2D plane. From my database, is it correct to say that the geometry type is represented in my database as a WKT text representation? With Point and SRID? Explaining better, is it correct to say that the entire the_geom column is represented in WKT format?
Or not, the WKT text representation is only for the POINT type, and does PostGIS extend this pattern by adding the SRID?


texts: Where can I download the Google Road network as shapefile?

texts: Is there any site from where I can get Google Road network as shapefile? As per my knowledge It is not distributed by google. So anybody can suggest how can I create this road network for a custom use?


texts: What was the rationale for the non-standard UTM zones near Norway?

texts:
UTM zones 30V, 31V, 32V, 31X, 33X, 35X, and 37X do not line up with the rest of the UTM zone grid.

Why was this allowed? According to Wikipedia, UTM was developed by the U.S. Army Corps of Engineers in the 1940s. What prompted them to have special zones for Norway and not for any other country?