texts: QGIS - qgis2web - IOError: [Errno 9] Bad file descriptor

texts:
I use QGIS 2.18.3. I am trying to use qgis2web plugin to send my map to the web (very simple map: 1 layer polygons, 1 layer points, 1 tile) but I get the following error:

Traceback (most recent call last):
File "C:/Users/pc/.qgis2/python/plugins\qgis2web\maindialog.py", line 162, in previewMap
MainDialog.previewOL3(self)
File "C:/Users/pc/.qgis2/python/plugins\qgis2web\maindialog.py", line 369, in previewOL3
cluster, params, utils.tempFolder())
File "C:/Users/pc/.qgis2/python/plugins\qgis2web\olwriter.py", line 302, in writeOL
print traceback.format_exc()
IOError: [Errno 9] Bad file descriptor


texts:
I already have the SNAP 7, and I can’t run sen2cor 2.08 directly from it. However, I can alternatively use the command prompt but it can only take the L1C images as input.
My problem now is that I want to eliminate clouds from my Sentinel-2 L2A image. My advisor wants me to do it with sen2cor. I know that I can always find new images cloud-free from the scihub_copernicus website, but he wants me to figure it out based on the L2A image. Any suggestions?


texts: ArcGIS Viewer for Flex - map won't display in preview tab

texts:
I'm new to ArcGIS Viewer for Flex and ArcGIS for Server, so far I have been able to get past the point where I get to choose my base map and operational layers from the ArcGIS for Server 10.2. I can see my map on the map extent tab, but when I go to the preview tab, all I see is information about the layers and other stuffs in the service directory. I am unable to see the map itself. I am using ArcGIS Viewer for Flex 3.6.

My setting page

web server base folder: - C:\inetpub\wwwroot\flexviewers
web server base folder url: - \http:/loc:6020/arcgis/rest/services/Jegol_web_mapping_trial/MapServer

My map page

I called the map from the same page as the web server base folder URL which is : - \http://loc:6020/arcgis/rest/services/Jegol_web_mapping_trial/MapServer

I don't seem to find where I made the mistake, can anyone help me out please


texts: Query polygons using ogrinfo but not generate a new layer

texts:
I want to list the polygons in my layer A that intersects other polygon B, using ogr.
The desired result is similar to "Select by location" tool in QGIS or ArcGIS Desktop, where I expect an index of table A, and NOT a new layer with the intersected geometry.

I'ld like to get the result like: 1, 4, 6 ... where the numbers represent the rows of A where B intersects, or the table A [ "Rows where A intersects B", ], and store in a file if possible.

Additionally I'll prefer to replace the B layer for a WKT string like
Polygon ((12 13, 11 12, ...)), just to avoid ask a new question.

Currently I'm using this expression using .vrt file:

ogrinfo vrt.vrt -dialect sqlite -sql "SELECT ST_Intersection(a.geometry, b.geometry) AS geometry FROM a a, b b" -geom=SUMMARY

And the results:

Layer name: SELECT
Geometry: Unknown (any)
Feature Count: 36
Layer SRS WKT:
PROJCS["Transverse_Mercator",
    GEOGCS["GCS_GRS 1980(IUGG, 1980)",
        DATUM["unknown",
            SPHEROID["GRS80",6378137,298.257222101]],
        PRIMEM["Greenwich",0],
        UNIT["Degree",0.017453292519943295]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",4.596200417],
    PARAMETER["central_meridian",-74.07750791700001],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",1000000],
    PARAMETER["false_northing",1000000],
    UNIT["Meter",1]]
Geometry Column = geometry
OGRFeature(SELECT):0
  MULTIPOLYGON : 3 geometries:
POLYGON : 168 points
POLYGON : 22972 points
POLYGON : 1990 points

OGRFeature(SELECT):1

OGRFeature(SELECT):2
....
OGRFeature(SELECT):17
  MULTIPOLYGON : 3 geometries:
POLYGON : 1251 points
POLYGON : 544 points
POLYGON : 15937 points


texts: How to remove false strip lines from UAV LiDAR data

texts:
I have a multiple return UAV LiDAR point cloud and it has false strip lines 0.5-1 m below ground. I think two different flight lines are not overlapping perfectly causing these strip lines, but not sure. Strip lines consist of many points, therefore I was unable to remove with lasnoise.

And lasground ground classification algorithms are picking up these noise points as ground points, as a result, correct ground points classified as non-ground points. During normalization of point cloud, height of correct ground points are calculated based on the falsely calculated ground points and shifting upwards which is causing me a problem.
Ground points:

And Vegetation points:

And  I created raster image from points below 1.3 m to show the problem clearly:

I have 100 1 ha plots and most of these plots have this issue. My question is how to get rid of these below ground strip lines?


texts: Leaflet MarkerCluster custom info window instead of spiderfy effect

texts:
So I have been experimenting with the Leaflet Marker Cluster and was wondering if there is a way to disable the spiderfy animation and change the event when cluster is clicked on? Instead the clicking of cluster would lead to an info window showing list of names(a list of each name of event within cluster, e.g. the cluster is 7 points stacked, when you click it displays an info window with the name of each event that occurred at the time) attributes instead.
It seems that setting spiderflyOnMaxZoom option to false seems to shut down the spiderfy animation.


texts: Missing some images of Landsat 8 Surface Reflectance data for some regions in google earth engine

texts:
I want to have the whole image collection of Landsat 8 Surface Reflectance data for my region (Path: 196, Row: 56) from 01.01.2014 to 31.01.2020. But GEE gives me only 66 images for this period while I can get 100 images on EarthExplorer (https://earthexplorer.usgs.gov/) for this period! The question is why GEE does not have the whole archive? I also checked for some other regions and I got the same problem but for the areas close to the equator there is no problem. is this a bug related to GEE?
It should be mentioned that I have not applied any filter as I need the whole archive regardless of clouds or etc.

here is my GEE script: https://code.earthengine.google.com/f954255b2565a6e4dfd807e4fa05af6b

It's appreciated if someone answers my question.


texts: Update base layer on toggle button change in OpenLayers

texts:
I have an html toggle button and try to toggle update between two layers (baselayer and earthlayer) when it's changed.

I am assigning a variable and change the value to the corresponding layer. Perhaps something is wrong and the layer is not getting updated.

var MapView = baseLayer;
$("#mapview").change(function () {//when toggle button is changed
  if ($(this).is(':checked'))
    MapView = earthlayer;
  else
    MapView = baseLayer;
});

var map = new ol.Map({
  target: 'map',
  layers: [MapView],
  view: new ol.View({
    center: ol.proj.fromLonLat([33.33386, 35.14710]),
    zoom: 14
  })
});


texts: ArcToolbox model run but warning output in Python

texts:
I have created an ArcToolbox and a Model in the toolbox. It uses a lot of Spatial Analyst tools and deals with raster images.

When I imported the toolbox in a Python script and ran the model from Python in Windows command prompt, I got a few warning messages as below but the model got run successfully.

<string>:1: SyntaxWarning: import * only allowed at module level

It was shown several times. I am unable to tell where this error comes from.

Any idea?


texts: Extract by mask: extract only those cells not covered by mask

texts:
I have two raster files of the same resolution: raster A and raster B. Raster A is the bigger raster and raster B covers some parts of raster A. I want to extract those grids from Raster A which are not covered by raster B.

I usually do using extract by mask but this time, it is the other way round. Can anyone tell me how can I do this.


texts: Finding routes in which more than one car passed (gps data)

texts:
I'm using ArcGIS 10.1. I have point GPS data of routes of several cars over a few days,  and I want to find only the routes where more than one or two different cars passed through (Each GPS points has the car id in it). How can I find a thing like that without going through the data manually ?

EDIT: The time interval (and the distance) between points differs - it moves from minutes to hours and from a few meters to tens of kilometers apart. I'm comparing routes right now by viewing the points by the eye (which isn't accurate or useful in any way). I thought about creating routes using network analyst, but as I have a lot of points I kind of hoped to fan out the less popular ways first, and calculating routes for only the popular ones...


texts:
We've been trying for a while to get the CartoDB tech stack working and we've made progress but getting stuck when trying to import data via the API. The logs show the POST working fine and we're getting this in the DB:

data_source            | /usr/local/src/cartodb/public/uploads/b24076067e8d1bd6cfd1/twitters.csv
data_type              | file
table_name             |
state                  | pending
success                | f
logger                 | 3525c925-f869-4349-bfd7-46447fea9345
updated_at             | 2015-04-07 04:21:51.628807+00
created_at             | 2015-04-07 04:21:51.62864+00
error_code             |
queue_id               |
tables_created_count   |
table_names            |
append                 | f
migrate_table          |
table_copy             |
from_query             |
synchronization_id     |
id                     | f03efe7d-1ce9-403b-b6f6-7da6147bb1a2
user_id                | f81a3066-7809-4bf7-a201-aebc948603da
table_id               |
service_name           | public_url
service_item_id        |
stats                  | {}
type_guessing          | t
quoted_fields_guessing | t
content_guessing       | f
server                 |
host                   |
upload_host            | ip-10-37-199-244
resque_ppid            |

But we're unclear what's supposed to happen next. The dashboard is showing a "Creating table" popup, seemingly based on the "pending" status, but nothing actually seems to be happening to import the data. I've found some code that seems related to kicking off import jobs (through resque?) but I'm having trouble making sense of it.

Any suggestions where to look, or what might be falling through the cracks?

Side note: it looks like the code eventually changes the status to "failed" after 60+ minutes of sitting in pending.


texts: I need to add a sentinel 2 image into Google Earth Engine Explorer graphical user interface

texts:
I am using the graphical user interface of Google Earth Engine (Explorer) to carry out classification of land cover using sentinel 2. My study site is somewhere in Kenya with an id:

S2A_OPER_MSI_L1C_TL_EPA__20150914T075319_20160317T153700_A001190_T37MDU_N02_01_01

However, the format required is:

COPERNICUS/S2/20151128T002653_20151128T102149_T56MNN

and so I attempted to truncate my id to:

COPERNICUS/S2/20150914T074233_20150914T074233_T37MDU

It is still not being accepted by Earth Engine graphical user interface.


texts: Is there a Reference Layer in the World Imagery Basemap?

texts:
How do I turn on the reference layer in the World Imagery Basemap?

It is not located on the Table of Contents.


texts:
I'm trying to reproject a 1/9 arc second DEM from EPSG:4269 to EPSG: 2264.  I also need to DEM to be in feet so I have been using the raster calculator to multiply the original DEM to feet.  When I assign a new projection to the layer and change the CRS of the project the measurement tool is way off.

before converting to feet and changing CRS

After.


texts: Creating Sentinel-2 NDVI time series chart in Google Earth Engine?

texts:
I am exploring Sentinel-2 time-series NDVI with Google Earth Engine.
In another post (here) I calculated and exported NDVI values from Sentinel-2.
Now I would like to create a NDVI time series chart, but when I run the code (see below) I get

  Computation timed out.

error in the console tab. I have no idea why.

Here is the code:

//Import GEE Feature Collection (Somaliland kml)
    var geometry = ee.FeatureCollection('ft:1bzF2pepsJHoGIXn-ojYYOkOBHWtfhqzS9QLuz_l3');

// Create image collection of S-2 imagery for the perdiod 2016-2018
var S2 = ee.ImageCollection('COPERNICUS/S2')

//filter start and end date
.filterDate('2018-09-01', '2018-10-31')

//filter according to drawn boundary
.filterBounds(geometry);

// Function to mask cloud from built-in quality band
// information on cloud
var maskcloud1 = function(image) {
var QA60 = image.select(['QA60']);
return image.updateMask(QA60.lt(1));
};

// Function to calculate and add an NDVI band
var addNDVI = function(image) {
return image.addBands(image.normalizedDifference(['B8', 'B4']));
};

// Add NDVI band to image collection
var S2 = S2.map(addNDVI);

// Extract NDVI band and create NDVI median composite image
var NDVI = S2.select(['nd']);
var NDVI = NDVI.median();

// Create palettes for display of NDVI
var ndvi_pal = ['#d73027', '#f46d43', '#fdae61', '#fee08b', '#d9ef8b',
'#a6d96a'];

// Create a time series chart.
var plotNDVI = ui.Chart.image.seriesByRegion(S2, geometry,ee.Reducer.mean(),'nd',500,'system:time_start', 'system:index')
              .setChartType('LineChart').setOptions({
                title: 'NDVI short-term time series',
                hAxis: {title: 'Date'},
                vAxis: {title: 'NDVI'}
});

// Display.
print(plotNDVI);

// Display NDVI results on map
Map.addLayer(NDVI.clip(geometry), {min:-0.5, max:0.9, palette: ndvi_pal}, 'NDVI');


texts:
I have two point layers, A and B. For every point in layer A, I want to count how many points from layer B lie within a 1-mile radius, and store that count in a new column of layer A's attribute table.

For context:

I am analyzing home values, and as part of this am learning to use GIS to add data about neighborhood attributes to the data about homes. The relevant layers here are Home_Parcels (a point layer where each point represents a home) and COBRA_2010 (where each point represents a crime).

The layers are large:
36,000 points in COBRA_2010
200,000 points in Home_Parcels

So if this would be too computationally intense, I can use a layer where points represent Census Blocks instead of Home_Parcels, which would bring me down to a couple thousand points for that layer.


texts:
I have 270 000 points on a map.
If you zoom out, the site hangs for a very long time and sometime the page just crash.

I use Cartaro with Drupal.
https://www.drupal.org/project/cartaro
Cartaro is a Drupal module, tha uses PostGis, GeoServer (drupal), GeoWebCache and Openlayers (for Drupal). Database is postgresql. In Ol i used Geoserver WMS for Layer.

This is my first attemt at all of the above, so im a very big noob. Hope the info provided is ok, or else, please ask again

Right now, i have the OL to zoom in automaticly on the user, so not so many points are loaded. I use bounding box, and the page loads instant. But the problem is, if the user zooms out a bit, the server dies or lag severly.

Im new to this, but stumbled over minscale and maxscale.
Could that be my solution ?
How and where do i insert minscale and maxscale ?


texts: Rounding column in attribute table using ArcGIS Field Calculator?

texts:
I want to use the Field Calculator in ArcMap to round an existing column to two decimals. Currently I have a column that is 6 decimal places long and would like to simply round it down to 2 decimals.

I had planned on using the Field Calculator (possibly using Python) to do this, but maybe there is an easier way?

The accepted answer is probably the easiest way to change a single field, however, here is how to do it with the field calculator for both python and VB.

VB:

round([column], 2)

Python:

round(!column!, 2)


texts: Extracting data from SNODAS file using QGIS

texts:
I'm trying to extract data from a raster file from "Snow Data Assimilation System (SNODAS)" and convert it to a MMQGIS (plugin) grid layer (Hexagon) with QGIS.
This raster file stores the snow depth of the whole of the United States.

I load my raster (tif)
I extract the data with a terrain analysis (Aspect)
I create my grid layer as a hexagon
I run Zonal statistics from my Aspect to my Grid
At this point my grid has "count", "sum" and "mean" attributes. I use the "sum" attribute.
I load it on TileMill and I style the "sum" attribute.

I would like to know if I have the correct process and the right information about snow depth on my Grid layer? Is "Sum" the right attribute to print?


texts:
I have a polyline shp file, which has 3 routes. I need to find the overlap segments of any two routes, and then merge all the overlap segments with the original 3 routes to get a new polyline shp file for the further analysis. Now my method is to use two iterations to obtain the results. The following picture shows the method I used for now. In the model, AllRoute.shp includes route0, route1 and route2.

However, I have to use the same method for 500 routes. Can anyone help me with adding another iteration to the existing model?


texts: Values change after re-projecting the NetCDF files to TIFF file

texts:
I am trying to re-project a NetCDF file to a uniform grid EPSG: 4326 and convert it to a TIFF file. However, after the GDAL re-projection the values in the file change, i.e. the values are not retained after the re-projection takes place. The commands that I am using are:

gdal_translate HDF5:"LC08_L1TP_029037_20190611_20190611_01_RT_L2GEN_2.nc"://geophysical_data/Rrs_443 -a_nodata -32767 output.tif

gdalwarp -t_srs EPSG:4326 -srcnodata -32767 output.tif reprojected.tif

The links to the files are:

NetCDF file: https://drive.google.com/open?id=1XKNOqf_-MIcgR-f_qBWP-92ZNP3pEsgj
output.tif file: https://drive.google.com/open?id=1nMcMIRiOVyiwjc8tHegzRyHxn2K9762o
reprojected.tif file: https://drive.google.com/open?id=1Uh_UfWkeiEDuwKvDp4qCTjyBh_Bra6cZ


texts: How to generate TTA mask for accuracy assement in eCognition?

texts:
I'm learning to use eCognition software recently. I have done the image segmentation for Worldview-2 image. Then, I have processed for classification for class names (e.g forest,built-up area, grass land)using Nearest Neighborhood.

Then I have tried to do accuracy assement. But failed. However, how can I generate accuracy assessment using TTA mask using eCognition?
Please help me..Thank you


texts: Making query layer with query parameters from generated values using ArcGIS ModelBuilder?

texts:
I want to join a relationship table (same database, not in SDE) to an ArcSDE feature class to display only needed features. This has to be done multiple times according to every row (value) of an other table (same database, not in SDE too).

The first step works, but I have to make a query layer first as it looks like I can't make joins with tables not in the SDE.

I do not want to create ~500 query layers (for the joins), so I'm trying to use ModelBuilder to create the query layers (and maybe also to make the joins in the future).

The figure below shows the model that I've created so far (at the moment I don't know if this approach would generally work/do its job).

I want to use the value catid to perform the query in the Make Query Layer (2).

How should/could I do that?

Here some pseudo things:

table a; -- contains some categories
table b; -- =m:n relationship table for "a" and "c"
feature class c;
table d_*; --contains all data from "b" related to a single category (row) of "a"

foreach row in a
{
 create table d_*;
}

In psql this would be kind of simple, but I need this done on the ArcGIS level as the join of the feature class relates on this.

I was also thinking about python, because using a programming language this would also be kind of simple -- but I have no idea about python and I thought I could avoid scripting.


texts: Importing data from Oracle Data Pump into empty Oracle ArcSDE Geodatabase?

texts:
I have an ArcSDE Geodatabase that contains empty layers created using ArcCatalog.

I would like to import data into those layers from a Data Pump of an Oracle database that has the same structure as the empty database.

How can I perform the spatial data import without getting problems linked to SDE?  Is it dangerous to use Oracle data import commands?  Should I use some special SDE command to perform the import correclty?


texts: Calculating length along multiple line segments in GeoPandas

texts: I would like to calculate the distance between two points. They are connected with multiple line segments. Each LineString only consisting out of two points. How do I connect these LineString segments between my start and my end point to one single LineSting (of which I can calculate the length).


texts: How to extract GeoJSON database column as geometry layer? (column as layer)

texts:
Have a one table of country with multiple categories for states. Table as below:

Country table:(converted shapefile to GeoJSON for use in Leaflet choropleth map)

So kind of table converted from shapefile to GeoJSON using QGIS for further procedure.

How can I convert these categorical columns into categorical layers for choropleth mapping?


texts: OpenLayers interaction is not working when initializing in map object

texts:
There is a problem to set interaction for the map object in OpenLayers. I've used this code but it's not working:

var osmLayer = new ol.layer.Tile({
    source: new ol.source.OSM()
});
var interaction = new ol.interaction.DragRotateAndZoom();
var center = ol.proj.transform([51.34, 35.65], 'EPSG:4326', 'EPSG:3857');

var view = new ol.View({
    center: center,
    zoom: 10
});

var map = new ol.Map({
    target: 'map',
    layers: [osmLayer],
    interactions: [interaction],
    view: view
});

But when I use this template the interaction is just fine. what's wrong with first approach?

var map = new ol.Map({
    target: 'map',
    layers: [osmLayer],
    view: view
});
map.addInteraction(interaction);


texts: centos 6.7 ArcGIS for Server 10.3.1: publishingtools will not start?

texts:
We're experimenting with ArcGIS for Server (10.3.1) on CentOS (6.7) and I've run into a wall; specifically this:

  Service failed to start
  System.PublishingTools.GPServer:
  LOCALHOST
  Could not create object 'esriGeoprocessing.GPServer'.

My analyst claims the sample service, "SampleWorldCities", worked previously, that he could connect to it AND preview the Geography in ArcCatalog. Unfortunately, I'm not entirely convinced that is true ..so consider that a possible red herring.

However the sample server HAS always been visible, of that I am certain..

Meanwhile, I've been going through the motions of adding PostGRESql, and today my analyst says the sample service is unavailable to him.

I suggested we try publishing a new service with a shapefile as an experiment, but ArcServer fails immediately. Anyway it seems System > PublishingTools is stopped and will not start.

Things I've tried:

stopping/restarting ArcGIS Server, e.g. service arcgisserver [stop|start]
rebooting the system entirely, e.g. reboot
logs have been vague and no help ..particularly when googling errors with "CentOS" among my search terms
stopping iptables entirely, just in case the documentation misled me and I needed to open additional ports, e.g.
service iptables stop
service ip6tables stop
double-checking the ArcServer Manager > System > PublishingTools > Edit Service > Pooling:..curiously, maximum instances shows 0 in the dashboard, even though minimum instances is set to 1 and maximum instances is set to 2. I tried changing the max to 3, saving, then restarting it, and it remained at 0 in the dashboard. I thought this was odd ..but perhaps residual of it failing to start in the first place.
double-checking ArcServer Manager > System > PublishingTools > Edit Service Security:Like iptables, I tried to side-step issues here by setting this to Public, available to everyone, then restarting it. Same.

Have I overlooked an obvious requirement necessary to enabling the PublishingTools.GPServer service? Alternatively, is there another troubleshooting step I might try to dig a little deeper?


texts:
I need to perform a postprocessing operation on an image produced by MapServer before it is handed off to the final user of that image.  I am looking for a way to hand a resulting image from MapServer directly to a postprocessing program, perhaps before it is written to disk.

I'm guessing something could be performed with MapScripts somehow with PHP instructing a different program to perform the postprocessing, but with the number of requests this will need to process, I'm hoping to keep the number of processes to a minumum.

Has anyone done something like this, or knows of a way to do something similar?


texts: Can I use leaflet to access a WMS service published as an OpenLayers javascript API?

texts:
I would like to add this layer to a site I'm working on that uses leaflet:

http://116.240.195.134/vicmapapi/UserGuide.jsp

The API itself is a javascript file based on OpenLayers: http://api.maps.vic.gov.au/vicmapapi/js/vicmapapi.js (I don't know OpenLayers at all, so I'm not really sure what I'm looking at).

Am I able to access the tiles using leaflet? This didn't work:

var map = L.map('map').setView([-37.3, 145.3], 5);

        var wms = L.tileLayer.wms("http://api.maps.vic.gov.au/vicmapapi/map/wms", {
        LAYERS: 'CARTOGRAPHICAL',
        FORMAT: 'image/png',
        FIRSTTITLE: false,
        transparent: true,
        WIDTH:512,
        HEIGHT:512,
        EXCEPTIONS:"application/vnd.ogc.se_inimage",
        SRS:"EPSG:3111"
        });
        wms.addTo(map);

In particular, it didn't work because the tiles it tried to fetch had URLs like this:

http://api.maps.vic.gov.au/vicmapapi/map/wms?service=WMS&request=GetMap&version=1.1.1&layers=&styles=&format=image/jpeg&transparent=true&height=256&width=256&LAYERS=CARTOGRAPHICAL&FORMAT=image/png&FIRSTTITLE=false&WIDTH=512&HEIGHT=512&EXCEPTIONS=application/vnd.ogc.se_inimage&SRS=EPSG:3111&srs=EPSG:3857&bbox=18785164.071364917,-3757032.814272983,20037508.342789244,-2504688.542848655

instead of:

http://api.maps.vic.gov.au/vicmapapi/map/wms?LAYERS=CARTOGRAPHICAL&amp;FORMAT=image%2Fpng&amp;FIRSTTILE=false&amp;SERVICE=WMS&amp;VERSION=1.1.1&amp;REQUEST=GetMap&amp;STYLES=&amp;EXCEPTIONS=application%2Fvnd.ogc.se_inimage&amp;SRS=EPSG%3A3111&amp;BBOX=2327867.7504022,2539132.2495979,2869735.5008043,3081000&amp;WIDTH=512&amp;HEIGHT=512&amp;cache_bbox=2296369.732799%2C2279293.001349%2C2931371.002802%2C2702627.181351

The server then returns error messages like "400: (x + width) is outside of Raster". If I increase WIDTH to 10000, then Apache returns 500 Internal Server error with a complaint about Java heap space. Oops.

So - am I barking up the wrong tree here, or is there a way to do this?


texts: Do zLAS files have to be unzipped to work in geoprocessing tools of ArcGIS Desktop?

texts:
I'm using ArcGIS Pro 2.4.3 and ArcMap 10.7.1 and I thought these versions can read zLAS files without having to unzip them.
When I put a zLAS file into a LAS Dataset (LASD) everything looks normal (I open the properties and see the classification codes as well as the number of points in each code. I can bring the LASD into a map and zoom to the level that the points become visible. I can filter so I'm just seeing ground points then change the symbology from elevation to classification and and watch the points change from blues and yellows to brown for ground).
However, when I use the Las Dataset to Raster tool and choose the elevation field to create a Digital Terrain Model (DTM)/bare-earth raster - the output is no data:( I know this because I use the identify tool and click on a pixel and the value is No Data. To further demonstrate, I choose the symbology for the DTM and choose to display no data cells as red and the entire raster turns red.
Do we have to use the Convert LAS geoprocessing tool and choose as compression type "None" to "unzip" the file and save it as a LAS file?
Because when I use that tool, then use the LAS file it creates, input that LAS file into a LAS Dataset, I can use the LAS Dataset as Raster tool and it creates the DTM raster.


texts: Does Cartodb have a function to change center location based on user location?

texts: I would like to have the map open up a the users current location. I only see documentation for setting the center on init, I would like it to be dynamic..


texts: Removing overlapping polygon parts in sequential order using ArcMap

texts:
I have 500 polygons in a feature class. There are a lot of overlapping polygons and I need to remove these overlapping areas according to three criteria.

What we do now: we have a column with the order in which they should be clipped (1 to 500), and we select each polygon in that order and clip it (editor clip -discard the area that intersects).

Is there a way to automate this?


texts:
I have this code for showing markers on map:

function addUserLayerCSV(){

layer1 = new OpenLayers.Layer.Vector("CSV", {
            projection: epsg6316,
            strategies: [new OpenLayers.Strategy.Fixed()],
            protocol: new OpenLayers.Protocol.HTTP({
                   url: "./../clientUploads/"+ FileURL,
                format: new OpenLayers.Format.Text()

            })
        }

and i have another js file where i convert commas in csv file into tabs, so Openlayers.Format.Text could do the job. Question is, how can i include that js file into this code?


texts: How to iterate through layers of an MXD?

texts:
I know model builder will not iterate through the layers in an MXD, which is kind of fundamentally stupid if you ask me, but I still need to do exactly that.  I am not a strong coder but i understand Python is how to handle this.

I have seen other questions addressing this, but they seem over-my-head in complexity (i have taken an 'intro to python' university course and did well in it, plus did all the ESRI tutorials, still doesn't help understand how this works).  I am hoping someone here can post some simple cut-n-paste python I can use to accomplish this.

I need to iterate through the ToC and export each layer as a shapefile to a specific folder.  The folder is always the same, but the ToC will change depending on the application, so I have different MXD's for each subset of layers that need to be exported, which is why i want to script/model this as then i can open any MXD and run a tool to get my outputs.


texts: Changing QGIS configs and plugins to another partition

texts:
I need to move all my QGIS configs and plugins to D: instead of C:
I installed QGIS in D: but QGIS is still using some folders in C:

I found this answer in StackExchange but I am unable to reproduce it since I don't understand how he managed to do that.

Could someone explain to me how to?

QGIS version: 2.18.18
Windows 10


texts: Cannot create proper URL for Geocoding API

texts:
Using Pentaho to send addresses to Mapquest for long/lat info.

This is the URL I am sending:

www.mapquestapi.com/geocoding/v1/address?key=MYKEY&inFormat=json&json={location: lookup_address}, options:{thumbMaps:false,maxResults:1}

The "lookup_address" field is the concatenated address for the API.

It keeps returning Invalid query.
Unable to get result from specified URL :

http://www.mapquestapi.com/geocoding/v1/address?key=MYKEY&json={location: {lookup_address}, options:{thumbMaps:false,maxResults:1}}&address=%70%20SouthBay%20Ave.%22%22Boston%22%22MA%22%2202118%22

If I put this into the URL it is fine.

www.mapquestapi.com/geocoding/v1/address?key=APPKEY&inFormat=json&json={location:%2270%20South%20Bay%20Ave.%20Boston,%20MA%2002118%22}

What is the proper syntax to send the string properly?


texts: Performance in calculating raster statistics in PostGIS

texts:
I'm trying to calculate raster statistics (min, max, mean) for each polygon in a vector layer using PostgreSQL/PostGIS.

This GIS.SE answer describes how to do this, by calculating the intersection between the polygon and the raster and then calculating a weighted average: https://gis.stackexchange.com/a/19858/12420

I'm using the following query (where dem is my raster, topo_area_su_region is my vector, and toid is a unique ID:

SELECT toid, Min((gv).val) As MinElevation, Max((gv).val) As MaxElevation, Sum(ST_Area((gv).geom) * (gv).val) / Sum(ST_Area((gv).geom)) as MeanElevation FROM (SELECT toid, ST_Intersection(rast, geom) AS gv FROM topo_area_su_region,dem WHERE ST_Intersects(rast, geom)) foo GROUP BY toid ORDER BY toid;

This works, but it's too slow. My vector layer has 2489k features, with each one taking around 90ms to process - it would take days to process the entire layer. The speed of the calculation doesn't seem to be significantly improved if I only calculate the min and max (which avoids the calls to ST_Area).

If I do a similar calculation using Python (GDAL, NumPy and PIL) I can significantly reduce the amount of time it takes to process the data, if instead of vectorizing the raster (using ST_Intersection) I rasterize the vector. See code here: https://gist.github.com/snorfalorpagus/7320167

I don't really need a weighted average - a "if it touches, it's in" approach is good enough - and I'm reasonably sure this is what is slowing things down.

Question: Is there any way to get PostGIS to behave like this? i.e. to return the values of all the cells from the raster that a polygon touches, rather than the exact intersection.

I'm very new to PostgreSQL/PostGIS, so maybe there is something else I'm not doing right. I'm running PostgreSQL 9.3.1 and PostGIS 2.1 on Windows 7 (2.9GHz i7, 8GB RAM) and have tweaked the database config as suggested here: http://postgis.net/workshops/postgis-intro/tuning.html


texts: More information from Summary statistics tool in ArcGIS

texts: I try to find the MAX value for every distinct element in a dataset (I done this using MAX in summary statistics) but I also want to have the id of the element corresponding of that MAX in my output table. Do you have any idea?


texts: Optional variable (with default value) in custom function

texts:
Is it possible to use an optional variable (value2) like

...
@qgsfunction(args='auto', group='Custom')
def my_test(value1, value2=True, feature, parent)
....

in a custom function? Using QGIS 3.8

...like


texts: how to overlay a polygon over SpatialPointsDataFrame and preserving the SPDF data?

texts:
I have a SpatialPointsDataFrame with some additional data. I would like to extract those points inside a polygon and at the same time, preserve SPDF object and its corresponding data.

So far I've had little luck and resorted to matching and merging through a common ID, but this works only becaues I have gridded data with individual IDS.

Here's a quick example, I'm looking for points inside the red square.

library(sp)
set.seed(357)
pts <- data.frame(x = rnorm(100), y = rnorm(100), var1 = runif(100), var2 = sample(letters, 100, replace = TRUE))
coordinates(pts) <- ~ x + y
class(pts)
plot(pts)
axis(1); axis(2)

ply <- matrix(c(-1,-1, 1,-1, 1,1, -1,1, -1,-1), ncol = 2, byrow = TRUE)
ply <- SpatialPolygons(list(Polygons(list(Polygon(ply)), ID = 1)))
ply <- SpatialPolygonsDataFrame(Sr = ply, data = data.frame(polyvar = 357))
plot(ply, add = TRUE, border = "red")

The most obvious approach would be to use over, but this returns the data from the polygon.

> over(pts, ply)
    polyvar
1        NA
2       357
3       357
4        NA
5       357
6       357


texts: Where can I find topology-related algorithms used in QGIS?

texts:
I have shapefile with polygons and I test topology - overlaps, duplicates and gaps. I tested it in QGIS and now I have to describe algorithms, which uses QGIS for this rules. I must write research algorithms.
I finded just general description:
http://docs.qgis.org/2.0/en/docs/gentle_gis_introduction/topology.html
and
http://qgis.org/en/docs/user_manual/plugins/plugins_topology_checker.html
but I can´t find algorithms.


texts: Slope between max and min values of a defined area within a DEM

texts:
Say I have a DEM raster.  I then also I have a set of polygons in which I would like to know the slope between the max and min DEM values of any given polygon.  If my thinking is correct, I could calculate the max slope as follows...

from within any given polygon:

(maxElev - minElev) / distance(maxElevPoint, minElevPoint)

I think the math is correct, but what tools within 3D Analyst could help me to automate this.  I can Python and know ArcPy, but I'm not sure where to start.  I found the Slope tool that deals with rasters within the 3D Analyst extension,  I could convert to raster but from what I can tell this gives the slope back for each cell compared to it's adjacent cells.

what is the best workflow for doing this type of analysis?


texts: Distance Between Two Points Same Value Different Feature Classes

texts: I have two point feature classes - one is the original data and one is the data after editing. I want to know how far each of the points moved in the editing process. How do I compare all of the points in the new data to the points in the old data that have the same ID# and display the distance between them? Also I have the standard license, so I don't have access to the Near or Generate Near Table tools.


texts: Selecting by attributes which have alphabet data using ArcGIS Field Calculator?

texts:
I only want to select values that are entirely alphabetic characters.

Select by attributes which I have alphabet data.

SELECT ONLY RED INDICATED VALUES


texts: Seeking tutorials online explaining how to draw squares (having specific length) in QGIS?

texts:
I need to produce in QGIS a square having a specific dimension (1km x 1km). First I created a new polygon layer, but when I start editing through the Advanced editing tool, the program crashes. I guess that there is something wrong in my procedure.

Beside that tool, does the 3.x version (I'm using Madeira) allow me to create such squares using another method?


texts: Is there a simple way to clip multiple layers and orthophotos?

texts:
(*disclaimer, I don't know any Python/scripting)  Ok, so I have clipped several layers in the past using Arcmap 10.+, but I am having a hard time wrapping my head around how hard it can be at times.

I used to own a photography/design business, so I'm used to manipulating crops via Photoshop. I think because of this, I have issues with simple clipping when composing a map.

For this current clip/crop, I would like to freehand a polygon around several connected streams and "crop" out everything around that polygon. Within the polygon, there are several layers having many polylines and points, and a Ortho_tile_tiff catalog layer aerial view. I believe it's this Ortho_tile geodatabase that's throwing me here. If someone could break it down for a semi_green GISer I would be much appreciative. Thanks!


texts: Simple data vs non simple data in ArcGIS

texts:
I want to enter some data in a feature class with ArcObjects. In my code I use EditSession for that, but in storing data an error returned that says I must use editsession. By researching I found that the feature class is non-simple data, and yesterday I find that this non-simple data option is because of creating a M-M relationship. At the end my goal is to enter data into this non-simple feature class with this code.

public string insertdata(string tablename, List<tableRow> fields, string primerykey)
    {
        IEnumDataset datasets = ws.get_Datasets(esriDatasetType.esriDTAny);
        IDataset dataset = null;
        ITable table = null;

            IFeatureWorkspace fws = ws as IFeatureWorkspace;
            IFeatureDataset fdataset = fws.OpenFeatureDataset("Electrical");
            datasets = fdataset.Subsets;
            do
            {
                dataset = datasets.Next();
                table = dataset as ITable;

                if (dataset == null) return "Table not exist";

            } while (dataset.Name != tablename);

        IWorkspaceEdit ews = ws as IWorkspaceEdit;
        ews.StartEditing(false);
        ews.StartEditOperation();

        IRow row = null;
            table.Fields.FindField(primerykey);
            string value = "123";  }
            }
            ICursor cursor = null;
            IQueryFilter2 queryFilter = new QueryFilterClass();
            queryFilter.WhereClause = primerykey + "='" + value + "'";
            cursor = table.Search(queryFilter, false);
            row = cursor.NextRow();

        foreach (tableRow field in fields)
        {
            Type type = field.type;
            int ind = row.Fields.FindField(field.name);
            if (ind == -1) return "Field not exist";
            row.Value[ind] = Convert.ChangeType(field.value, type);

        }

        row.Store();

        ews.StopEditOperation();
        ews.StopEditing(true);

        return "Row inserted succesfull";
    }

This method worked on all datasets but in one of them an error returned that says: "Objects in this class cannot be updated outside an edit session."

I want to know what the difference between "simple data" and "non-simple data" is in an ArcGIS geodatabase? And how can I convert non-simple to simple data?


texts: Understanding the unit of area for standard deviational ellipse

texts:
I am using the aspace package to calculate the standard deviational ellipse of a set of points. Test shapefile with 5 points zipped here.

library(maptools)
library(sp)
library(aspace)
library(rgdal)
shpfile <- readShapePoints("testSHP.shp")
coords <- shpfile@coords
calc_sde(points=coords)
plot_sde(plotnew = TRUE)
r.SDE$Area.sde

The area is 2824712729.

I think I am getting tripped up by the coordinate system:

coords
#  coords.x1 coords.x2
#0   3667152   8861968
#1   3713150   8878606
#2   3714374   8843128
#3   3683545   8834809
#4   3690396   8854383

The shapefile is set to +proj=longlat +datum=WGS84 +no_defs.

What are the units for this area?
How can I project the coordinates in decimal degrees?


texts: MMQGIS geocoding location error in QGIS 2.18

texts:
I´m getting a locational error during the process of geocoding a csv file (with 3 records) with MMQGIS in QGIS 2.18.

The error is that the geocoding result are in the same locational place.

I already made the geocoding with Google and OpenStreeMap servers, having the same error. Any ideas...?

I´m using a single csv file (UTF-8 encoded) with 5 records (3 in Comodoro Rivadavia city, 1 in Puerto Madryn city and 1 in Trelew city). Sorry but I don´t know how to add a screen shot here. It´s seems that the geocoding is well done in a City scale, but not so good with the addresses detail (in the city where we have 3 addresses, the 3 of them are in the same spot.

CSV detail: Direccion;Ciudad;Provincia;Pais Belgrano 574;Comodoro Rivadavia;Chubut;Argentina Sarmiento 230;Comodoro Rivadavia;Chubut;Argentina Belgrano 874;Comodoro Rivadavia;Chubut;Argentina Libertad 147;Puerto Madryn;Chubut;Argentina Laprida 147;Trelew;Chubut;Argentina


texts: Creating boundaries for Land-use classification

texts:
I have done the land use classification for sentinel imagery. Now I want to create boundaries, within that classified image, enclosing different classified section.
For example boundaries for water, forest etc.
If a category has more than one portion inside the image, then it can have more than one enclosing boundaries. Like, the forest can have more than one enclosing boundaries depending on where foreset has been detected inside the image.

How can I perform this operation? (in QGIS)

For example - I have shapefile of the complete image given below. I want to create shapefiles for each unique colour in the image.


texts: looking for datum to georeference maps of 19th century Maine

texts: I am trying to georeference three (1840, 1844, and 1877) northern Maine maps. I really am not able to find which datum to use. Can anyone point me to a list of historic and regional datums?


texts: Batch map creation, iterated over attributes of a layer, in QGIS

texts:
I am using QGIS 2.0.1 under Mac OS-X 10.8.2. I have a layer of counties in New Mexico to which I have joined a csv file that includes presence / absence of plant species by county. I've set up the layer style to show presence (coded "1") as grey and absence (coded "0") as white. So I can go into "Layer Properties", go to "Style", select the appropriate plant from the "Column" pull-down list, and I have a map that is nicely filled to show counties in which the plant is present.

The question is--how do I create an output image for each species showing its distribution?

At present, what I would have to do is the following:
1) Select a species from the pull-down list in "Layer Properties":

yielding a map that looks like this:

2) Enter "Print Composer", get all that set up, export as an image:

3) Switch back to the main QGIS window, open up "Layer Properties", and select the next species:

4) Switch back to the Print Composer, and export the next image:

Now, that all works, but given that I have several thousand species to deal with, this is not a remotely efficient way of approaching the issue. So I am trying to figure out how to automate map creation so that I can create a set of images incremented over all the attributes in this layer. So far, I've come up empty. The "Atlas Generation" plugin that's been integrated into the Print Composer seems to solve essentially the same issue, but so far as I can tell it can only increment across rows (features) and has no ability to increment across columns (attributes). Any suggestions?


texts:
I am trying to import CSV files into PostGIS. Following this post, I have created tables before. I found other suggestions saying that I can run the copy command.

If I run this command:

COPY table FROM '/Users/macbook/file.csv' DELIMITERS ',' CSV HEADER;

it didn't copy the table at all. It says that "table" is not recognized.

I tried this:

COPY moulding
(Borough,Block,Lot,CD,CT2010,CB2010,SchoolDist,Council,ZipCode,FireComp,PolicePrct,Address,ZoneDist1,ZoneDist2,ZoneDist3,ZoneDist4,Overlay1,Overlay2,SPDist1,SPDist2,LtdHeight,AllZoning1,AllZoning2,SplitZone,BldgClass,LandUse,Easements,OwnerType,OwnerName,LotArea,BldgArea,ComArea,ResArea,OfficeArea,RetailArea,GarageArea,StrgeArea,FactryArea,OtherArea,AreaSource,NumBldgs,NumFloors,UnitsRes,UnitsTotal,LotFront,LotDepth,BldgFront,BldgDepth,Ext,ProxCode,IrrLotCode,LotType,BsmtCode,AssessLand,AssessTot,ExemptLand,ExemptTot,YearBuilt,BuiltCode,YearAlter1,YearAlter2,HistDist,Landmark,BuiltFAR,ResidFAR,CommFAR,FacilFAR,BoroCode,BBL,CondoNo,Tract2010,XCoord,YCoord,ZoneMap,ZMCode,Sanborn,TaxMap,EDesigNum,APPBBL,APPDate,PLUTOMapID,Version)
FROM
'/Users/macbook/file.csv'
DELIMITERS
','
CSV HEADER;

but didn't work either.

An example of such data set can be downloaded from this link:

Should I create a model and then execute it?


texts: Adding a new and removing an old marker every time the user click on the map

texts:
On Leaflet, I am drawing a marker everytime the user click on the map like this:

function getClickedLanLon(e) {

    var lat,
        lon,
        zoom;

    lat = e.latlng.lat;
    lon = e.latlng.lng;
    zoom = map.getZoom();

    marker2 = new L.Marker(new L.LatLng(lat, lon));
    map.addLayer(marker2);
}

map.on('click', getClickedLanLon);

but the map get populated with so many markers. Instead, I want to delete the old marker before adding the new one. So, I add the following if condition in the beginning of the function:

if (!map.hasLayer(marker2)) {
   map.removeLayer(marker2);
}

But I get the following error:

  Uncaught TypeError: Cannot read property '_leaflet_id' of undefined

My question is how to add a new and remove an old marker every time the user click on the map?


texts:
How do I resize & replace a layer that is displaced & too big?

The purple background above is actually a layer representing a state. I have tried using "Define Projection" and "Project" tool, but I wasn't successful with it.

I am really new to this software.


texts: How to traverse all Landsat 8 pixels in each MODIS pixel space on GEE

texts: I have two scenes of Landsat 8 and MODIS images in a designated area. I want to traverse all Landsat 8 pixels within the space of each MODIS pixel. How should I deal with this?


texts: Simple way to modify GeoWebCache Capabilities?

texts:
I have been working with WMTS OGC Standard in Geoserver and now i'm worried about the metadata of my service. The capabilities document is managed by GeoWebCache :

http://localhost:8081/geoserver/gwc/service/wmts?request=GetCapabilities&service=WMTS&version=1.3.0

and there's no visible way to modify some elements such as title, abstract, keywords,... in Geoserver.

I try to find some info by looking for through the internet but no luck. Is there a way to modify the capabilities generated by gwc?


texts:
I am trying to execute in mapcalc a simple band ratio exercise.
I have set new map as "Band 12"
and entered this equation into mapcalc:

Band12== Band1@MApset/Band2@Mapset

I hit "Run" but I keep getting this answer: Invalid Map .

What am I doing wrong?


texts: How do I save a polygon once I have created it in editor?

texts: I am attempting to create a new polygon. I go to ArcCatolog>New>Shapefile and make sure Feature type is set to polygon. I create the polygon and if I click save edits or stop editing, it disappears. What am I doing wrong?


texts:
I am totally new to either PostGIS or PostgreSQL. I am trying to install PostGIS on my computer, which is running with MacOS Big Sur version 11.2.3
I have installed PostgreSQL 13 with the EDB.
Then I tried to install PostGIS 3.0 with the application Stack Builder, but I get an error message. I understand that it is a permission problem that often happens with Big Sur. I click on cancel the installation, and then I get this message :

"Please excute sudo /tmp/rZxm7Y/postgis-pg13-3.0.2-1-osx.app/Contents/MacOS/osx-x86-6"

I open the Terminal and run it but again I get an error message

"Command not found"

I have deleted "postgis-3.1.1.tar.gz" after a first attempt to clear everything, but now I don't know where it should be (Applications ? Library ? User ?).
I have also installed PostgreSQL and PostGIS with homebrew at some point, and it worked, but I've uninstalled it since, I don't feel confident enough to work without a GUI.
I have installed Postgres with success, but I am afraid it will conflict with my previously installed version of PostgreSQL. I've tried to uninstall it, but I can't find the uninstaller, even when I run find . -name 'Postgre*' -type d in the Terminal.


texts: Python Band.SetNoDataValue(), sets all bands at once

texts:
Can someone tell me what I am doing wrong?

I want to set bands 1 and three to 255, and band 2 to 0. The code below only should only set band 2, but gdalinfo shows that all bands get set to 0.  Ubuntu 14.04, gdal 1.11.0.

--------------------------------------
import sys
import os.path
import gdal

input = sys.argv[1]

dataset = gdal.Open( input, gdal.GA_Update )
if dataset is None:
    print 'Unable to open', input, 'for writing'
    sys.exit(1)
band = dataset.GetRasterBand(2)
band.SetNoDataValue(0.0)
band.FlushCache()
band = None
dataset = None

----------------------------------------------

$ gdalinfo test.tif
...
Upper Right (-11546254.150, 4881966.253)
Lower Right (-11546254.150, 4244288.152)
Center      (-11950987.074, 4563127.202)
Band 1 Block=15107x1 Type=Byte, ColorInterp=Red
  NoData Value=0
Band 2 Block=15107x1 Type=Byte, ColorInterp=Green
  NoData Value=0
Band 3 Block=15107x1 Type=Byte, ColorInterp=Blue
  NoData Value=0


texts:
I am looking for a method for convert orthophotomap (.adf files) to ERDAS .ecw file.

Have you got any ideas ?

Or maybe you have got any other idea how i can reduce a size of orthophotomap .adf files?
Now i have got 16 GB file I want to compress it to something like 3-4GB- is it even possible?

I read that .ecw is highly compressed for raster formats...


texts: Failed to load ArcGIS JS API polar satellite map

texts:
I am trying to draw satellite map with polar (Arctic) projection using ArcGIS JS API. Firstly I load default (Mercator) satellite map and it works fine. Then I change map spatialReference to polar the following way:

require([
"esri/Map",
"esri/views/MapView",
"esri/geometry/Polygon",
"esri/geometry/SpatialReference",
"esri/layers/FeatureLayer"
],
function (Map, MapView) {
    var map = new Map({
        basemap: "satellite"
    });

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

    var view = new MapView({
        container: "map",
        map: map,
        spatialReference: viewSpatialReference,
        center: centerPoint,
        scale: 2000
    });
});

When I try to load this map in browser, I get blank map and two errors in console

  [esri.views.2d.layers.TileLayerView2D] #resolve() Failed to resolve layer view (layer title: 'World Imagery', id: 'satellite-base-layer') message: "The spatial reference of this layer does not meet the requirements of the view"

  dojo.js:310 [esri.views.LayerViewFactory] Failed to create view for layer 'World Imagery, id:satellite-base-layer' of type 'tile'. message: "The spatial reference of this layer does not meet the requirements of the view"

Does anyone have experience working with ArcGIS JS for Arctic region? My original idea was to redraw basic satellite map using ArcGIS; I hope it is possible.


texts: Creating grid in QGIS 'anchored' to specific point

texts:
I need to create grids using QGIS with specific intervals 'anchored' to, or beginning at, a specific UTM coordinate. This is for a survey sampling procedure. For example: a point was determined - UTM 0561111, 6111111, as the starting point and I must create a grid with 92m spacing extending north and east for the extent of the polygon. The intersection of those gridlines will be unbiased sample points.
I can't seem to 'anchor' the grid to the starting point using "Create Grid" or MMQGIS "Create Grid Layer". I'd like to be able to do this with QGIS 3.4. Here's a sketch of what is needed.


texts: Displaying images in QGIS - Actions - unable to run command

texts:
I was following the video below:
https://www.youtube.com/watch?v=aW7CI7SiqQ4&ab_channel=WanjohiKibui
In order to prepare the action in my QGIS project. I want to have my images popped-up instantly.
Unfortunately, I am getting an error like this:

My action properties look as follows:

The image directory is related to the project path.
I also tried the command:
   eog [%Image%]

and
   display [%Image%]

as per the thread here:
Using Hyperlink URL in attribute table of QGIS?
the result the same everywhere
What causes the problem? Is there something with the hotlink?
Issues with hotlinking in QGIS


texts:
I have problem with ran tool from ArcToolbox. When I run any tool, I get error:

  Content within this application coming from the website listed is
  being blocked by Internet Explorer Enhanced Security Configuration.

I work on Windows Serwer 2012 R2 and use ArcGis for Desktop 10.3.1

I try:

reinstall ArcGis for Desktop 10.3.1
to activate ActiveX (http://support.esri.com/technical-article/000011895)
add http://localhost to trusted sites (http://www.makeuseof.com/tag/configure-trusted-sites-internet-explorer-group-policy/)

But problem still exists.
How to fix this?


texts: Finding row that has shapefile containing lat/long point in PostGIS geometry column

texts:
I am fairly new to PostGIS and I'm trying to return a row in my table that has the shapefile containing a lat/long point I pass in via st_contains. Below is a query I wrote myself that brings back 0 rows with no errors as seen below.

SELECT *
FROM es_zones
WHERE ST_Contains(geom, ST_GeomFromText('POINT(-73.952545 40.774576)'))=true;

Having found a query similar of that to which I am trying to accomplish here, below is my last query ran.

SELECT *
FROM es_zones
WHERE ST_Contains(es_zones.geom,ST_transform(ST_GeometryFromText('POINT(-73.952545 40.774576)',26918),26918))=true
AND ST_Distance(es_zones.geom,ST_transform(ST_GeometryFromText('POINT(-73.952545 40.774576)',21918),21918))=false;

which also brings back 0 rows as seen below

----+-------------+-------+--------+--------------+-------------+--------------+-------------+---------------+---------------+-----------+-----+----------+-------+--------

gid | schooldist | boro | initials | creat_date | edit_date | zoned_dist | boro_num | shape_leng | shape_area | remarks | dbn | esid_no | label | geom

(0 rows)

Really not sure where I am going wrong.

I was able to get QGIS set up and connected to my PostGIS data. I don't understand what the 997487, 221261 coordinates are.

How do they correspond to the lat/lon points I got back from my GeoCoder Gem?

I passed them in as my point but got the following error;

  ERROR:  Operation on mixed SRID geometries

Verifying the SRID shows

SELECT ST_SRID(geom) FROM es_zones LIMIT 1;
 st_srid
---------
  926918
 (1 row)

So it sounds like something may have happened to my spacial_ref_sys table, or perhaps I didn't get it configured correctly to begin with. A Google search for my last error ERROR:  GetProj4StringSPI: Cannot find SRID (4326) in spatial_ref_sys returned this StackOverflow question. Which points to running spatial_ref_sys.sql to regenerate the rows in that table.

I found some elaboration on how to go about executing spatial_ref_sys.sql in the PostGIS documentation here, but I'm confused why it says it can't find SRID 4326 when looking in PGAdmin I can see it as shown below!


texts:
I am searching for a tool to generate vector-tiles that I can then display using the Mapbox Map-View.

We currently have a postgres database with all of our data. What is the best tool to generate vector mbtile tiles?
What would be awesome it a tool where we can (ideally with a GUI) define what kind of data is embedded into the tiles at which zoom level.

If we could gather the information from our postgres database would also be cool, so we can pack all the information we need into the vector tiles.


texts: Copy feature class to feature class in loop and a list of fields in modelbuilder

texts:
I have a several file geodatabase with seven features in it (polygon).
I want to copy them to a new feature in a other empty file geodatabase (instead of deleting the fields...had same problems with that tool)

What I want:
I just need 20 out of 60 fields and in the new feature in a new order.
I tried to set up this in Modelbuilder, but there a several problems which I don't get.
I will use the iterator for features thats clear. No problem.
Then feature to feature tool. Up to this point is clear to me.
I create a variable with my necessary field names in the order I want them. clear to me.
But what is not clear:

How to say take this list every time for each feature for the field map?
How to reorder them? like variable-list has it?

Is it possible with model builder or do I need a script for that?
I have around three million polygons in each feature...

I see my english expression was not the best. I am german ;-)
I will simple it a bit easier.

One Filegeodatabase contents: 7 features(polygons) and around 3 million attributes and 60 fields.
And I have 4 of this filegeodatabases...
Each field in each feature is same by name.
It is a result of many intersection steps.
But for the end product, I just need 20 of this fields in my desired order.
Yes, I agree to you, I can use the tool, but so I have to click 40 times of each to have my new feature but i cannot move the order in the "fieldmap"-window.

The other option what I could do is, to copy the feature first into the new filegeodatabase and delete the unwanted fields, but at the end it is same clicking time which i not have, so it must be working automated or?
In the delete field tool, it brings check boxes, I can select all or unselect all and search through the fields where are my wanted and so i could get it.
I can do it for one. And i can not recall the tool from the resultslist, because then he reads the list of the feature again, and i can start choosing again...
is this cleared explained?
On searching the net, i didnt find any useful solution. The example-code of delete-tool in the arcgis-help resource i could not get to run...


texts: Inserting Timeline (similar to Legend) on map using ArcGIS for Desktop?

texts:
is there any way to insert a timeline in ArcGIS 10.1? I mean timeline more as a kind of legend showing at which time (or after how many minutes) point features occured (not animated or dynamic!).

I have plotted a number of whale blows in different colours, depending on the time they were detected. However, it doesn't show how much time there was between blows and how long the entire sighting was.

Any ideas?


texts: Automating ArcMap topology summary report

texts: I am trying to automate a topology check process in ArcMap. As the last step, I hope to generate the error report automatically rather than clicking the 'generate summary' and 'Export to file' buttons manually. Is it feasible to automate with arcpy...?


texts: Calculating area under each polygon in QGIS

texts: I would like to know how to calculate the area of the vector layer (red in the photo) under each polygon (white square in the photo). The red layer has the area (km2) in its attributes. But what I want to know is how much (area) of the pink layer is within the square grid. Both my layers are in meters already.


texts:
I am new to scripting.

I have around 4500 small shapefiles to be intersected with one big shapefile. The primary key is to make 4500 small intersected shapefiles.

How to do it correctly by using ArcPy?

I tried using batch processing, but I quit after few tries, as it is lot of work. I wrote some kind of "script".

#system modules
import arcpy, os, datetime
from arcpy import env
from arcpy.sa import *

arcpy.env.workspace = r"D:\test\split\los"
arcpy.env.overwriteOutput = True
print "start"
infeature = ["big.shp", "small1.shp"]
join_attributes = "ALL"
cluster_tolerance = "-1 Unknown"
output_type = "INPUT"

for i in infeature:
    outfeature = r"D:\test\split\los_" + i
    arcpy.Intersect_analysis(infeature, outfeature, join_attributes, cluster_tolerance, output_type)

print "end"


texts: Managing Data Driven Pages with ArcObjects

texts:
I can't find any resources on managing Data Driven Pages with ArcObjects. E.g.: is it possible to use ArcObjects (.NET) for the following tasks:

Check if the current map document has enabled Data Driven Pages?
Loop through all pages and get the name, number, extent and scale of each page?
Enable DPP and set the index layer?
Disable DPP?

Any hints would be appreciated.


texts: Arcpy get point feature at start and endpoint of Line feature

texts:
I am having trouble getting a feature that is at the start and end point of a line. I can get the start and end point of the line but I am not sure how to grab the Point feature that is at the start and end point of a line. I need to get the point feature that is at the start and end point of a line and then grab some data from the point feature and update the Line feature data. How can I get a point feature that is at the start or end point of a Line and then get that feature data?

    with arcpy.da.SearchCursor("GravityMains",["Material","SHAPE@"]) as curs:
        for row in curs:
            if row[1]:
                startptx = row[1].firstPoint.X
                startpty = row[1].firstPoint.Y
                arcpy.SelectLayerByLocation_management("Manhole","INTERSECT",arcpy.PointGeometry(arcpy.Point(startptx,startpty)))

             else:
                 print "NULL"

I am using ArcMap 10.3.1


texts: How to calculate the distance to the nearest water source?

texts: I have a polygon hydro layer of my study area and am trying to calculate the distance to nearest water bodies from my site. Since my site are 2km*2km grids I am not sure how to calculate it? Whether to calculate from centroids of the grids or convert the layer to raster and somehow calculate mean distance of raster cells within the grid. Would be a great help if somebody could provide steps :)


texts: Checking maximum number of features for shapefile to be opened in ArcMap?

texts:
I have a shapefile that contains 10,608,268 polyline features.

What is the maximum limit for records that a shapefile can contains?


texts: Iterating through geodatabase and assigning field value based on spatial location using ArcPy - Part 2

texts:
This currently loops through my dataset and calculates the installation ID field from my dataset with the installation code from my Bragg Area feature that are within the polygon of each row. It works on the first feature in my dataset but on other features, in that dataset, it only calculates the features that fall into the last row of the Bragg feature class. The picture below shows how it looped through the monitor stations twice. It should have also shown the same for a couple of other features.

How do I get it to calculate every feature in my dataset within every row in the Bragg area as it does with the first feature?

import arcpy

# Set overwrite option
arcpy.env.overwriteOutput = True

# Define your workspace. This should be your SDE database with the feature class
arcpy.env.workspace = arcpy.GetParameterAsText(0)
fcZone = arcpy.GetParameterAsText(1)
aFieldID = arcpy.GetParameterAsText(2)
zFieldID = arcpy.GetParameterAsText(3)
flZone = arcpy.GetParameterAsText(4)
flAsset = arcpy.GetParameterAsText(5)
# This can be done once for all iterations
arcpy.MakeFeatureLayer_management(fcZone, flZone)

for fc in arcpy.ListFeatureClasses():
    fcs = fc
    fields = arcpy.ListFields(fcs, aFieldID)
    for field in fields:
        if field.name == aFieldID:
            # KEEPING ORIGINAL LAYERS (NOT CREATING ADDITIONAL "JOINED" LAYER)
            # Create FeatureLayers
            arcpy.MakeFeatureLayer_management(fc, flAsset)
            # Create a search cursor for the Area
            rows = arcpy.SearchCursor(flZone)
            for row in rows:
                # What you'll do is select each Area one at a time, and then select all the Features in that Area and calculate the field
                # NOTE: If you are using not using shapefiles, then you'll have to change the FID in the line below to OBJECTID (or similar)
                arcpy.SelectLayerByAttribute_management(flZone, "NEW_SELECTION", "\"OBJECTID\" = " + str(row.getValue("OBJECTID")))
                arcpy.SelectLayerByLocation_management(flAsset, "WITHIN", flZone, "", "NEW_SELECTION")
                arcpy.CalculateField_management(flAsset, aFieldID, "'{0}'".format(str(row.getValue(zFieldID))), "PYTHON_9.3", "")
                arcpy.AddMessage(fcs)


texts: Intersection of geometries as GeoSPARQL result

texts:
I'm querying two geosparql conform datasets at endpoint http://onto.fel.cvut.cz:7200/sparql trying to get intersectred geometry as result using GeoSPARQL language.

PREFIX rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#>
PREFIX geo: <http://www.opengis.net/ont/geosparql#>
PREFIX geof: <http://www.opengis.net/def/function/geosparql/>
PREFIX databaseTableFVU: <http://onto.fel.cvut.cz/ontologies/town-plan/databaseTableFVU/>
PREFIX databaseTableSSVU: <http://onto.fel.cvut.cz/ontologies/town-plan/databaseTableSSVU/>
PREFIX geosparql: <http://www.opengis.net/ont/geosparql#>
PREFIX common: <http://onto.fel.cvut.cz/ontologies/town-plan/common/>
PREFIX xsd: <http://www.w3.org/2001/XMLSchema#>
SELECT distinct ?var1 ?var2
WHERE  {
       {
        ?var1 databaseTableFVU:kodfp1_a ?val1.
        ?var1 geosparql:hasGeometry ?geometry1.
        ?geometry1 geosparql:asWKT ?WKT1.
        FILTER(?val1 = "OB" || ?val1 = "OV")
       }
       {
        ?var2 databaseTableSSVU:kod ?val2.
        ?var2 geosparql:hasGeometry ?geometry2.
        ?geometry2 geosparql:asWKT ?WKT2.
        FILTER(?val2 != "BD" && ?val2 != "BQ" && ?val2 != "BRV" && ?val2 != "BRR")
      }
FILTER(geof:sfIntersects(?WKT1, ?WKT2)).
}

I know that this query returns only features, that intersects, but I have no idea how to actually get intersection as geometry representation.


texts: Add comma to number labels, multiple records in attribute table

texts:
A QGIS question:

How can I insert commas in 4 digit numbers in an attribute table that has multiple records, with some numbers that are only in the hundreds? It's for an animation so multiple labels are being displayed.

I tried using the format_number option under the string menu in field calculator, but it only works when I enter a specific number in the expression dialogue, so i imagine taht this tool only works for attribute tables with one record.

I don't know much about coding so hopefully there's an option that doesn't need that, or a simple way to describe how to do the coding.

I don;t think this is the same as the topic "Forcing decimal numbers to use comma in ArcGIS Field Calculator? " because i still want to keep the decimal places.


texts: Joining 17 million points to one dissolved polygon layer taking over 24 hours to run!

texts:
I am trying to count how many points, out of 17 million, fall inside a polygon (see image below)

which I created from dissolving various buffers and features. The polygon layer is thus one dissolved feature. Both the points and the polygon are feature classes saved inside a file geodatabase. I am also using ArcGIS 10.3 on a 64 bit desktop.

To count the number of points I first added a column called 'Count'and populated all fields with 1. I then right clicked my polygon layer, clicked 'Join and Relate', and then 'Join'. I then chose the join option: 'Join data from another layer based on spatial location'. The parameters I entered are below:

. As you can see, I am summing all fields. Thus, I will be able to sum the 'Count' column to determine how many point features lie inside the polygon feature as shown above.

I started this process 23 hours ago, and so far only half of the 17 million points are processed. I am also concurrently running two other similar processes with copies of the 17 million points and two separate polygon layers (both point/polygon feature classes are saved in separate geodatabases).

I was able to successfully do a similar join with the 17 million points and a fishnet with a grid size of 1 mile by 1 mile covering the state of Ohio. That process only took me 2.5 hours yesterday.

Why is it taking so long to sum all points inside a dissolved polygon feature class? Is running such a process 3 times at the same time adding to the slow speed? Is there a way to speed up this process in ArcGIS 10.3?

I am the only GIS Analyst at my company and am at a loss in understanding why it is taking so long.


texts: Spectral Profile Chart Using Landsat 5 TOA Earth Engine

texts:
I'm attempting to create a spectral profile chart for 4 features in my feature collection that represent vegetation. I'm not really sure how to move forward, this is one of my first few times coding for a class and am having difficulty creating my code for Landsat 5 C1 T1 TOA. I've been mostly crafting based on code given to me by my professor, but it is designed to be used for Landsat 8. I've looked for answers online, but I'm not particularly skilled in deciphering long code language and most questions are regarding Landsat 8.
https://code.earthengine.google.com/d6fbedd85ee315948e322b1d64e302ca
var landsat = ee.ImageCollection("LANDSAT/LT05/C01/T1_TOA")
var town = ee.Image(landsat.filterDate('1984-07-01', '1984-07-31').filterBounds(town).sort('CLOUD_COVER').first());
//What we did there was declare a variable, called park. We then called Image from our ‘landsat’ import and filtered by date, then cloud cover and took the first image.
var trueColor = {
  bands: ['B3', 'B2', 'B1'],
  min: 0,
  max: 0.45
};
//Add Image
Map.addLayer(town, trueColor, 'True-color image');
var falseColor = {
  bands: ['B4', 'B3', 'B2'],
  min: 0,
  max: 0.45
};
Map.addLayer(town, falseColor, 'False-color image');
//TOA
var toaImage = ee.Image('LANDSAT/LC08/C01/T1_TOA/LC08_220015_20170802').select('B([1-7])');
print('TOA', toaImage);
Map.addLayer(toaImage, {bands: ['B4', 'B3', 'B2'], min: 0, max: 0.45}, 'toa');
var points = ee.FeatureCollection([
 ee.Feature(greenpea, {'label': 'Green Pea'}),
 ee.Feature(potato, {'label': 'Potato'}),
 ee.Feature(sorghum, {'label': 'Sorghum'}),
 ee.Feature(marsh, {'label': 'Marsh'}),
]);
Map.addLayer(points);
var wavelengths = [0.44, 0.48, 0.56, 0.65, 0.86, 1.61, 2.2];
var options = {
 title: 'Landsat 5 TOA spectrum of crops in Hajdůbőszőrmėny, Hungary',
 hAxis: {title: 'Wavelength (micrometers)'},
 vAxis: {title: 'Reflectance'},
 lineWidth: 1,
 pointSize: 4,
 series: {
 0: {color: '00FF00'}, //Marsh
 1: {color: '0000FF'}, //Green Pea
 2: {color: 'FF0000'}, //Potato
 3: {color: 'f0af07'}, //Sorghum
}};
var chart = ui.Chart.image.regions(
 toaImage, points, ee.Reducer.mean(), 30, 'label', wavelengths)
 .setChartType('ScatterChart')
 .setOptions(options);
print(chart);


texts: leaflet-routing-machine add waypoints, routing summary

texts:
I'm playing around with leaflet-routing-machine (https://github.com/perliedman/leaflet-routing-machine) and have 2 questions about it's usage:

1)

I have the following code now

routeControl = L.Routing.control({
     waypoints: [
                 L.latLng(48.133205,11.565225),
                 L.latLng(48.131830,11.588227)
             ],
    useZoomParameter: false,
    autoRoute: false,
}).addTo(map);

which works fine, but what when i want to start without any waypoints? When i leave out the waypoints from the constructor and try to add it then with

routeControl.spliceWaypoints(0, 0, e.latlng);

routing is not working anylonger (i know, i need to start routing manually, because I set autoRoute to false). But if I start with those initial points and add further waypoints, routing is working.

2) I'm trying to get totalDistance, which is part of the Interface IRouteSummary.

According to the API, I need to access the IRoute Interface, but i didn't find out how to get IRoute from L.Routing.Control

I found this post (https://stackoverflow.com/questions/26570980/getting-distance-between-2-markers) which uses an array of routes[], but i also don't know how to get this array?! It also refers to this Tutorial http://www.liedman.net/leaflet-routing-machine/tutorials/interaction/ but here i don't find any information about how to access this array as well.


texts: Where are TIGER's "Consolidated Cities"?

texts:
The TIGER reference docs make numerous references to Consolidated Cities. It's also one of the checkboxes on TIGERWeb.

But there are only 5 files in the Consolidate cities directory and available for download (http://www2.census.gov/geo/tiger/TIGER2011/CONCITY/). And I can't find any in the Places file (i.e. no lines have MTFCC='G4120').

So can someone please explain how to find which "Consolidated City" (if any) any given (lat,long) belongs too in the US ?

Thanks


texts: Convert multi KML and KMZ into shapefile using ArcPy?

texts:
I try to convert with arcpy (python 2.7.8) 50 KML and KMZ files, that i downloaded from google earth into shapefiles. All files located in many sub-Folders. All Sub Folders are located in one large directory. When i run this code:

import arcpy,fnmatch,os

rootPath = r"C:\Project\gis\layers"
pattern = '*.km*'

for root, dirs, files in os.walk(rootPath):
    for filename in fnmatch.filter(files, pattern):
        print(os.path.join(root, filename))
        arcpy.KMLToLayer_conversion( r"C:\Project\gis\layers" ,r'C:\Project\gis')

i get an error:

C:\Project\gis\layers\ESD.kml
ExecuteError: Failed to execute. Parameters are not valid.
ERROR 000814: Invalid file type
Failed to execute (KMLToLayer).

I red:

Python Script to convert a mass of KMLs to a single fGDB

http://resources.arcgis.com/en/help/main/10.2/index.html#/KML_To_Layer/00120000004w000000/ but didn't understand my code mistake


texts: MapInfo equivalent of SQL 'AS' statement

texts: I'm trying to execute a SQL query in MapInfo, and want to change some of the column names to something more readable. The SQL command 'AS' doesn't seem to work - is there a MapInfo equivalent?


texts: Specifying geoid in ArcGIS Raster Project tool?

texts:
Background:

ArcGIS Desktop 10.6, Spatial Analyst

I have two bare-earth DEMs, derived from lidar, acquired on different dates.  I need to subtract one DEM from the other as part of a stream erosion study.

The DEMs are in different horizontal coordinate systems, so I will project one raster to match the other prior to any analysis.

Additionally, they have different geoids (03 vs. 12B), although they have the same vertical coordinate system (NAVD88).

According to the first answer to Are LiDAR datasets with same vertical datum, but based on different geoids comparable?, I need to correct for the different geoids or else I will have systematic deviations in the elevation values, adding error to any resulting analysis.

The Problem:

The ArcGIS Project Raster tool does not have the option to specify the output geoid (in my case, 03 is required).  Below is a screenshot of the tools's Z-Coordinate System tab.

How/where is the output geoid specified?


texts:
I have a large xyz file. it takes me long time to load it to Qgis. I wanted to know that how I can crop the file to my area of interest.

I tried to save the file using extent option but still have problem and Qgis crashes during the saving process.


texts: How does one use leaflet layergroups with CartoDB layers?

texts:
I'm trying to use leaflet layergroups with a cartodb layer, ala: http://leafletjs.com/reference.html#layergroup but haven't been able to successfully do it.

This is what I have so far:

    var map = new L.Map('map', {
      zoomControl: true,
      center: [41.460, -81.67],
      zoom: 11
    });

    map.locate({setView: true, maxZoom: 18});

    var background = L.tileLayer('http://{s}.tile.stamen.com/toner-lite/{z}/{x}/{y}.png', {
      attribution: 'Stamen',
      maxZoom: 15,
      errorTileUrl: 'http://fashiontech.files.wordpress.com/2008/08/plaid_b_single.jpg?w=256'
    })
     .setZIndex(0),

    cdb = cartodb.createLayer(map, 'http://cleveland-metroparks.cartodb.com/api/v2/viz/05af5a7c-668d-11e3-8914-cb7f69afe66f/viz.json')

     .on('done', function(layer) {

      layer.setInteraction(true);

      layer.on('featureOver', function(e, pos, latlng, data) {
        cartodb.log.log(e, pos, latlng, data);
      });

      layer.on('error', function(err) {
        cartodb.log.log('error: ' + err);
      });
    }).on('error', function() {
      cartodb.log.log("some error occurred");
    }),

    labels = L.tileLayer('http://{s}.tile.stamen.com/toner-labels/{z}/{x}/{y}.png', {
      attribution: 'Stamen',
      maxZoom: 18
    })
     .setZIndex(600);

    var existing = L.layerGroup([background, cdb]);
    var baseLayers = {
        "Existing": existing
    };

    L.control.layers(baseLayers).addTo(map);

This answer hints that it's possible, but my js-fu is really lacking:
How does Leaflet layer ordering work with CartoDB layers?


texts: Problem on changing EPSG:3376 to EPSG:4326 using Proj4js

texts:
I been using the proj4js to transform the projection1 (EPSG::3376) to projection2 (EPSG::4326) with the data location of [66333, 170175] the outcome is supposed to be [110.281820983;1.5390299426] where I test it on the https://mygeodata.cloud/cs2cs/ and test it on https://epsg.io but when it transforms it turns to [115.59856224058234, 5.5388802885425905].

Can I know what are the problem with it?

Here is my code

    var firstProjection     = '+proj=omerc +lat_0=4 +lonc=115 +alpha=53.31580995 +k=0.99984 +x_0=0 +y_0=0 +no_uoff +gamma=53.13010236111111 +ellps=GRS80 +units=m +no_defs';
    var secondProjection    = '+proj=longlat +datum=WGS84 +no_defs';
    var point               = proj4(firstProjection,secondProjection,[66333,170175]);


texts: Result of 'Rasterize' shows only "0" data

texts:
I have point shapefile "TAEBAEK".
And i want convert it to raster.

Rasterize works well.
but result show only "0" data.

I thought result shows "0" and "1" data.
(the cell, containing the point, show "1". And the others show "0")

help me~


texts: Adding OpenStreetMap basemap from Python script in ArcGIS?

texts: Maybe it's too simple, but: Is there a way to add the OpenStreetMap basemap from Python-Script in ArcMap?


texts: Create a One-to-Many Spatial View with ArcGIS 10.1

texts:
We have two datasets: parcels and owners.  This is a many-to-many relationship, where a parcel can have multiple owners and there may be multiple properties on a single parcel (also with multiple owners: condos or strata property)

We would like to publish this data to ArcGIS Server, and issue queries on it, so that we can get a single result for a single condo owner.  Currently, when mapping the condo owners back to the parcel, we get all other condo owners, which gets confusing if the result was obtained by a text search.

In ArcGIS 10 and earlier, it's possible to create a spatial view that generates a unique Object ID field for repeated geometries.  This lets us publish a parcels layer with stacked parcels for owners.

Is this possible with ArcGIS 10.1?  The online documentation implies that ObjectID and Geometry must be obtained from the same feature class in a spatial view.

The underlying database is SQL Server.  We do not want to convert the data storage type.


texts: Spatial join based on a categorising attribute

texts:
I have two shapefiles that I want to spatial join (closest match), however, in the second shapefile I have records that are split by attribute X.
I want the closest to records from each unique X record to shapefile nodes. i.e. the closest match usually pick the closest record from shapefile 2 to shapefile 1 nodes, but I want the closest from each X attribute value to each of shapefile 1 nodes.

An example would be:
Shapefile 1: different schools around the city.
Shapefile 2: Kids around the city Attribute X: Age category (1,2,3).
I want the closest kid from each age category to every school.

How could this be done?


texts: Displaying forest landcover using Google Earth Engine

texts:
My question concerns code from Displaying landcover classes using Google Earth Engine.
I tried to apply the code in the answer to the new V20 data set, but I doesn't seem to work the same way. If I change the version it doesn't visualize the data and everything seems to be masked.
To try and solve the problem and because I'm interested in the 2018 data I changed the code to:
//Load Corine Image 2018
var corine = ee.Image('COPERNICUS/CORINE/V20/100m/2018');
print('CORINE:', corine);
Map.addLayer(corine, {}, 'CORINE');

// Get location of forest (=23-25)
//print(corine.propertyNames());
//var lc_value = corine.get('landcover_class_values');
//print(lc_value);

//Mask all layers that are not forest
var forest = corine.updateMask(corine.gte(23).and(corine.lte(25)))
print('Forest:', forest);
Map.addLayer(forest, {}, 'Forest');

but I still can't visualize the data.
Can you tell where the problem is?
I'm new to earth engine and still learning the syntax.


texts: Creating map booklets or flyers (not mapbooks)

texts:
I'm wondering if there is existing extensions or add-ons for ArcGIS 10.1 where I can create more sophisticated printable maps. I'm aware of datadriven pages and mapbooks, but I'm looking for something that's more like a booklet or a folded pamphlet  that would be given away at say a tourism booth. I'd like to have separate pages for a legend and map (or different "panels"). Ideally, this would be integrated into ArcGIS so that as points change the legend would update as well. I'm aware I could do this within ArcGIS with some clever measuring - but what if I wanted to print double sided all in one go?

So far I've found this:

MapLogic Layout Manager

As another option I would not be opposed to exporting a map as a PNG and inputting it into another program to make booklets. Are there any programs that are more map centric?


texts: How to save geodatabase attachments to disk with updated exif data

texts:
I have a folder of photos without xy information. Using a spreadsheet and lots of patience, I attached the photos to their corresponding GPS point as geodatabase attachments. I now need to produce a kmz file of these photos. Normally I would add the folder to Picasa and export to Google Earth, but since the original photos are not geotagged I cannot do that.

After looking at this question, it would appear that gaining access to the temporary files created when previewed in ArcMap would solve the problem. However, this would entail clicking on each point with the HTML popup on so that a temporary jpg is generated. Upon inspecting the properties of the temporary jpg, I saw that the coordinates of the associated point is not written into the EXIF data, so even if I were to save these temporary files, they still would not have location information.

Is there a way to geotag these attachments with the coordinates of their attached points, and save these files back to disk? I would like to do this in Python, but the part about getting the attached jpg files to generate without physically clicking on each point is a concern.

Edit I found a post on the ESRI forum which extracts the attachments to disk exactly as requested. I will add in the solution to extracting the coordinates and geotagging the photos, and post it here.


texts: I want re-project my Esri Arcinfo ASCII (.asc) raster file to Albers Equal Area Projection Which is currently in WGS84 and output it as .asc file. But gdalwarp cannot output as .asc file but gdaltranslate can. Is there anyway that i can use gdalwarp output as gdaltranslate input.


texts: Generate list of legend's RGB colors and values in ArcGIS

texts:
I am writing a node/typescript app that analyzes raster data.  It uses leaflet and esri-leaflet quite a bit (modified to work server-side).  For a given map bounds, it can download a raster image from ArcGIS server, and for a given latlng, it can tell you a pixel's rgba value from the raster image.
What I want to do is then map the RGBA value to its meaning according to the map's legend. Take this layer for example:

There are an absolute ton of values available in the legend.  I don't want to have to guess at the rgba values, or use some kind of color-getter app to do it one by one.  I need a legend that has the actual RGB(A) values of these colors listed, so I can then map my latlng value back to a meaningful value.
Most of this code runs in a node server, but to give you something to play with, I reacreated a playground in a codesandbox. Let's say I'm working with layer 30 of the LANDFIRE dataset, which is vegetation condition class. I created a map which fetches the raster image based on the map's bounds and zoom, writes the image to a canvas, and then reads the RGB values on the fly:
Codesandbox Playground
Here's the app without the code.  It acts a little wonky sometimes to don't be afraid to reload, and give it a second on load, zoom, and pan events.
As you mouseover the map, you get the rgb values of the pixels instantly, because the whole dataset is downloaded in one swoop as an image.  For this particular raster data, arcgis online has this legend.  I want to be able to map the pixels back to the actual values represented by the colors, i.e. instead of seeing the rgb values in my little box, I want to see the ground cover type, tree cover %, or whatever the legend contains.
How can I get the legend's color values as a list of RGBA: Layer Value pairs programatically?


texts:
When I zoom too much on a raster the layer does not appear anymore.

How to make it appear anyway?


texts: Unique List of Image Dates from Image Collection in Google Earth Engine

texts:
I would like to create a list or array object containing a unique set of dates represented in an image collection for a given region of interest. I have a function that extracts and converts dates for each image to a number, and stores these as a list. However, there are often multiple images (tiles) for the same date, meaning the date list has many duplicates. How can I filter out replicate dates to get a unique date list?

See code here: https://code.earthengine.google.com/20de3ef240fd96d1b25f0b52b9dce5f3

Map.addLayer(roi,{})
Map.centerObject(roi,7);

//Dates of Interest
var start = ee.Date("2014-10-01");
var finish = ee.Date("2018-05-01");

///--------------------- Landsat Collection ---------------------------------------///
var landsat = ee.ImageCollection('LANDSAT/LC08/C01/T1_TOA')
.filterDate(start, finish)
.filterBounds(roi);

// Year-Month-Day Extract function
function ymdList(imgcol){
    var iter_func = function(image, newlist){
        var date = ee.Number.parse(image.date().format("YYYYMMdd"));
        newlist = ee.List(newlist);
        return ee.List(newlist.add(date).sort())
    };
    return imgcol.iterate(iter_func, ee.List([]));
}
print(ymdList(landsat))


texts: If keep resolution = False in gdalogr:cliprasterbymasklayer

texts:
What does it mean if i put keep_resolution = False in gdalogr:cliprasterbymasklayer?
In the tool itself in QGIS I could choose then a resolution which I want to have for the output, but here via scripting you can juts choose True or false, so either it keeps the reolution from input file or not. But if not, what res. does it take?


texts: Collecting street data to populate a graph stucture for routing

texts:
I'm implementing a basic digraph with weighted edges for street routing in a small geo area (600mx600m). How and where can I get the required data to populate such a structure? I've been reading about OpenStreetMap, but am unsure as to how I can export the required data into my application.

The kind of data I would need is - all the streets (nodes) in the specified area, their connecting streets (subsequent nodes), and the distances between them.


texts:
The formula for NDBI for landsat 8 is (NIR-swir)/(NIR+swir).
I am confused if it is (Band5-band6)/(Band5+band6) OR (Band5-band7)/(Band5+band7)?

(For Landsat 5TM, NDBI formula is (Band5-Band4)/(Band5+Band4), and wavelength of band 5 is between 1.55-1.75.
So, by refering the wavelength, I am guessing, below formula to be correct for Landsat 8.

(Band5-band6)/(Band5+band6), but I want to be sure.
Please help me clear my confusion.


texts: How to insert raster in PostGIS from java JDBC

texts:
I have been looking for a solution to store dynamically generated raster data into PostGIS from java and PostGIS JDBC.

I am working with geotools toolkit but i am unable to find any documentation that fit my needs.

Here is what I've tried:

I am currently importing com.vividsolutions.jts.geom.GeometryCollection using java jdbc and preparedStatement like so:

PreparedStatement p =
   conn.prepareStatement("INSERT INTO myGeometryTable (geom) VALUES (?)");

p.setObject(1, (PGgeometry) convertToPGgeometry(geometryCollection));
p.executeQuery();

So I though it would be the same for raster, but there is no PGraster type. My raster object is a GridCoverageLayer that implement RasterLayer.

I got a simple table created with this sql command

CREATE TABLE myRasterTable ("rid" serial PRIMARY KEY,"rast" raster);
PreparedStatement p =
   conn.prepareStatement("INSERT INTO myRasterTable (rast) VALUES (?)");

-- PGraster DOES NOT EXIST
p.setObject(1, (PGraster) myRasterType);
p.executeQuery();

Almost every solution to store raster in postgis involve using raster2pgsql. raster2pgsql create a sql file that serialize a raster image (geotiff in my case). The output is like so:

INSERT INTO "public"."test" ("rast") VALUES ('01000001000000F ... I think it is Hex format '::raster);
CREATE INDEX "test_rast_gist" ON "public"."test" USING gist (st_convexhull("rast"));

I don't know if I can do the same using java


texts: Importing CSV table to geopackage with PyQgis using CSV file name

texts:
I'm trying to import some CSV with the same name of the file in a geopackage using PyQGIS before opening a project, but what I'm using doesn't work

I'm using this code:

from qgis.core import *
import glob
import os

def openProject():
    ruta = QgsProject.instance().absolutePath()
    tmp = ruta + '/tmp'
    try:
        for csv in glob.glob(tmp + '/*.csv'):
            rutacsv = 'file:///' + csv + '?delimiter=;'
            tabla = QgsVectorLayer(rutacsv, os.path.basename(csv), 'delimitedtext')

            opt = QgsVectorFileWriter.SaveVectorOptions()
            opt.actionOnExistingFile = QgsVectorFileWriter.CreateOrOverwriteLayer
            QgsVectorFileWriter.writeAsVectorFormat(tabla, ruta + '/cdau.gpkg|layername' + os.path.basename(csv).split('.csv', 1)[0], opt)

I tried also this code:

from qgis.core import *
import glob
import os

def openProject():
    ruta = QgsProject.instance().absolutePath()
    tmp = ruta + '/tmp'
    try:
        for csv in glob.glob(tmp + '/*.csv'):
            rutacsv = 'file:///' + csv + '?delimiter=;'
            tabla = QgsVectorLayer(rutacsv, os.path.basename(csv), 'delimitedtext')

            opt = QgsVectorFileWriter.SaveVectorOptions()
            opt.actionOnExistingFile = QgsVectorFileWriter.CreateOrOverwriteLayer
            opt.LayerName = os.path.basename(csv).split('.csv', 1)[0]
            QgsVectorFileWriter.writeAsVectorFormat(tabla, ruta + '/cdau.gpkg', opt)

Although if I use this code, the table is stored but with the name of the geopackage (cdau).

from qgis.core import *
import glob
import os

def openProject():
    ruta = QgsProject.instance().absolutePath()
    tmp = ruta + '/tmp'
    try:
        for csv in glob.glob(tmp + '/*.csv'):
            rutacsv = 'file:///' + csv + '?delimiter=;'
            tabla = QgsVectorLayer(rutacsv, os.path.basename(csv), 'delimitedtext')

            opt = QgsVectorFileWriter.SaveVectorOptions()
            opt.actionOnExistingFile = QgsVectorFileWriter.CreateOrOverwriteLayer
            QgsVectorFileWriter.writeAsVectorFormat(tabla, ruta + '/cdau.gpkg', opt)

But I'm not able to see where the mistake is.


texts: Georeferencing vector PDF/SVG/DXF linework using QGIS?

texts:
I want to  create a small web application to visualize election results on a per poll site basis for my city. The result is going to be very similar to what has been done for Berlin here. Hopefully this helps in understanding what I'm going for.

The election raw data on the scale I need is published and I know how to visualize this on the grounds of a shapefile with separate polygons lining out the election districts. My problem is - as you could probably guess - that something like this doesn't already exist.

What actually is published though, is a PDF file lining out the borders of the districts I'm interested in:
http://muenster.de/stadt/stadtplanung/pdf/a3_stimmbezirk.pdf

I'm not keen on redrawing the linework myself, so I'm looking for a way to make a shapefile out of this. Luckily, the PDf includes the linework as separate vector objects, I extracted these into a SVG file:
https://www.dropbox.com/s/bn7698yrdh5tdqj/a3_stimmbezirk_ungrouped.svg

I already tried converting the SVG to a DXF and importing it into qgis 2.0.1. This actually works, but (as expected) places the vectors in the atlantic ocean and I can't find a way to edit it and fix this.

What can I do now to georeference these vectors with open source tools and finally produce a shapefile?


texts:
I am trying to compute the area between two longitudes (using the angle between them) but am having a hard time. I can do that more easily between latitudes but longitude reveals to be tricky.

Anyone has a good reference formula for that one?


texts: Pulling Data off PNG map scan using ArcGIS Desktop?

texts:
I have scanned a map with several trade routes on it and georeferenced it on top of a hydrology (river system) map of the U.S. The scanned map is a PNG file.

How can I pull off the trade routes?

End Goal: story map showing the evolution of the Rocky Mountain Fur Companys trade network.


texts: Raster calculator in Process Toolbox QGIS

texts:
I am not able to find the Raster Calculator in my process toolbox anymore, even though I worked with the Raster Calculator through the toolbox in the same version of QGIS (2.8) one month ago.

It wouldn't matter much because I still can run the raster calculator in the raster menu, However I have a problem in the Graphical Modeler because I need to use this algorithm and it is not available in the GDAL/OGR group.

Any ideas?


texts: Openlayers dynamic filtering/styling of geojson layer

texts:
I'm working with Openlayers and overlaying a geojson file representing species type over OSM. I have an HTML form that allows users to select a species type and I then need to filter the geojson layer based on the selection. I then want the map to redraw, displaying only the filtered features.

I've tried CQL and messed around with styling but I'm not getting very far. Does anyone have any experience of doing this?

Thanks in advance.

James


texts: How to iterate through feature class to produce multiple point density files in ArcGIS?

texts:
I'm trying to create separate point density layers for each feature within a feature class owning to differing density characteristics within each feature.

I would like to be able to accomplish this without having to create a new feature class for each feature before computing point density. Here is my current model:

Edit: At the moment my model returns an error

How can I create new point density feature classes for each feature independent of each other either in Model Builder or ArcPy?


texts: ArcGIS JavaScript 3.28 extent not working

texts:
I'm trying to get an extent to work on a web map app but keep getting an error when trying to pass the extent parameter:

<!DOCTYPE html>

    <meta http-equiv="Content-Type" content="text/html; charset=utf-8">
    <meta name="viewport" content="initial-scale=1, maximum-scale=1,user-scalable=no"/>
    <title>Simple Map</title>
    <link rel="stylesheet" href="https://js.arcgis.com/3.28/esri/css/esri.css">

The error returned is:

a.equals is not a function

And is happening right at map = new Map("map", {


texts: Comprehensive Intro to GIS Data Structures

texts:
I am looking for a comprehensive introduction to design and analysis of spatial data structures.

If possible, please only link to online resources unless the offline resources is especially of note.


texts: How to obtain MGRS 100Km square id for every such square on world map?

texts:
I'm trying to compile a list of MGRS 100km square identifiers for every such square (eg 4QFJ for Honolulu Hawaii); excluding the poles (85 to -85 latitudes). If possible it'd also be great to know which of these 100km squares intersect land.

Is there any way to compute this or is there a recommended data set to look at?