texts:
In my current project I'm updating my app from ArcGIS runtime 10.2 to 100.7. In 10.2 we add UI elements to the collection of overlays and have no issues zooming, but in 100.7 whenever I add a UI element to the collection of overlays there is a dead zone wherever the control is. Is there a way to work around this?

Below is the code to replicate the issue, I've tried the same code using the 10.2 ArcGIS Nuget package and the zooming works

MainWindow.xaml

<Window x:Class="EsriTestApp.MainWindow"
    xmlns="http://schemas.microsoft.com/winfx/2006/xaml/presentation"
    xmlns:x="http://schemas.microsoft.com/winfx/2006/xaml"
    xmlns:d="http://schemas.microsoft.com/expression/blend/2008"
    xmlns:mc="http://schemas.openxmlformats.org/markup-compatibility/2006"
    xmlns:esri="http://schemas.esri.com/arcgis/runtime/2013"
    xmlns:local="clr-namespace:EsriTestApp"
    mc:Ignorable="d"
    Title="MainWindow" Height="450" Width="800">
<Window.Resources>
    <Style TargetType="Button">
        <Setter Property="Padding" Value="5" />
        <Setter Property="Margin" Value="5,5,0,0" />
    </Style>
</Window.Resources>
<Grid Name="MainGrid">
    <Grid.ColumnDefinitions>
        <ColumnDefinition Width="Auto"/>
        <ColumnDefinition Width="*"/>
    </Grid.ColumnDefinitions>
    <Grid.RowDefinitions>
        <RowDefinition Height="Auto"/>
        <RowDefinition Height="*"/>
    </Grid.RowDefinitions>
    <esri:MapView x:Name="MyMapView"  Grid.Column="0" Grid.ColumnSpan="2" Grid.Row="0" Grid.RowSpan="2">
        <esri:Map>
            <esri:ArcGISTiledLayer x:Name="Basemap"
                                   Source="https://services.arcgisonline.com/ArcGIS/rest/services/World_Topo_Map/MapServer"/>
        </esri:Map>
    </esri:MapView>
    <Button Content="Create Rectangle"
            Click="OnButtonClick"
            Grid.Column="0"
            Grid.Row="0"/>
</Grid>

MainWindow.xaml.cs

private void OnButtonClick(object sender, RoutedEventArgs e)
{
    Random rng = new Random();
    MapPoint point = new MapPoint(rng.Next(1, 100), rng.Next(1, 100), SpatialReferences.Wgs84);
    var rect = new BigRectangle();
    GeoView.SetViewOverlayAnchor(rect, point);
    MyMapView.Overlays.Items.Add(rect);
}

BigRectangle.xaml

<UserControl x:Class="EsriTestApp.BigRectangle"
         xmlns="http://schemas.microsoft.com/winfx/2006/xaml/presentation"
         xmlns:x="http://schemas.microsoft.com/winfx/2006/xaml"
         xmlns:mc="http://schemas.openxmlformats.org/markup-compatibility/2006"
         xmlns:d="http://schemas.microsoft.com/expression/blend/2008"
         xmlns:local="clr-namespace:EsriTestApp"
         mc:Ignorable="d"
         d:DesignHeight="450" d:DesignWidth="800">
<Grid>
        <Rectangle Width="100"
                   Height="50"
                   Fill="FloralWhite"
                   Stroke="Black"/>
</Grid>

I'm not quite sure what to do here. I've tried a few different things like attempting to reroute the MouseWheel event on the BigRectangle.cs class the MapView.MouseWheel event but all of the code to send events to that event's handler seems to be internal (MapView.OnMouseWheel() etc). The only way I've been able to get around this is to create an event handler for the MouseWheel event on BigRectangle.cs and call the MapView SetViewpointScale method to zoom the map, but I can't quite figure out the scale I need to accomplish this.


texts: GeoServer vector tile preview fails to load SourceMap

texts:
Very new to GeoServer and web mapping in general. My end goal is to publish data from PostGIS through GeoServer and OpenLayers. I have PostgreSQL/PostGIS and GeoServer on localhost.
When I try to preview data (either my own or the data that is already available upon installing GeoServer) as cached Mapbox vector tiles (pbf), I get the error:

DevTools failed to load SourceMap: Could not load content for http://localhost:8080/geoserver/gwc/rest/web/openlayers3/ol.js.map: HTTP error: status code 404, net::ERR_HTTP_RESPONSE_CODE_FAILURE

I interpret this as that it is supposed to show some sort of base map, but fails? How do I solve this?


texts: Evaluating elevation change between multiple LiDAR datasets

texts:
I am hoping to get some insight from others who are familiar with comparing multiple LiDAR surveys to one another. I have two LiDAR surveys with similar extents that I would like compare for elevation change. I have the point cloud .las files for each of them, so I am able to manipulate the original data in any way that is necessary.

So far, I have created 3 ft. resolution (1/2 the sampling density) rasters for each of the surveys using the same (IDW) interpolation method. I was hoping that by creating a grid of the same size for each of the surveys that I would get rid of any artifacts in the subtraction product. This has not been the case thus far as there is grid-like striping in the subtraction product. These artifacts result in about a +/- 3 feet error in the product. Maybe with a 3 feet grid spacing this is the best that I can hope for??

Does anyone have a comparison method that works well for them?


texts: Getting Error 000539 from Raster Calculator?

texts:
I am new to Arc and 10.2 but can not seem to figure out why this won't work.

Executing: RasterCalculator Con(("SLOPE_SLO.tif"  >=  .015) & ("SLOPE_SLO.tif" <= .15,1,0)) D:\GEOG_480\Final2\RastCalSlope
Start Time: Tue Jul 29 22:24:21 2014
Con((Raster(r"SLOPE_SLO.tif")  >=  .015) & (Raster(r"SLOPE_SLO.tif") <= .15,1,0))
ERROR 000539: Error running expression: rcexec()
Traceback (most recent call last):
  File "<expression>", line 1, in <module>
  File "<string>", line 5, in rcexec
  File "c:\program files (x86)\arcgis\desktop10.2\arcpy\arcpy\sa\Functions.py", line 3620, in BooleanAnd
    in_raster_or_constant2)
  File "c:\program files (x86)\arcgis\desktop10.2\arcpy\arcpy\sa\Utils.py", line 47, in swapper
    result = wrapper(*args, **kwargs)
  File "c:\program files (x86)\arcgis\desktop10.2\arcpy\arcpy\sa\Functions.py", line 3617, in Wrapper
    ["BooleanAnd", in_raster_or_constant1, in_raster_or_constant2])
TypeError: expected a raster or layer name

Failed to execute (RasterCalculator).
Failed at Tue Jul 29 22:24:21 2014 (Elapsed Time: 0.02 seconds)


texts:
Most of my data is geocoding in cartodb except a handful.  I've verified the location exists using google maps but it's still not coding it.  For example, here's a location it's unable to code:

3641 S TACOMA AVE
TULSA, OK  74107
US

I've attempted to code using city & postal code with no success.


texts: Client Side Reprojected Tiles shown at the wrong place

texts:
Description

I have to use a WMTS-Service in EPSG:31254 and show it on an OpenLayers 5 Map in EPSG:3857.
The WMTS-Service to use is this one.

My Solution

I use proj4 and configure "EPSG:31254" according to epsg.io/31254

proj4.defs("EPSG:31254", "+proj=tmerc +lat_0=0 +lon_0=10.33333333333333 +k=1 +x_0=0 +y_0=-5000000 +ellps=bessel +towgs84=577.326,90.129,463.919,5.137,1.474,5.297,2.4232 +units=m +no_defs");

Then i configure OpenLayers to use proj4 and create an OpenLayers Map in EPSG:3857 (i know, its the default projection) with a default background at a given center and zoom level.

Then i add the TileLayer in EPSG:31254:

read WMTS Capabilities
build Options from Capabilities
build WMTS-Source with Options
add TileLayer with WMTS-Source

Problem

An Image is shown, but it's shown at the wrong place (e.g. the Italian Peninsula should overlap!)

Source Code

I created an example showing the problem over at codesandbox.io
What is wrong with my solution?


texts: Trimble coordinates in field not accurate

texts: I am using a Trimble Juno 3E and when I load the points taken in the field on GIS, it is off by a metre or so from the actual object (It doesnt fall on the exact object location). Any help to improve the accuracy please?


texts: Open source option to ArcGIS tool MosaicToNewRaster_management in python?

texts: Is there any open source option in python to mosaic rasters to one raster? The equivalent of the MosaicToNewRaster_management in ArcGIS.


texts: Calculating frost frequency in Google Earth Engine?

texts:
I want to work with frost and several metrics

First, I want to map the number of nights that have temperatures below 0 in a year. This code seems to do this.

Second, I want to calculate the number of days between the first and last night with temperatures below zero.

How can I get the frequency of frost for the month of October since 2002 or better still in a decidical way? Sum of values below 0 ° C divided amount of values without clouds.

Code: https://code.earthengine.google.com/20e6accd99765fca9fef98371e02c366


texts:
I want to fill the NA values using the average values of the nearest neighbors:

r <- raster(matrix(1:16, nrow=8, ncol=8))
r[r==12] <- NA


texts: Splitting discontiguous vector feature in QGIS

texts:
I'm a QGIS novice. I'm trying to correct a map of Europe in QGIS, and it's incorrectly assigned Corsica to Germany. I can't seem to split Corsica off of Germany because when I use "split features" to draw a line around Germany or Corsica, it says "No feature split done" because (I think) the line didn't actually cross a feature.

How do I split them?


texts: How to calculate the building coverage ratio (BCR) if buildings overlap with two or more land parcels?

texts:
I use the following SQL query to calculate the BCR of land parcels in PostGIS:

CREATE VIEW public.v_bcr AS SELECT
    land_parcels.gid,
    ST_Collect(buildings.geom)::geometry(Multipolygon, /*SRID*/) AS geom,
    ST_Area(ST_Collect(buildings.geom))/(ST_Area(land_parcels.geom)/100) AS bcr
FROM land_parcels, buildings
WHERE ST_Within(ST_Snap(buildings.geom, land_parcels.geom, 1e-5), ST_Snap(land_parcels.geom, buildings.geom, 1e-5))
GROUP BY land_parcels.gid;

This query considers two possible cases:

Land parcels contain more than one building.
Vertices of buildings are located at segments of land parcels.

But sometimes buildings overlap with two or more land parcels. Is there any way to split these buildings in the query posted above?


texts:
I've been writing a script that clips large datasets based on a selected feature. It works fine when I run it in a test workspace where only related layers are active and running. However, when I run the script in my organisation's main workspace where alot of unrelated layers are active I get the error "Runtime error : list index out of range". This has baffled me as I have the same layers active in this workspace as I have in the test workspace yet I am getting this error, could it be something to do with the dataframe? or something else? Anyway, here is the script (apologies about the length):

import arcpy
from arcpy import env
from arcpy import mapping

env.workspace = "C:\MyArcGIS\Clip_tool\uFMfSW\uFMfSW_Banded_Vector_SE00.gdb"

arcpy.env.overwriteOutput = True

edrn_node = "C:\\MyArcGIS\\EDRN_Shp\\EDRN_NODE.shp"
ID = "EDRN_043066"
query = '"LCC_DRN_ID"' + " = '" + str(ID) + "'"
polygon_extent102 = "C:\\MyArcGIS\\Clip_tool\\poly_extent102"
depth_1in100 = "uFMfSW_SE00_DEPTH_1in100_BV"
depth1in100 = "C:\\MyArcGIS\\Clip_tool\\uFMfSW\\uFMfSW_Banded_Vector_SE00.gdb\\dep1in100"
depth_1in30 = "uFMfSW_SE00_DEPTH_1in30_BV"
depth1in30 = "C:\\MyArcGIS\\Clip_tool\\uFMfSW\\uFMfSW_Banded_Vector_SE00.gdb\\dep1in30"
depth_1in1000 = "uFMfSW_SE00_DEPTH_1in1000_BV"
depth1in1000 = "C:\\MyArcGIS\\Clip_tool\\uFMfSW\\uFMfSW_Banded_Vector_SE00.gdb\\dep1in1000"
extent_1in30 = "uFMfSW_SE00_EXTENT_1in30_BV"
extent1in30 = "C:\\MyArcGIS\\Clip_tool\\uFMfSW\\uFMfSW_Banded_Vector_SE00.gdb\\ext1in30"
extent_1in100 = "uFMfSW_SE00_EXTENT_1in100_BV"
extent1in100 = "C:\\MyArcGIS\\Clip_tool\\uFMfSW\\uFMfSW_Banded_Vector_SE00.gdb\\ext1in100"
extent_1in1000 = "uFMfSW_SE00_EXTENT_1in100_BV"
extent1in1000 = "C:\\MyArcGIS\\Clip_tool\\uFMfSW\\uFMfSW_Banded_Vector_SE00.gdb\\ext1in100"

mxd = arcpy.mapping.MapDocument("CURRENT")

arcpy.SelectLayerByAttribute_management ("edrn_node", "NEW_SELECTION", query)

df = arcpy.mapping.ListDataFrames(mxd, "Layers") [0]
df.zoomToSelectedFeatures()
df.scale = 2000

dataframe = mapping.ListDataFrames(mxd, "*")[0]
frameExtent = dataframe.extent
XMAX = frameExtent.XMax
XMIN = frameExtent.XMin
YMAX = frameExtent.YMax
YMIN = frameExtent.YMin
pnt1 = arcpy.Point(XMIN, YMIN)
pnt2 = arcpy.Point(XMIN, YMAX)
pnt3 = arcpy.Point(XMAX, YMAX)
pnt4 = arcpy.Point(XMAX, YMIN)
array = arcpy.Array()
array.add(pnt1)
array.add(pnt2)
array.add(pnt3)
array.add(pnt4)
array.add(pnt1)
polygon = arcpy.Polygon(array)

arcpy.CopyFeatures_management(polygon, "polygon_Extent102")

arcpy.Clip_analysis(depth_1in30, "polygon_Extent102", "Depth1in30")
arcpy.Clip_analysis(depth_1in100, "polygon_Extent102", "Depth1in100")
arcpy.Clip_analysis(depth_1in1000, "polygon_Extent102", "Depth1in1000")
arcpy.Clip_analysis(extent_1in30, "polygon_Extent102", "extent1in30")
arcpy.Clip_analysis(extent_1in100, "polygon_Extent102", "extent1in100")
arcpy.Clip_analysis(extent_1in1000, "polygon_Extent102", "extent1in1000")

arcpy.Delete_management("polygon_Extent102")

I've also noticed when I convert the above to a script tool that I have to reconnect the geodatabase files before I get the expected results.


texts: Discrepancy in values between Summed Line Lengths and Clip in QGIS

texts:
I have

a polygon vector layer containing four polygons, corresponding to four Protected Areas in a state in India.
a line vector corresponding to all the rivers in that state.

Desired output: I want to find out the proportion of the total river length flowing through Protected Areas.

Results so far:

I tried using 'Sum Line Lengths' under Vector Analysis to get at this, and then used the 'Basic Statistics for Numeric Tools' to get the summed length of all the rivers falling within the Protected Areas, by putting the river layer in the 'Lines' layer and the Protected Area layer in the 'Polygons' layer, in the 'Sum Line Lengths' dialog box. The value I got was 368.62 km.

However, the new layer looks the same as the polygon layer of Protected Areas in the workspace, so I decided to use 'Clip' as well, where the river layer was the Input layer and the  Protected Area layer was the Clip layer.

I then ran a '$length' operation on the field calculator to get the length of all the river segmens, and then used the 'Basic Statistics for Numeric Tools' to get the summed lengths. Now I'm getting a value of 552.33 km.

I'm a little confused with this, as I was assuming I should get the same value for both. Does the 'Sum Line Lengths' work differently than a straight up 'Clip'?


texts: Error using PST(Place Syntax Tool) in QGIS

texts:
I'm working with the PST(Place Syntax Tool) that I want to create segment map (Geometry is LineString) in QGIS.

I'm having the problems with these errors below

An error has occurred while executing Python code:

Exception: QGISModel: Column '' not found!
Traceback (most recent call last):
  File "/Users/chenyuru/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/pst/ui/pages/progresspage.py", line 64, in start
    task.run(dlgt)
  File "/Users/chenyuru/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/pst/analyses/createsegmentmapanalysis.py", line 74, in run
    extra_column_type = model.columnType(props['in_network'], copy_column_name)
  File "/Users/chenyuru/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/pst/model/model.py", line 67, in columnType
    field = layer.fields()[self._safeFieldIndex(layer, column)]
  File "/Users/chenyuru/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/pst/model/model.py", line 338, in _safeFieldIndex
    raise Exception("QGISModel: Column '%s' not found!" % name)
Exception: QGISModel: Column '' not found!


texts: How to build a web site or mobile app with a circular routing capability?

texts:
Hoping this is an easy one...

I want to build a web/mobile app for a college project which allows the user to enter a start point and a distance to travel e.g. distance to travel = 6km and the site would then plot a 6km circular route which brings the user back to where they started without (if possible) backtracking on any of the route already travelled.

I've been looking at OpenStreetMap and OpenLayers or Leaflet.js.

I'm not sure if this explains it very well but I suppose the app might work by calculating a 3km "end" point with a 1.5km waypoint on the east and a 1.5km waypoint on the west. Prefer to use open source technology.

I'm new to this so appreciate advice.


texts:
I am fairly new to actively work with databases. I wanted to try out some stuff on my home computer where I have installed pgAdmin and PostgreSQL previously. I can connect to my localhost in pgAdmin but I can't do the same in QGIS. I just get an error message saying:

"Connection failed -  Check settings and try again"

I am new in this game and don't really know where to look for solutions. My database connection info is standard settings:

Name: postgres
Host:Localhost
Port:5432
Database:test

Any ideas what could be wrong?


texts: PostGIS: Having trouble getting a value from a raster

texts:
I'm working on a subsetted grib file that I retrieved using this URL:
http://nomads.ncep.noaa.gov/cgi-bin/filter_gfs_0p25.pl?file=gfs.t12z.pgrb2.0p25.f001&lev_925_mb=on&var_RH=on&subregion=&leftlon=-126&rightlon=-66&toplat=50&bottomlat=24&dir=%2Fgfs.2017020912

First I opened the grib file in zyGrib to confirm that the file has the right data, i.e., relative humidity is appearing for lower 48 US states at an altitude of 925hPa. So far so good.

Next, I try loading the grib file into postgres:

raster2pgsql -s 4326 gfs.t12z.pgrb2.0p25.f001 | psql

Then, I confirm the new raster exists, has one band, and SRID=4326:

SELECT ST_Height(rast) As height,
ST_Width(rast) As width,
ST_SRID(rast) As srid,
ST_NumBands(rast) As num_bands,
ST_BandPixelType(rast,1) As btype
FROM public."gfs.t12z.pgrb2.0p25";

Next, I try querying a point located in the USA:

SELECT ST_Value(rast, ST_SetSRID(ST_Point(-87.9080059, 41.9773201), 4326))
FROM public."gfs.t12z.pgrb2.0p25";

I'm getting the following warning:

NOTICE:  Attempting to get pixel value with out of range raster coordinates: (-1288, 32)

And a blank value is returned in the resultset.

So next, I try gdalwarp as suggested in this answer: https://gis.stackexchange.com/a/33524/12291

gdalwarp -t_srs EPSG:4326 gfs.t12z.pgrb2.0p25.f001 warped.grib2

Then import:

raster2pgsql -s 4326 warped.grib2 | psql

Same as before, confirm new raster exists, has one band, and SRID=4326. Results match last attempt.

SELECT ST_Height(rast) As height,
ST_Width(rast) As width,
ST_SRID(rast) As srid,
ST_NumBands(rast) As num_bands,
ST_BandPixelType(rast,1) As btype
FROM public.warped;

Querying a point located in the USA from the new warped raster:

SELECT ST_Value(rast, ST_SetSRID(ST_Point(-87.9080059, 41.9773201), 4326))
FROM public.warped;

Gives the same result ... warning message and blank value in resultset :-(

NOTICE:  Attempting to get pixel value with out of range raster coordinates: (-1288, 32)

What am I doing wrong here?

EDIT: here's the gdalinfo for this grib file:

gdalinfo gfs.t12z.pgrb2.0p25.f001
Driver: GRIB/GRIdded Binary (.grb)
Files: gfs.t12z.pgrb2.0p25.f001
Size is 241, 105
Coordinate System is:
GEOGCS["Coordinate System imported from GRIB file",
    DATUM["unknown",
        SPHEROID["Sphere",6371229,0]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433]]
Origin = (233.875000000000000,50.125000000000000)
Pixel Size = (0.250000000000000,-0.250000000000000)
Corner Coordinates:
Upper Left  (     233.875,      50.125) (233d52'30.00"E, 50d 7'30.00"N)
Lower Left  (     233.875,      23.875) (233d52'30.00"E, 23d52'30.00"N)
Upper Right (     294.125,      50.125) (294d 7'30.00"E, 50d 7'30.00"N)
Lower Right (     294.125,      23.875) (294d 7'30.00"E, 23d52'30.00"N)
Center      (     264.000,      37.000) (264d 0' 0.00"E, 37d 0' 0.00"N)
Band 1 Block=241x1 Type=Float64, ColorInterp=Undefined
  Description = 92500[Pa] ISBL="Isobaric surface"
  Metadata:
    GRIB_COMMENT=Relative Humidity [%]
    GRIB_ELEMENT=RH
    GRIB_FORECAST_SECONDS=3600 sec
    GRIB_PDS_PDTN=0
    GRIB_PDS_TEMPLATE_NUMBERS=1 1 2 0 96 0 0 0 1 0 0 0 1 100 0 0 1 105 84 255 0 0 0 0 0
    GRIB_REF_TIME=  1486641600 sec UTC
    GRIB_SHORT_NAME=92500-ISBL
    GRIB_UNIT=[%]
    GRIB_VALID_TIME=  1486645200 sec UTC

EDIT #2: Here is the information zyGrib shows for this file. Somehow it is generating the correct area for this file despite missing datum in the source file.

File : gfs.t12z.pgrb2.0p25-7.f001
Size : 31811 bytes

Weather center 7 - Model 96 - Grid 0
1 dates:
    from Thu 2017-02-09 13:00 UTC
    to Thu 2017-02-09 13:00 UTC

Avalaible data :
* Relative humidity: 925hpa

Grid : 241x105=25305 points : 0.25°x0.25°
Area :24.00°N 126.00°W  ->  50.00°N 066.00°W

Reference date: Thu 2017-02-09 13:00 UTC


texts:
I am trying to use proj4js to convert between WGS84 (4326) to NAD27 BLM 14N (ftUS) (32064).  I am using this definition for EPSG 2674

`proj=tmerc +lat_0=0 +lon_0=-99 +k=0.9996 +x_0=500000.001016002 +towgs84=-8,161,176 +y_0=0 +ellps=clrk66 +to_meter=0.3048006096012192 +no_defs "`

I have also tried:

proj=tmerc +lat_0=0 +lon_0=-99 +k=0.9996 +x_0=500000.001016002 +datum=NAD27 +y_0=0 +ellps=clrk66 +to_meter=0.3048006096012192 +no_defs "

I got the transformation Parameters from here and have gone through these three similar SO posts one, two, three. My proj4js transformation remains about 100ft off. I'm validating the result using both corpscon and blue marble.  I can accept the limitations of a 7 parameter transformation but can't find suitable values for the Gulf of Mexico.  3 parameters are not tight enough.  How do I form the proj definition to get this to work.

The coord pair I am using to test is 25.340631, -96.647033 the expected result for EPSG 26714 is 9201257.09usft,2417471.02usft.


texts:
I have created a map to display 2011 census data using .shp files of output areas in England and Wales and joined this to a csv file which has the census data for each output area. I then display he data using graduated colours.

This works very well and displays maps showing the data for each output area.

Now I would like to also be able to display how many people who have attribute X within the output area, but also (separately) how many there are within a radius of N metres of each output area (let's say how many people there are over 65 within 1 mile of Output Area G (or the centre of the Output area if that is simpler)

Is there a way to do this?


texts: Received corrupt data from Overpass (incomplete polygon) with Overpass API Python wrapper

texts:
I'd like to collect forest data on overpass turbo. But my code isn't working properly.
Here is the first version:

import overpass

api = overpass.API()

res = api.get("""
 area[name="Charente-Maritime"]->.searchArea;
 (
 way[landuse=forest](area.searchArea);
 relation[landuse=forest](area.searchArea);;
 way[natural=wood](area.searchArea);
 relation[natural=wood](area.searchArea);
 );
 (._;>;);
 out body;
""")

Here's the error message:

Traceback (most recent call last):
  File "D:/3_PROJETS/00_TD/TD7-analyse_spatiale/02_TRAVAIL/forest.py", line 16, in <module>
    """)
  File "C:\python35\lib\site-packages\overpass\api.py", line 130, in get
    return self._as_geojson(response["elements"])
  File "C:\python35\lib\site-packages\overpass\api.py", line 214, in _as_geojson
    raise UnknownOverpassError("Received corrupt data from Overpass (incomplete polygon).")
overpass.errors.UnknownOverpassError: Received corrupt data from Overpass (incomplete polygon).

I've tried many different syntaxes but none of them have worked so far.


texts: Is there an active equivalent to Yahoo GeoPlanet available for commercial use?

texts:
It looks like GeoPlanet is dead (Is GeoPlanet dead?).

Is there a suitable alternative that covers Parents (preferably Ancestors) and Neighbours? If there's a name for this type of data, that'd be helpful.

Spatial data would be a plus, but I can get by without it for now.

This project may come to nothing, but in case I am ever able to charge for the service I need to make sure I don't shoot myself in the foot with a restrictive licence.


texts: Network Analyst provides different values for Total length and Shape length

texts:
I'm using the Network Analyst extension in ArcMap 10.2 to calculate the shortest routes between ~50000 pairs of stops. As the impedance is set to "Length (Meters)" the resulting Routes Featureclass has a field for "Total_Length" and "Shape_Length". What is confusing me, is that the values in these two fields are not the same, they differ by +15m to -0.5m. At first I assumed that this had something to do with the distance of the stops to the network, but looking more closely at some of the routes, the distance to the network doesn't match the difference in length. I've also re-run my analysis with a 0m offset of the Stops from the network, but the length differences stayed the same.

Has anybody else come across this before or does anybody have an idea of why this is happening?


texts: PostGIS multiple conditional statement on spatial join

texts:
I want to make a multiple conditioned spatial query in PostGIS 2.0/PostgreSQL 9.3., OS Windows 7.
There are 2 tables, cadastre (Polygon) and roads (Linestring).

Query: Select the sum (one number) of the length of all roads which are longer than 500 m, and which intersect cadastral parcels with area greater than 2000 m2.

This code produces duplicate sums for all roads >500 intersecting parcels >2000, which is obvious if the sums are grouped by parcel id.

SELECT sum(ST_Length(r.geom)) AS total_length
FROM roads AS r
JOIN cadastre AS c
ON ST_Intersects(r.geom,c.geom)
WHERE ST_Length(r.geom)>500 AND ST_Area(c.geom)>2000;

Duplicate values of road lengths:

So, how to avoid the multiplication of values, and get the real aggregated number?


texts: Export to OpenLayers plugin doesn't work with raster and WMS

texts:
I am trying to export a QGIS_2.4 project to OpenLayers using the "Export to OpenLayers" plugin. I cannot export a raster layer. I am not getting any error message, but the raster layer doesn't appear in the preview window.

I have followed all rules in the plugin web page. Export to OpenLayers

Raster layers: They will be exported to JPEG files and reprojected into the default projection of the OL3 map to create (as EPSG:3857). The symbology is not used, so the layer must be a rendered image already. You can use the QGIS Save as rendered image functionality to get such an image from your data and symbology.

Any idea of what I am missing?


texts: GEE error code while exporting Sentinel-2 cloud masked image

texts:
I have filtered a Sentinel-2 image and created a cloud mask for it based on its QA60 band. It works fine and displays with this mask applied. However, when I try to export it to a drive, I get the following error message:

Error: Image.clipToBoundsAndScale, argument 'input': Invalid type.
Expected type: Image. Actual type: ImageCollection.

Why is this happening and how can I fix it? My code is below.
    //define the cloud mask layer
function maskS2clouds(image) {
  var qa = image.select('QA60');
  var cloudBitMask = 1 <<10;
  var cirrusBitMask = 1 <<11;

//layer called mask = the cloudy bits
  var mask = qa.bitwiseAnd(cloudBitMask).eq(0)
  .and(qa.bitwiseAnd(cirrusBitMask).eq(0));

//'image' is a layer with the mask applied to it
  return image.updateMask(mask).divide(10000);
}

//image collection is filtered
var dataset = ee.ImageCollection("COPERNICUS/S2")
  .filterBounds(geometry)
  .filterDate('2015-08-26', '2015-08-28')
  .map(maskS2clouds);

//image with mask applied
Map.addLayer(dataset);

//exported to drive
Export.image.toDrive({
    image:dataset,
    description: '2015_CLOUDFREE',
    folder: 'users/emilywest2',
    region: geometry.geometry,
    scale: 10,
});


texts: Seeking web service that returns elevation profile between two points and can process 5000 requests per second

texts: I need a web service that provides elevation profile data between two coordinates with request rates at a minimum of 2000 requests per second but preferably at 5000 per second. The distance between the coordinates is on average 50 kilometers, and I need perhaps 10 or 20 elevation sample points between the two coordinates. What are the logistical and financial means to make this possible?


texts: Assigning categories to a single band raster in Google Earth Engine?

texts:
I'm trying to assign bins to a raster containing values between -1 and 1 in GEE. The bins would be 0.41<x<0.60, 0.61<x<0.80 for example. Is there a function for this?

The only function I can find is image.remap which works for individual values rather than ranges.


texts: What are the units of the QGIS roughness tool results?

texts: When I estimate Roughness, using QGIS [Raster > Analysis > DEM (Terrain models)], what is the unit of the values I get?


texts: How to extend a (river water level) raster, transfering its values to the added area?

texts:
Exploring the raster tools in QGIS and related software, I found raster buffer functions to surround existing rasters. The values of the raster buffer are restricted to be distances from the interior raster or constant however. I would rather like to derive the values of the added cells from the boundary of the interior raster...

Data: 1. river flood level raster; 2. polygon at the side of 1.; 3. bended river axis as polyline; 4. DEM.

Aim: I'd like to laterally extend the river flood level raster, which is restricted between dikes, constraining the flood plain. Some dikes are to be removed and additional flood retention area (delineated by polygon 2.) will be flooded.

It's not about modelling the flow or change of water level. - I just want to transfer the given flood levels at the boundary of the raster to the neighboring area, though keeping the overal water level inclination along the main flow direction of the river.

The laterally extended flood water level will be used to approximate inundation depths by subtracting DEM values in a second step. I don't need very precise water level results, because the study is about sedimentation and soil contamination mainly.

Any idea how to laterally extend the (flood level) raster values by use of QGIS & friends?


texts:
I have secured my layers on Geoserver. Users can not access my layers without authenticated username and password.

Now, I created a proxy script with PHP. Proxy script sends request with username password.

http://myproxydomain.com/proxy.php is my proxy address.

And I am using this proxy on my OpenLayers application to access WMS and WFS:

  OpenLayers.ProxyHost="http://myproxydomain.com/proxy.php?url=";

  var myWMSLayer = OpenLayers.Class(OpenLayers.Layer.WMS, {

      getURL: function (bounds) {
         var url = OpenLayers.Layer.WMS.prototype.getURL.call(this, bounds);
         if (OpenLayers.ProxyHost && OpenLayers.String.startsWith(url, "http")) {
             url = OpenLayers.ProxyHost + encodeURIComponent(url);
         }
         return url;
      }
   });

As a result, everyone can use my services that knows my proxy url. My application is in Javascript and you can see my proxy url from browser view code. And you can add your my proxy to your own Openlayers application and access my layers.

How can I solve this problem?


texts: Aggregating geometry features on geometry simplification using PostGIS?

texts:
PostGIS has ST_Simplify functionality which gets one multiline as input and returns one output multiline where the general shape of lines are retained while eliminating unnecessary detail.

Instead of one multiline, I have a table with a LineString geometry column (street segment) and an integer column (traffic velocity of that street). I can create a single multiline out of that LineString column and simplify it, but I need to know traffic velocity of simplified geometries too. Is there any tool to aggregate data of integer column into simplified lines?


texts: Binding IActiveViewEvents::ItemReordered to method in my ATL class?

texts:
How I can wire this event to my method?

I have an ATL class and I need to know when a layer added to map object so that I can update my list of layers in this class.My class implements IExtension and I've put this code in Startup method of this interface:

//Some code to obtain a reference to Map object
__hook(&IActiveViewEvents::ItemReordered,pMap/pointer to map/,&MyClass::MyEventHandler);

and in Shutdown method of IExtension I unhook it,but I receive compile error says the source and handler method not compatible while both signatures,arguments and names are exactly same.

Is it true a all this method?any simpler method for this wiring?or any reference on net?


texts: can i use layername instead of layerid in identifyParams using arcgis api javascript 3.17?

texts:
Is it possible to use layer name instaed of layerid in identifyparams.
However, over time, the layer order will change because I'll be adding new layers. Therefore, I want to refer to the layers by name, instead of id.can anyone suggest me please (am using arcgis api javascript api 3.17  example
https://developers.arcgis.com/javascript/3/sandbox/sandbox.html?sample=find_popup)


texts:
In QGIS there are "Advanced Editing Tools", which are able to measure distances and angles of the feature being digitized as you are drawing it.

Does ArcGIS Pro also have this capability?

I have a Basic level license. I haven't been able to find any sort of advanced editing options where measurements are provided, or angles given for that matter.

I would like to use ArcGIS Pro as I can easily work with CAD (.dwg) data, which is quite important for me, as well as having access to ESRI's basemaps versus Quick Map Services in QGIS as it being used for commercial purposes.

Image of the "Advanced Digitizing Tool". While digitizing a feature (i.e. when you start a line or polygon) "d" for distance in metres starts at 0 and updates as you move your cursor/next point. This allows for precise dimensions to be created. Similarly it can be set to snap to various angles while calculating the distance. If ArcGIS Pro has this (fairly basic) functionality it would be fantastic.

]1


texts: Cropping Oblique Photos to Geographic Feature

texts:
Let me preface this by stating that I know I'm in the GIS neighborhood and I'm asking a code heavy question, but my hiccup here isn't the programming - it's the GIS and associated math.

I am writing some C# code that will crop an oblique photograph to a 1500 pixel square surrounding a geographic point and I am having trouble with finding the right formula to translate geographic coordinates to pixel coordinates with the following known variables:

Feature geographic shape and centroid (State Plane)
Photo geographic centroid (State Plane)
Photo geographic footprint (Polygon shapefile [State Plane])
Aircraft position at time of photo (State Plane)
Pixel coordinates of photo centroid
Aircraft altitude in feet
Focal length of lens (50mm)

The oblique photo metadata is stored in shapefiles and served to my application via ArcGIS server

Here is the code I have so far

        double originX = xg_0;  //Geographic centroid of photo (State Plane) (X)
        double originY = yg_0;  //Geographic centroid of photo (State Plane) (Y)
        double x0 = -1 * (originX);
        double y0 = -1 * (originY);

        aircraftX = aircraftX + x0;  //Aircraft position in state plane adjusted to (0,0) geographic origin at centroid
        aircraftY = aircraftY + y0;  //Aircraft position in state plane adjusted to (0,0) geographic origin at centroid

        double pixelCenterX = xp_0;
        double pixelCenterY = Math.Abs(yp_0);

        double featureX = xg_1 + x0;  //Feature point (or centroid) in state plane adjusted to (0,0) geographic origin at centroid
        double featureY = yg_1 + y0;  //Feature point (or centroid) in state plane adjusted to (0,0) geographic origin at centroid

        double pixelX = 0.0;
        double pixelY = 0.0;

        double groundDistanceToCenter = Math.Sqrt((Math.Pow((featureX - 0), 2)) + (Math.Pow((featureY - 0), 2)));  //Distance in feet from photo centroid to feature
        double groundDistanceToAircraft = Math.Sqrt((Math.Pow((featureX - aircraftX), 2)) + (Math.Pow((featureY - aircraftY), 2))); //Distance in feet from aircraft to feature
        double height = 5000.0; //aircraft altitud in feet
        double imageAngle = Math.Atan(height / groundDistanceToAircraft) * 180 / Math.PI; //viewing angle from feature to aircraft

        double u = featureX - (featureY * (1 / Math.Tan(imageAngle)));
        double v = featureY * (1 / (Math.Sin(imageAngle)));
        double imageStatePlaneX = u + (v * Math.Cos(imageAngle));
        double imageStatePlaneY = v * Math.Sin(imageAngle);
        pixelCoordinates = new Point((int)(pixelCenterX - imageStatePlaneX), (int)(yp_0 - imageStatePlaneY) + 4912);  //4912 is pixel height of image


texts:
I'm somewhat confused about regression models. I decided to create a GWR model for a project using 35 polygons of health data. However I find now that my model has flaws and I cannot create a GWR model probably because the dataset is too small.

When creating a OLS are we meant to remove all variables that are related? if I'm modelling obesity, I would think that income, exercise are factors that relate to it so what purpose does the OLS model serve? I created a model after removing auto correlated variables but I'm kind of confused as to what its meant to be showing.


texts: Calculating difference between bearing values in degrees

texts:
Suppose I have 3 geographic points, A, B and C.

I can calculate the bearing of B from A (in degrees) using this SQL snippet:

degrees(ST_Azimuth(geogA,geogB)) as brg_a_b

...and the same of C from B:

degrees(ST_Azimuth(geogB,geogC)) as brg_b_c

Now, I'm wanting to check if A, B and C are close to being in a straight line, let's assume I want the B -> C bearing to be within say, 2 degrees different to the A -> B bearing to consider them in a straight line.

Azimuth is measured clockwise from north, so North = 0º, East = 90º, South = 180º, West is 270º.

Say the A-B bearing is 1º, how do I calculate if B-C is within, say 2º (either way) - IE between the values of 359º and 3º?

I need to construct a SQL "WHERE" clause which effectively does this:
WHERE brg_b_c within 2º of brg_a_b (but copes with the 360-0 crossover)

I'm sure this is really simple but am having a mental block working it out.

Thanks for any pointers (or suggestions of achieving the same thing in a different way ?)


texts: Stop QGIS from printing true in python console when updating records

texts:
How can I stop QGIS from printing True twice in python console for each updated record? This is slowing down the calculations.

Example:

layer = iface.activeLayer()
with edit(layer):
    for feature in layer.getFeatures():
        if feature['Randint'] % 2 == 0:
            new_value =  'Even'
        else:
            new_value = 'Odd'
        feature.setAttribute(feature.fieldNameIndex('Parity'), new_value)
        layer.updateFeature(feature)

True*10000:


texts: How to georeference an areal photograph given known coordinates of the four corners and center of the image

texts:
First, let me say that I am new to this field.  I have an existing image/video viewing application.  My primary goal is to implement an efficient method to determine the latitude and longitude of any pixel in an image or video frame.  I know the latitude and longitude of the four corners and the center of the image/frame.  I also know the projection system that the coordinates are based in.

The images in question are not maps.  They are areal photographs.  I cannot make assumptions about the viewing angle or anything like that.

I'm not looking to use any existing application or tool.  I need to implement the necessary functionality in my own application.  I am open to using a third party library if that's the best way, but my application is proprietary.  So licensing is a concern.

I've been looking into this for a while now.  I've read up on world files and using them for rectification.  But I don't know how to create one for my image.  It seems like the tutorials that I have seen make assumptions for simplicity sake that I cannot make.

I'm basically looking for a description of the necessary math to perform the required analysis and/or a good library I could use that is compatible with proprietary software.  My existing application is written in c++ and c#.

I apologize if any of terminology is incorrect.  Any advice would be appreciated.


texts: Layers with the same Geographic Coordinate system don't overlap

texts:
I am a GIS newbie. I'm having a problem that, from the looks of things on this website, is a fairly common problem. I have 3 layers, all in the same GCS that do not line up. I have read all the the strings that others have posted on this website regarding this kind of problem, but I still haven't been able to figure out how to fix my particular issue. I am just a social studies student, trying to do a project. I don't have GIS training, so I feel very lost at sea, with all the complex technical terminology.

Here are the details:

I am working with data for the City of Chicago, USA. I have a street map for Chicago (chicagostreets1999, derived, I think from the US Census), a map of census tracts for the year 2000 (tracts2000, from the US census), and a zoning map (Zoning_nov2012, provided online by the Chicagy City Hall).

Census tract boundaries run along streets. The tract boundaries in tracts2000 line up nicely with the street map from chicagostreets1999. Zoning boundaries also follow streets, but the zoning map does not line up properly with these two other files. When I compare the zoning map to the street file, its lines are about 175 feet (53 meters) south, and 148 feet (45 meters) east of the relevant lines on the street map.

I looked at the "source" information on each file in arcGIS. Here is the info on their Geographic Coordination Systems.

Tracts 2000*****
Data Type:  Shapefile Feature Class
Geometry Type:  Polygon

Geographic Coordinate System:   GCS_North_American_1983
Datum:  D_North_American_1983
Prime Meridian:     Greenwich
Angular Unit:   Degree
![Misalignment, zoom in][1\

*Chicago streets 1999****

Data Type:  Shapefile Feature Class
Geometry Type:  Line

Geographic Coordinate System:   GCS_North_American_1983
Datum:  D_North_American_1983
Prime Meridian:     Greenwich
Angular Unit:   Degree

**Zoning_nov2012*******************

Data Type:  Shapefile Feature Class
Geometry Type:  Polygon

Projected Coordinate System:

NAD_1983_StatePlane_Illinois_East_FIPS_1201_Feet
Projection: Transverse_Mercator
False_Easting:  984250.00000000
False_Northing: 0.00000000
Central_Meridian:   -88.33333333
Scale_Factor:   0.99997500
Latitude_Of_Origin: 36.66666667
Linear Unit:    Foot_US

Geographic Coordinate System:   GCS_North_American_1983
Datum:  D_North_American_1983
Prime Meridian:     Greenwich
Angular Unit:   Degree

I tried setting the coordinate system for all three files by double clicking the "Layer" meta-category in the data view window, and then selecting "Coordinate System" within the data frame properties window. The coordinate system was already defined as:

*Layers****
NAD_1983_StatePlane_Illinois_East_FIPS_1201_Feet
Projection: Transverse_Mercator
False_Easting: 984250.000000
False_Northing: 0.000000
Central_Meridian: -88.333333
Scale_Factor: 0.999975
Latitude_Of_Origin: 36.666667
Linear Unit: Foot_US

GCS_North_American_1983
Datum: D_North_American_1983

This seemed to be the same as all three layers, so I did not change it.

I asked the GIS librarian at my university what the problem might be and what, if anything could be done, to fix the misalignment. He said that he did not know how to fix it, and that the problem was probably due to the City Hall GIS people scanning in their zoning map "at a different scale" than the Census people did. I don't know whether that is true, but I really would like to get these files to align properly with each other. It will make the process of analysis a lot easier (ie not horribly visually confusing), and also more accurate.

I am happy to share the files in question with anyone who might be able to help. The Zoning data, for instance, is freely available here on the City of Chicago's City Hall government website:

https://data.cityofchicago.org/Community-Economic-Development/Boundaries-Zoning-Districts/p8va-airx

I can email the street file to anyone who is willing to help.

Thanks so much for your advice!


texts: sdequery.exe: When using functions as column names (-C parameter), how do you escape the comma for functions with parameters?

texts:
Using sdequery, I'm specifically trying to return sys_context('USERENV', 'SESSIONID'). user and sysdate  work fine. Functions with single parameters work fine.

But as soon as I use a function with multiple parameters, I get:   Error: Stream not found (-107). The SDEINTERCEPT shows the query columns as 2 columns (split at the comma). Makes sense since the  comma is usually a column delimiter, but how do you escape commas that aren't delimiters?

I've tried preceding the comma with a backslash.  Double commas. Wrapping the entire thing in double quotes. Wrapping the comma in various quote combinations.

UPDATE
USERENV('SESSIONID') works and gives me what I need, but I'd still like to know if/how you can escape the comma for future reference.


texts: I'm working in QGIS 2.8.6 with layers with legends with many entries.Is it possible to mark or unmark elements of the legend with a single click?


texts: Adding points along line at varying distances in ArcMap?

texts:
I have a shapefile for a single (3 km) line generated from a start and end location of a trail. I would like to create a new shapefile with points along the line that are not of an equal distances apart. I have all the point labels with the distances in an Excel speadsheet.

I am familiar with "Contruct Points" in the ArcMap editor, but this does not allow for varying distances. For example, I would like the first point at the beginning of the line, the next 47 meters from the starting point of the line, the next point 46 meters from the next point point, the next 48 meters, etc.

How can I do this in ArcMap 10.5


texts: What kind of geo information is stored in a GeoTIFF

texts:
When given a geotiff, what kind of information is normally stored in it as part of the geo information?

Is it the position of the 4 corners?
Is it the position of each pixel?
Is it a mesh of geo points?

(Latest OpenCV include gdal support and the option to Open geotiff files and I was wondering what kind of information i will be able to find in it? as I have a goal of writing a small program that can tile a DEM model up in XYZ tiles using OpenCV)


texts: Accuracy assessment (confusion matrix) of image classification from Landsat imagery in SAGA GIS

texts:
I have classified Landsat images in two different classes, namely pervious and impervious, in SAGA GIS, using training sample in polygon shapes. This training sample has only one column of attribute table comprising of class value of each polygon (value = 1 for impervious and 2 for impervious). The total of polygon shapes is 65.
I created a new test sample for evaluating the accuracy of the classification result. The test sample also has attribute value with one column specifying the class value of each polygon. The total of polygon shapes is 45.
I tried to do the evaluation using Confusion Matrix (Polygons/Grid) module from SAGA GIS.
I specified the grid as the result of the classification and the polygons are from the test samples.
The result of the confusion matrix is shown in this figure:

I am struggling to understand this confusion matrix.
There are around 46000 features are assessed in the evaluation shown in the table.
I am wondering what are these features?
What method I should use to evaluate my classification? Should I just manually check my test samples instead whether it s correctly classified or not?
PS: This is the grid system of the classified Landsat image, and the classification is done using random forest classifier.


texts: How to save the shapefile with points and labels together?

texts:
I'm using QGIS and trying to save a point layer with labels. When I save that point layer and then try to open it separately. But now there are only the points no any labels which I fix earlier.

Why is that?

Before saving the point layer.

After saving the point layer and opened it.

What can I do in here to avoid this situation?


texts: How to install PostgreSQL and PostGIS on Ubuntu 13.10?

texts:
I have upgraded my OS from Ubuntu 12.10 to 13.10. I had a clean installation of 13.10 version.   Can someone please suggest which version of PostgreSQL and PostGIS can be installed on Ubuntu 13.10 version and how to install it?

I have looked at the steps given in How to install PostgreSQL 9.2 with PostGIS 2.0 on Ubuntu 11.10 (or higher)?? But the three repos are not adding after update.

Hence that answer is not helping me.


texts: CartoDB API to make a table public after an import?

texts:
I am currently using the import API to import my layer data into CartoDB from another application.

My question is whether there is a way to make this table PUBLIC using an API instead of the UI?


texts: PyQGIS map rendering export shows no map labels?

texts:
I'm trying to learn a little bit of Python in QGIS with the following tutorial: QGIS Maprendering.

Now the tutorial works, with minor adjustments.
In QGIS I have four layers, and one them is showing labels (country names), but when I export to PNG or PDF the labels of the layer won't show in the output file.

My code:

from qgis.PyQt.QtGui import *
from qgis.PyQt.QtWidgets import *
import processing
from qgis.PyQt.QtCore import *

img = QImage(QSize(1200, 1000), QImage.Format_ARGB32_Premultiplied)
color = QColor(250, 250, 250)
img.fill(color.rgb())
p = QPainter()
p.begin(img)

p.setRenderHint(QPainter.Antialiasing)
render = QgsMapRenderer()

layers = iface.mapCanvas().layers()
print layers

layers_ids = []
for layer in layers:
    print layer.id()
    layers_ids.append(layer.id())

render.setLayerSet(layers_ids)
rect = QgsRectangle((iface.mapCanvas().extent().xMinimum()),(iface.mapCanvas().extent().yMinimum()),(iface.mapCanvas().extent().xMaximum()),(iface.mapCanvas().extent().yMaximum()))
rect.scale(1.1)
render.setExtent(rect)

render.setOutputSize(img.size(), img.logicalDpiX())
render.render(p)

p.end()

img.save("map2.png","png")

Is there anyone who can help me with this problem?


texts: Transfer a cartodb legacy map to carto builder and keep the same url publish link

texts:
I have several cartodb legacy map layers that are already embedded in a separate website, is there a way to update the map to the new carto builder without changing the embed link?

Essentially, take the an old map in the old interface and update it to the new interface and functions without changing the link address.

I have tried adding the same dataset to a new map but that creates a new link and I have tried exporting the map and uploading it in a new map, which gives me the new functionality, but also gives a new link to the map.


texts:
I created a TIN raster for an height values. The color ramp is:

but i want to delete all the zero numbers and stay with only integers like this( i do it manually in the layer properties):

How can i change the values in automatic way?


texts: Convert shp 2 sqlite with correct character encoding

texts:
I want to convert Shape data with the QGIS browser in the SQLite database . Unfortunately, the umlauts then appear with a mark (see screen shot ) . I do not want to change the UTF-8 encoding for each Layer .
Can someone help me ?


texts:
How do I connect my WFS service to CARTO - dashboard/datasets?
http://services.cuzk.cz/wfs/inspire-cp-wfs.asp?


texts: ArcMap not recognizing first column of CSV?

texts:
ArcMap is not recognizing the first column of my csv. Field name "AppID" and the values are all in the format of wxyz-#### (letters-numbers). I've checked the files for erroneous carriage returns. If I duplicate the first column in my csv before bringing it in to ArcMap, it comes in the way I wanted it to originally, but still doesn't recognize the very first column.

I should clarify that when I say the field is not recognized I don't mean it doesn't show up in arcmap, only that when the table is opened the second column is the first shown column of the attribute table and you have to manually scroll left to see the first column. I am able to perform joins with all other fields in the table, except for the one that is not recognized. When I join using that field, nothing is returned in the attribute table after the join. Exporting to dbf worked, but this is not an ideal solution since it makes my workflow longer


texts: Best practice for opening Sentinel 2 data with Python

texts:
I'm just trying to learn how to process Sentinel 2 data with Python. All of the tutorials seem to use the data in GeoTIFF format. So I opened the data with SNAP and tried to export as a GeoTIFF but it says that I would need to resample to do that due to the differing pixel sizes. I don't think there should be too much of a problem doing that as long as I'm resampling down the the smallest pixel size (10m) but I'm not sure I'm going about this the right way.
What is the best practice for taking raw Sentinel (2) data and opening it with Python/rasterio?
Edit: Level 2A, in case that matters


texts: Why does QgsSpatialIndex.nearestNeighbor(point, 1) sometimes return 2 lines?

texts:
I'm trying to get the closest line to a point using PyQGIS in a plugin. The active code file is this: http://pastebin.com/TDzmsKVM

When I request the closest line (only 1)

nearestIds = self.spIndex.nearestNeighbor(self.feat2.geometry().centroid().asPoint(),1)

I sometimes get a list with two lines and the first one is not the closest but as far as I can see the second one is. Why is this?

Does QgsSpatialIndex use only boundingboxes for lines?
How do I get the real closest line always, not only most of the time?


texts: gdal_merge with python for PDS4 does not generate xml

texts:
When I use gdal_merge to merge two LOLA PDS3 DEMs* into a PDS4 DEM, I do not get a PDS4 xml file. It creates the img file just fine but it does not create the xml file. I am using Python 3.7.6 in Windows with the gdal tool installed from osgeo4w. I have also attempted to installed from QGIS. This is what I typed in the command line:
python gdal_merge.py -of PDS4 -o merged.img SLDEM2015_256_0N_60N_240_360_FLOAT.LBL SLDEM2015_256_60S_0S_240_360_FLOAT.LBL

and the output is:
0...10...20...30...40...50...60...70...80...90...100 - done.
ERROR 1: Cannot find pds4_template.xml and TEMPLATE creation option not specified

I know that I will need an xml file to open this merged DEM; I have tried opening the img file directly but without metadata it doesn't work. Is there a way to resolve this issue so gdal_merge.py can generate the xml file automatically?
*Data files and associated PDS3 metadata I am using are downloadable here:
http://imbrium.mit.edu/EXTRAS/SLDEM2015/TILES/FLOAT_IMG


texts: Why does the python editor "Idle" reference two Python installations?

texts:
I have installed python 2.6 with ArcGIS 10.0 but also Python 2.7 for use with gdal.

I have a problem at the moment with the editor "IDLE". When I right-click on a .py script, the idle editor references python 2.7 as shown on the first line here:

  Python 2.7.6 (default, Nov 10 2013, 19:24:24) [MSC v.1500 64 bit (AMD64)] on win32
  Type "copyright", "credits" or "license()" for more information.

Subsequently, all scripts fail with the following error:

  ImportError: DLL load failed: %1 ist keine zulässige Win32-Anwendung.

However, if I type "idle" in my Windows 7 search field, and start idle independently without a script, the editor starts in this mode:

  Python 2.6.5 (r265:79096, Mar 19 2010, 21:48:26) [MSC v.1500 32 bit (Intel)] on win32

So how can I make idle reference Python 2.6 when right-clicking on the script?


texts: Create automatic path layer for deriving digital elevation model using GoogleEarth and QGIS

texts:
I recently found this post to create a digital elevation model (DEM) with Google Earth and QGIS. First it requires to manually create a path for sampling elevation data:
http://monde-geospatial.com/create-contours-and-dem-using-google-earth-and-qgis-2-10/
While this works fine for smaller study regions, it its much of work for a large area to get enough sampling points, which are regularly distributed over the whole area (e.g. a point every 30m). Is there a way to create such a layer automatically and regularly (e.g. by predefining point distance)?


texts: Feature is visible in a layer but geometry is None in QGIS using Python

texts:
I am using the feature's geometry in a layer to get the centroid point of the polygon. One particular feature is visible in a TOC layer and also able to select a feature. But, while reading that feature's geometry it shows "NONE". And also, while checking the nodes using "Node Vertex", it shows all the vertex of that polygon in vertex editor table. Kindly suggest us to get the geometry of that feature. Here I have attached a sample code and a error screenshot.
layer=self.iface.activeLayer()
for ft in layer.selectedFeatures():
    print ft.geometry().centroid().asPoint()


texts: Precision when working with spatial data in GIS

texts:
Is there a basic or introductory study that examines and compares the precision of spatial data when working

with varying precison of the data input, like having 1, 2, ... decimal places?
with the varying implementations of floating points (float, double)?
with data near the equator in comparison to data near the poles?
with geographical distances computed with tunnel distance, great circle distance, vincenty, bowring, lambert?

All work I found so far stated that these are sources of error but don't give  exact error boundaries one can expect.


texts: Symbolizing 3d polyline by z-value using ArcMap

texts:
I have a polyline with z-values.  I need to symbolize the polyline by z-value.  There doesn't appear to be any way to do this.
How can I do this using ArcMap 10.2?
A color ramp symbology is the sort of thing I'm looking for. The individual line segments are too long of effectively symbolize by some summary statistic like Zmin/max/avg/etc. I need the symbology to be continuous or nearly so.


texts:
I've got a CSV file with 29 rows and 2146 columns including coordinates and a number of reflectance for each wavelength starting from 350 up to 2500.
I imported the CSV file in QGIS without any problems. But when I want to write my data as a new shapefile, I got the following error:

  (OGR-Fehler: Can't create field W2392 in Shape DBF file, reason unknown.)

Is there field number limitation of attribute table for exporting it as new shapefile using QGIS?


texts: Get the measure M along a CircularStringM in PostGIS without first converting to a LineStringM

texts:
I've been struggling to get PostGIS 2.5 to run any linear referencing functions on geometry containing curves (CircularStringM, CompoundCurveM for example).  Is it possible to get the measure M along a CircularStringM in PostgreSQL/PostGIS without first converting it to a LineStringM?  LineStringM isn't accurate enough to represent the curve for the specific application I'm working on.

Looking to get the closest PointM along a CompoundCurveM or CircularStringM from a user-added point some distance away.  Too bad ST_ClosestPoint doesn't support PointM and Curves.

Willing to settle for a python solution that works nicely with PostGIS.


texts: Adding vector layer in QGIS with standalone PyCharm script?

texts:
I've been struggling with a standalone script in PyCharm, trying to add a vector layer. I've set up PyCharm with Python3.7 to work with QGIS3 using: http://spatialgalaxy.net/2014/10/09/a-quick-guide-to-getting-started-with-pyqgis-on-windows/#ide-example. I've tried tweaking paths in the bat file and script itself based on various posts to no avail.

At this point I am at a loss. I am new to both QGIS and Python.

Below I have provided all code and output. For reference, I am on a Windows 8 laptop.

Script:

import sys
import os

sys.path.append("C:/OSGeo4W64/apps/qgis/python/")
sys.path.append("C:/OSGeo4W64/apps/qgis/python/plugins")

# set path because apparently bat file didn't?
os.environ['QT_QPA_PLATFORM_PLUGIN_PATH'] = 'C:/OSGeo4W64/apps/Qt5/plugins'

from qgis.core import *
from qgis.gui import *
from qgis.utils import *
from PyQt5.QtCore import *
from PyQt5.QtGui import *

print("Import successful")

# Initialize QGIS Application
qgs = QgsApplication([], False)
QgsApplication.setPrefixPath("C:/OSGeo4W64/apps/qgis", True)
QgsApplication.initQgis()

print("QGIS initialized")

pathto = "C:/Users/USER/Documents/Research/bufftest.shp"
layer = QgsVectorLayer(pathto, "bufftest", "ogr")
if not layer.isValid():
    print("Layer fail!")

# Valid file?
print("File Path points to a file: ", os.path.isfile(pathto))

# Valid layer?
print("Output of layer.isValid(): ", layer.isValid())

qgs.exitQgis()

Bat file:

@echo off
SET OSGEO4W_ROOT=C:\OSGeo4W64
call "%OSGEO4W_ROOT%"\bin\o4w_env.bat
call "%OSGEO4W_ROOT%"\bin\qt5_env.bat
call "%OSGEO4W_ROOT%"\bin\py3_env.bat
@echo off
path %PATH%;%OSGEO4W_ROOT%\apps\qgis\bin
path %PATH%;C:\OSGeo4W64\apps\Qt5\bin
path %PATH%;C:\OSGeo4W64\apps\Python37\Scripts

set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\qgis\python
set PYTHONHOME=%OSGEO4W_ROOT%\apps\Python37
set QT_PLUGIN_PATH=%OSGEO4W_ROOT%\apps\qgis\qtplugins;%OSGEO4W_ROOT%\apps\qt5\plugins

start "PyCharm aware of Quantum GIS" /B "C:\Program Files\JetBrains\PyCharm Community Edition 2018.3.1\bin\pycharm64.exe"

Debug output:

  Import successful
  QGIS initialized
  Layer fail!
  File Path points to a file:  True
  Output of layer.isValid():  False
  sys:1: ResourceWarning: unclosed family=AddressFamily.AF_INET, type=SocketKind.SOCK_STREAM, proto=0, laddr=('127.0.0.1', 57568), raddr=('127.0.0.1', 57566)>

  Process finished with exit code 0


texts: Getting inset map selected and magnified in ArcPy?

texts:
I want to create a mapbook(PDF) of all states in USA, and also create an inset map by using a second data frame. I got my mapbook to run, but then I can't get the inset map in the second data frame selected and zoomed in.  can i use a "select by attribute" loop? or any other way?

# Mapbook

#Imports the operating system module
import os

#Sets the environment
import arcpy

#Allows us to re-run the code multiple times, overwriting previous material
arcpy.env.overwriteOutput=True

#Defines the folder path
folder_path=r"E:\Geog173Programming\Final Project\LabData\tracts (1)"

#Allows to access the arcmap document that I created and saved
mxd=arcpy.mapping.MapDocument(folder_path+"\Mapbook.mxd")

#Allows it to access the data frame
df=arcpy.mapping.ListDataFrames(mxd,"Map")[0]

#Creates a final PDF document and titles it
finalPDF_frame=folder_path+"\Mapbook.pdf"
final_PDF=arcpy.mapping.PDFDocumentCreate(finalPDF_frame)

#Defines the temporary pdf file
tmpPDF=folder_path+r'\tmp.pdf'

#Reads that we it will put information in the new empty text box that
#I created
titlePageinfo=arcpy.mapping.ListLayoutElements(mxd,"TEXT_ELEMENT")[0]

#Gives the cover page a title, and the authors name
titlePageinfo.text="Population Weighted Centroid of Each State" + "\r" + "By: XXX"

#Assigns it to my temporary pdf
arcpy.mapping.ExportToPDF(mxd,folder_path+ r'\tmp.pdf')
#And then appends it to the final pdf
final_PDF.appendPages(folder_path+r'\tmp.pdf')

#Selecting the north arrow,lengend, and scale bar
#and issues a postion of x and y coordinates
North_Arrow=arcpy.mapping.ListLayoutElements(mxd,"MAPSURROUND_ELEMENT","North Arrow")[0]
North_Arrow=elementPositionX=10.1713
North_Arrow=elementPositionY=7.8034

Legend=arcpy.mapping.ListLayoutElements(mxd,"LEGEND_ELEMENT","Legend")[0]
Legend=elementPositionX=0.5284
Legend=elementPostitionY=0.5784

Scale=arcpy.mapping.ListLayoutElements(mxd,"MAPSURROUND_ELEMENT","Alternating Scale Bar")[0]
Scale=elementPositionX=7.3671
Scale=elementPositionY=0.5532

#Labels the state
lyrlist=arcpy.mapping.ListLayers(mxd,"State",df)
for lyr in lyrlist:
    if lyr.name=="State":
        lyr.showLabels=True

#Defines variable. Searches in the states shp.
inputRows=arcpy.SearchCursor(folder_path+'/states.shp')

#Defines the tile variable. Reads the same text box that
#created earlier
titleinfo=arcpy.mapping.ListLayoutElements(mxd,"TEXT_ELEMENT")[0]

df2=arcpy.mapping.ListDataFrames(mxd,"Inset")[1]

df2.elementPositionX=8.6642
df2.elementPositionY=0.4446
df2.elementHeight=2.0162
df2.elementWidth=2.8323

#lyrlist2=arcpy.mapping.ListLayers(mxd,"states",df2)
#for lyr in lyrlist2:
    #field="STATE_NAME"
    #cursor=arcpy.SearchCursor(folder_path+'/states.shp')
    #row=cursor.next()
    #while row:
        #feature2=(row.getValue(field))
        #arcpy.SelectLayerByAttribute_management("states","NEW_SELECTION","STATE_NAME='Montana'")
        #row=cursor.next()
        #print feature2

#Searches through each row in the state shp.Prints the state name and state FIPS
#for each page
for row in inputRows:
    titleinfo.text="State Name:"+str(row.STATE_NAME)+"\n"+"State FIPS:"+\
                    str(row.STATE_FIPS)
    #Positions the the title on each page at a particular x and y coordinate
    titleinfo.elementPositionX=9.9883
    titleinfo.elementPositionY=8.2111
    #Pulling out descripitve information and getting value from the state shp.
    feature=row.getValue(arcpy.Describe(folder_path+"\states.shp").ShapeFieldName)
    #Extent of that particualr feature. Zooms in ot the lake.
    df.extent=feature.extent
    #Finally appends it to the temporary pdf, which then appends to the final pdf
    arcpy.mapping.ExportToPDF(mxd,folder_path+r'\tmp.pdf')
    final_PDF.appendPages(folder_path+r'\tmp.pdf')

    #for row in inputRows:
        #arcpy.mapping.ExportToPDF(mxd,folder_path+r'\tmp.pdf')
        #final_PDF.appendPages(folder_path+r'\tmp.pdf')

    #Saves and closed the final PDF document
final_PDF.saveAndClose()

print "Done"


texts:
I am using google earth engine toolbox to estimate land surface temp for landsat 8 imagery. I have been able to run all the necessary calculation without error but when I add the band as a layer all that is rendered is a solid blue box. I am having trouble determining if it is something that I am doing or an issue with the library itself. I have looked through the code and the calculations seem to be correct.

var geet = require('users/elacerda/geet:geet');

var image2 = ee.Image(
  l82.filterBounds(point)
    .filterDate('2018-07-01', '2018-7-31')
    .sort('CLOUD_COVER')
    .first()
);
var new_toa_radiance = geet.toa_radiance(image2, 10);
var brightness_temp_img = geet.brightness_temp_l8c(new_toa_radiance, true);
var l8_ndvi = geet.ndvi_l8(brightness_temp_img);
var img_pv = geet.prop_veg(l8_ndvi);
var lse = geet.surface_emissivity(img_pv);
var surfTemp_img = geet.surface_temperature_oli(lse);
var im = surfTemp_img.select("LST");
print(im);
Map.addLayer(im, {min: 7000, max: 8000, palette: ['blue', 'green', 'red']});


texts:
Is there a 'friendly' way to create dot-density maps, where a polygon attribute data defines the number of dots (e.g. 1 dot = 100 ) using QGIS 1.7+? I say friendly because after a couple of hours researching all the alternatives seem to require programming. Maybe I am not looking in the right place. I am already aware of the fTools 'Random Points' feature and it is not what I am looking for.

That is precisely what I am seeking. However every time I attempted to run this function. I received the following:

An error has occured while executing Python code:

Traceback (most recent call last):
File "/usr/share/qgis/python/plugins/fTools/tools/doRandPoints.py", line 118, in accept
self.randomize(inLayer, outPath, minimum, design, value)
File "/usr/share/qgis/python/plugins/fTools/tools/doRandPoints.py", line 210, in randomize
else: points = self.loopThruPolygons(inLayer, value, design)
File "/usr/share/qgis/python/plugins/fTools/tools/doRandPoints.py", line 256, in loopThruPolygons
sExt = sGeom.boundingBox()
AttributeError: 'NoneType' object has no attribute 'boundingBox'

Python version:
2.7.1+ (r271:86832, Apr 11 2011, 18:37:12)
[GCC 4.5.2]

QGIS version:
1.8.0-Trunk Trunk, exported


texts:
Within ModelBuilder I've got an orginal file, out of this file a few locations are extracted by mask and the data has been changed afterward. Now I am trying to put this data back by first erasing the mask locations in the orginal file, and then merging this with the changed data.

This is supposed to paste it back into place, and the geometry does...but the data gets created as an empty output. I have no idea what is going wrong!

I've added my model as extra information.

My goal with this is to then delete all the features within the orginal file and append the features from this file. Therefor updating the orginal file without overwriting or deleting it.

Running: Arcgis 10.3.1


texts: Using a proxy for unsecure Mapquest layers in Leaflet L.tileLayer call

texts:
I have a web mapping client that needs to be secure and use only HTTPS requests for layers. I would like to include some basemaps such as Mapquest OSM and OpenSeaMap which do not (as far as I can see) have a HTTPS urls.

I have a proxy setup but I can't see any option with L.tileLayer to apply it to the layer.  Is there a way to call a Leaflet tileLayer through a proxy?

If I use L.esri.tiledMapLayer to call the service I get different errors presumably because these are not ArcServer published services so this is not an option.


texts: How to calculate the area in kernel density?

texts:
I have estimated and reclassify my kernel densities and I would to calculate the area of each value. How could I do it?

My Projected Coordinate System is WGS84 World Mercator. My pixel size is 10000 x 10000. Linear unit: Meter (1).


texts: What should I include in the GetCapabilities of a WFS to make it transactional

texts:
I have created a (reversed engineered) WFS that I can load in QGIS. However I now want to make my WFS transactional.
QGIS however does not display the editing of my current WFS as an option, I guess this is because I need to indicate that my WFS is editable in the GetCapabilities.
I have been searching for examples on how I should indicate this in the GetCapabilities, but I wasn't able to find an example. Does anyone know what I should place in the GetCapabilities XML file to indicate that the WFS supports transactions? Or maybe give an example XML of a WFS that supports transactions?


texts: NA values and Focal Statistics for a Neighborhood of Cells

texts:
I have a bathymetry raster, and I am trying to use it to calculate a parameter called the "contour index". This parameter is calculated using the max and min values from a moving window of cells (let's say, 9x9). It's defined as
"(max depth - min depth)/max depth * 100". Once calculated, the values should range from 0-100.

To try to calculate the contour index, I first used ArcMap-->Spatial Analyst-->Neighbourhood-->Focal Statistics to create a new layer for the min water depth and max water depth in each 9x9 window. I then tried to use the "raster calculator" to calculate the contour index value for each cell.

However, the resulting layer is definitely wrong, with values ranging from -70,000 - 180,000. I honestly am at a loss for what I am doing wrong. Could it be related to NA values in my raster (eg. coastlines and edges of the raster)?


texts: How to make Leaflet maps accessible and WCAG2 compliant

texts:
I have followed several threads on this topic but have not been able to
fix the following issues:

How to prevent Leaflet map tile images names from being read out by screen readers?
How to make vector layers (like Circles and Rectangles) tabbable?
How to make tile images not selectable - is this related to #1?

I have tried with current Leaflet distribution 0.7.5 and 1.0 Beta2 - without any success.  Here are some  related discussions:

Leaflet issue: accessibility improvements

StackOverflow: how to make Leaflet tile layers accessible

Does anyone have a link to share - to a Leaflet based site/page that is WCAG2 compliant?


texts: Returning TopoJSON from server compressed

texts:
I am using TopoJSON Server for my geo Spatial project.
By executing the geo2topo (version 3.0) conversion from GeoJSON to TopoJSON on the server side back to the caller of an endpoint that normally returns GeoJSON.
However, I noticed that some TopoJSON files are the same size as the GeoJSON.
If I should be compressing the responses, then what compress method (gzip?) should I use?


texts: Should Python Add-In installed in map which is packaged be part of *.mpk?

texts:
Using ArcGIS for Desktop 10.1 SP1, I just created a map package (.mpk) from a map (.mxd) which had a Python Add-In (.esriaddin) installed.

When I unpacked the .mpk the .esriaddin was neither installed nor present in the folder that was created in C:\Users\PolyGeo\Documents\ArcGIS\Packages.

Intuitively, I would expect it to be at least present in the package, so am wondering if it that capability has been omitted by accident or design?


texts:
I have edited the app_config.yml
dropbox:
    app_key:              'myappkey'
    app_secret:           'myappsecret'
    callback_url:   'https://urlofmyinstance/user/admin/api/v1/imports/service/dropbox/oauth_callback/'

but in config there is also:
dropbox_api_key: "" but i cannot find on dropbox site this key!
how to have a oauth_callback for all user and not only for 1 user?
i restarted the services of carto but i still cannot see the dropbox feature


texts: Distance to nearest hub using raster data QGIS

texts:
I'm trying to figure out if in QGIS 3.10, there is a way to obtain the distance of each raster pixel from a line feature (a vector)?

I'm quite familiar with the "distance from nearest hub" using vector data but I'm a bit lost when it comes to raster data. The idea would be to obtain the distance of each pixel from a line representing the border of Vietnam. Coordinate set at UTM48N for Vietnam.


texts: Mapping methodology and/or clustering technique needed

texts:
I have a scenario in which I have two service centers (A & B) in an area. The capacity of these service centers are different. Service Center A can serve twice as many as Service Center B.

They both are servicing in the same area. The issue is I want exclusive servicing. Meaning the locations Service Center A services cannot be the same as those Service Center B services. Priority would be distance from centers, however because of capacity restraints I understand Service Center A may have to go to locations technically closer to Service Center B.

Is there a particular type of analysis I should begin to explore?

I feel like K-Means almost took me there, but those results generated my service centers as seeds and my serviceable locations (observations) in terms of proximity only and clustered from there.


texts: How to copy the contents of an attribute from one vector layer to another one that touches it or intersects it?

texts:
I have 2 layers that intersect, called layer A and layer B.

Layer A has an attribute (X) that I want to copy wholesale (for the whole layer) to layer B, so that when object A1 intersects object B1, attribute X1 from A1 is copied to B1.

I am trying to do this in QGIS.  I used to do this in Mapinfo with an UPDATE query... then I switched to Mac.  I am new to QGIS.


texts: How to reverse the drawing order of markers?

texts:
I saw this example which demonstrates how ordering can be controlled with markers added using styling options (graphicZIndex) applied to a vector layer. Due to the limitations in the combination of event listeners that I can apply using this method, I am forced to use the markers layer instead.

There are 100 points with associated data that I need to present both on the map (as markers) and listed below it (as detailed data). The data is sorted with the most important data appearing first. It makes sense that the corresponding marker is also on top of the heap. But by looping through the array, and calling the addMarkers function while appending data to the bottom of the map, the most important marker would be buried under the heap.

Is there a way to reverse the ordering of the markers without having to loop through the data twice?


texts: measuring a distance on QGIS 2.2. in meters not degrees

texts:
Im using QGIS 2.2.0 and when I want to measure the distance on a map, the units are in degrees instead of meters. I have tried this :

Right click on layer -> "Save as..." -> CRS: Selected CRS (WGS84)

but nothing has happened.

Does anyone have any ideas?


texts: How to find closest point on multi line String from a point in PostGIS

texts:
I want to find out the nearest point on my road network which is a multiline string from a point with lat lon value. I have tried the following codes but it giving a error message that my road network is not a line.

SELECT ST_AsText(ST_Line_Interpolate_Point(foo.the_line, ST_Line_Locate_Point(foo.the_line, ST_GeomFromText('POINT(90.402 23.752)'))))
FROM (SELECT ST_GeomFromText(geom) As the_line from roads) As foo;

Where geom is my geometry column with multiLineString
roads is the table name
Long 90.402 and Lat 23.752
Projected Coordinate system UTM_Zone_46N


texts:
I'm trying to upgrade my Openlayers 2.12 map to Openlayers 3 to take advantage of the fantastic transition effects on transparent layers (Something OL2 and Leaflet can't do attractively). This and I don't want to lag behind upgrading my sites when OL3 is officially released.
In my current site (OL2.12) I use check boxes in a simple HTML menu to toggle layer visibility. I push each layer to an array (I hope I'm right in thinking OL3 now automatically generates an array for the layers called 'layers') and each checkbox calls this function (checkboxes are given a value that is representative of their layers array number):

function layerswitch(evt){
    layers[evt.value].setVisibility(evt.checked);
}

In OL3 this no longer works, and I can't find any examples or documentation that details how to set layer visibility.


texts: Create furthest point in multilineString

texts:
I am very new to PostGIS. I have a table with sections (multilineString) drawed on a map and a table with observations (point).

Users can make observations and I have to automatically place the point on the section.

My first try was to use ST_Centroid to place the observations in the middle of the section and it technically works.

Unfortunately, it is not usable later. I print all observations on a map, and they all have same positions, so the top observation only is visible and clickable.

I would like to create the furthest new point from existings points that is still in the multilineString (So it will naturally distribute on the multilineString)

So for exemple:

1st point         : Center
2nd and 3rd point : begin or end of
multilineString (0% or 100% lenght of the section)
4th and 5th      : 25% or 75% length of the section
etc.

I don't need to find the 1st observation, only the followings.

I saw a lot of resources saying that we should use ST_Distance DESC but it's not working as my points are new and not already in the database.

Could someone help me, please?


texts: How do I install ArcGIS Developer help in Visual Studio Express 2008?

texts: I would like to install ArcGIS 10 Developer help in VS Express (2008) so that it works with the F1 key. Is it possible? How?


texts: Using Update Cursor gives IndexError: list assignment index out of range?

texts:
I am trying to use an update cursor and receiving an error. I'm not very familiar with using cursors, so please bear with me. I am trying to write a 1 or 0 to a field, based on a condition.
The error I am receiving is:

  IndexError: list assignment index out of range.

I have googled around and have seen solutions, but I'm still not understanding why its throwing an error.

#Define Environments
rasters = r"C:\Data.gdb"
points = r"C:\Data.gdb\Points"

#Set workspace
arcpy.env.workspace = rasters
arcpy.env.overwriteOutput = True
# Get a list of the rasters in the workspace
rasters = arcpy.ListRasters()
# Loop through the list of rasters
for inRaster in rasters:
    print str("Processing " + inRaster)
    rastername = str(inRaster)
    fieldname = str(inRaster + "_perf_flag")
    # Add surface information to Well FC
    arcpy.AddSurfaceInformation_3d(points,inRaster , "Z", "LINEAR")
    print "Formation Depth Added"
    # Create positive values from SS values
    arcpy.CalculateField_management(points, "Z", '!Z!*-1', "PYTHON")
    print "Depth Values Made Positive"
    # Add new field name for perf flag
    arcpy.AddField_management(points, fieldname, "LONG", 9, "", "", fieldname, "NULLABLE", "NON_REQUIRED")
    print "Flag Field Added"
    # Raname Z field as formation grid name
    arcpy.AlterField_management(points, 'Z', rastername, rastername)
    print "Field Altered"
    # Calculate if well was perfed in the formation - Flag 1 as yes - Flag 0 as no
    with arcpy.da.UpdateCursor(points, ["Base", "Top",rastername, fieldname]) as cursor:
        for row in cursor:
            if (row[3] >= row[0] and row[3] <= row[1]):
                row[4] = 1
            if not (row[3] >= row[0] and row[3] <= row[1]):
                row[4] = 0
                cursor.updateRow(row)
        print "Processing perf flags for " + rastername


texts: Clipping: Image Offset when clipping raster in both Python and QGIS

texts:
I am trying to clip two raster with the same co-ordinates using Python. But even though they are from the same co-ordinates they are properly not align with each other. Its the same problem with QGIS 3.16.6 and also I tried with QGIS 2.8.7. The cell size for the image is 10,10 and 0.04,0.04 for the other one.Is it because of the mismatch of the cell size or something I am doing wrong. I can't understand what causing it how to fix.

import rasterio
from rasterio.plot import show
from rasterio.plot import show_hist
from rasterio.mask import mask
from shapely import geometry
from shapely.geometry import box
import geopandas as gpd
from fiona.crs import from_epsg
import pycrs
import pandas as pd
import json
from shapely import speedups
speedups.disable()
df = pd.read_csv("coordinates.txt")

fp = r"resample_UAV.tif"
data = rasterio.open(fp)

i = 1
while True:
    minx = df.iloc[i][0]
    miny = df.iloc[i][1]
    maxx = df.iloc[i+1][0]
    maxy = df.iloc[i+1][1]
    bbox = box(minx, miny, maxx, maxy)
    geo = gpd.GeoDataFrame({'geometry': bbox}, index=[0], crs=from_epsg(32614))
    geo = geo.to_crs(crs=data.crs.data)
    def getFeatures(gdf):
        return [json.loads(gdf.to_json())['features'][0]['geometry']]
    coords = getFeatures(geo)
    out_img, out_transform = mask(dataset=data, shapes=coords, crop=True)
    out_meta = data.meta.copy()
    crs = from_epsg(32614)
    out_meta.update({"driver": "GTiff",
                    "height": out_img.shape[1],
                    "width": out_img.shape[2],
                    "transform": out_transform,
                    "crs": crs}
                            )
    out_tif = r"{}.tif".format(i)

    with rasterio.open(out_tif, "w", **out_meta) as dest:
            dest.write(out_img)
    i += 1
    if i > (len(df)-2):
        break

I also tired with GDAl_translate but no improvement
Geotransform
For the images with cell size 0.04
Dimensions
X: 3454 Y: 3454 Bands: 4
Origin
690789,1.8704e+06
Pixel Size
0.04343,-0.04343
No Data Value
0
Data Type
Byte - Eight bit unsigned integer
Pyramid overviews
Layer Spatial Reference System
+proj=utm +zone=14 +datum=WGS84 +units=m +no_defs
Layer Extent (layer original source projection)
690788.9751400001114234,1870251.9970800001174212 :
690938.9823600001400337,1870402.0043000001460314

for the image with cell size 10
Dimensions
X: 15 Y: 15 Bands: 4
Origin
690785,1.87041e+06
Pixel Size
10,-10
No Data Value
0
Data Type
Int16 - Sixteen bit signed integer
Pyramid overviews
Layer Spatial Reference System
+proj=utm +zone=14 +datum=WGS84 +units=m +no_defs
Layer Extent (layer original source projection)
690784.9840000000549480,1870258.6740900001022965 : 690934.9840000000549480,1870408.6740900001022965


texts: How to use model-builder iterator with a custom bit of arbitrary code?

texts:
How do I integrate a in-line variable from modelbuilder into a python script? Specifically, I'm using an iterator, and I'm trying to include the 'name' variable in a python script. I believe the problem is on step four, when I add it to the model as a a script (as shown here).

Script code I'm trying to include looks like this:

pop_rank = 1
rows = arcpy.UpdateCursor ("path_to_dataset", "", "", "", "pop A")
for row in rows:
    row.setValue("RANK", pop_rank)
    pop_rank += 1
    rows.updateRow(row)

Specifically, I'd like to make the "path_to_dataset" into the 'Name' variable from the iterator.

I can pull the path using the 'parse path' model builder only tool, but I'm not sure how to feed that into the Python Script.


texts: Creating shapefile from Excel file in QGIS?

texts:
I am using QGIS 2.6.1 and have worked in an Excel file database (which was actually converted from a Shapefile). I have managed to convert it as CSV file and load it into QGIS (displays data fine) but I have trouble to convert it to a shapefile.

I can view it but I can´t add or alter data in this file. I think I have to save it as a database.

Any clues how to convert the xlsl-data to a database and finally to a shapefile?


texts: Transforming CRS of Spatial Polygon Data Frame?

texts:
The code below shows the outline of Canada obtained from GDAM as a Spatial Polygons Dataframe with a CRS of EPSG:4326. I wish to transform this to EPSG:3857 (Mercator). Yet, I obtain the following error.

Why?

Error in UseMethod("st_transform") :
  no applicable method for 'st_transform' applied to an object of class "c('SpatialPolygonsDataFrame', 'SpatialPolygons', 'Spatial', 'SpatialVector')"

==========

library(mapdata)
library(sf)

# Get Canada as spatial polygons data frame
# GADM data is +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
# CRS("+init=epsg:4326")
canada <- getData("GADM", country = "CAN", level = 1)

# Transform to Mercator projection *** ERROR ***
canada_t <- st_transform(canada, "+init=epsg:3857")


texts: How can I achieve Geocoding with the JavaScript OpenLayers API, in the same fashion as the Google Maps API?


texts:
Could somebody help me?
How to open InfoTemplate programmatically???
Using of FeatureLayer.selectFeatures(...) not helps.
Here is example:
http://jsfiddle.net/gfV2K/
When I click on the map then single feature selected and InfoTemplate is opened.
But when I click some feature in the grid then single feature selected and InfoTemplate isn't opened!
So how to open InfoTemplate in this case? May be there is a way to simulate this behavior (for example by calling map.click() event)?


texts: What is needed to develop ArcGIS application in C#?

texts:
I want to develop a piece of software that uses the ArcGIS Engine. This application needs to import user data (GeoTIFF, etc) and allow users to put points and polygons on a map. I would think this would be VERY straightforward.

I am unable to find sufficient detail on how to begin this process and was hoping someone here could help me.

There are APIs, there is the ArcGIS Engine, there is ArcObjects, also SDKS. I don't know which I need to develop my application. I would love to make a Windows Form application in C#.

What would I need?

Do I need ArcGIS Engine?

I am working on computers that have ArcGIS desktop installed.

How can I know what else is installed (Engine, Objects, etc...)?

How can I get ArcEngine if that is what I need?


texts:
I am trying to download protected area data for the Democratic Republic of the Congo (ISO3 code = COD) using the wdpar R package.
library(wdpar)
library(dplyr)
library(ggmap)
library(sf)

cod_raw_pa_data <- wdpa_fetch("COD")

I keep getting the following error message.
Error in wdpa_data[[1]] : subscript out of bounds

How do I fix this?


texts: Check if two layers are identical in Google Earth Engine - Is a visual inspection useful?

texts:
I made a copy of an Image that was shared with me as an asset in Google Earth Engine.
I expect the new copy to be identical to the image that was shared. To check that nothing was changed during export, I made a new image that is the difference between my two images, and I expect the result to be an image with uniform zero value.
var original = ee.Image("..."); //shared layer
var copy = ee.Image("...");

var difference = original.subtract(copy)

However, if I add difference to the map, it shows a pattern if I take a look from far away:

Inspecting the value of pixels shows values other than zero.
However, the pattern fades away if I zoom in, and I only find zeros if I inspect pixel values. I've added some transparency to the layer in the capture, but you can see that it's a plain grey.

Still, if I print the properties of my difference layer, it shows -255 and 255 as min and max values.
I think this is something about how GEE works with pyramids or pixel representation at different zoom scales.
Indeed, the original and copy layers also look a bit different from the distance, but they look the same from close distance:
original

copy

I'll zoom in to the island where the arrow points
original

copy

Can I trust that the two layers are identical in their values and their difference is zero for every pixel?


texts: How to avoid random transparent pixels in georeferenced grayscale raster images?

texts:
I'm in the process of georeferencing some old black and white aerial pictures from 1990 in QGIS 2.4. So far, this is business as usual, no problems.

Of course I get the usual black triangles, because I'm reprojecting them to EPSG:25833. So I checked "use 0 for transparency when needed". But now my georeferenced pictures are spotted with occasional transparent pixels, because of course there are legitimate black pixels all over my old aerial pictures.

How do I avoid this? Is there any workaround available? Maybe change the no-data color value to some better value? And most important, how can I do this for ~50 files as economically as possible?


texts: arcgisscripting.ExecuteError: ERROR 001270: Consolidating the data failed

texts:
I am working on PuTTY and trying to create service definition file in server. My Python script is creating .sddraft but I am receiving 'Consolidating the data failed' error when I run stage service script.
I already checked the following link and fixed my issues.
 http://resources.arcgis.com/en/help/main/10.1/index.html#//00vp0000002r001270

import arcpy

wrkspc = '/home/gis/arcgis/server/arcpy/arcpy/idtech/Example/'
mapDoc = arcpy.mapping.MapDocument(wrkspc + 'test.mxd')

service = 'test'
sddraft = wrkspc + service + '.sddraft'
sd = wrkspc + service + '.sd'

analysis = arcpy.mapping.CreateMapSDDraft(mapDoc, sddraft , service, 'MY_HOSTED_SERVICES')
arcpy.StageService_server(sddraft, sd)


texts: Searching/adding temperature values to CSV file with locations and dates in Google Earth Engine

texts:
I am brand new to GEE and JavaScript.
My goal: To add a column in a CSV file with locations and dates with corresponding temperatures
My problem: I am inexperienced with coding and have no history with JavaScript or GEE
I have a very simple CSV file that is structured as follows:

eventDate
Elevation

2015/02/15
300

2013/08/04
420

and the geometry has been correctly read into GEE (I'm able to draw the points as a layer).
My specific issue is this: I want to create a loop or function that will allow me to add another column, "MaxTemp", and find corresponding maximum temperature data (for the moment from GRIDMET, but I honestly do not really care) for the locations/dates of my events.
This is how far I've gotten on my own:
//converting my data table of Mountain Chickadee sightings into a FeatureCollection
var MtnChickadee_points = ee.FeatureCollection(table);

//This was the only way I could limit the data after filterBounds() not working for some reason
var clipToCol = function(image){
  return image.clip(geometry);
}

//my attempt at grabbing the 'max temperature' band in the date range (2001-2016) that I care about
//seems to have worked
var dataset = ee.ImageCollection('IDAHO_EPSCOR/GRIDMET')
                    .select('tmmx')   // select 'max temperature' band only
                    .filterDate('2001-01-01', '2016-12-31')
                    .map(clipToCol);

//This was just for me to confirm I had imported the data table correctly
///Map.addLayer(MtnChickadee_points.draw({color: '006600', pointRadius: 5}), {}, 'Mountain Chickadees in Montana');

//using the visualization parameters from the example code of the dataset
var maximumTemperatureVis = {
  min: 273.15,
  max: 314.0,
  palette: ['d8d8d8', '4addff', '5affa3', 'f2ff89', 'ff725c'],
};

//This is where I'm stuck. Here's some psuedo code for what I'm trying to do
var MtnChickadee_points$MaxTemps = dataset @ MtnChickadee_points(Lat/Long, eventDate);

Now there's all sorts of things that I would need to further do. Likely there are many eventDates that do not have a corresponding temperature data point and so on, so I'd need to make a buffer. But for the moment I am completely stuck at the most simple part of my problem: to add a new column with corresponding temperatures to my events.


texts: How to fix "Error, query must return columns 'id', 'source', 'target' and 'cost'"?

texts:
I have one spatial table with gid, FNODE_, TNODE_ and LENGTH column in my Postgresql/PostGIS database. I am using pgRouting. Problem is, when i try with shortest path query it shows following error. But it should not be. I double check my geographical database there is a connection between 505793 and 496514 nodes. But, whenever i query it shows following error.

SELECT * FROM shortest_path('
                SELECT "gid" as id,
                         "FNODE_"::integer,
                         "TNODE_"::integer,
                         "LENGTH"::double precision as cost
                        FROM roads_dist',
                505793, 496514, false, false);

ERROR:  Error, query must return columns 'id', 'source', 'target' and 'cost'

********** Error **********

ERROR: Error, query must return columns 'id', 'source', 'target' and 'cost'
SQL state: XX000

Thanks, in advance.


texts: Calculating inverse polynomial transforms for pixel sampling when map georeferencing

texts:
In map georeferencing, one typically chooses a set of ground control points that represent pixel locations in the map image along with their real world geographic coordinates, which are then used to estimate a polynomial forward transformation of 1st/affine, 2nd, or 3rd order, for instance using least squares. The forward transform is then used to calculate the bounds of the geographic coordinate system.
Next, one needs the equivalent backwards transformation for determining which pixel coordinate in the original map image to "sample" for each of the pixels in the output warped/resampled map image.
My question is what is the correct or common approach used to calculate the backwards/inverse of polynomial transforms when map georeferencing? How is it done in common GIS software/libraries such as gdalwarp? The forward transform is simple enough using a least squares approach, but I can't figure out how to get the inverse transformation, especially for 2nd and 3d order polynomials.
My solution so far has been to simply estimate a separate backwards transformation using least squares to predict the image coordinates based on the geographic coordinates of the GCPs (instead of the other way around). The problem with this however, is that reestimating a new transformation does not create an exact inverse of the original forward transform, i.e. forward and then backward transforming an image pixel coordinate does not return the same image pixel that I started with, it ends up slightly offset. For the image transformation to be accurate the backwards sampling needs to be the exact inverse of the forward transform, so you can go back and forth without losing information.
Is there another mathematical or programmatic approach to get the exact inverse of such polynomial map transforms, involving for instance matrix algebra or numpy? I'm not looking to outsource the problem to external libraries, I want to implement the nuts and bolts, preferably in numpy.
I've dived under the hood of GDAL and gdalwarp but so far I've been unable to find the source code for calculating the inverse transform, so any pointers here would be useful as well.