texts: Using function defined in function editor of Advanced Python Field Calculator?

texts:
I defined a new function SetPopulation using the function editor.  I have used it in the field calculator, and it works as desired.

I have a model (deliberately very simple for testing purposes) in which I wish to use this new SetPopulation function in the Advanced Python field calculator to calculate the values of a new attribute created by the Advanced Python field calculator.

Why does the the Advanced Python field calculator algorithm in the model not recognize SetPopulation as a legitimate expression?


texts: Options to export a polygon defined in ArcGIS SDK for .NET WPF into a KML file?

texts:
Using ArcGIS SDK for .NET WPF.

I have a ArcGIS SDK for WPF user interface that allows users to create polygons on a map.

What are my options to export the ploygon shape as a KML or SHP file?

I saw examples showing how to read and display KML files. But I am looking to export my plygons to KML (or SHP) files

If the ArcGIS SDk cannot directly export to KML, I am fine with running another tool or SDk to export in memory plygon data into an KML file.


texts: Adding independent panel in ArcMap layout view?

texts:
In ArcMap 9.3, I want to make a map in which I start from a wider location (for example, a country) and then, I give some zooms until one specific area.

I have all the data I need: country administrative boundaries, state/province boundaries and city boundaries (in Brazil and Paraguay).

I already got to make my first panel, but how do I plot an independent "frame/panel" into my ArcMap layout view?

I added a screenshot below, where I drew in red (paintbrush) the region I would like to zoom and plot in a separate panel.


texts: Exporting raster with transparency from QGIS

texts:
As a follow on from my earlier question on reprojecting raster images - is there a way to export a reprojected PNG file (with transparency) from QGIS, with an associated world file?
I can now successfully reproject a TIF, but I'd like to keep the transparencies of the PNG.
I'd also like to have better control over the export resolution, which appears poor from using the Export > Save As > GeoTIFF option, even with 'rendered image' selected.
Edit: I know I could export a world file with transparency using the print composer; but I want to export the raster in a different CRS, so it needs to be a layer export option, I think!


texts: Load OSM data into PostgreSQL using Python script

texts:
I want to load OSM data into PostgreSQL database using Python script. When I am trying this Python script, it doesn't load OSM data into database. Can anyone guide me? I know how to load data in PostgreSQL using osmosis but currently I am looking for something using Python.

import psycopg2
from osgeo import ogr

# connect to the database
connection = psycopg2.connect(user="postgres",
                              password="password",
                              host="localhost",
                              database="example")
# create cursor
cursor = connection.cursor()
cursor.execute("DROP TABLE IF EXISTS trial1")
cursor.execute("CREATE TABLE trial1 (id SERIAL PRIMARY KEY, geom Geometry)")
cursor.execute("CREATE INDEX trial1_index ON trial1 USING GIST(geom)")
print("Successfully created ")
connection.commit()

# define OSMfile path
osm = ogr.Open("pedestrian.osm")
layer = osm.GetLayer(1)

#  delete the existing contents of the table
cursor.execute("DELETE FROM trial1")
print(str(layer.GetFeatureCount()))
for i in range(layer.GetFeatureCount()):
    feature = layer.GetFeature(i)
    # Get feature geometry
    geometry = feature.GetGeometryRef()
    wkt = geometry.ExportToWkt()
    # Insert data into database,
    cursor.execute("INSERT INTO trial1 (geom) VALUES (ST_GeomFromText(" + "'" + wkt + "', 4326))")
    print("Data inserted successfully")
connection.commit()


texts: How can i spatially add text to feature in argis 10?

texts: I have alot of text that write with drawing tools and same feature . i want to find a method to add texts to features spatially like spatial join for feature to feature . is there any method to add these texts?


texts: Converting Excel Sheet to Shapefile without XY Coordinates?

texts:
I am currently trying to build a model that will join two excel tables (without xy coordinates) to one shapefile, which will save some time on a very long process that I have to repeatedly do each day at work. I know that there is a way to manually do it if I right-click the shapefile in the table of contents and add a join from that menu; however, there seems to be no way of doing it in the ModelBuilder.

I am guessing it is a similar situation to the two different merge options in the toolbox and in the editor drop down menu, where they share the same name yet function differently. Therefore, I was wondering if there was a way to make a conversion of an excel table to a shapefile in model builder without coordinates? I know that there is a new tool in ArcGIS 10.2 that would probably help me with this, but the updates are not easily accessed at my work (I am currently using 10.0), so are there other ways of doing this?

This may seem like a simple question, but things are a lot more difficult using model builder due to the fact that I need the process to work automatically each time I use the tool and cannot think of a way that would add xy coordinates and delete them each time I use the tool while not messing with the structure of the original input shapefile. What are my options?


texts: Prevent pixel doubling for a mobile display - OpenLayers

texts:
I'm programming an OpenLayers (4.6.5) map with a layerswitcher in order to switch between my different layers. On mobile, the display of my layer is blurry, OpenLayers is doubling pixels and I don't want that.

I also tested Leaflet (1.3.1) and the same layer with the same parameters (in the URL) has no blurring problems.

In OpenLayers I tested this code :

var hiDPI = ol.has.DEVICE_PIXEL_RATIO > 1;
var layerHdpi = hiDPI ? 'gny_bright_mobile' : 'gny_bright';
var tilePixelRatio = hiDPI ? 2 : 1;

var wms_gnybright_auto = new ol.source.TileWMS({
    url: '...',
    params: {'LAYERS': layerHdpi, 'FORMAT': 'image/png'},
});

I found no way to prevent pixel doubling and the blurry problem.

I tested an another way supposed to fix my problem but it's not :

var map = new ol.Map({
    layers: [new ol.layer.Group({
                title: 'Fonds de carte',
                layers: [
                    new ol.layer.Tile({
                     title: 'gny_bright_auto',
                     source: wms_gnybright_auto,
                     type: 'base',
                     tilePixelRatio: pixelRatio
                    }),
                    //...

Captures :

Leaflet version without blurry problem :

OpenLayers version with blurry problem :

Do you have any idea ?


texts: Measure distance between points along a line in QGIS

texts: I have a road network and different points connected to the roads. I would like to measure the distance between several points along the roads. Basically I want to do the same like the road graph plugin, but I want to measure not only the distance between two points and also would like to get the results in an attribute table.


texts: Trying to start Collector for ArcGIS, can't publish service definition because GIS is on virtual machine with no internet

texts: I'm trying to set up Collector for ArcGIS. We are storing data in databases using SDE. The GIS is stored on a VM with no internet connection. I am trying to figure out how to publish the service definition without internet.. I have saved it as a .SD file and published to our default ArcServer Manager site? localhost:6080/ etc..  It publishes fine on the site but I do not know how to get it to show in Collector app on mobile. It doesn't appear. We have an organizational account and i've created a map but the features are not editable.. i assume it's because they weren't published as feature services. Is there any way around this?


texts: Adding attributes to line feature according to adjacent polygon

texts:
I'm trying to add attributes to a field in a polylines attribute table according to directly adjacent polygons....

In the pic, the purple line has an empty field labelled 'FG' and I'd like to fill it with whatever is in the 'FVG' field in the adjacent polygon that shares the boundary ie. 'Rf'.. 'Mit'...'Con'. The purple polyline is split into segments at each polygon boundary.

This should be quite simple but can't for the life of me figure it out. Is it spatial join I'm after? Haven't had much experience with that


texts: Problem with QGIS 2.2 joins between shapefile and CSV

texts:
Didn't see question like this asked, so...

I've joined a CSV file to a shapefile and everything seems OK.  I've modified attribute form so that fields from joined CSV appears in there.

So far so good- it's there. But when I reopen the project and identify some feature, instead of joined fields there's "Error:Field does not exist in datasource" even though I see it when opening the attribute table.

If I open the layers properties and then press OK, the attributes are back in my form.

Is there any way to avoid this annoyance?


texts:
My immediate geographic area is lacking an organized local GIS Users Group. I belong to a group that meets a couple times a year, but its about 1.5 hour drive. I'm thinking about starting up a Users Group for our County, we have about 50 users just in my organization that I think could benefit, not to mention the other users from different organizations.

Problem is, I've never done anything like this and don't exactly know where to begin. Does anybody have tips or pointers for how to start a group or how to make a GIS group successful?


texts: Sum values of a field based on unique values in another field?

texts:
I have a layer made of a field with repeted string values ("building") and a field with float values ("flatssurface"). I'd like to write a script to do some arithmetic operations (like sum, max and min, etc.) on all the values in column B which are associated to same values in column A but it doesn't work.

Dataset:

BUILDING       FLATSURF    SUM (wanted result)
    1           40.0       115 (i.e. 40+45+30 m^2)
    1           45.0       115
    1           30.0       115
    2          200.0       300 (i.e. 200+100 m^2)
    2          100.0       300
    3           60.0       140 (i.e. 60+80 m^2)
    3           80.0       140

Code trial:

building=layer.dataProvider().fieldNameIndex('BUILDING')
flatssurface=layer.dataProvider().fieldNameIndex('FLATSURF')

uniquevalues=layer.uniqueValues(building,limit=10000)
for uv in uniquevalues:
    feat=layer.getFeatures()
    for f in feat:
      tot=sum([flatssurface])
      print tot

however it returns only the index of the field flats' surface (f[1]).
For simplicity, I neglected the part on adding the field but I just tried to print the results, in order to check if it worked.
If possible, I don't need answers related to GroupStats (Calculating sum of parts of column based on another column using QGIS Field Calculator?) or SQL or statistics (https://gis.stackexchange.com/editing-help) but I am looking only for solution related to Python.


texts: How to call the method of DescribeFeatureType in OWSLib

texts:
I want to execute DescribeFeatureType operation of WFS services. However, I can’t find the method to fulfil this operation in OWSLib. Can anyone kind enough to tell me how to do this?

The following is my code snippet:

from owslib.wms import  WebMapService
from owslib.wfs import  WebFeatureService as wfsServericelayer

wfslayer=wfsServericelayer('http://***:8080/geoserver/test/ows')
for i in wfslayer.operations:
    print i.name

Gives:

GetCapabilities
DescribeFeatureType
GetFeature

But how to execute the DescribeFeatureType method?

operationT=wfslayer.getOperationByName('DescribeFeatureType')


texts: Splitting KMZ file into two files using QGIS?

texts:
I have a KMZ file that contains about 10,500 features.
- My mapping application only allows 10,000 (a google api limit)
- So how do I split the KMZ file into two files?
- Example, file_A = 5,000 and file_B = 5,500

In QGIS v3, I can open the file, then open the Attribute Table:
- I see the 10,500 lines..
- I can select the first 5,000 lines..
- But how do I save these to a new KMZ file?

Or maybe delete the last 5,500 and save to a new file_A..?

Surely there is a way to quickly/easily break a KMZ into two files?


texts: Converting *.adf files to *.xyz for XPSWMM?

texts:
i am trying to convert adf files to xyz files in order to import them to XPSWMM.  i also need the elevation values.  however, when i use QGIS to translate the raster files (adf files) to xyz these errors show up,

  Traceback (most recent call last):   File
  "C:/Users/Mark.Valiente/.qgis2/python/plugins\BLNExporter\BLNExp.py",
  line 208, in EnablePlugin
      if (layer.geometryType() == QGis.Polygon) or (layer.geometryType() == QGis.Line) or (layer.geometryType() == QGis.Point): AttributeError: 'QgsRasterLayer' object has no attribute 'geometryType'

and it says,

  Warning 6: A dataset opened by GDALOpenShared should have the same
  filename (H:\DEM\Tiles_1kmx1km\e492161301001\w001001.adf) and
  description (H:\DEM\Tiles_1kmx1km\e492161301001)

when the conversion has finished only a solid black color was displayed on the work area.

here are the names and formats in every folder of the directory,

dblbnd.adf
hdr.adf
metadata.xml
prj.adf
sta.adf
w001001.adf
w001001x.adf


texts:
is there any means of telling QGIS to use a different Python installation than the default, system-wide one?
E.g. use some virtualenv Python /myvirtualenv/bin/python instead of /usr/bin/python?

Background: I'm developing a custom plugin that uses some modules which are only installed in a certain virtualenv, but not system-wide.

The virtualenv is allowed to access system-wide packages, so the qgis core libs should still be accessible when using that virtualenv.

The only solution that I see right now is to install the required modules in the system-wide Python. However, I'd like to keep them in the virtualenv only, if possible.


texts:
I am trying to figure out how to insert Geometry Array in Geometry [] column in PostGIS.

CREATE TABLE public.array_data_demo
(
id bigint PRIMARY KEY,
txt_arr text[] COLLATE pg_catalog."default",
bigint_arr bigint[],
jsonb_arr jsonb[],
geom_arr geometry[]

)

--------------------------------
INSERT INTO public.array_data_demo
(id, txt_arr, bigint_arr, jsonb_arr, geom_arr)
VALUES
(3, '{txt1,txt2}', '{1,2,5,3}',  (array['{"sender":"pablo","body":"they are
on to us"}'::json, '{"sender":"arthur"}'::json]),
{ST_GeomFromText('POINT(23.3 19.234)',4326)}
);

SELECT * FROM  public.array_data_demo ;

Gives the following error:

ERROR:  syntax error at or near "{"
LINE 6:     {ST_GeomFromText('POINT(23.3 19.234)',4326)}
            ^
SQL state: 42601
Character: 243


texts: Delete observations that fall outside of shapefile

texts:
I have data that locates places by xy coordinate.

I have plotted them on a shapefile.

As you can from the sample map I have uploaded, many of the xy coordinates are incorrect and locate places outside of the shapefile, and therefore need to be removed from my dataset.

Some of the incorrect coordinates have a y coordinate = 0 (the vertical line of points), and are easily removed, by dropping observations where y is = 0.

However, there are many other coordinates (not all represented by this map), that are incorrect but are not equal to 0 and therefore not easily identified.

It would be very timely to identify 10 years of measurement errors individually.

So, I was thinking that an easy way to remove these points is to delete the observations situated outside of the shapefile/island.

I am still an ArcGis beginner, so I am not confident yet on how to manipulate the data in such a way.

Do you know if it is possible, and how to go about it?


texts:
I have just installed the OpenGeoServer Suite from Boundless and when I go to open the dashboard @ http://localhost:8080/dashboard I get a 404 error not found.

edit: just realized that apache is installed on http://localhost:8080/, I have tried to uninstall it but no luck


texts: Cannot populate a Ext GridPanel through GeoExt WMSCapabilitiesStore

texts:
I've been working through a few tutorial/starter map samples on the OpenGeo and GeoExt websites and have come upon a problem I cannot solve. I'm trying to load a list of layers from my GeoServer into a GridPanel using GeoExt.data.WMSCapabilitiesStore. No matter what I try, I cannot get anything to show in the GridPanel. I thought I had cross domain issues at first, but the problem persists even when the webpage and the GeoServer data are on the same server.

Here's a link to my very barren looking web map: http://24.105.210.45/OpenGeoGeoExtMap1.html.

Here's the JavaScript that's supposed to access the wms and load the layer list:

    Ext.onReady(function() {
var grid = new Ext.grid.GridPanel({
    title: "Available WMS Layers",
    region: "north",
    height: 150,
    viewConfig: {forceFit: true},
    store: new GeoExt.data.WMSCapabilitiesStore({
        url: "http://24.105.210.45:8080/geoserver/wms?service=wms&version=1.1.1&request=GetCapabilities",
        autoLoad: true
        }),
        columns: [
            {header: "Name", dataIndex: "name", sortable: true},
            {header: "Title", dataIndex: "title", sortable: true},
            {header: "Abstract", dataIndex: "abstract"}
            ],

Any help with this, or just ideas of where to look next would be greatly appreciated.


texts: Create a Line FeatureClass from IPointCollection

texts: I want to create a new line layer from a point layer. I store the points in an IPointCollection. I want to use that IPointCollection to create a line IFeatureClass. I know how to create a single Feature and how to put it in a layer, but this doesn't work.


texts: Calculating distance between two feature layers in specific direction using ArcMap?

texts:
I'm working on a project in which I'd like to calculate the distance from specific locations within a city (points) to the nearest body of water (polygons) following the path of the wind (direction).

I am using ArcMap 10.5 Desktop and the only related solution I can find is the Near tool in the Proximity toolset.

However, while this tool can output the direction of the nearest feature, it does not allow you to specify the direction in which you would like to calculate distance.

Does anyone know of a tool or workaround that I could use to do this?


texts:
I am getting the error while I am trying to stack individual raster bands into one image.

  ValueError: expected array of dim 2

import os
from osgeo import gdal
import numpy as np
os.chdir(r'C:\osgeopy-data\Landsat\Washington')
band1_fn = 'p047r027_7t20000730_z10_nn101.tif'
band2_fn = 'p047r027_7t20000730_z10_nn201.tif'
band3_fn = 'p047r027_7t20000730_z10_nn301.tif'

in_ds = gdal.Open(band1_fn)

in_band = in_ds.GetRasterBand(1)

gtiff_driver = gdal.GetDriverByName('GTiff')
out_ds = gtiff_driver.Create('nat_color.tif',
      in_band.XSize, in_band.YSize, 3, in_band.DataType)
out_ds.SetProjection(in_ds.GetProjection())
out_ds.SetGeoTransform(in_ds.GetGeoTransform())

in_data = in_band.ReadAsArray()
out_band = out_ds.GetRasterBand(3)

out_band.WriteArray(in_data)
in_ds = gdal.Open(band2_fn)
out_band = out_ds.GetRasterBand(2)

out_band.WriteArray(in_ds.ReadAsArray())

out_ds.GetRasterBand(1).WriteArray(
gdal.Open(band3_fn).ReadAsArray())

out_ds.FlushCache()
for i in range(1, 4):
    out_ds.GetRasterBand(i).ComputeStatistics(False)

out_ds.BuildOverviews('average', [2, 4, 8, 16, 32])
del out_ds


texts:
How is't possible to change the geomarker in another symbol?
For examle 'a bike'

http://openlayers.org/en/latest/examples/feature-move-animation.html


texts:
I'm trying to create a simple polygon in sf and select only points within that polygon. What am I doing wrong here?
library(sf)
library(concaveman)
library(ggplot2)

foo.df <- data.frame("long"=c(136,137,137,136),"lat"=c(36,36,37,37))
foo.sf <- st_as_sf(foo.df, coords = c("long","lat"))
poly <- concaveman(foo.sf) ## in case points are out of order
point.df <- data.frame("long"=c(136.2,136.5,137.5),"lat"=c(36.5,36.5,36.5))
point.sf <- st_as_sf(point.df, coords = c("long","lat"))

good_points <- st_join(point.sf,poly,join=st_within)

The st_join function doesn't seem to do anything
ggplot() +
  geom_sf(data = poly) +
  geom_sf(data= good_points)

The problem isn't with the concaveman package
good_points <- st_join(point.sf,foo.sf,join=st_within)

ggplot() +
  geom_sf(data = poly) +
  geom_sf(data= good_points)

And this attempt to create a polygon throws an error
another_poly <- st_polygon(list(as.matrix(foo.df)))
good_points <- st_join(point.sf,another_poly,join=st_within)

What am I missing?


texts: Combine rasters in PostGIS with algebra operations

texts:
I have 2 raster tables with one band. What I want to do is clip the rasters with a polygon and multiply pixel values with different values and combine those multiplied rasters into a new one. The main goal is to analyze DEM data for the specified areas. I can do those things with separate queries which you can find below. But I was wondering if it is possible to achieve all those things in one query with ST_MapAlgebra
--create tables for clipping rasters

CREATE TABLE schema.t1(
    rid integer primary key,
    rast raster,
    new text
);
CREATE TABLE schema.t2(
    rid integer primary key,
    rast raster,
    new text
);

--clip the data with attribute-based polygon and insert into table

INSERT INTO schema.t1(rast, new, rid)
SELECT (ST_Clip(rast, sf.geom, TRUE)) as new, rast, rid
FROM schema.rast1
JOIN schema.polygon sf ON ST_Intersects(rast1.rast,sf.geom)

WHERE sf.name_1 = 'Van'
;

--clip the data and insert into a different table for other algebra operation

INSERT INTO schema.t2(rast, new, rid)
SELECT (ST_Clip(rast1.rast,sf.geom, TRUE)) as new, rast, rid
FROM schema.rast
JOIN schema.polygon sf ON ST_Intersects(rast1.rast,sf.geom)
WHERE sf.name_1 = 'Van'
;

-- now update the pixels with multiplication

UPDATE schema.t1
    set rast = ST_MapAlgebra(rast, 1,'32BF'::text, '[rast.val] * 3'::text, -1);

UPDATE schema.t2
    set rast = ST_MapAlgebra(rast, 1,'32BF'::text, '[rast.val] * 5'::text, -1);

-- now create another table and insert multiplied tables
CREATE TABLE schema.t3(
  rid SERIAL primary key, rast raster
);

INSERT INTO schema.t3(rast)
SELECT st_mapalgebra(b1.rast,1, b2.rast, 1,'([rast2] + [rast1.val])')
AS rast

FROM schema.t1 AS b1 cross join schema.t2  AS b2;


texts: Error in converting shapefiles into spatstat window

texts:
I wish to import shapefiles of the Seoul Capital Area to use as a window in spatstat. The SCA is the combination of the Seoul, Incheon, and Gyeonggi-do level-1 administrative areas of South Korea.
I download shapefiles for South Korea from GADM (unless there is a better source for Korean provinces). The zip file contains KOR_adm1.shp for level-1 administrative areas of South Korea.
library(rgeos) #required by maptools
library(maptools)
library(spatstat)

sca <- readShapeSpatial("KOR_adm1.shp") #upload Korean provinces
sca <- sca[sca@data$NAME_1 %in% c("Seoul", "Incheon", "Gyeonggi-do"),] #subset required 3 areas
plot(sca) #plot reveals an anomaly in the bottom right corner

sca <- as.ppp(sca) #convert 3 areas to spatstat ppp format returns an error

Error in as.ppp.default(sca) : Can't interpret X as a point pattern

seoul <- as.ppp(sca[sca@data$NAME_1=="Seoul",])  #convert 1 area to spatstat ppp format returns the same error

Error in as.ppp.default(sca[sca@data$NAME_1 == "Seoul", ]) :    Can't
interpret X as a point pattern

My question is how to convert Seoul Capital Area shapefiles into a single window for spatstat and secondarily whether it's possible to preserve the boundaries for the 3 sub-areas within that window. This vignette advises (p2) to "[multiply] the longitude coordinates by the cosine of the latitude of the centre of the region." Is this also recommended in this case?


texts: Why are web maps able to display features from a REST endpoint that is disabled?

texts:
Let's say I have web service X: maps.example.com/servername/rest/services/folder/ServiceName/MapServer/0

If I navigate to X in my browser's search bar, I get the following error:

{
    "Error": "Services Directory has been disabled",
    "Code": 403
}

However, when I add X into an ESRI web map using "Add data from the Web", everything displays fine. According to this ESRI documentation, one can disable a services directory. By my reasoning, this should prevent users from accessing or querying the data in the map or through an HTML page.

Why is it possible to add a service from a disabled directory to a web map? Is it a weird quirk with this particular X service, or does this have broader context?


texts: QGIS Multi-Basin Delineation (Pour-Point) Hydrology

texts:
I'm working on a QGIS plugin and struggling to find some existing QGIS processing functionally I could leverage. I need to replicate multi-basin delineation using pre-defined pour points. Specifically, I would like to delineate basins (or contributing area) for each cell representing a river segment. The ESRI watershed tool was used before for the function, referencing a input flow direction grid, and Pour_Points (a raster with each cell being a unique ID).

previous ArcPy implementation:

gp.Watershed_sa(fdirnd, Pour_Points, wshed_g, "")

Graphic below shows needed output. Blue line is a vector river segment. The blue line is converted to raster and matched cell size to a flow Direction raster. Each cell of the river represents a Pour point, in which a basin is delineated.


texts:
I want to realize "Portable gdal2tiles.py" in both Windows and Mac environment.
"Portable" means, run script without installing something.

To realize this, I believe, we have to have such components:
1. Python runtime which is portable
2. libproj4 which is portable
3. libgdal/libosr which are portable

Are there such components? Is this possible?


texts:
I'm a forestry technician and I've one question dealing with forest road planning... in order to calculate areas that can be cover by a specific machine (with a given work distance) I do as follow:

"proximity" routine to calculate distance from road network;
ratio between proximity output and a raster containing the cosine of slope angle in order to create a "orographic distance" raster,
reclass of the output in served areas and not served areas, for example considering tractor equipped with forest winch served areas are all zones with distance from road till 80 m (work distance), distances greater than 80 m define instead not served areas...

Now I need to define two different work distances, one for uphill and one for downhill extraction... is there a way to do that directly with reclass or other QGIS tools?

At this moment I create two different rasters for uphill and downhill and manually clip them according to slope morphology.


texts: Sequenced path using pgrouting (going through the POI of the user)

texts:
I'm using the pgr_astar function of Postgres and I want to personalize the optimal path proposed according to the preferences of the user.
For example if the pg_astar propose the red path, I want to make this function propose the orange path if I specify that the user wants to go through the blue point.


texts: Is it possible for ArcMap to treat data types (numbers as text) correctly when importing from a text file?

texts:
I am creating a text/csv file (using a python script). The file contains data to import into ArcMap. I have several fields that look like numbers but have leading zeros ("0012") etc. that I want to preserve. I thought I could force this on ArcMap by adding a dummy line at the beginning of the file that has obvious text in it (like an "X").

In the screenshot below, I have an X in several of the fields at the beginning (including ACCTID in the first column) for the file but ArcMap has made them into integers anyway. To clarify, I'm just trying to add the file into ArcMap directly using "Add Data".

Excel seems to treat them similarly. I do not want to open this in Excel and fiddle with the formatting of over 100 fields if I can avoid it and it seems preposterous that I would have to do that. I found similar questions but nothing exactly like this. Some folks had the reverse issue (wanting numbers and getting text) and others were working with Excel files rather than plain text.

The code I'm using to write the CSV (if that matters) looks like this ...

    with open(filehandle, 'wb') as csvfile:
    data_writer = csv.writer(csvfile, delimiter = ",")
    if printHeaders:
        csvfile.write(', '.join(these_headers) + '\n')
    if dummy:
        csvfile.write(dummy + '\n')

further down I have ...

if linecount < linelimit and printline:
    data_writer.writerow(theDict[acct])
    linecount += 1

I am actually running a search cursor through this file to use to run against an update cursor using a key (the ACCTID).


texts: Is it possible to convert GeoJSON to KML format using geotools. I am trying to figure out this using Java but I can get results only using node js.


texts: Display a color ramp legend (of a vector layer) in the qgis map composer

texts: I would like to have a color ramp instead of 10 categories of color when displaying my vector layer in the map composer of QGIS. Is it possible ?


texts: Printing only selected features in report QGIS 3?

texts:
My layer has more than 2900 elements, and I would like to release a map with 1 or 10 or 50 selected elements but never 2900 in the time.

I can't select them other than visually (no column to identify).
Do you have an idea to output a report with only the selected elements?

I have tried the different solutions already seen on this site but nothing that allows me to release a report without saving many unnecessary layers.


texts: NAD27 data distributed in NAD83 by Canadian government

texts:
Natural Resources Canada (NRCAN) produced 1:50000 topo maps for all of Canada according to the Natural Topographic System grid configuration. The maps are now considered a legacy product, but GeoTIFFs are still available via FTP (e.g. here).
Text on the GeoTIFF indicates that the map's datum is NAD27, however gdalinfo reports that GeoTIFF's SRS is based on NAD83. The XML file that ships with the GeoTIFF includes the following:

The source of this raster file (paper map) is geocoded according to
the North American Datum of 1927 (NAD27).  The resulting raster file
is georeferenced to the North American Datum of 1983 (NAD 83)
according to the Universal Transverse Mercator (UTM) projection

If I load the GeoTIFF into QGIS in a new project it will default to using NAD83 UTM 9N as the project SRS. If I click on a grid line intersection that the map says is 0600000,6070000, QGIS tells me its coordinate is roughly 0600000,6070200. If I change the project SRS to NAD27 UTM 9N and click at the same point I get a much more appropriate coordinate.
What is the explanation for NRCAN distributing NAD27 data with a NAD83 SRS, and even acknowledging in the XML documentation that the wrong datum is used in the GeoTIFF? Is this likely to be a mistake, or is there something I'm missing? I know that NAD27 / NAD83 can be a complicated problem, and this is a product from the federal government, so I'm guessing it has to be my mistake.


texts:
I followed this tutorial downloading all the resources
https://docs.geoserver.org/latest/en/user/tutorials/imagemosaic_timeseries/imagemosaic_timeseries.html

I was able to reproduce everything, but I have an issue with colors.
Everything I put inside ImageMosaic become GREY.

So instead of visualizing https://docs.geoserver.org/latest/en/user/_images/step_3.png
with blue colors

I see everything GREY.


texts: TIFF file can only be opened if downloaded through browser

texts:
I'm downloading a TIFF file from https://floodmap.modaps.eosdis.nasa.gov/getTile.php?location=060W030S&day=120&year=2021&product=3
When I download it in my browser I am able to successfully download it and open it. When I use a simple download.file command in R it downloads a file that is the right size but can't be opened.
When I run tiffdump on the two versions I see why, which is that the browser downloaded file has all the correct tags while the R downloaded one only has a subset of the tags and is missing some important ones like StripOffsets.
Here is an example:
download.file("https://floodmap.modaps.eosdis.nasa.gov/Products/060W030S/2021/MSW_2021120_060W030S_A14x3D3OT.tif", "./MSW_2021120_060W030S_A14x3D3OT.tif")


texts: qgis 2.10.1. mac - geotiff singleband pseudocolor

texts:
i've installed mac version from kyngchaos including all required packages ...
... and just noticed the following:

while windows version behaves like expected when trying to render an ordinary tif via "singleband pseudocolor" (see the image)

the same procedure on the current mac-release returns:

interesting detail: pressing "apply" in mac version doesnt change anything and if you press "OK" and re-open layer properties the "render type"-setting has  swiched back to "singleband grey".

even more interesting: "singleband grey" is faulty too:
none of the options "color gradient", "min", "max", "contrast enhancement" does affect the way the tif is displayed,
neither do the "load min/max values" settings !

anyone any ideas?


texts: 'Generate Near Table' and 'XY To Line' Not Working As Expected

texts:
I have two point layers (layer-1 and layer-2).

I used 'Generate Near Table' to identify connections between points in layer-1 and layer-2. I used a search radius of 5000m and maximum number of closest of 10.

I then used the 'XY To Line' tool to create connection lines between the connections identified in the first step.

Issue: some points in layer-1 are connecting to less than 10 points in layer-2 although there are more layer-2 points within the 5000m radius.
In the attached image the three green lines are generated by ArcGIS and orange broken-lines are not generated but are within the 5000m radius.
What could be the reason and how do I correct this?
Software: ArcGIS Pro 2.6.0


texts: Using raster calculator for subtracting landslide mass from raster topo?

texts:
I have a topographic raster data of an area prior to a landslide event. Unfortunately, I don't have a topo-data after the event. So, in order to simulate the landslide event and its run-out distance, I need to subtract the landslide source volume from the main topography raster file.

I have the shape of the landslide source area as a shapefile, and I even have a raster file clipped in the shape of the source file. I have tried to subtract this file from the main file with a uniform depth (let's say 2m) but the result is only on the clipped area.

Is it generally possible to do such calculation (subtraction of selected area) on raster files and, if yes, then how?


texts: QGIS: Export multiple points layers to Excel

texts: I am trying to export multiple points layers to separate Excel files. I could do it manually but since I have 500+ points layers it would take too long. Is there any Python code to solve this?


texts:
I used Torque to animate points over time in CartoDB

I have a map of polygons and I want to animate the changing colors.

Over five years the intensity of the type of incident I'm mapping changes. I can make a slideshow that shows a visualization of each year in an individual choropleth map, but I'd love to be able to animate it in torque instead.

Is this even possible?


texts: Using Google Earth server data with custom GIS app

texts:
I have implemented my own GIS-like app that is able to show tiled data (Bing & Esri tiles for example) and georeferenced images. However, customer wants to use the app in an offline environment with their Google Earth server.

Is this possible?

I tried to find a way to retrieve data from such server, but all the examples are using Google-based  solutions for client side. My application is written from scratch in C++/OpenGL.


texts: How to label vector features programmatically?

texts:
I am using python for QGIS 2.0.1 and I have a vector layer of points representing the positions of state labels. I want to display the 'state abbreviation' attribute on the corresponding point for each feature.

Does anyone have an example using QgsPalLabeling to do this? Any other suggestions are also appreciated.

This is a fairly simple process using the new Labeling tab of the GUI, but I have not found any documentation to make this work programmatically using pyQGIS.


texts: How to write to Sql server from MapBasic and insert rows

texts:
I have created a table in Sql server 2014 with the columns: Id (int), Geom (geometry), RouteId (int), Desc (varchar(MAX))
I have succeeded in making the table mappable through MapInfo
I have checked .TAB file for the table, that it has IsReadOnly = false

Stil when i try to execute an insert query on the table I get the error messsage that the table is read only. Anyone knows how to fix this?

Also I cannot figure out how to write a correct insert query I have tried doing that on a copied version of the table code is as follows

Dim rId As Integer
Dim oPline as Object
rId = 20
Dim t As String
t = "test"
Create pline into variable oPline 2 (1,2)(1,3)

Insert into TestGeom1 (Geom, RouteId, Desc) VALUES (oPline, rId, t)

Writing like above I get the error messsage field Geom does not exist in table TestGeom1
If i replace Geom with obj or just remove it completely i get the error Expression does not evaluate to a column or a table name.
So please who can tell me what is the correct syntax for this?


texts: PyQGIS: Importing own modules results in ModuleNotFoundError or KeyError

texts:
I've installed QGis 3 on Mac using osgeo4mac and wrote a script in the Processing Script Editor. The script works fine but when I try to divide it into several modules, I get an error when importing them.

import module_name
or
from module_name import function_name

results in ModuleNotFoundError.

Relative import using import .module_name results in KeyError: "'__name__' not in globals"

All scripts (main.py, module_name.py etc.) are in the same folder. Any ideas what I'm doing wrong?


texts:
Is there any software that can be used to view point clouds in bpf format? My preference is for free software that has low-learning curve for quick visualizations.

Currently, I convert bpf to las and use a viewer like lasview or FugroViewer, but that's a little inconvenient.

There is a discussion on lidar software at Working with LiDAR data using other than Esri software? but it doesn't cover bpf format specifically. I've looked through some of the packages mentioned for visualization capabilities and bpf support but nothing pops out. A google search of some combination of search terms like "bpf point cloud viewer" hasn't been fruitful either.

Edit: The binary point file (bpf) format specification can be found at https://nsgreg.nga.mil/doc/view?i=4202 or in the source document https://nsgreg.nga.mil/NSGDOC/files/doc/Document/NGA.SIG.0020_1.0_BPF.zip


texts:
I am trying to make maps and histograms for Landsat EVI. My code (extract below and full thing attached) works but the results seem very strange. The histogram outputs are bimodal, despite the exact same operation with Sentinel-2 imagery returning a bell curve. Random patches appear very unhealthy which should not be the case. I am wondering if it is to do with the quality of my images despite being cloud-filtered, but don't know how to check this or change it.  The total values that come out of my histogram are also quite different when I run the full code so perhaps this has something to do with it?
My question therefore is how do I correct these EVI values and what has gone wrong?
link to full code

var forest = BROADLEAF.merge(MIXED).filterBounds(geometry4);

// Load a raw Landsat collection for a single year.
var collection8 = ic2
    .filterDate('2018-06-01', '2018-08-31')
    .filterBounds(geometry4);

var c2018 = ee.Algorithms.Landsat.simpleComposite({
  collection: collection8,
  percentile: 0,
  cloudScoreRange: 0,
  asFloat: true
});

// Display the composites.
Map.addLayer(c2018.clip(geometry4));
var evi2018 = c2018.expression(
  '2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))', {
    'NIR': c2018.select('B5'),
    'RED': c2018.select('B4'),
    'BLUE': c2018.select('B2')
});

var evi2018real = evi2018.clip(forest);

var histogram5 = ui.Chart.image.histogram({
  image: evi2018real,
  region: geometry4,
  scale: 30,
  maxPixels: 1e13})
    .setSeriesNames(['18']);

print(histogram5);
Map.addLayer(evi2018real, VIS, "2018");


texts: Why simplify() returns unbounded polygon? Google Earth Engine

texts:
Why does simplify() return an unbounded geometry when the maxError parameter is above a certain threshold? I have a feature from RESOLVE ecoregions, when I specify maxError of 4 km it returns a nice simplified geometry, but when I up the maxError to 5 km it returns an unbounded geometry.
simplify(4000) vs simplify(5000):

Code and link below:
var ecoregions = ee.FeatureCollection("RESOLVE/ECOREGIONS/2017");
var geometry = ee.Geometry.Point([-151.06, 61.67]);

var simpError1 = 4000;
var simpError2 = 5000;

var myRegion = ee.Feature(ecoregions.filterBounds(geometry).first());

Map.addLayer(myRegion,{},'an ecoregion', false);
Map.addLayer(myRegion.simplify(simpError1),{},'low error simplify');
Map.addLayer(myRegion.simplify(simpError2),{},'high error simplify');

https://code.earthengine.google.com/7636a43eef5b2054e483bcc1ed7f12bd


texts: How can I add a total (sum) at the bottom of this table? I would like to sum the vines and Acerage columns. Using Arcmap.


texts:
I currently have a trigger defined in my database which executes when my edits in ArcMap are saved using the AFTER INSERT ON method.

Is there anyway to execute this trigger without having to select Save Edits in ArcMap? I have tried both BEFORE INSERT and UPDATE as well.


texts: Aerial photogrammetry book with all the math

texts:
Can anyone recommend a book that hits all the math for the general situation of doing serious photogrammetry calculations on aerial photographs?

Particularly an inexpensive book or one I can get very cheaply used, as my budget is quite constrained.

I am a mathematician. I am working on an aerial photogrammetry project which includes making corrections as finely as possible. I need a book that covers the general situation with all the necessary mathematics to do actual computations, particularly serious mathematical discussions of atmospheric refraction correction and radial lens correction, .... simple stuff like triangulating a point I can do myself.

I want a serious book because after significant searching I have found that free or online resources on the topic fall into two camps:

1) scholarly articles of current research that involve specialized methods that are not what I am after - I want the general big picture described precisely, not an article that takes 2 hours to read only to discover that it is only useful if you used their special calibration technique.

2) Books and summaries that that contain equations that are thrown in, but are not expected to be used and probably weren't used by the author, e.g. contradictory equations, undefined constants, variables that are defined but don't show up in the equations, etc.... These tend to look like they have the math, but were clearly not written with the idea that math would need to actually need to be used by the reader. For instance I have two free sources with radial correction equations, but they contradict each other, reversing the undistorted and the distorted x and y. I don't know which is really correct.

I have found that there are a number of photogrammetry books that have a sprinkling of formulas to seem more technical, but not what you would really need.


texts: How to set radius about 100km by buffer? (QGIS 3.0.1)

texts:
I am trying to create a buffer with radius about 100km (from the previous version of QGIS 2.12 it should be equal buffer distance 150000 with EPSG3857 (Pseudo Mercator) - fixed distance buffer.
I tried to use the same for buffer in the latest version of QGIS but its creating one big buffer around all the way or buffers which are less then 2km in radius. Could you please give me an advice what number should I set up please and in which format?
The other way which I am trying is used to fixed distance buffer as well but when start to run the analyse the the following error apears Error:

WKT: unknown projection (Mercator_Auxiliary_Sphere)


texts: How to write dynamic SQL for a function that has dynamic SQL parameter?

texts:
I am trying to execute a string SQL query in function i.e.

Create or Replace Function Route(edge_table text)
Returns table(sequence integer) as $$
Begin
Return Query Execute 'with r as(
SELECT seq, id1 AS node, id2 AS edge, cost FROM pgr_dijkstra( ' || '''SELECT id, source, target, cost FROM ' || edge_table || ' ''' || ', 7, 12, false, false );
)
select seq from r';
End;
$$ language plpgsql;

When I run function gets created, but when I call this function it says verify your query is not correct.

Actually dijkstra takes an string sql.like that,

SELECT seq, id1 AS node, id2 AS edge, cost
FROM pgr_dijkstra( **'SELECT id, source, target, cost FROM edge_table'**, 7, 12, false, false );

I want table name to be dynamic, how can I achieve this?


texts: Calculating minimum distance between points using PostGIS

texts:
I have a table (tab_Cities with cities which own more than 50000 inhabitants. The table also contains lat, long and the_geom as columns. I can already see the points via QGIS.

Now I want to write a postgresql query which gives me the minimum distance of the points but I only now how to get the distance of two points. But now I have around 60 points and I want to have the smallest distance of ALL points.

How can I do this?

I found some similar questions but they were using two tables. I am using one table only.


texts: Create a mapping between polyline and raster in ArcObjects

texts:
In my case I have two sources of data: a vector shapefile with a polyline (expressing the path) and a raster with some integer metadata values available for each cell.

What I would like to do is to query the metadata values from the raster cells for each segment of the polyline and then to create another column for this data in my shapefile, where this data would be finally stored. In other words, I want to perform a mapping between two different datasources. How can I perform such actions using ArcObjects?

I don't really need a code sample but rather a strategy of using available tools and functions.


texts: Using relative path in QGIS 3.0 Processing Modeler?

texts:
Using QGIS 3.0, we are working on creating a tool that will allow people to check that their snapping is correct with several linear features stacked on each other. A co-worker of mine figured out how to do it. However, when you first set it up, you have to open the graphical modeler and set the style manually. What she currently has is a path to her own desktop to the folder that has all the possible styles. I am trying to find a way to make it universally, find the CURRENT users profile, THEN search for the file path.

Currently the path is set to something like this:

C:/Users/jwagner/Desktop/Staging Folder/Layer Styles/CableVertex_Incorrect.qml

I tried editing it to use something like this:

Environ("USERPROFILE") & "\Desktop\V1.2\Layer Styles\CableVertex_Incorrect.qml"

We have also tried using "%," and "~" to no avail. The end goal is to have this path already set when we give people the model file. That way we aren't having to personally show every single person on our project how to set styles. I know from personal experience that it will lead to problems later on down the line when trusting everyone else to remember how to properly do it on their own.


texts: Conversion of (lat, lon, height) from GPS to another CRS

texts:
I have latitudes, longitudes and ellipsoidal heights derived from GPS readings, referenced to the ITRF08 system with vertical datum (CGVD2013) - http://webapp.geod.nrcan.gc.ca/geod/hmod/hmod.php

How can I convert the data to another reference system and another vertical datum?

I might need to convert the lat lon height data at every sampling interval in a RINEX file to another reference system and datum. The amount of data may be large and hence I doubt a web-based conversion platform would be suitable.


texts: Create dissolved buffer from multi-geometry (Union by shared attribute and spatial intersection)

texts:
I've to create dissolved buffers from multi-point input features. In the example below, the input table contains 4 features. Feature #2 consists of two point geometries. After creating a buffer, I get 4 polygon geometries:

Is there a way to group the result? The buffers of the points #1 and #2 are dissolved and should be a single multi-polygon feature (a).

What I've done so far:

-- collect all buffers to a single multi-polygon feature
-- dissolve overlapping polygon geometries
CREATE TABLE public.pg_multibuffer AS SELECT
    row_number() over() AS gid,
    sub_qry.*
FROM (SELECT
    ST_Union(ST_Buffer(geom, 1000, 8))::geometry(MultiPolygon, /*SRID*/) AS geom
FROM
public.multipoints)
AS sub_qry;

EDIT:

-- create sample geometries

CREATE TABLE public.multipoints (
gid serial NOT NULL,
geom geometry(MultiPoint, 31256),
CONSTRAINT multipoints_pkey PRIMARY KEY (gid)
);

CREATE INDEX sidx_multipoints_geom
ON public.multipoints
USING gist
(geom);

INSERT INTO public.multipoints (gid, geom) VALUES
(1, ST_SetSRID(ST_GeomFromText('MultiPoint(12370 361685)'), 31256)),
(2, ST_SetSRID(ST_GeomFromText('MultiPoint(13520 360880, 19325 364350)'), 31256)),
(3, ST_SetSRID(ST_GeomFromText('MultiPoint(11785 367775)'), 31256)),
(4, ST_SetSRID(ST_GeomFromText('MultiPoint(19525 356305)'), 31256));


texts: Applying a two step floating catchment with boundary limitations

texts:
The literature states that gravity models and the two step floating catchment area (which in itself is a version of the gravity model) method for analysing the accessibility of healthcare services by populations is considered the best. There are many variations taking into account distance decay functions, travel times and supply specific weightings.

The equation (summarised here by Kanuganti et al. (2016))

My questions applies to all of them and is something I can't quite work out from the literature or the mathematics.

Question

If you have an island for example then this is not a problem but if you are calculating the accessibility of a county that is surrounded by other counties (which themselves have supply and population locations) then they will have an influence on the accessibility of populations near the boundary. Is it true then that you must assess all counties in the analysis otherwise your results on the periphery will be distorted?

References

Kanuganti, S., Sarkar, A.K. and Singh, A.P., 2016. Quantifying accessibility to health care using two-step floating catchment area method (2SFCA): A case study in Rajasthan. Transportation Research Procedia, 17, pp.391-399.


texts: Changing index fails to show raster layer above vector layer in OpenLayers?

texts:
I have 6 different layers. Some are layers with markers and some other with raster data. In some point, when I zoom in the map, I create a vector layer. Upon adding the layer on the map I want a specific WMS layer to appear above the vector layer.
That's why I do this:

map.setLayerIndex(wms_layer_network, map.layers.length+100);

Although the layer indeed change position and comes to the top, on my map I don't see this change. The vector layer keeps being on top of the WMS.
This is how I check the layer array and their position (assuming the last value of the array is the one which will appear on top):

                    var mLayers = map.layers;

                for(var a = 0; a < mLayers.length; a++ ){
                    alert(mLayers[a].name)
                };

Does anyone know why this happens?


texts: GeoWebCache throwing an error 400: All backends (1) failed. Last request: when trying to send a WMS service request

texts:
We have a stand alone GeoWebCache (version 1.18) that is connected to a GeoServer (version 2.16.0). The data source for the layer in GeoServer is from PostGIS. Recently we are experiencing an issue where in if we send a WMS request to GeoWebCache, it throws an error.
400: All backends (1) failed. Last request:'http://*servername*:8080/geoserver/wms?'

Here is an example of the service request:
https://*servername*/geowebcache/service/wms?service=WMS&version=2.0&request=GetMap&layers=*layername*&bbox=-12503874.838186786,4960457.388858138,-12494090.898563793,4970241.328481132&width=256&height=256&srs=EPSG:3857&format=application/vnd.mapbox-vector-tile

If we try to request in EPSG:4326, the error is different:
400: X,Y values for the tile index were calculated to be {769.9999997724444, 1486.4759872853333} which had to be rounded to {770, 1486} and exceeds the threshold of 10%. Perhaps the client is using the wrong origin ?

Upon further investigating the issue we found out that when GeoWebCache sends the request to GeoServer, it does not use the same min/max bounding box as defined for that layer in GeoSever.
For example, if we publish a service with data source from PostGIS on GeoServer with min/max bounding box of -136.097427368164,28.8160457611084,-70.327880859375,58.095703125, the GeoWebCache has the following properties in the WMS GetCapabilities XML for that service:
<BoundingBox SRS="EPSG:3857" minx="-2.003750834E7" miny="0.0" maxx="0.0" maxy="2.003750834E7"/>
<BoundingBox SRS="EPSG:4326" minx="-180.0" miny="-90.0" maxx="0.0" maxy="90.0"/>

Notice that the maxx and miny in the above properties is 0.0. Due to this the bounding box of GeoServer and GeoWebCache misalign and we get the errors.
I have checked the gridsets in GeoServer for EPSG:3857 and EPSG:4326 and the min/max bounding box and resolutions look correct.
Is there a way or any setting that we need to change so that GeoWebCache picks up the correct bounds? Or is there some other way around to fix this issue?


texts: Merging layers via Python OGR ExecuteSQL

texts:
I have two tables with a common ID field that I would like to merge.

Is it possible to do this through python OGR? I am attempting to use the ExecuteSQL() command but am getting an error that table2 doesn't exist.

Here is my code:

connection = ogr.Open(connString)
connection2 = ogr.Open(connString2)

layer1 = connection.GetLayer()
layer2 = connection2.GetLayer()

newlayer = connection.ExecuteSQL("SELECT table1.field1, table1.field2, table2.field2 FROM table1 INNER JOIN table2 ON table1.field1 = table2.field1",  dialect = "SQLITE")

Is this something that can be done through OGR? Or should I be looking into something else?


texts: Action to directly export selected features as GPX

texts:
I would like to directly export selected tracks of a line layer as gpx. The common way is:

Select the line-layer
select the tracks for export
right-click on the layer -> save selected features as...
choose gpx-format, a filename and EPSG:4326 as CRS
hit OK-Button

But this takes too much time. So I'm looking for a fast solution. I know, the Layer-Actions could be perfect for this task!
How do I export selected features as gpx to a fixed file path with QGIS actions?


texts: Reducing/merging vector features after QGIS Vectorize?

texts:
i have a high res raster layer. Each pixel has an attribute value between 0 - 10000. I would like to vectorize the layer and finally get a geojson object. But atm i get for each pixel an own vector object and finally a huge geojson file. I think it's because of the different values.

How can I reduce the number of generated vector elements?


texts: No data in the attribute table when creating buffer in QGIS

texts:
I wanted to create a buffer around a polygon using QGIS. This buffer zone is something I will use to estimate the distances between patches of habitat but every time I do so there's no data in the attribute table of the buffer, see an image below.
Is there a way to extract the buffer zone data in QGIS?


texts: Query ShapeFile to get informations from coordinates

texts:
I have a shapefile that represents all the Hydrological Basins in the world. You can download the file from here:

http://www.fao.org/geonetwork/srv/en/resources.get?id=38047&fname=Major_hydrological_basins.zip&access=private

In that file, there's the whole world divided into areas (all the Hydrological Basins in the world). I need to create a function that query the shapefile with latitude and longitude, and gets back the values of that point, that are:
MAJ_BAS, MAJ_NAME and MAJ_AREA

I need to do this server side, without a GUI. To me seems simple, when I do in my QGIS I just click on a point and get back the informations I need, so I thought "there'll be a function that you'll give latitude and longitude, and the function will answer with the datas of that point". But I've been looking for 2 days and still get no solution about this problem.

I'm quite new to QGIS.


texts: What projection does the global climate region map from Wikipedia use?

texts:
I would like to use this global climate region map from Wikipedia as a "basemap" (by plotting over the image in R).

Is it possible to determine the projection of this map?


texts: Counting and storing events contained within polygon using ArcPy?

texts:
I have a fairly large project and I would like to automate some tasks using Python.  Problem is, I don’t have any experience with Python for ArcGIS.  And I rarely use ArcGIS until something big like this comes up, and usually I am under some time constraint (smh)  I am looking for some guidance (such as python syntax to use, methodology, etc.) to help me achieve the following, using the screenshots as a reference…

The red dots denote multiple events that overlap each other (contain multiple records) that fall within the polygon buffer points (labeled A and B, respectively):

The events are part of a shapefile specific to a road corridor, so that each corridor shapefile may contain several event points, representing motor vehicle crashes at an intersection. Each event contains a field “YEAR” that has either 2007, 2008, 2009, or 2010 in the field.:

The polygon buffer shapefiles have a separate field for each year listed above. Most corridors have about 6-10 rows, representing checkpoints along the corridor.

What I wish to accomplish is to:

Count all the events that fall within in the polygon points.
Populate the polygon fields for each year with the corresponding
counts of the “YEAR” field in the event shapefiles to the “CNT_20XX”
fields in the polygon shapefiles
Example: Polygon A has a total
of 10 events within it.  The counts for each year are 2007= 3, 2008=
1,  2009= 2, 2010= 4; so the Polygon for point A would have its
corresponding “CNT_20XX” fields populated with the sum of the counts
for each year in the events shapefile.

So...can this be accomplished in python?  Is there a non-python way to accomplish the same thing?  I have 16 corridors.


texts:
I have a Trimble Yuma tablet running Windows 7.  I need to build an application that allows viewing of points on a map and editing of a geodatabase (feature layers and tables). I'm hoping to use ArcGIS Engine SDK.

How do I leverage the GPS of the Trimble within my ArcGIS Engine app?  Anybody done this with a Trimble before?  Any code samples out there?


texts: GDAL translate in python - where do I find how to convert the command line arguments to use in python scripts?

texts:
I'm trying to use gdal.translate() to convert a raster from unit16 to int16 and I am trying to figure out what arguments to use.

I think in the command prompt version I would use -ot Int16 but that isn't working in python?

The documentation on the gdal website is for the command prompt and not to use in python. How do I know how to convert these argmunets so they work in python? I've come across this problem a few times while using gdal but was lucky enough this time to find an example of what I wanted to do. This time not so much.

Here is what I have so far. It then gives the warning: "Warning 6: driver GTiff does not support creation option outputType"

co = [
        'outputType=Int16'
        'noData=-10000'
    ]

# save as a new raster that is 16bit SIGNED to allow for -10000 value
print("copying to new raster" + newimage)
gdal.Translate(newimage, image, creationOptions=co)


texts: Arcpy to PostgreSQL using psycopg2 with MULTIPOLYGON using WKT

texts:
I am trying to create a MULTIPOLGYON in PostgreSQL from WKT format using arcpy.da module and psycopg2. I have had success creating point and multi point shapes this way however I cannot seem to get how to format the code for a MULTIPOLGYON.

I have read over the documentation http://postgis.net/docs/manual-1.4/ST_GeomFromText.html but it does not give any examples of inserting a variable holding the WKT geometry.

this is my part of the code I am having trouble with

cur = conn.cursor()  #setting up connection cursor
cur.execute('''drop table polygon''')
cur.execute('''CREATE TABLE polygon
            (pid serial primary key,
             Name TEXT,
             geom geometry);''')
conn.commit()
for row in arcpy.da.SearchCursor(shp, ["name","SHAPE@"]):
    name =  row[0]
    geome = row[1].WKT
    cur.execute('''insert into polygon(Name, geom)values(%s, ST_GeomFromText('MULTIPOLYGON'(%s)))''', (name,geome))
    conn.commit()

this is the error

traceback (most recent call last):
  File "<module1>", line 20, in <module>
ProgrammingError: syntax error at or near "("
LINE 1: ...values('Minnesota', ST_GeomFromText('MULTIPOLYGON'('MULTIPOL...

I have tried taking out the 'MULTIPOLYGON' but it still has not worked! what else can I do!


texts: Saving pyqgis QgsGeometryAnalyzer.dissolve output in memory?

texts:
How to save pyqgis QgsGeometryAnalyzer.dissolve output in memory (temporal layer)?

Usage example:

from qgis.analysis import *
shp_path = "C:\shape_file.shp"
QgsGeometryAnalyzer().dissolve(layer, shp_path, onlySelectedFeatures=False, uniqueIdField=-1)


texts:
I periodically use a nice basemap served by the Utah state government in QGIS. Recently they upgraded their server from ArcGIS Server 9.3 to ArcGIS Server 10.1. Ever since then, I can't seem to get the WMS basemap to load in QGIS in the same way. Previously I would get nice raster tiles of Utah. Now I just get plain county, municipal, and road identification layers without the pretty images. All of WMS basemaps I've tried from their server seem to work fine in ArcGIS, though.

Am I missing something? The URL to the basemap I generally use is below.

http://mapserv.utah.gov/arcgis/services/BaseMaps/Vector/MapServer/WMSServer?request=GetCapabilities&service=WMS


texts:
I would like to update a field with new values. With the code below it assigns the value '2' to all the rows and not just if the field == 2. But by 'print i' and 'print "it is 2"' everything looks correct. Any ideas?

    field = "field"
    field2List = []
    rows = arcpy.SearchCursor(input_layer)
    for my_row in rows:
        field2List.append(my_row.getValue(field))

    cursor = arcpy.UpdateCursor(input_layer)
    for row in cursor:
        for i in field2List:
            print i
            if i == 2:
                print "it is 2"
                row.setValue(my_field,2)
                cursor.updateRow(row)
            else:
                print "not 2"


texts:
How to read netCDF raster data in GRASS or QGIS?

I have a file with this characteristics:

Command: gdalinfo C:\Users\XXXX\gisdata\Test\mb.grd

Output:

Driver: netCDF/Network Common Data Format
Files: C:\Users\XXXX\gisdata\Test\mb.grd
Size is 412, 433
Coordinate System is `'
Origin = (-87.389553008819831,11.262298195006528)
Pixel Size = (0.000778017639680,-0.000768428196547)
Metadata:
  NC_GLOBAL#Conventions=COARDS, CF-1.5
  NC_GLOBAL#description=
    Projection: Geographic
    Grid created by mbgrid
    MB-system Version 5.5.2318
    Run by <gery> on <gery> at <Wed Nov  1 19:56:46 2017>
  NC_GLOBAL#GMT_version=5.2.1 (r15220) [64-bit] [MP]
  NC_GLOBAL#history=mbgrid
  NC_GLOBAL#title=Topography Grid
  x#actual_range={-87.38916399999999,-87.06939875009137}
  x#long_name=Longitude
  y#actual_range={10.929953,11.26191398090825}
  y#long_name=Latitude
  z#actual_range={-4783.08984375,-760.27880859375}
  z#long_name=Topography (m)
  z#_FillValue=1.#QNAN
Corner Coordinates:
Upper Left  ( -87.3895530,  11.2622982)
Lower Left  ( -87.3895530,  10.9295688)
Upper Right ( -87.0690097,  11.2622982)
Lower Right ( -87.0690097,  10.9295688)
Center      ( -87.2292814,  11.0959335)
Band 1 Block=412x1 Type=Float32, ColorInterp=Undefined
  NoData Value=nan
  Metadata:
    actual_range={-4783.08984375,-760.27880859375}
    long_name=Topography (m)
    NETCDF_VARNAME=z
    _FillValue=1.#QNAN

My first attemp to read it was with the standard gdal procedure

command:
r.in.gdal -o -k input="C:\Users\XXXX\gisdata\Test\mb.grd" outpu\   |
 |    t="mb" memory=300 offset=0 num_digits=0

but it looks to have just zeros, not the z#actual_range={-4783.08984375,-760.27880859375}

Here the r.info map=mb

|----------------------------------------------------------------------------|
 |                                                                            |
 |   Type of Map:  raster               Number of Categories: 0               |
 |   Data Type:    FCELL                                                      |
 |   Rows:         433                                                        |
 |   Columns:      412                                                        |
 |   Total Cells:  178396                                                     |
 |        Projection: Latitude-Longitude                                      |
 |            N: 11:15:44.273502N    S: 10:55:46.447629N   Res: 0:00:02.76634 |
 |            E: 87:04:08.435069W    W: 87:23:22.390832W   Res: 0:00:02.80086 |
 |   Range of data:    min = 0  max = 0                                       |
 |                                                                            |
 |   Data Description:                                                        |
 |    generated by r.in.gdal

There is a link with the file https://ufile.io/3i2ix and this belongs to this question How to restore a scarp surface grid to calculate material loss

System info:

GRASS version: 7.2.0
GRASS SVN revision: r70156
Build date: 2016-12-29
Build platform: x86_64-w64-mingw32
GDAL: 2.1.2
PROJ.4: 4.9.3
GEOS: 3.5.0
SQLite: 3.14.1
Python: 2.7.5
wxPython: 2.8.12.1
Platform: Windows-7-6.1.7601-SP1 (OSGeo4W)


texts: Batch travel distance calculation in QGIS

texts:
How to calculate the walking distance between postcodes and destination in QGIS?

I have the following shapefiles:

1 layer containing the destination;
1 layer containing 1,000 postcodes.

I've tried Directions to Points in the ORS Tool set but it does not handle multi-point layers. Is there another way to calculate the walking distance?


texts:
Is it safe to delete the contents of the ArcGIS Desktop Staging folder to free up disk space?
This would be the example file location: C:\Users\user\AppData\Local\ESRI\Desktop10.8\Staging


texts:
In QGIS 2.18, I have a shapefile of lines whose attributes are the lengths (in Km) of the lines. When displaying the labels of this attribute they appear repeated, as seen below:

How do I get the information to appear only once?


texts: Minimum bounding geometry outputs in meters

texts:
Which CRS is best to use on QGIS to get the minimum bounding geography outputs to present in meters?

I need the length and width measurements of > 1000 polygons, but can figure how to get the outputs in meters.


texts: How can I find the median column in a database?

texts:
I have data on period of construction for dwellings in dissemination areas. I have transferred this data to my study areas and would like to determine the median period of construction for each study area. The only problem is that the columns' information is the number of dwellings and I do not want the median of these, but the title of the column to populate the row in a new field (or something else indicating it, preferably the period in text format, but it's not the end of the world if it's just a number indicating the period).

I am attaching a picture of the attribute table with the relevant fields highlighted. There is a MEDIAN field but the data I am using was created by someone else and poorly documented, so I am not certain as to whether the calculation was already conducted or not.

(The seven columns represent seven non-overlapping time periods ordered chronologically.  The [Median] field appears to index the period at which total construction was half-complete; that is, it records the median time.)


texts: Why is my georeferenced image not visible in ArcGIS?

texts: I georeferenced a tif in ArcGIS 10, copied the folder into another PC opened the PGDB with 10.2.2 version, everything is loaded but i can't see the image itself. However the zoom to layer takes me exactly where it has to be. I figured it should work if i just add again the georeferenced image, how can i do that in ArcGIS?


texts:
Bit of a strange question but hope this is OK to ask here.

Has anyone heard of a '1-dimensional' projection of the world map - that is mapping all the points on the globe to a single line?

I was thinking of doing such a thing - trying to keep cities that are 'close' on the globe 'close' on the line.

Before I do this, I wondered what the state of the art might be in this area?


texts: Add elevation to contour vectors using sequential numbering by drawing a line

texts: I have several contour lines and I want to add the elevation as attribute. I found a tool in ArcGis Pro that would do exactly what I need, select the contours with a line, specify the starting value and increment with a value. The problem is that I need to increment with 0.5 and I am allowed to use only integers. Any ideas on how to add the values?


texts: using REST "delete features" in arcgis api for javascript

texts:
I am developping a web application that uses a feature service, what i want to do is save the data in the feature server, and then delete all the features at the start of the application.
I am able to save data in the feature service, my problem is how to delete the features in the start of the application.
I found a solution to that suggesting that I use the operation: "Delete Features" of the REST API, but I think I am making a mistake in the syntax, here is my code sample:

var theurl= "***/arcgis/rest/services/edt/FeatureServer/0/deleteFeatures"
var request = esri.request({

url: theurl,
// Service parameters if required, sent with URL as key/value pairs
content: {
f: "json",
where: "1=1"}

});

I don't know what I am missing can you please help


texts:
I would like to connect a PostGIS-PostgreSQL database created in pgAdmin to GeoNetwork.
Which files in the webInf under the geonetwork folder I am supposed to modify and how?


texts:
I have a big table that contains a set of osm_nodes_keys, like the following working example:

 a <- c(5993651, 27886377, 104537102)
 b <- c("a", "b", "c")
 df <- data.frame(osm_nodes_keys = a, myid= b)

I would like to be able to make a map in R, starting from the osm_nodes_keys column. How do I do a batch geocode?


texts: How to Unregister Replicas that have an invalid data-set information using ARCGIS REST APIs Version(10.2.1)

texts:
I am writing an application which will UnRegister replicas that are currently available. I managed to unregester new replicas which have the latest Server Url address stored in the DATASET_INFO field, but some of the replicas have an old address stored from a year ago, and those do not get unregistered, I get an exception.

The way I get the Replicas is by querying the gdb_items table, I gather the UUIDs of each replica and the value stored in DATASETINFO which is basically the url address that the replica is on.

The replicas that contain the new address of the server get unregistered just fine, all the replicas that has old address wont get unregistered because the address provided for the replica is not responding since it no longer exists. How can I Unregister those replicas ?

Here is a part of the code where I Unregister them

  public string UnregisterReplica(string user,string password,string url,string token,string replicaUUId) {

        const string unregUrl= "unRegisterReplica";
        if(!url.EndsWith("/")) {
            url += "/" + unregUrl;
        } else
            url += unregUrl;

        var values = new NameValueCollection
{
   { "username", user},
   { "password", password },
   {"replicaID",replicaUUId},
   {"token",token},
   {"f","json"},

};
        using(var wc = new WebClient()) {
            try {
                var result = wc.UploadValues(url,"POST",values);
                var response = Encoding.Default.GetString(result);
                return response;

            } catch(Exception ex) {
               LoggingManager.LogException("Exception Occured during request made on Url:" + url,ex);
                return null;
            }
        }
    }

The exception I get when trying to unregister is as follows "Unable to connect to the remote server" "{"A connection attempt failed because the connected party did not properly respond after a period of time, or established connection failed because connected host has failed to respond 133.100.10.03:8424"}"

Edits:-
So I looked through the sde tables to see if i can make up the relationships between the tables that hold the replicas, unfortunately I ended up at xml tables that I can't figure out. What I am trying to do now, is if the url is incorrect, I should just manually erase those replica records from the sde tables, but I am not sure which tables will take care of this problem.


texts: How to create a new 'gis' database in PostGIS?

texts:
I'd like to create a new database in PostGIS, so I can load stuff into it while the current database is being used. According to the docs

  Some packaged distributions of PostGIS (in particular the Win32
  installers for PostGIS >= 1.1.5) load the PostGIS functions into a
  template database called template_postgis. If the template_postgis
  database exists in your PostgreSQL installation then it is possible
  for users and/or applications to create spatially-enabled databases
  using a single command.

In my case this appears not to be so:

$ createdb -T template_postgis my_spatial_db
createdb: database creation failed: ERROR:  template database "template_postgis" does not exist

In the past I have messed around with copying the primary gis database, then deleting the contents of all the tables. There must be a better way. What do you do if you accidentally drop it?


texts: GDAL: C++ API to Convert Text file of Coordinate points to Raster Pixel line?

texts:
I am using GDAL to "Grow" areas form starting seed points.  I have a text file listng the coordinates of my seed locations.  The X,Y coordinates are in the coordinate system of a raster image.  I want to burn these point locations into the raster image.  I don't kow how to convert my point coordinates to the proper raster image pixel/line locations.  I have the geotransformation of the image.

I see GDAL classes in the API tutorial that return transformed coordinates from image row/col locations but not the reverse.  I am famialiar with the utility that does this, but I don't see a C++ programmers API to do this. Does anyone know how to do this?  This capability is  fundamental to what I'm trying to do.

Thanks in advance for any help you can offer.


texts: Making output shapefile taking fields of input, and adding new field populated with list using PyQGIS?

texts:
I found a couple of good answers like here for example, but I do not manage to adapt them and make them work...

What i am trying to do: take the input, create new output with same fields, add a new field with values calculated from another field.
For the moment, I managed to get the new values for the new field in a list:

from collections import Counter as counter

layer = qgis.utils.iface.activeLayer()
features = layer.getFeatures()

# find the index of the 'relict_id' field.
idx = layer.fieldNameIndex('RELICT_ID')

relict_id = []
for feat in features:
    attrs = feat.attributes()
    relict_id.append(attrs[idx])

relict_id_new = [x for x in relict_id]

c = counter(relict_id)
idx = 0
while idx in range(len(relict_id)):
    rel_id = relict_id[idx]
    if c[rel_id] != 1:
        for x in range(c[rel_id]):
            relict_id_new[idx+x] = '{}_{}'.format(relict_id[idx],x)
        idx += c[rel_id]
    else:
        idx +=1

I would like now to adapt this answer, but all I get is a bunch of errors. I tried to get rid off Processings since it does not seem to work for me (QGIS 2.18.2). I am not understanding pretty well the implication of the processes, I do not know if it is the good way to do it and what are the problems.

Here is the code for creating the new shapefile and new field:

##testes=group
##Ratio =name
##layer=vector
##output_file=output vector

from PyQt4.QtCore import QVariant, QFileInfo
from qgis.core import QGis, QgsFeature, QgsField, QgsFields, QgsProject
import processing

layer= iface.activeLayer().dataProvider().dataSourceUri()
layer = path[:path.rfind('|')]

mylayer = QgsVectorLayer(layer, "inventaris_points_in_roi_single_new", "ogr")
if not mylayer.isValid():
    print("Layer failed to load!")

vectorLayer = processing.getObject(mylayer)
provider = vectorLayer.dataProvider()

New_layer= processing.getObject(layer)
fields = QgsFields()
fields = provider.fields()
fields.append(QgsField('AREA', QVariant.Double,'',20,3))
fields.append(QgsField('PERIMETER', QVariant.Double,'',20,3))
fields.append(QgsField('RATIO', QVariant.Double,'',20,3))

writer = processing.VectorWriter(output_file, None, fields.toList(),
                              QGis.WKBMultiPolygon, New_layer.crs())

features = processing.features(vectorLayer)
for feat in features:
    writer.addFeature(feat)

del writer

And the error :

Traceback (most recent call last):
  File "<input>", line 1, in <module>
  File "c:/users/hocquetfp/appdata/local/temp/tmp_oajpe.py", line 18, in <module>
    vectorLayer = processing.getObject(mylayer)
  File "C:/Program Files/QGIS 2.18/apps/qgis/./python/plugins\processing\tools\dataobjects.py", line 216, in getObject
    ret = getObjectFromUri(uriorname)
  File "C:/Program Files/QGIS 2.18/apps/qgis/./python/plugins\processing\tools\dataobjects.py", line 244, in getObjectFromUri
    if normalizeLayerSource(layer.source()) == normalizeLayerSource(uri):
  File "C:/Program Files/QGIS 2.18/apps/qgis/./python/plugins\processing\tools\dataobjects.py", line 222, in normalizeLayerSource
    source = source.replace('\\', '/')
AttributeError: 'QgsVectorLayer' object has no attribute 'replace'


texts: Attribute Error when trying to run a script with Zonal Statistics

texts:
When I run my script I get the following attribute error:

Traceback (most recent call last):
HouseGradeBasementElev.tbx#HouseGradeBasementElev.py", line 31, in
<module>
AttributeError: 'module' object has no attribute 'ZonalStatistics_sa'

Failed to execute (HouseGradeBasementElev).

Here's my code:

import arcpy
import getpass

arcpy.CheckOutExtension("Spatial")
from arcpy.sa import *

from arcpy import env
arcpy.env.overwriteOutput = True

username = getpass.getuser()
GDB      = "C:\\Users\\"+username+"\\Documents\\ArcGIS\\Default.gdb\\"
arcpy.env.workspace = GDB

buildings   = arcpy.GetParameterAsText(0)
rd_cxn_line = arcpy.GetParameterAsText(1)
front_line   = arcpy.GetParameterAsText(2)

city_dem = "\\\\HRN-A001\\TUPW\\Infrastructure Services\\Infrastructure
Management\\Water Resources\\Wet Weather Database\\Major
Proxy\\Results\\Areas\\Ottawa East\\LOS.gdb\\DEM_1m"

arcpy.AddField_management(buildings, "Front_Grade", "FLOAT")
arcpy.AddField_management(buildings, "Basement_Elevation", "FLOAT")

#Selecting front building face
arcpy.FeatureToLine_management(buildings, "blds_lines")
arcpy.SplitLine_management("blds_lines", "blds_splt_lines")
arcpy.MakeFeatureLayer_management("blds_splt_lines",
     "blds_splt_lines_layer")
arcpy.SelectLayerByLocation_management("blds_splt_lines_layer", "INTERSECT",
     rd_cxn_line, "1 Meter")
arcpy.Select_analysis ("blds_splt_lines_layer", front_line)

#Front building face elevation
arcpy.ZonalStatistics_sa(front_line, "FID_Ottawa_East_Buildings", city_dem,
     "MEAN", "NO DATA", "bldg_fl_zonal_a")
arcpy.EucAllocation_sa(buildings, {}, {}, 1, "fieldID", {}, {},
     "bldg_euc_allo")
arcpy.ZonalStatistics_sa("bldg_euc_allo", "Value", "bldg_fl_zonal", "MEAN",
     "NO DATA", "bldg_euc_zonal_b")
arcpy.ZonalStatistics_sa(buildings, "fieldID", "bldg_euc_zonal", "MEAN", "NO
     DATA", "bldg_zonal_c")

To be clear, the attribute error is occurring on the first instance of the zonal statistics tool being used.


texts:
I have a CSV (value,lat,lon) that I am trying to create a raster grid from.

CSV File Sample:

temp,lat,lon
24.1,40.503,-69.248
21.7,38.461,-74.703
19.2,40.694,-72.048
18.5,41.443,-70.187
22.5,40.251,-73.164
14.3,44.287,-67.307
13.9,44.106,-68.109
20,40.369,-73.703

R Script:

library(raster)
library(rasterVis)

csv_path <- "file.csv"
csv_file <- paste(csv_path, sep="")
temp <- read.csv(csv_file) # takes a while
str(temp)
temp_coords <- cbind(temp$lon, temp$lat)
temp_pts <- SpatialPointsDataFrame(coords=temp_coords,
data=data.frame(temp$temp))
names(temp_pts) <- "temp"

cell_size <- 0.75
lon_min <- -125.0; lon_max <- -65.0; lat_min <- 25.5; lat_max <- 50.5
ncols <- ((lon_max - lon_min)/cell_size)+1; nrows <- ((lat_max - lat_min)/cell_size)+1
gridtemp <- raster(nrows=nrows, ncols=ncols, xmn=lon_min, xmx=lon_max, ymn=lat_min, ymx=lat_max, res=cell_size, crs="+proj=longlat +datum=WGS84")
gridtemp

us_fire <- rasterize(temp_coords, gridtemp, fun=mean)
us_fire

plot(us_fire, col=brewer.pal(9,"PuBuGn"), sub="CSV2RASTER")

Instead of displaying the correct values from the CSV file, it appears that the rasterizing converted the numbers to bogus values:

I need to create a raster with the correct values (column 1) from the CSV.


texts: How to check shape length before running code in Python/ArcPy?

texts:
I am using a code to convert a polygon into evenly spaced parallel lines

for point in range(0,numPoints+1):
        shortPoint1 = currentLines[firstLong].positionAlongLine(lineSpaceNum*point)
        shortPoint2 = currentLines[firstLong + 2].positionAlongLine(currentLines[firstLong + 2].length - (lineSpaceNum*point))
        parallel = arcpy.Polyline(arcpy.Array([shortPoint1.centroid,shortPoint2.centroid]), SR)

but the problem is that sometimes the first line it will produce has a shape_length of 0 which causes other scripts to break. If the object is square it will have a shape_length > 0 but if it is an unusual polygon it wont.


texts: Unable to run R in QGIS 2.6 - rgdal unavailable for R 3.1

texts:
I'm attempting to use R within QGIS2.6 on OSX Mavericks but cannot get R to run because rgdal is not available as binary for R3.1.2 (see QGIS log output from R execution console below).

I have installed rgdal on my system R installation using an alternative to install.package but QGIS doesn't seem to be seeing this. Does QGIS use an different R/rgdal location to my main system versions?

R version 3.1.2 (2014-10-31) -- "Pumpkin Helmet"
Copyright (C) 2014 The R Foundation for Statistical Computing
Platform: x86_64-apple-darwin13.4.0 (64-bit)

...some packages downloaded here...

The downloaded binary packages are in
    /var/folders/q4/r2k75kq956l54c4xtcxkscl40000gn/T//RtmpOSE4WR/downloaded_packages
tryCatch(find.package("spatstat"), error=function(e) install.packages("spatstat", dependencies=TRUE))
[1] "/Library/Frameworks/R.framework/Versions/3.1/Resources/library/spatstat"
tryCatch(find.package("rgdal"), error=function(e) install.packages("rgdal", dependencies=TRUE))
package 'rgdal' is available as a source package but not as a binary
Warning message:
package 'rgdal' is not available (for R version 3.1.2)
tryCatch(find.package("raster"), error=function(e) install.packages("raster", dependencies=TRUE))
[1] "/Library/Frameworks/R.framework/Versions/3.1/Resources/library/raster"
library("raster")
Loading required package: sp
library("rgdal")
Error in library("rgdal") : there is no package called 'rgdal'
Execution halted

added whole error log output January 7, 2015 23:07:23
R execution console output

R version 3.1.2 (2014-10-31) -- "Pumpkin Helmet"
Copyright (C) 2014 The R Foundation for Statistical Computing
Platform: x86_64-apple-darwin13.4.0 (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

options("repos"="http://cran.at.r-project.org/")
tryCatch(find.package("rgdal"), error=function(e) install.packages("rgdal", dependencies=TRUE))

package 'rgdal' is available as a source package but not as a binary

Warning message:
package 'rgdal' is not available (for R version 3.1.2)
tryCatch(find.package("raster"), error=function(e) install.packages("raster", dependencies=TRUE))
[1] "/Library/Frameworks/R.framework/Versions/3.1/Resources/library/raster"
library("raster")
Loading required package: sp
library("rgdal")
Error in library("rgdal") : there is no package called 'rgdal'
Execution halted


texts: Cropping TIFF in Python using GDAL returns NoneType object

texts:
I tried to clip a TIFF in Python, but I always get a NoneType object instead of a gdal.Dataset.
The most promising solution I've come about so far is this:
t = os.path.join(dir_omk(plot)) # gives path to the tif file
tif = gdal.Open(t)
tif_crop = gdal.Translate(destName = "/vsimem/clip.tif", srcDS = tif,
                                  projWin = [extent[0], extent[2],
                                  extent[1], extent[3]])
tif = None
tif_crop = gdal.Open(r"/vsimem/clip.tif")

In the first part, I opened a tif file as gdal.Dataset. Then I want to crop it to a certain extent. The "extent[i]-part" returns the following:
[637791.2537999959, 7364288.8728, 637909.596199995, 7364406.1047]

which are the x and y limits to which I want to clip the tif.
I read I could use "/vsimem/clip.tif" in order to save the clipped raster to a temporary file (I'll only use it in the script and don't want a copy on the hard disk).
I tried the above shown code and similar solutions with small variations (e.g., omiting tif_crop =). Since I am totally new to Python, I don't have any clue what to try further...


texts: Get the distance from a point on a polyline to polyline's from point

texts:
I have a polyline and a point on the polyline. I need to find the distance from the from point of the polyline to that point.

ICurve.QueryPoint method does the opposite thing. It gives a point based on given distance. In ICurve.QueryPointAndDistance method, specified point is a outside point and it gives shortest distance. I could not find any method that gives the distance I require. So I have tried something different. I took the polyline and split at specific point using IPolycurve.SplitAtPoint Method. I wanted to get my desired part's length from splitted part. But I guess it was a wrong procedure. I found a geometry collection and those geometry were not type of esriGeometryPolyline, they were type of esriGeometryPath. And the length of the path is not similar of my desired part's length.
 Can anybody give me some pointer about how to get the length?


texts:
I can load a PostGIS layer on QGIS' map canvas but "Toggle Editing" button comes "inactive". How can I load vector layer as editable?

The code given below loads PostGIS layer on map canvas & "Toggle Editing" button comes "active" but once I run the script below in Python console; when I press "Toggle Editing" button, QGIS gives me an info message: "Start editing failed: Provider cannot be opened for editing":

sql = "(SELECT * FROM tablename)"
uri = QgsDataSourceURI()
uri.setConnection("localhost", "5432", "database", "postgres", "password")
uri.setDataSource("", sql, "geom", "", "gid")
vlayer = QgsVectorLayer(uri.uri(), "abc", "postgres")
QgsMapLayerRegistry.instance().addMapLayer(vlayer)

# should the layer be set as "active" ???

vlayer.dataProvider()
vlayer.dataProvider().capabilities()
vlayer.startEditing()
iface.actionToggleEditing().setEnabled(True)

Can someone please give me a hint on how to load an editable layer programmatically?


texts: Is there a correlation between the absolute errors of nearby GPS receivers?

texts:
Suppose I have two different GPS receivers located x meters apart and their exact location is know. Then, would there be a correlation between the absolute errors of nearby GPS receivers? By absolute error I mean difference in distance and angle between the the GPS reading and the exact location.

If so, how to calculate this correlation? how much does it decrease with x?

Background: some one has told me that since the two receivers are close to each other, they would receive same signal with same error, thus they would most likely drift from the exact location at the same direction and error distance (assuming two identical receivers).


texts:
I am trying to use open source project named SciGRID to construct a power network. The project requires power tag data to be filtered (using osmosis) and loaded into a postgres database using osm2pgsql. I have currently managed to get osm2pgsql version 0.92 installed in my system. However the tables created by the new version (0.92) is not in line with tables required by SciGRID. They specify the table osm_planet_nodes table (created by osm2pgsql and .style file when loading filtered power data)  to contain a tag column, which it doesnt (contains only id, lat and lon). According to the following post,
Column "tags" is missing in PostgreSQL table created with osm2pgsql
the new version of osm2pgsql is likely creating some additional tables with potentially these information (separated into 3 or 4 columns). However my main goal is to get tables created by osm2pgsql to work with SciGRID, which require the tags column to be present in the planet_osm_nodes table.
I am therefore wondering if there is a way get osm2pgsql version 0.82 (requested by SciGRID) installed in my system. I am using ubuntu 18.04 and have osm2pgsql 0.92 installed.


texts: Web mapping application with time slider and raster data

texts:
I need to create a web mapping application that can show raster layer of historical maps with a time slider.

Basically, something like this: http://goo.gl/GYHOX

I have different maps of different years of the same are, that can be showed on a time series and moving on with a time slider. The interface can be simple, the important thing is the time slider.

I have programming experience only in HTML/CSS, little of PHP.

Any suggestions of which tools to use? I had a look at TileMill but the problem is the time-slider that I have no clue on how to build it.

Thanks in advance for who will help me.


texts: Problem in showing marker and Waypoints own OSRM server

texts: We have installed the OSRM step by step as provided in the link below benjaminberhault.com//post/2018/12/08/… but we are able to see the map but not the waypoints and markers.


texts: How to let the user select a spatial reference using ArcObjects?

texts:
I need to build a tool where the user selects a coordinate system for a given CAD file. What is the easiest way to do this with ArcObjects?

Obviously, I could populate a list with the hundreds of projection constants (i.e. esrSRProjCSType), but what user is going to want to scroll through such a huge list? Must be a better way... any tips?

Update:

Thanks everybody. I went with the SpatialReferenceDialog approach below. Now, anybody know how to use NewGeoTransformationDialog to apply datum trans?

private void btn_getCs_Click(object sender, EventArgs e)
    {

        string StatusMsg = "";
        bool Status = true;
        SpatialReferenceDialog oSpatialReferenceDialog = null;
        ISpatialReference oSpatialReference = null;
        ISpatialReferenceFactory2 oSpatialReferenceFactory2 = null;
        IGeographicCoordinateSystem oGCS = null;

        try
        {
            oSpatialReferenceFactory2 = (SpatialReferenceEnvironmentClass)new SpatialReferenceEnvironmentClass();
            oGCS = oSpatialReferenceFactory2.CreateGeographicCoordinateSystem(esriSRGeoCSType.esriSRGeoCS_WGS1984.GetHashCode());
            oSpatialReferenceDialog = new SpatialReferenceDialogClass();
            oSpatialReference = oSpatialReferenceDialog.DoModalEdit(oGCS, false, false, false, false, true, this.Handle.ToInt32());

            txtBx_srDwg.Text = oSpatialReference.Name;
        }
        catch (System.Runtime.InteropServices.COMException CE)
        {
            StatusMsg = CE.ErrorCode.ToString() + ": " + CE.Message;
            Status = false;
        }
        catch (Exception E)
        {
            StatusMsg = E.ToString();
            Status = false;
        }
        finally
        {
        }
}


texts: Data not showing when project is reopened

texts: I am a very new QGIS user and I thought I had the hang of it, until I saved a project then when I reopened it to continue work, nothing would show up. All the layers are still there and boxes checked to "show all layers", but none are appearing in the map area of the program. I have tried everything and done numerous searches to figure it out and no luck.


texts:
I want to model a 3D scene, showing layers floating above a flat map. Currently, I am able to show the floating layers, but I cannot add the bottommost map layer. I just want to add a "basemap", like ArcMap allows, but none of the files seem to be compatible. Just saving a basemap as a lyr file in ArcMap and importing it in ArcScene does not work, and raster images obtained from the Internet do not line up properly when added to the scene.

I know that the online basemaps that are available in ArcMap do not work in ArcScene. I was wondering if there was any way to add a simple, flat map to a scene in ArcScene?

The illustration below illustrates what I want as my final product:


texts: Turning off Feature Classes in ArcGIS Pro using ArcPy?

texts:
I'm struggling to create a python script that will simply turn off all feature classes in the table of contents within my ArcGIS Pro project. I think part of my problem is the syntax change from ArcMap to Pro. I know that ESRI has a web page dedicated to that very issue (http://pro.arcgis.com/en/pro-app/arcpy/mapping/migratingfrom10xarcpymapping.htm), but the diagram on that page just seems to confuse me more.

Here is my attempt at creating a script to to turn off feature classes in ArcGIS Pro.

import arcpy
p = arcpy.mp.ArcGISProject("Current")
for lyrx in arcpy.mp.ListLayers(p):
    print lyrx.name('')
    lyrx.visible = False

When I run the code above I receive the following message

File "<string>", line 4
  print lyrx.name('')
         ^
SyntaxError: invalid syntax


texts: Rotating geometry to specific angle using pyqgis?

texts:
I do not know how to rotate a geometry using pyqgis.

Is there any sample/example code available to show me how?


texts: Extract a list areas and indices of multi-part polygons (SpatialPolygonsDataframe) R

texts:
I'd like to extract a list of polygons from a SpatialPolygonsDataFrame that include more than one area (so a polygon with a hole for example, or just 2 areas in one polygon entry in general). There was a question answered previously which was very close to what I want but it wasn't testing for whether or not more than one area was listed, and it did not preserve the indices for which polygon had more than one area listed within it.

In the answer posted by @cengel, there are ways to list all the areas of the polygon, to get a list of polygons in the file, and to get the areas of a prespecified polygon, which is close but not what my specific question is:

-Example poly from cited question (without holes or multiple areas)

mp <- readWKT("MULTIPOLYGON (((30 20, 45 40, 10 40, 30 20)),((15 5, 40 10, 10 20, 5 10, 15 5)))")

-total area

sapply(slot(mp, "polygons"), slot, "area")

-get list of individual polys

p <- lapply(mp@polygons , slot , "Polygons")

-areas of individual polygons

lapply(p[[1]], function(x) slot(x, "area"))

My structure looks more like this:

And the polygons with multiple areas like this:

There are two areas even though for some reason they don't show up in the list. I can plot them individually by using:

plot(poly@polygons[[112]]@Polygons[[2]]@coords,type = "l", asp=1)

(this is showing that I'm plotting the 2nd of 2 areas)

So basically, I'd like a list of which polygons have areas with more than one entry, but also which polygons those are (the index).


texts: Slow performance viewing attribute table of PostgreSQL/Postgis views in QGIS 2.14

texts:
A Postgis view was created with a query that executes in about one second.

It is a parcel layer composed of about 20,000 polygons (the original layer loads without any problems).

When loading this view in QGIS 2.14 the geometry loads quickly enough (perhaps a few seconds longer than a regular table), but when viewing attribute table, QGIS fails to respond, freezes for a few seconds, and then shows the table. However, it is not possible to scroll through the table as it freezes or simply does not respond to mouse clicks.

At this point I will have to abandon the use of views since the performance is terrible, unless I'm doing something wrong on the QGIS or DB side.

Is there a way to improve the performance of PostGIS views in QGIS 2.14?


texts: Updating ArcGIS Online feature service Layer Error with FME?

texts:
I am having issues updating a feature service with Safe FME. I looked up the example model and as far as I can tell I am following it. It seems to be running most of the model and throws an error as it starts trying to update the Feature Service. I am not sure what the issue is. I set the Feature Layer to allow edits. I have the reader for the ArcGIS Feature Layer and the Reader for the shapefile source, also the FeatureMerger to add the FID Field to the update table.
LOG: http://pastebin.com/U1Hu06GS