texts: How to export elevated contours to a DXF file in QGIS

texts:
Having tried Project -> Export DXF as well as right clicking my vector layer and using the save by dxf feature I have found that the contours do not show in AutoCAD to be elevated from 0

I have checked the attributes of the contours vector layer and find and attribute set containing ID,ELEV and going further I have tried the command:

ogr2ogr -f "DXF" contours-out.dxf contours.shp -zfield ELEV

which I have run successfully although I found that I had to download some version of gdal tools and add the path to my system to the various gdal subfolders I am finding that still the contours do not elevate from 0 and all are still drawn to a 0 elevation in AutoCAD.

I also received at some point an error in my console about the creation of an additional attribute ID which is not supported reported by the ogr2ogr command.

Is there a way to export contours to DXF with elevations?

Before I posted this question I read:
Exporting QGIS contour attributes in DXF file?
and a few others like it.


texts: Writing massive kernel files for focal statistics in ArcGIS Desktop

texts:
I would like to run focal statistics on a raster and apply a weighted kernel file.
The issue is that I have a dense raster (and I would like to keep it this way for a reason) on the other hand I am interested in a vast focal area (eg. 1000 x 1000 px). I am about to write a kernel file to specify my weights and it seems completely crazy to type 10000000 digits into a notepad text file.

Am I missing something?

Is there any function or tool to create massive kernel txt files?

Is there a smarter way to get the same distance-weighted focal statistics result?


texts: Mapnik/Tilemill: Align or orient a Marker along a Line, such as a Gate

texts:
In Mapnik (actually Tilemill), using OSM data, I'm trying to draw a gate across a way, such that the gate is perpendicular to the road.  But the gate always gets drawn with the same orientation, unless I manually rotate it to match a specific road.

I've been using a Marker, placed along a line, both svg, arrow and ellipse.

Is there some way when querying the barrier points from the database, that the orientation of the line can be pulled at the same time and stored as metadata with the point, so that it can be used to then perform a unique rotation for each point?


texts: Why is osgeo4w (w QGIS 3.6.2) asking for admin rights?

texts:
On my work computer, I get prompted for an elevated login when I try to run osgeo4w-setup-x86_64.exe on Windows 7

Sounds like I should be able to install without admin rights. As of 5/8/19 are they or are they not required?


texts:
I have taken a raster image and coverted it to a vector file.  It is an Ordnance Survey map which shows field boundaries etc.

As the vector file has been created from the pixelated raster file the lines are very uneven.

Is there a way to tell QGIS to staighten the lines by taking the mid points?

Is that clear. Here are two images.

and


texts: Nominatim disk space Needed for whole planet-osm

texts:
How much disk space needed for installing planet-osm Nominatim file ?

After trying to install nominatim with planet-latest.osm.pbf I got this error:

Done 15249 in 44 @ 346.568176 per second - Rank 26 ETA (seconds): 80273.023438
index_placex: UPDATE failed: ERROR:  could not extend file "base/32849/1731142": No space left on device
HINT:  Check free disk space.

I have a SSD hard-disk with 1.2T

df -h /var/lib/postgresql/9.3/main
Filesystem                                              Size  Used Avail Use% Mounted on
/dev/disk/by-uuid/5d9be0b7-27a3-4085-837f-a96ae3c2400b  1.1T  1.1T     0 100% /

Is there temporary tables that I can removed in order to make some free place ?


texts: Does line feature have one or two dimensions?

texts:
In the coursera class fundamentals of gis it's said that lines are one-dimensional feature and polygons are two-dimensional.
Is this correct? if the line has a curvature two dimensions are needed to represent it right?


texts: Why the street label is not showing in a feature in QGIS 2.14.2?

texts:
Why this street label is not visible even with zoom?

It has 'name' in it...

Here_qgis2.14.2_file


texts:
I have two polygons created at st_concavehull(st_collect(geom), 0.99)
and I wonder how to extract some points from union area... example below

I have no idea how to achieve this but I think it's possible.

Postgres 9.5 / PostGIS 2.2.1


texts:
I have encountered a problem with the raster>extraction>clipper function in QGIS when a mask layer is used.

In the image below the diagram on the left shows part of the original raster.

The diagram in the middle shows the original with one of the red areas overlain by a polygon.

The diagram on the right shows the raster created by the clipper function using a mask.

Notice that the pixels in the clipped raster (on the right) are not in their correct location - offset downwards.

There seems to be an offset across all of the output file.  The offset is approx 130 m.
Does anyone have any ideas as to why this happens?  It does not happen when the extent option is used instead of mask.

I am using QGIS 2.8.1 (same happens in a previous version - 2.4).  The mask layer and raster layer are in the same CRS.


texts: Fuzzy Clip Layer to maintain feature geometry

texts:
I have a number of large datasets that I want to clip to a specific area to make more manageable and carry out further processing. I have a polygon outline of my AOI, however I don't want to hard cut the features to this, rather select and then export all the features that 'intersect' it, therefore retaining the original datas geomtetry and attributes - in the same way that a select by location does.
(I realise that I could individually export the selected features but I'm hoping there is a tool or process that will allow me to batch multiple layers to export at the same time as I need to repeat this process for multiple locations)
So my question:
How can I 'fuzzy' clip data rather than 'hard' clipping it?
In the example below I want to export the features that intersect the blue polygon (my AOI) including the features highlighted by the red boxes that cross the boundary of the AOI, but without cutting or editing their geometry.


texts: Importing Rasters to ArcSDE 9.3.1 using sderaster

texts:
I am trying to import a raster to 11g using the sderaster command.  It is a standalone jpg with a world file associated, and pyramids.
Am I right in assuming I cannot import a JPG, I can only convert to a jpg on the import process i.e. store it as jpg, but not import a jpg?
I am currently lost as importing the jpg file I get an "Error: Invalid input file specified" error.
I know I can do this via sdeexport - sdeimport, but I believe the correct process should be the sderaster -o import process
I am sure we were able to do this with ArcSDE 8.x


texts: How can I get the coordinates of a grid cell?

texts: How can I get the real coordinates (i.e. X and Y) of a given grid cell when I know its row and column position?


texts:
We have many datasets, each containing about 150'000 polygons. We would like to dissolve each dataset down to about 2500 polygons, based on a given attribute.

We've successfully done this using the ArcGIS Dissolve tool, but we would like to move towards an open-source solution. The problem is that every open source solution we have tried so far is incredibly slow. It takes ArcGIS about one to two minutes to dissolve one dataset. With other methods, it takes at least one hour, at which point I stop it since my guess is it has crashed completely due to memory issues.

These are the methods we have tried so far:

QGIS Dissolve tool (NOT the Saga or GRASS tools). I use QGIS 2.18,
would QGIS3 make a big difference?
GeoPandas/Shapely.
ogr2ogr (as described here http://darrencope.com/2015/01/23/dissolve-shapefiles-using-ogr/)
We would gladly try other tools such as GRASS, but we require this to be able to run on a machine for which we do not have admin rights. This means no installations or sudo commands.

Are there are there other tools or methods we can try? Are there some "tricks" to try out? For example, would dissolving in tiles and then dissolving the entire shapefile be more efficient?


texts: Display error when uploading images to google earth engine

texts:
I am adding square images (1024x1024 pixel, left image) from the xbd dataset (a disaster dataset comprising of satellite images) to the google earth engine. When displaying them they are shown as a rectangle with uneven aspect ratio. One can see that the image does not align with the underlying map because the roads (especially the needle curve in the lower part of the image) does not align with the streets surrounding them.(middle image)
Additionally when displaying many images of the dataset empty spaces between images appear and they overlap. (right image)
I checked the metadata of the image with gdal info and the pixel size is equal in both x and y. I await that therefore the image should be displayed as a square but it is not. (You can see the output of gdalinfo below)
Is there an explenation to this behavior?

gdalinfo output:
Driver: GTiff/GeoTIFF
Files: woolsey-fire_00000009_post_disaster.tif
Size is 1024, 1024
     Coordinate System is:
        GEOGCS["WGS 84",
            DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
        AUTHORITY["EPSG","7030"]],
    AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]]
Origin = (-118.812958604038045,34.105219763502319)
Pixel Size = (0.000004577876532,-0.000004577876532)
Metadata:
   AREA_OR_POINT=Area
   Image Structure Metadata:
   INTERLEAVE=PIXEL
   Corner Coordinates:
   Upper Left  (-118.8129586,  34.1052198) (118d48'46.65"W, 34d 6'18.79"N)
   Lower Left  (-118.8129586,  34.1005320) (118d48'46.65"W, 34d 6' 1.92"N)
   Upper Right (-118.8082709,  34.1052198) (118d48'29.78"W, 34d 6'18.79"N)
   Lower Right (-118.8082709,  34.1005320) (118d48'29.78"W, 34d 6' 1.92"N)
   Center      (-118.8106147,  34.1028759) (118d48'38.21"W, 34d 6'10.35"N)
   Band 1 Block=1024x1 Type=Int16, ColorInterp=Gray
   NoData Value=-99
   Band 2 Block=1024x1 Type=Int16, ColorInterp=Undefined
   NoData Value=-99
   Band 3 Block=1024x1 Type=Int16, ColorInterp=Undefined
   NoData Value=-99


texts: Calculating "relative slope" per Mora-Vahrson method in ArcGIS

texts:
I'm a longtime lurker who finally decided to register on the forum (in no way driven by my need to have a certain question answered, of course...).

I was contacted by a colleague asking for help calcuating an input for a landslide model he's working on. Specifically, the model he's following (the Mora-Vahrson paper/model) asks for the slope input to be not in degrees or percents, but "relative slope." This is obtained by using the formula:

Sr= [(ElevMAX-ElevMIN)/(Area)],

where the the min and max elevations are the lowest and highest points within a bounded area (most commonly 1 square km); these area searches/calculations are carried out across the desired study region without overlap.

I've tried multiple ways to calculate this using my existing skillset (limited; I'm an archaeologist with a background in GIS who was asked to help out...) and program (ArcGIS). I'm using a 90m resolution SRTM DEM for the area and have currently settled on using the Block Stats tool to calculate the MIN (and the MAX) for for a specified area 1000 x 1000m, and then entering those resultant rasters into the Raster Calculator as part of the Mora-Vahrson equation), and dividing by 1 (as the area of the Block Stat calculation, expressed in square km).

What I would like to gauge from more experienced GIS users is whether or not the Block Stats process I described seems to meet the requirements for the model's relative slope component. If anyone else has had experience calculating relative slope or a similar process, I'd be very interested to hear about it. Thanks.


texts: Maintaining Form Designs with one Geopackage to Multiple Different Projects

texts:
Long time reader, first time poster.
I work for an Archaeological consulting company in Australia and one of the major parts of our job is to go out and survey for aboriginal heritage sites before constructions can be approved. The requirements for recording are standardised from the government, so this makes creating forms quite easy.
I have been trying to move our documentation to a digitised method using QGIS, Mergin and Input. Which works great for what we need it for.
To do this, I am hoping to create a master form as a GeoPackage with all of the forms required for each of the site types we are likely to encounter in the field (Artefacts, Scarred Trees, Shell Middens etc). Then transfer this Geopackage to  the QGIS project for each separate job. I have created some nice looking forms that function well in Input.
The only problem is, I have to use a master project to keep the formatting and value map lists intact. If i copy and add the Geopackage with all the survey layers in it, all the form editing i did gets destroyed and i am only left with the fields again. The forms have quite a few entries and rather large value maps, so it is not something I can do at the start of each job.
Is there a way to ensure that the formatting for these forms gets carried over into new projects within QGIS, so i don't have to enter in all the value map information and ranges every time i start a new job?
Thanks so much for your help


texts: PyQt4 pushbutton not connecting for my QGIS plugin

texts:
I'm trying to display information in a text browser with the push of a button in my QGIS plugin. The .ui is created using PyQt4 (button's object name is calcConso). No error message, but the button isn't connecting. The button simply displays values returned by my stats function, but the function itself isn't actually called anywhere. If it should be, I don't see where yet.

Any ideas what's missing?

def initGui(self):
# connect button to calculation function
    QObject.connect(self.dlg.ui.calcConso, SIGNAL("clicked()"), self.buttonClicked)

def buttonClicked(self):
    global stats, values
    self.dlg.clearTextBrowser()
    self.dlg.setTextBrowser( str(stats, values) )

def stats(self):
    global stats, values
    ...
    return stats, values

in a separate dialog class i've defined clearTextBrowser() and setTextBrowser(), which is imported.

def setTextBrowser(self, output):
    self.ui.txtFeedback.setText(output)

def clearTextBrowser(self):
    self.ui.txtFeedback.clear()


texts: ArcGIS10 Highest point via polygon features

texts:
I am using ArcGIS 10.0 with 3D Analyst, Spatial Analyst and Military Analyst.

I feel like I am halfway to solving a problem but do not seem able to completely solve it.

I have a polygon shapefile with about 100 polygons and each one has a unique ID.
I want to overlay that on a raster and then for each polygon extract the highest value.
I used the Military Analyst>Highest Point and also Spatial Analyst>Zonal Statistics and they both work.
However, the main thing I need in the output is what polygon it is linked to.
It seems like the zonal statistics work on zones but the id of the zone from the polygon is not transferred to the output.

So how can I get the max values back to the polygon shapefile to update a field?

Update: The solution I found was Zonal Statistics To Table in the Spatial Analyst Toolbox. I just selected ALL statistics.


texts:
I would like to write a cost functiont that takes into account the previous edge, rather than only the current one, to calculate for instance:

  Route with the smallest number of turns

For that I would have to check for turns every time, ie. check the difference in angle between my current edge and the previously visited one.

I have seen similar things for [pgr_trsp][1] but these seem to totally restrict certain given turns rather than returning an aggregate on the number of turns.

Is there an easy way I am not seeing for calculating this?


texts: How to join the attributes table from one input layer to another layer with unique id using python

texts:
I have one input layer with the attibute values and i wanted to join the attribute table of the input layer to the another layer with unique values using python.
gis.stackexchange.com/questions/89706/how-to-merge-shape-files-with-attributes-from-a-python-script-in-qgis
This link gives the information about about how to merge the features and attributes manually. but i want to join the atttribute tables. How can i do that using python?


texts: Change the symbology of a set of points in the legend

texts:
I use QGIS 3.01.

I have these set of points, where they have the values "yellow, red and green". I was wondering how I can change the symbology in the QGIS layout editor to only show three symbols, with yellow, red and green as label? I am trying to edit the legend, not the actual map items.


texts: Calculating max and average distance from point to points in ModelBuilder?

texts:
I have two point layers. One contains nest locations and the ID of the individual using that nest. The second point layer contains movement points for all individuals in the area.

I would like my ModelBuilder model to generate a table with a row for each nest showing the MAX and AVERAGE distance from the nest to the movement point which are from the individual ID using the nest.

The images below are a mockup of the basic idea of what the input looks like and what I'd like the results to contain.

So far, and with some help on this board, I've used Point Distance to output a table which has the distance from each nest point to all other points. Then I've been playing with Summary Statistics to try to generate the output table. The problem is that the Point Distance tool, I believe, is creating distance from the nest to ALL points and not just the ones that match the ID. I need to put a filter into the SurveyPoints at the beginning to iterate through the nest IDs. But ideally I want this all to wrap into a single final table.

My untried idea is to try to do a iterate through SurveyPoints using a filter based on the nest id (I don't know how to do this). Then somehow append the results of the Point Distance or the Summary Statistics to get a final table.

---UPDATE---
I tried the approach as suggested which looked like it was going to work. However, I'm getting a weird snag. Let me know if this should be a new question.

My Iterate Field Values is setup just like the suggestion with the Iterator using the nests UnitID field which is a string.

When I run the model (after validating) it gives me an error because the Value = 1 for some reason. Actually if I click on the Value circle before I run the model the greyed-out value in the Value box is '1'. None of the UnitID fields are 1.

So when I run the model I get the error below. "Nests_meter" is the true name of the "Nests" labelled circle in my model.

ERROR 000358: Invalid expression UnitID = 1
An invalid SQL statement was used.
An invalid SQL statement was used. [Nests_meter]
Failed to execute (Make Feature Layer).

The weird thing is that the first time I added the iterator, it did contain a valid value. But after that it just starts with a value of 1. I recreated the model from scratch but still get the '1'.


texts: How to exclude outliers from raster using ArcGIS for Desktop?

texts:
I calculated terrain curvature from a DEM. The DEM I use is interpolated from point elevation data. Generally, the DEM seems to be good, but some artifacts from interpolation exist. Those artifacts cause extreme values for curvature. Around 90% of the curvature values are between -25 and 25, but the other 10% are between -500 and 500.

I can use histogram equalize to properly visualize the curvature raster. For further analysis however I'd like to replace this extreme values.

So my question is: How can I exclude this extreme values from the raster?
I imagine something like this would work (pseudo-code):

m = median("curvature")
IF("curvature.value") > m * 10 THEN ("curvature.value") = m * 10
IF("curvature.value") < m * 10 THEN ("curvature.value") = m * 10

What tool would I use for this? I tried Raster Calculator, but I couldn't get the syntax right.

Software: ArcGIS 10 with spatial and 3D analyst


texts: Installing QGIS on Mac OS Sierra 10.12.5. "The GDAL 2.1 framework is required."

texts:
I use the version on web("http://www.kyngchaos.com/software/qgis").
Whether I come to install QGIS 2.18.9-1 or QGIS 2.14.15-1,  I click on the QuantumGIS.pkg and it opens the install dialogue box but with a warning "The GDAL 2.1 framework is required."

Could anyone give me a clue how to properly install QGIS on  Mac?
The version is Mac OS Sierra 10.12.5.


texts: Using Spyder (Anaconda) environment to execute GRASS Python script in Windows 10

texts:
How can I set my spyder environment in order to use GRASS Python script in Spyder for Windows 10? I need to use some Python script from GRASS GIS. I am really new in Python and manipulate with scripts.
I follow the following setup to set my spyder environment but I always got error messages (No module named grass.script)

@echo off
rem set OSGEO4W_ROOT=c:\OSGEO4W64
rem set PATH=%OSGEO4W_ROOT%\bin;%PATH%
call c:\osgeo4w64\bin\o4w_env.bat
call %OSGEO4W_ROOT%\apps\grass\grass-7.6.0\etc\env.bat
set GDAL_DRIVER_PATH=%OSGEO4W_ROOT%\bin\gdalplugins\1.9
set PATH=%PATH%;%OSGEO4W_ROOT%\apps\qgis\bin
set PATH=%PATH%;%OSGEO4W_ROOT%\apps\grass\grass-7.6.0\lib
set PATH=%PATH%;%OSGEO4W_ROOT%\apps\Python27\Scripts\
set GISBASE=C:\OSGeo4W64\apps\grass\grass-7.6.0
set GISDBASE=C:\Users\Wael Hassan\Documents\grassdata
set GISRC=C:\Users\Wael Hassan\AppData\Roaming\GRASS7
set GRASS_SH=C:\OSGeo4W\apps\msys\bin\sh.exe
set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\qgis\python;
set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\Python27\Lib\site-packages
set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\grass\grass-7.6.0\etc\python
set QGISPATH=%OSGEO4W_ROOT%\apps\qgis
start c:\OSGeo4W64\apps\Python27\Scripts\spyder.bat


texts: PostGIS SQL Query Filter by Bounding Box

texts:
I have the following query on a PostGIS Database

SELECT row_to_json(fc)
FROM (
  SELECT 'FeatureCollection' As type,
         array_to_json(array_agg(f)) As features
  FROM (
    SELECT 'Feature' As type,
      ST_AsGeoJSON(lg.geom)::json As geometry,
      row_to_json((id, name)) As properties
    FROM lines As lg) As f
  ) As fc;

which returns a GEOJSON Object of all features in the database.

How do i filter the results by a Bounding Box already in the query?


texts: Can HTML annotations be autogenerated for each location in QGIS 2.4

texts: Can HTML Annotation Boxes be autogenerated on the whole map, linked to locations, and always have the same html file loaded for all of them?


texts: How to find Bounding Box given 4 Lat/Long and Bearing?

texts:
I need to find a bounding box given 4 lat/long points and a bearing (as seen in the example picture). I always know which two points are lined up by the bearing (1 and 2 in the example), so I will always know the length of the bounding box.  The width however is arbitrary, with the points being anywhere along the lines (3 and 4 in the example).

My first thought is that I'll have to calculate the angles between the points (1 & 3, 1 & 4, 2 & 3, 2 & 4) and then use a series of "law of cosine" equations to calculate the corner points.  Is there a simpler way?  Would that even work?


texts: Import only planet_osm_polygon osm2pgsql

texts:
My purpose its to only use queries with the ST_Contains function, so for me to have only the planet_osm_polygon its enough i think. How can i tune the osm2pgsl cli command or the style file to import only this kind of polygon data without POIs, roads, etc...?
I dont know what tables are strictly necessary or which ones can be avoided.


texts:
I'm new to pyqgis/python and programming in general..  and i am working on my first small project. I want to create plugin so that I can add layers from my PostGIS database to a QGIS project.

I suppose I need to create a list with my layers and the connection parameters to QGIS and populate the combobox with the list. I am not sure how to go about it, are there someone willing to give me some pointers?

Best regards


texts: Getting node information of a way in OSM2PGSQL data

texts:
I have a feeling this is not possible, but I am looking to get information on the nodes which are part of a line in the planet_osm_line table.

The issue which I am trying to handle is I am trying to map out train routes, which I can identify through the "route=train" tag in the planet_osm_line table.

However, I want to identify the train stops. Nodes which indicate a train stop are marked "railway=station" tag in the planet_osm_point table.

I'm looking for a way to associate a train stop, or nodes in general, in the planet_osm_point table to the train route, or any way in general, in the planet_osm_line table.

As a final note, I am working through SQLAlchemy, but I have a feeling that is not important to this particular problem.


texts: I've seen that gdalwarp has the -multi option to use multiple threads. But I don't see that option for gdal_translate. Does gdal_translate support multi thread?


texts: Automated QGIS Atlas with layer styled according to parameter

texts:
I have a unique layer (a map of a country divided in regions) which I'd like to to style differently according to a variable.
I would like to color country regions according to a variable that changes from a map to another.

Here is the sample table that I have:

I would like to automatically create 3 maps with regions that are colored:

Map 1: regions colored according to car registrations in 2010
Map 2: regions colored according to car registrations in 2011
Map 3: regions colored according to car registrations in 2013
etc..

I say automatically because I have to create hundreds of maps of the same country, colored each time differently.

I'm currently running QGIS 2.16.


texts: "save as" a vector layer with new source files but retain formatting

texts: I often need to duplicate a vector layer, creating new data source (shapefile) files so that I can then delete some features in the duplicate layer without deleting them in the original data. I need to retain all of the formatting settings, etc. as I do this. I have been using "save as" (after right-clicking on the layer to be duplicated), but then need to go through and format the symbol, labels, etc.. Is it possible to do this while retaining the formatting settings?


texts:
Is there any map projection of the earth that has been invented by a woman?

Context: I am starting a math outreach program and one of the activities I am proposing involves playing with different maps projections. I also want to give a historical perspective on the different maps projections and the related debates stemming from these.


texts: Converting individual x,y coordinates into individual polygons

texts: I am currently struggling to convert x,y coordinates into individual polygons. Feature polygon and points to line is not what I am trying to do. I have sampled vegetative plots on the landscape and have recorded the X,y coordinate for the bottom left corner of the sampling polygon. I have then taken the bearings from the baseline and first transect to orientate the polygon and I have the length measurements for the baseline and transects of the polygon. I have all the information to produce and orientate a polygon from each GPS coordinate, but can not seem to find a function in arcGIS to accomplish this task. Obviously I could create feature polygons for each x,y coordinate and define each length and angle, but I have a several hundred x,y coordinates to convert to polygons. I have added a table of the information I want to use, and a photo to visually depict what I want to do.


texts: Upslope and downslope distance or area using ArcGIS Desktop?

texts:
Is there a possibility to determine the surface area above the road and the surface of the terrain below the road or between two or more roads? I want to pull the trees down the hill and up the hill to the roads. The price of skidding trees is different so I would like to know what the surface (ha) is above and below the road?

Here is an example, orange is uphill and red is downhill.
I think this is solved with a raster calculator

Here is the method, but I do not know how to solve it in ArcGIS Desktop


texts: How to update a bulk selection of fields using field calculator in QGIS 2.8.1

texts: I have an attribute table showing vegetation communities.  The table has three columns ("ID", "Veg Comm" and "Area_ha".  Instead of having to copy and paste text after each polygon created into each field (i.e. coastal wattle scrub), is there a way to highlight a selection and mass copy?


texts:
import arcpy
import sys
import os

#arcpy.env.workspace = "C:\\Test\\addresspoint"
#arcpy.env.overwriteOutput = True

input_database = "C:\\Test\\addresspoint\\connection\\Test.sde"
out_layer_name = "addPointtbv"

Query = """ SELECT S_ID, S_Key, Address, City, State, ZipCode, Latitude, Longitude FROM [ST].[dbo].[ST] """
Unique_Identifier_Field_s_ = ["S_ID", "S_Key"]
Shape_Type = "POINT"
SRID = "3857"
Coordinate_System = """PROJCS['WGS_1984_Web_Mercator_Auxiliary_Sphere',GEOGCS['GCS_WGS_1984',DATUM['D_WGS_1984',SPHEROID['WGS_1984',6378137.0,298.257223563]],PRIMEM['Greenwich',0.0],
UNIT['Degree',0.0174532925199433]],PROJECTION['Mercator_Auxiliary_Sphere'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',0.0],
PARAMETER['Standard_Parallel_1',0.0],PARAMETER['Auxiliary_Sphere_Type',0.0],UNIT['Meter',1.0]];-20037700 -30241100 10000;-100000 10000;-100000 10000;0.001;0.001;0.001;IsHighPrecision"""

# Process: Make Query Layer (Make Query Layer)
add_pnt = arcpy.MakeQueryLayer_management(input_database, out_layer_name, Query, Unique_Identifier_Field_s_,Shape_Type, SRID, Coordinate_System)
print(out_layer_name + " query layer is created")

X_Field="Latitude"
Y_Field="Longitude"
Spatial_Reference = Coordinate_System

# Process: Make XY Event Layer (Make XY Event Layer)
in_layer = out_layer_name
out_layer = "xyeventlayer"

arcpy.MakeXYEventLayer_management(in_layer, X_Field, Y_Field, out_layer, Spatial_Reference)
print(out_layer + " is saved")

Error:

  arcgisscripting.ExecuteError: ERROR 000223: Cannot extract event table
  properties

  Failed to execute (MakeXYEventLayer).

It does work in ArcGIS Pro based python 3.6.9 (64 bit) anaconda, It is not working at all in 2.7.5 and 2.7.16 which desktop based installation and it is 32 bit.


texts: Getting attribute list ONLY of particular published layer (WFS or WMS)

texts:
I do NOT want to use GetFeatureInfo as it is implemented using map click. I just want to get the columns list.

I tried

http://example.com/geoserver/wfs?
  service=wfs&
  version=2.0.0&
  request=DescribeFeatureType

but it did not worked for me.


texts: Syntax error in parameter formula of field calculator funtion

texts:
I'm using a QGIS 2.8.1 in windows.
Have some doubts with syntax when use string in python console. I want to use the field calculator funtion like:

processing.runalg('qgis:fieldcalculator', input, field_name, type_field, width, precision, Create new field, formula, output)

In the formula parameter I want use the condition:

case when "Municipio" is NULL then 'PENICHE'

I know that PENICHE is a string, so have to be between '' (plicas) but all the expression should be too between '' (plicas) like:

' case when "Municipio" is NULL then 'PENICHE' '

The problem is that there happened a syntax error.

How can I correct the expression?

The code I´m using is:

processing.runalg('qgis:fieldcalculator', input, "Mun", 2, 27, 0, True, ' case when "Municipio" is NULL then 'PENICHE' else "Municipio" end ', output)

I used formula = "case when \"Municipio\" is NULL then 'PENICHE' end" know its working!


texts: ArcMap toolbars always stacking up instead of remaining where placed when opening MXD on multi-monitor configurations?

texts:
In ArcMap 9.3, we had issues with our toolbars always stacking up instead of remaining where they were placed.  This was mainly an issue on dual monitor workstations and we would occasionally see the issue on a single monitor workstation.  When we reported the issue to Esri, they said it was a known bug: NIM003161.

I won't have an opportunity to install ArcGIS 10 on a dual monitor workstation for a couple of months and I was wondering what the user experience is when opening saved MXDs on dual (or triple) monitor workstations.

Are toolbars staying put?

I looked up bug #NIM003161 and its status is on-hold and last modified on March 25, 2010 and the version fixed is still blank.  That part worries me a little, but I'm hoping Esri just hasn't gotten around to updating that bug.


texts: Lost the properties in ImageCollection after Mosaicking a Image Collection in Google Earth Engine

texts:
As per the question from Mosaicking a Image Collection by Date (day) in Google Earth Engine

I tried the answering code from Rodrigo E. Principe.
That code is very good. I can mosaic all the images (tiles) all my area.

But my problem is that I lost all the properties in each image list after Mosaicking.

In order to analysis with Time series, I need to have metadata('system:time_start')
but I can not find this Metadata in in each image list from the new ImageCollection(mosaicking).

So, Please advise how to transfer all properties to new ImageCollection after mosaicking

The Mosaic code is modified from 'Rodrigo E. Principe' by the following code:

var ROI =
    ee.Geometry.Polygon(
        [[[100.14445909748613, 14.689307337967634],
          [100.14445909748613, 14.04009719681829],
          [100.98491319904863, 14.04009719681829],
          [100.98491319904863, 14.689307337967634]]], null, false);

Map.centerObject(ROI, 10);

var start = ee.Date('2017-09-03');
var finish = ee.Date('2017-10-30');

var S1 = ee.ImageCollection('COPERNICUS/S1_GRD')
                     .filterDate(start, finish)
                     .filterBounds(ROI)
                     .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VH'))
                     .filter(ee.Filter.eq('orbitProperties_pass', 'DESCENDING'))
                     .filter(ee.Filter.eq('instrumentMode', 'IW'));

print('S1',S1)

/////////////////////Mosaicking a Image Collection by Date (day) ////////////////

// Difference in days between start and finish
var diff = finish.difference(start, 'day');

// Make a list of all dates
var range = ee.List.sequence(0, diff.subtract(1)).map(function(day){return start.advance(day,'day')});

// Funtion for iteraton over the range of dates
var day_mosaics = function(date, newlist) {
  // Cast
  date = ee.Date(date);
  newlist = ee.List(newlist);

  // Filter collection between date and the next day
  var filtered = S1.filterDate(date, date.advance(1,'day'));

  // Make the mosaic
  var image_m = ee.Image(filtered.mosaic());

  // Add the mosaic to a list only if the collection has images
  return ee.List(ee.Algorithms.If(filtered.size(), newlist.add(image_m), newlist))
}

// Iterate over the range to make a new list, and then cast the list to an imagecollection
var S1_mosaic = ee.ImageCollection(ee.List(range.iterate(day_mosaics,ee.List([]))));

print('S1_mosaic',S1_mosaic)

// Adds a time band to the image.Original from (https://developers.google.com/earth-engine/reducers_regression)
var createTimeBand = function(image) {
  return image.addBands(image.metadata('system:time_start').divide(1e18));
};

var linearFit = S1_mosaic.select(['system:time_start', 'VH'])
  .reduce(ee.Reducer.linearFit());

print(linearFit)

Code link : https://code.earthengine.google.com/dbd55ea2354984025b0fc10874f23ea9


texts:
I would like to acquire the coordinates for a specific corner of a Section.  The Section is Section 12 located as part of/within Township 8 South, Range 2 West.  This section is located in the southwest portion of Mobile County, Alabama.  This Section 12 is shown on the City Map of Mobile (http://maps.cityofmobile.org/citymap).  It is not part of the PLSS data base.

I have ask for this information from the Mobile City/County GIS office.  I think the information given me was 'measured' from the map view itself.  This is not acceptable.

I believe this information (coordinates) could be got from a shapefile in some way.  I have as yet been unable to figure out how to do such. The corner of this section is the beginning point of a survey description of some property my brother has.  We are not surveyors but are interested in gaining more knowledge in our pursuit of understanding the survey itself.

Can you help?


texts: Using OpenLayers 3 or 2 to add layers from Google Maps and Mapserver?

texts:
I am beginner in OpenLayers and am starting with the www.openlayers.org documents, which use OpenLayers 3 to teach the basics of OpenLayers. Now, I would like to add layers from Google Maps and Mapserver, but I can't find that in the documents.

I have found many site related to my questions, but all of them use different syntax that I guess are related to OpenLayers 2.

for a beginner would you recommend ol3 or ol2 ?
if you recommend ol3, where I can find what's new in ol3 syntax?


texts: Advantages of OSM over Google Maps for custom map in app

texts:
I wanted to map my city on a private map. I have read about maps and I found that I will need tiles. Looking at available service, could not find one that allows me to host them in my server. Which means I have to comply with their own terms.

So I was left wondering what is the advantage then, of using OSM over Google Maps?

Am new to GIS so it might be no brainer question. In that case I apologize

UPDATE

So here is what I want basically (Just to give an idea since the idea itself is huge). I take OSM/GM data into my database. I add different stuffs in different locations local to that City and correct them to better accuracy where it is not. Add bus stops, hotels, homes, et al.

Then I make simple app to guide someone to go to place they want. That is the basic of the idea.


texts:
I have two points on a EPSG:4326 layer. Measuring the distance between these two points using the Measure line tool gives me 10'000 meters, which is correct.

I convert the layer with a EPSG:3857 CRS using Save as in the layers panel. I open this new layer and measure the distance between the two points again and the Measure line tool gives me 13'000 meters.

Why? I don't understand. It means that I can not use the buffer tool on my EPSG:3857 layer because it gives me incorrect values and I have to use meters units.

Same issue with geopandas. Here are two screenshots to illustrate the problem.

My two points on Google Map (8,92km)

My two points on my Jupyter notebook, with a buffer equals two the distance displayed on Google Map (8920 meters)

The distances don't match! What am I doing wrong?


texts:
I'm using QGIS 2.18 Las Palmas, and I'm trying to build up a Distance Matrix of a series of points for a Mantel test. Specifically, I have a series of points and I want to calculate the distance of each of them to all the other ones. Now, when I use the 'Distance Matrix' tool in QGIS it gives me a three column table with all the distances, in this way:

InputID,C,254   TargetID,C,254  Distance,C,254
1               1               0
1               2               16.2788206
1               3               19.92485885
1               4               36.138622
1               5               69.23149572

This kind of output file contains all the information I need, but in the wrong format. What I'm looking for is a triangular distance matrix, as I need it as an input file for a Mantel test, like this one:

    1           2           3           4           5

1       0
2       16.2788206  0
3       19.92485885 23.84674    0
4       36.138622   29.4846735  55.68871    0
5       69.23149572 45.22356    66.254128   78.6542212  0

If it is possible to do this with QGIS, how do I do it?


texts: Speeding up test debugging when loading large module (ArcPy) as part of program?

texts:
I'm using PyScripter to do some python geocoding automation and I'm coming across the situation where I am waiting several tens of seconds for a large module to be imported.

This is very inconvenient since it makes editing and debugging much slower than could be.

I'm hoping there would be a workaround where this large module could be saved in memory so that after the first run, consecutive re-runs wouldn't need to reload the module.

Does such a workaround exist?

I found this Importing arcpy on application launch to speed up subsequent arcpy methods.  However it doesn't seem to directly address what I'm experiencing, or perhaps I do not understand the answer.

If someone could clarify the answer or suggest a workaround that would be great!

The code in question is this:

import arcpy

#use arcpy stuff


texts: Opening an modis hdf file in Envi and Calibration

texts:
I have downloaded an modis hdf file.it is l1 and its name is  MOD021KM ..I want to convert it to radiance.so I performed this steps in envi

Open>EOS>MODIS
selected modis hdf file.then I tried to run Radiometric Calibration but when I select the image it says image must have gain and offset for each band.I tried to opent image from open>generic format>hdf but there is not any hdf format there.

Anyway.I think envi does not recognize gains and offsets.I have checked metadata of opened file and the table of gains and offsets is empty.

I am using Envi 5.1

Can you please help me to fix this problem?


texts: Clear selection - argument as a layer name or layer object?

texts:
I am writing a script tool that takes multiple layers, selects every item from each layer, zooms to it and exports the image. After the last row in layer, the selection is cleared and new layer comes in a loop.

mxd = mapping.MapDocument("Current")
df = arcpy.mapping.ListDataFrames(mxd)[0]
arcpy.env.overwriteOutput = True

for lyr in arcpy.mapping.ListLayers(mxd):
    if lyr.name == "that_layer":
        oid_fieldname = arcpy.Describe(lyr).featureClass.OIDFieldName

        with arcpy.da.SearchCursor(lyr.name, ["OID@","adr_for"]) as cursor:
            for row in cursor:
                expression = """{} = {}""".format(oid_fieldname, row[0])
                arcpy.SelectLayerByAttribute_management (lyr.name, "NEW_SELECTION", expression)
                df.zoomToSelectedFeatures()
                arcpy.mapping.ExportToPNG(mxd, out_folder + "\\" + lyr.name + "\\" + row[1][8:] + "#"+ str(row[0]) +'.png', resolution = 300)

        arcpy.SelectLayerByAttribute_management(lyr.name, "CLEAR_SELECTION")

    #here go few elifs

Everything goes well, however one of the layers responds only to

arcpy.SelectLayerByAttribute_management(lyr.name, "CLEAR_SELECTION")

while given layer name (as in above example), in opposite to the rest where first argument given is a layer object:

arcpy.SelectLayerByAttribute_management(lyr, "CLEAR_SELECTION")

I cannot find any differences between layers - all contain polygon feature classes stored in the same personal geodatabase. All are in the Table of Content. I don't want to risk that in other project, magic "something" will be different and lyr.name will not work.

In the ESRI example layer object is used, but this particular layer does not respond for that.

What causes that behavior and what is the proper way to refer to the layer in this case?


texts:
I have a LASCatalog that is producing some funky outputs so I ran a lascheck() on it. The output tells me that catalog has inconsistent number of VLR, but what does this mean? I know that variable length records are part of the las file structure, but what problems could come from having an inconsistent number of them in my catalog?

Output from lascheck():

> lascheck(chunk2)

 Checking headers consistency
 - Checking file version consistency... ✓
 - Checking scale consistency... ✓
 - Checking offset consistency...
    ⚠ Inconsistent offsets
 - Checking point type consistency... ✓
 - Checking VLR consistency...
    ✗ Inconsistent number of VLR
 - Checking CRS consistency... ✓
 Checking the headers
 - Checking scale factor validity... ✓
 - Checking Point Data Format ID validity... ✓
 Checking preprocessing already done
 - Checking negative outliers... ✓
 - Checking normalization... no
 Checking the geometry
 - Checking overlapping tiles... ✓
 - Checking point indexation... no

Some tiles have the following VLRs:

Variable length records:
   Variable length record 1 of 3
       Description: by LAStools of rapidlasso GmbH
       Tags:
          Key 1024 value 1
          Key 3072 value 26918
          Key 3076 value 9001
          Key 4099 value 9001
   Variable length record 2 of 3
       Description: NIIRS10 Tile Index
   Variable length record 3 of 3
       Description: NIIRS10 Timestamp
       Extra Bytes Description:

While others only have:

Variable length records:
   Variable length record 1 of 3
       Description: by LAStools of rapidlasso GmbH
       Tags:
          Key 1024 value 1
          Key 3072 value 26918
          Key 3076 value 9001
          Key 4099 value 9001

LASCatalog before processing:

LASCatalog after normalization showing strange results:

Figure from spplot(chunk2_ctg, "Number.of.variable.length.records"):


texts: Adding building footprints to DEM in ArcGIS Desktop?

texts:
I'm currently working on a project wherein I am trying to assess the feasibility of creating a DSM from a DEM by adding building footprints to it (LiDAR manipulation is something I'd like to avoid for now). I have a DEM in ESRI Grid format, and another ESRI Grid raster that I created from a building footprint polygon shapefile (set the z-value to 20m for all buildings, then used Polygon to Raster conversion). Both rasters have the same projection and cell size. I would like to take the 20m values for the buildings, and add them on top of the elevation values present in the DEM (vs. simply replacing them). I will be using the resulting DSM for flow path modeling, so I need a continuous raster surface.

I thought to use Raster Calculator to combine the two, but I'm unsure if I am using the correct syntax. Simply trying to add DEM + Footprint Raster didn't work, and doesn't even create a new raster.

This was my syntax

"BldFootClp" + "CityDEMClip"

This is what my layers look like.


texts:
After several intents to reinstall everything (even qgis3.8.2 zanzibar), and verify the path (several times) for the lasstools plugin, still cannot use it, as consequence of the error message the system cannot find the specified path. I did doble check the path in my laptop, but still nothing works. Here is the error notice from the console:

QGIS version: 3.8.2-Zanzibar
QGIS code revision: 4470baa1a3
Qt version: 5.11.2
GDAL version: 2.4.1
GEOS version: 3.7.2-CAPI-1.11.0 b55d2125
PROJ version: Rel. 5.2.0, September 15th, 2018
Processing algorithm…
Algorithm 'lasground' starting…
Input parameters:
{ 'ADDITIONAL_OPTIONS' : '', 'BY_FLIGHTLINE' : False, 'CPU64' : False, 'GRANULARITY' : 1, 'GUI' : False, 'HORIZONTAL_FEET' : False, 'IGNORE_CLASS1' : 0, 'INPUT_LASLAZ' : 'C:\\lidar\\lazlid\\700515NW.las', 'NO_BULGE' : False, 'OUTPUT_LASLAZ' : 'C:/lidar/lazlid/70015.las', 'TERRAIN' : 4, 'VERBOSE' : True, 'VERTICAL_FEET' : False }

LAStools command line
C:/LAStools\bin\lasground -v -i "C:\lidar\lazlid\700515NW.las" -city -o "C:/lidar/lazlid/70015.las"
LAStools console output
The system cannot find the path specified.

Execution completed in 0.07 seconds
Results:
{'': None}

Any guess?


texts: What can I use to read MapInfo files in Java?

texts:
I want to read MapInfo files using Java.
 Before I used geotools for this but now I have problems with the GDAL library in Win7x64. Can i find geotools or gdal only for Java? To it not usede gdal's dlls only Java classes.
 I found the MITAB library but its only for C++.
Does anyone know something that can help me?

UPADTE

Another question: i read geotools quickstart and dont se any words about Gdal. But when i try use geotools i get error:

Exception in thread "main" java.lang.NoClassDefFoundError: org/gdal/ogr/ogr
at org.geotools.data.ogr.OGRDataStore.<clinit>(OGRDataStore.java:169)

When i add gdal.jre in build path all works fine.
Whats wrong?


texts: What are the advantages/disadvantages of populating a PostGIS database with Raster information?

texts:
...as opposed to keeping the rasters as locally stored GEOTIFFs?

I have a large project that will require making available approximetely 50,000, 10mb GEOTIFFs as WMS to an OpenLayers based web map. Only 1 WMS will be displayed at any given time so performance is a non-issue. Considering MapServer and GeoServer for this application.

I would like to know what are the advantages and disadvantages of storing this type of data. Will updating a PostGIS database when the raster data changes be more difficult, rather than replacing a single (or series of) raster files. WMS load times drastically improved with PostGIS?


texts: Hillshading all tiles by default in QGIS?

texts:
I've downloaded a zip of img files from the TN Lidar page, but see no useable tilemap to know which tiles I actually want to look at (the whole zip was over 4GB amd 600 tiles). I'm only wanting to look at terrain details, so just need to hillshade the tiles. I only seem to be able to choose one at a time and manually check off hillshade.

Is there a way to set styling for all layers I have loaded at once?


texts: Algorithm for clustering points into pairs (exactly two) by proximity without re-use

texts:
I have a reasonably large list of centroids that I want to cluster into groups of two by proximity (minimizing proximity).

I've explored k-means, which does cluster them by proximity, but the count of members in each group varies.  With k-means you set a number of clusters, not a number of members in each cluster.

The nearest-neighbor problem solves this issue for two items from the set, but not against the entire data set.

K-nearest neighbors seems to break them into groups of N, but it appears to allow for points to be reused.  In my scenario there can be no overlap.

Is there a particular algorithm, or suite of algorithms designed to address this?  I'm pretty handy when I know what I'm working against, but I don't have a good sense of how to approach the problem.

To add more about the context and what we're trying to solve:

The points represent a number of sites throughout the USA.  Each of these sites is a competitor (supply).  Independently, we've aggregated demand (from census data, etc).  We want to average the nearest pairs so that we can use the aggregated supply when calculating our supply/demand indexes for a given spatial extent (defined by the demand polygons).

We need to use at least two points so that individual data from a given site is obscured.  This is a licensing/privacy requirement.  We would otherwise analyze every point individually.  We don't want to use more than two, because that further obscures the data.  By using two, we adhere to licensing requirements, while minimizing the effect of averaging across a cluster.


texts: How can I Enable Google Earth 3D features in QGIS?

texts: I want to ask if there is any plugin or any other way where I can use Google Earth features in QGIS. When I'm saying Google Earth features I mean to be able to see roads, towns, etc and to be able to use the 3D option (Street view tool) I have already installed the OpenLayers Plugin but I did not find something to do it.


texts: How to work with OGR and GDAL in GP Service in ArcGIS for Server 10.1

texts:
I have got a lot troubles with GDAL and OGR in the context of ArcGIS Server on 64 bit server machine.
I can run script that use GDAL in ArcGIS as Tool and in the debug mode in PyScripter, then publish script as GP Service. I know that ArcGIS for Dekstop use 32 bit version of Python, and Server use 64 bit. I can see both version in Python27 folder. So in the desktop mode (32 bit) everything is Ok. But when I use submit job by REST I get an ERROR

  esriJobMessageTypeError:
  Traceback (most recent call last):
  File "C:\arcgisserver\directories\arcgissystem\arcgisinput\project\Analize.GPServer\extracted\v101\Analize.py",
  line 14, in from osgeo import gdal File "C:\Python27\ArcGISx6410.1\lib\site-packages\osgeo__init__.py",
  line 21, in _gdal = swig_import_helper() File "C:\Python27\ArcGISx6410.1\lib\site-packages\osgeo__init__.py",
  line 17, in swig_import_helper _mod = imp.load_module('_gdal', fp, pathname, description)
  ImportError: DLL load failed: Nie można odnaleźć określonego modułu.

I know that the reason is GDAL library. How to configure environment on ArcGIS Server machine to force to use 32 bit GDAL in the case of ArcGIS for Desktop and 64 bit for ArcGIS Server?

I installed:
GDAL for Desktop 32 bit according this solution How to run OGR in ArcGIS Scripts?
GDAL for ArcGIS 64 bit the same packages but for 64 bit

When I change PATH to new GDAL 64 bit script in ArcGIS and Service does not work

If PATH is set to GDAL 32 bit script works, but I still get error by REST. It seems that ArcGIS Server does not see the GDAL library on the server machine.


texts: Why are These KML Polygons Not Filled With Assigned Colour in Google Earth?

texts:
In ArcGIS Desktop 10.2.2 I converted a Feature Class into a KML using "Layer to KML (Conversion)".

The polygons in that KML look like this in Google Earth.

Why are (most of) the polygons rendering "oddly", that is, lacking a solid fill colour, and often mostly entirely absent?

@MichaelStimson provided the solution (to run Repair Geometry tool).  See results below.


texts: Oracle - find N neighbors (points) of a point within a radio

texts:
I have this table, there it goes a sample output to show that it is a point

select SDO_UTIL.TO_WKTGEOMETRY(geo)
from places
where id=1
_______________________________________________________

POINT (-88.706148 32.374127)

Given a point (i.e. latitude: 32.378097, longitude:-88.709356) I need to get the N neighbors that are around that point in a radio distance of 1 kilometer (0.6 miles)

So far I have this

SELECT id, description, latitude, longitude, SDO_NN_DISTANCE(1) dist
FROM places P
WHERE sdo_nn (geoloc,  MDSYS.SDO_GEOMETRY(1, 8307, MDSYS.SDO_POINT_TYPE( -88.709356 , 32.378097, NULL), NULL, NULL), 'distance=1 unit=KM', 1) = 'TRUE'

But, that query does not return points from places that are near, like 400 meters

Is there something I'm missing?


texts:
I was trying to use the lines to polygons tool in qGIS and it failed and gave me the following error:

  TypeError: setExteriorRing(self, ring: QgsCurve): argument 1 has
  unexpected type 'QgsPolygon

It seems to do with QgsPolygons being somewhere they shouldn't be(?) Can I just remove the offending polygons somehow? How can I fix this?


texts: downloading waypoint layer from garmin gps 60CSx to QGIS 1.8.0Lisboa

texts:
I want to download waypoint from my Garmin GPS but it responded

  the layer GPS.gpx?type=waypoint is not a valid layer and can not be
  added to the map

whereas my colleagues do the same without a problem.


texts: In JOSM, adjusting one line to follow another

texts: I would like to adjust a rudimentary drawn area (neighborhood) so that it uses the same nodes as the nearby streets. Doing it by hand seems very combersome - I use "a" or a nearby plus symbol to create a new node on the way, then I select it, and shift+select nearby street node, and hit "m" (merge).  How do I select "bad" line, and make it auto-join all the nearby nodes along the way, possibly with clicks or a keystroke.  I know about "follow", but that only works for a new drawing.


texts:
I would like to be able to script some things that are not exposed via arcgisscripting or ArcPy.

How do I access ArcObjects from Python?


texts: Geocoder's autocomplete only shows 26px in Leaflet?

texts:
I've got a map on GitHub (I have removed everything that isn't needed). The problem is that the geocoding widget doesn't show the search results completely:

I guess this happens due to a custom zoom bar control that includes a Zoom Home function which controls the width to 26px. Removing the.leaflet-bar a part leads to smaller icons but to fully result display.

Does anyone have a hint on how to have big icons and full size result display?


texts: When topological editing how can I automatically snap lines to lines close to them?

texts:
Using ArcMap 10.3, advanced license, I have a data set of road lines and another data set of road surface conditions, and I need to apply the road conditions to the road line as an attribute where it overlaps.  The problem is the lines do not overlap.

I have built a topology for it and have it identifying where one line does not get covered by the other line, but unfortunately that is most of the town.  I would therefore like to use automatic error correction and be able to go through the roads and let ArcMap do the spatial adjustment to align the conditions to where the road center line is.  I gave the center line a topology importance on 1 so it would not move and the road condition 5 so it would move.

However, when adding the topology and verifying it, then I choose the error correction tool and right-click on the error in the list, or on the error line, but I do not have a "fix error" type option.  Where do I find the trigger to have ArcMap move the one line onto the other?

I have also tried to use the Integrate geoprocess but it keeps failing with little info as to why.  I did remove the data from the topology before trying that.

Picture of data is below.


texts:
I want to toggle the Leaflet layer control when the user clicks in a button.

var layerControl = L.control.layers(
  {
    '': OSM,
    '': Esri_WorldImagery
  }, null, {collapsed: false}
).addTo(map);

<button onclick="toggleLayerControl()">toggle control</button>

function toggleLayerControl(){
   if(layerControl.display=="block"){
      layerControl.display = "none";
   }else{
      layerControl.display = "block";
   }
}

What's the simplest way of doing that?


texts: ArcPy can I delete the 9th column if I don't know its field name?

texts: Is there a way to delete a field based on its number, rather than field name? For example, if the field name is FID_myshp I would do arcpy.DeleteField_management(newshp, ['FID_myshp']) but the name of the field to delete changes, and I would just rather delete the 9th column because I am sure it will always be the 9th column.


texts: Earth Engine: issues with NAIP image collection

texts:
I am trying to visualize the NAIP data available by direct import in Earth Engine. I am unable to use the simple print method: Collection query aborted after accumulating over 5000 elements.

I suppose the problem comes from the fact that the NAIP is an image collection containing images with different bands over time, but not sure how to solve this, and how to visualize the details/information on the discussion?

var NAIP = ee.ImageCollection('USDA/NAIP/DOQQ')

// problem 1:
print(NAIP, "NAIP")


texts: What command do I use to do map algebra on rasters from the QGIS python console?


texts:
I have one image that I filtered it in MATLAB and finally, I saved it in tiff format from MATLAB. My image is RGB. it has three bands;
When I open this image in QGIS, and then go to ‘new print composer’, I do not know why I cannot add coordinate to it. I change ‘CRT’ but still did not work or some black area cover my image.
Do you have any idea?


texts: TauDEM tools in QGIS 2.10.1 Pisa (Ubuntu 15) are not in the Processing Toolbox

texts:
My TauDEM tools have disappeared in my latest version of QGIS.

Any ideas how to get them back? They used to be in the processing Tab.

Using
QGIS 2.10.1
TauDEM
Ubuntu 14 stable

My other tools are there (Gdal, grass, etc)

This appears to be the instructions but I am not really linux strong enough to follow them.

For example, Number 1. Uninstall openmpi-common (causes conflicts) and install mpich2.


texts:
What are the alternative vector tiles pipeline tools (producing -> rendering) to mapbox-gl as it is no longer open-source starting from its version 2.0.0?

Source: https://github.com/mapbox/mapbox-gl-js/blob/main/CHANGELOG.md#200


texts: NumPyArrayToFeatureClass gives RuntimeError: create table

texts:
I've tried to change polygon feautre  to point by using tools:

FeatureClassToNumPyArray
NumPyArrayToFeatureClass.
My code as below:

    import arcpy
    import numpy as np
    in_file='C:/ios/OneDrive/siatki_obliczeniowe/siatki.gdb/siatka_005'
    a = arcpy.da.FeatureClassToNumPyArray(in_file, ["*"], explode_to_points=False)
    out_file='C:/ios/OneDrive/siatki_obliczeniowe/siatki.gdb/xxxx3'

    arcpy.da.NumPyArrayToFeatureClass(a, out_file,["shape"])

I got an error:

Traceback (most recent call last):
  File "C:/ios/OneDrive/python/feature_to_point.py", line 7, in <module>
    arcpy.da.NumPyArrayToFeatureClass(a, out_file,["shape"])
RuntimeError: create table

I had changed out_file to shp file and it works but it is not solution for me. I would like to be able write to .gdb.
Any idea how to solve this problem?

I modified my code as below:

import arcpy
import numpy as np
in_file='C:/ios/OneDrive/siatki_obliczeniowe/siatki.gdb/strefy_miejskie_siatka_01'
flds = ['OID@', 'SHAPE@X', 'SHAPE@Y']
other_flds = [i.name for i in arcpy.ListFields(in_file) if i.type not in ('OID', 'Geometry')]
flds.extend(other_flds)
a = arcpy.da.FeatureClassToNumPyArray(in_file, flds, explode_to_points=False)
out_file='C:/ios/OneDrive/siatki_obliczeniowe/siatki.gdb/axxxyz9'
arcpy.da.NumPyArrayToFeatureClass(a, out_file,['SHAPE@X', 'SHAPE@Y'])

Unfortunately I got the same error again.


texts:
I am attempting to join 2 layers from the same ArcGIS Online service together, using ESRI's limited sample for MapImageLayer.  I would like to join the first [ID: 0] and second [ID: 1] layers in my service together using the "FLOORID" field, which both layers have.

var floors = new MapImageLayer({
    url: baseURL, //my service URL
    title: "All Layers",
    sublayers: [{
            title: "InteriorSpace",
            id: 1,
            opacity: 0.75,
            source: {
                type: "data-layer",
                dataSource: {
                    type: "join-table",
                    leftTableSource: {
                        type: "map-layer",
                        mapLayerId: 0
                    },
                    rightTableSource: {
                        type: "data-layer",
                        dataSource: {
                            type: "map-layer",
                            mapLayerId: 1
                        }
                    }
                },
                leftTableKey: "FLOORID",
                rightTableKey: "FLOORID",
                joinType: "left-inner-join"
            }
        }]
});

The console doesn't log any errors, but I do not get any tiles displayed.  What am I missing?


texts: Build static raster tile overview cache from WMS server

texts:
Often WMS servers will not give you anything but a transparent image at lower zooms. When using WMS servers on the web this can often lead to a poor user experience.
I would like to build a local cache of the missing lower level zooms of the remote WMS server, so that we can also display map tiles on the web when zoomed out.  I do not control the remote server so we can assume it doesn't have any special features.
What is an easy way to pull down the lowest resolution available from the WMS and then building a pyramid of the missing lower zoom level raster tiles as static files (preferably .mbtiles)?


texts: Saga:processing resample extent in python editor QGIS

texts:
I would like to resample many geotiff layers to another geotiff layer's extent and pixel size. For this I would like to use the QGIS python console (editor). When I run my script it raises the following error:

File "/Applications/QGIS.app/Contents/MacOS/../Resources/python/plugins/processing/core/parameters.py", line 165, in setValue
    tokens = text.split(',')
AttributeError: 'list' object has no attribute 'split'

When I run the saga resampling from the processing toolbox

enter code here

import processing
input    ="/Users/rutgerhofste/GitHub/AlgeriaTest01/Data/Raster/CHIRPSMewina/CHIRPSMewinaY2012M01.tif"
keep_type = True
target = 0
scale_up_method =  0
scale_down_method = 0

user_xmin = -17.0999999594
user_xmax = 35.7000007799
user_ymin = 14.6999995737
user_ymax = 37.

output_extent = [user_xmin,user_xmax,user_ymin,user_ymax]

user_size = 0.050000
user_grid = "/Users/rutgerhofste/Desktop/werkmap/testgrid.tif"

processing.runalg('saga:resampling', input, keep_type, target, scale_up_method, scale_down_method, output_extent, user_size, user_grid)

When I run the algorithm from the processing toolbox everything works fine


texts:
I am new in the area. I am working with GRASS in R. Now I am trying to load an vector corine land cover data and my before calculated catchment areas of the DEM and the created Stream network. I understand that my error tells me that my input paramters are multiple values, but I am not sure if I need to change a setting or if the code is wrong or need to be changed. I can't find an problem solving answer as I do not understand the GRASS key support on this subject. Anyone understands my question an can help me?

Problem in R with GRASS: The input layer has multiple values can I change the code to read all of the values or do I need to create a layer with all the information before and if so can someone please help me with that

for(i in 1:length(locID))
+ print(i)
[1] 1
[1] 2
[1] 3
[1] 4
[1] 5
[1] 6
[1] 7
[1] 8
[1] 9
[1] 10
[1] 11
[1] 12
[1] 13
[1] 14
[1] 15
[1] 16
[1] 17
[1] 18
[1] 19
[1] 20
[1] 21
[1] 22
> cname <- paste0("sites_catchm_", locID[1:22])
> if(cname %in% rast){
+     if(!cname %in% vect){
+
+       # change raster to vector
+       execGRASS("r.to.vect", flags = c("overwrite","v"),
+                 parameters = list(
+                   input = cname,
+                   output = paste0(cname,"_v"),
+                   type = "area"
+                 ))
+       # calculate area
+       execGRASS("v.db.addcolumn",
+                 parameters = list(
+                   map = paste0(cname,"_v"),
+                   columns = "area double precision"
+                 ))
+       execGRASS("v.to.db",flags = "quiet",
+                 parameters = list(
+                   map = paste0(cname,"_v"),
+                   option = "area",
+                   columns = "area"
+                 ))
+       tarea <- sum(as.numeric(execGRASS("v.db.select",flags = "quiet",
+                                         parameters = list(
+                                           map = paste0(cname,"_v"),
+                                           columns = "area"
+                                         ), intern = T)[-1]))
+     }
+
+     # intersect
+     execGRASS("v.overlay", flags = c("overwrite","quiet"),
+               parameters = list(
+                 ainput = paste0("clc_"),
+                 binput = paste0(cname,"_v"),
+                 operator = "and",
+                 output = "intersect_clc",
+                 olayer = "1,0,0"
+               ))
+     # calculate areas from the features
+     execGRASS("v.db.addcolumn",
+               parameters = list(
+                 map = "intersect_clc",
+                 columns = "area double precision"
+               ))
+     execGRASS("v.to.db",flags = "quiet",
+               parameters = list(
+                 map = "intersect_clc",
+                 option = "area",
+                 columns = "area"
+               ))
+     # read data from table
+     a <- execGRASS("v.db.select",flags = "quiet",
+                    parameters = list(
+                      map = "intersect_clc",
+                      columns = paste0("CODE_12",",area")
+                    ), intern = T)[-1]
+     a <- do.call(rbind,strsplit(a,split = '\\|'))
+     a <- data.frame(a,  stringsAsFactors = F)
+     a <- apply(a,2,as.numeric)
+     if(class(a) == "matrix"){
+       a <- round(tapply(a[,2],a[,1], sum),2)
+     } else {
+       nam <- a[1]
+       a <- a[2]
+       names(a) <- nam
+     }
+   }
Error in doGRASS(cmd, flags = flags, ..., parameters = parameters, echoCmd = echoCmd,  :
  Parameter <input> has multiple values
In addition: Warning messages:
1: In if (cname %in% rast) { :
  the condition has length > 1 and only the first element will be used
2: In if (!cname %in% vect) { :
  the condition has length > 1 and only the first element will be used


texts: Troubleshooting failed ArcGIS Server Geoprocesses?

texts:
There are a few geoprocesses on a recently restored ArcGIS Server 10.0 SP5 that are behaving strangely. If we turn on 'show messages' and try to run the process from the rest interface, the process will report 'failed due to crash or termination of the server object'. This happens at different points in the script and sometimes it will even return messages after the job has officially failed. Every now and then, the final message is 'the job has succeeded' even though the 'Job status' field has been set to 'failed' half way through.

Some information about the setup:
Data and server directories are accessed via UNC path
Authentication is handled via domain accounts
Geoprocesses are backed by Python scripts

I've been searching around for file access issues in both the data store and server output directories. Any advice on best practices to narrow down what may be going on?


texts: DMS grid when SRC is not in decimal degrees using QGIS print composer

texts:
I have a QGIS project where the layers are in projected SRC, and when I go to the print composer to draw the grid, the units are in millons of meters, and I want them in decimal degrees.

Is it possible to draw the grid´s coordinates in DMS or DD without reproject the layers of the project?


texts: ArcGIS 10.0/10.1 - arcpy - automate projection tool with different datum or projection

texts:
I am getting stuck with a task that is slightly different from the posts I found so far on this same website.

Assumption: all layers have a defined coordinate system

Problem: I would like to automate a python script (using arcGIS 10.1 but 10.0 is fine as well) that projects a list of vector layers (some of which with different coordinate systems) to a common coordinate system, given by a target vector layer. The main issue I would need help with is how to automate something this:

1) check if the coordinate system of the current vector layer is the same as the target layer

2) if not, check if the datum is the same

3) if it's the same datum, then project the current layer to the target layer without any datum transformation

4) if, instead, they do not have the same datum, pick a pre-defined transformation and include that in the projection function

Any help would be much appreciated! thanks!

Here is my sample code:

import arcpy
from arcpy import env

arcpy.gp.overwriteOutput = True
env.workspace = 'C:\\Temp'
Target_lyr = 'water_Broward.shp'

try:
 desc = arcpy.Describe(Target_lyr)
 TargetSR = desc.spatialReference

 for inputFC in arcpy.ListFeatureClasses('*'):

  if inputFC != Target_lyr:
   inputFC_SR = arcpy.Describe(inputFC).spatialReference

  if inputFC_SR.Name == TargetSR.Name:
   continue
  else:
   #check if current layer and target layer are the same type
   #thus if they are both 'projected' or 'geographic'
   if inputFC_SR.type == TargetSR.type:

    #HELP STARTING FROM HERE...
    #Check if both have same datum..

    #if they ARE, then
    #outFC = inputFC + '_prj.shp'
    #arcpy.Project_management(inputFC, outFC, TargetCS)

    #if datum is NOT the same use
    #transformation = for example "WGS_1972_To_WGS_1984_1"
    #arcpy.Project_management(inputFC, outFC, TargetCS,transformation)

   #if current layer and target layer are NOT same kind
   else:
    #HELP
    #repeat as above same steps as above to check the datum

except:
 print arcpy.GetMessages()


texts:
I have a map showing administrative boundaries (available here).  I've loaded it into R and have attempted to run gTouches (from rgeos) to determine which areas are adjacent to which, but receive an error:

library(rgeos)
library(rgdal)

> Districts <- readOGR(... )
OGR data source with driver: ESRI Shapefile
with 9235 features and 16 fields
Feature type: wkbPolygon with 2 dimensions

> adj <- gTouches(Districts)
Error in RGEOSBinPredFunc(spgeom1, spgeom2, byid, "rgeos_touches") :
  TopologyException: side location conflict at 34.353274826787519 -15.739538287828307

After Googling that error I found a lot of advice about "bad geometries."  Following that I took the file into QGIS (following advice here from Faunalia's "Bad, bad polygon" document) and ArcMap (using Python code here) and run the geometry checking/repair tools. Neither program found an error.

I've zoomed my map into the offending coordinates but don't see anything that appears wrong--it's a location where two lines come together near another vertex.

How do I repair the geometry, ignore that particular region, or some other R function that might work better than gTouches?


texts: Are there any general rules of thumb for selecting a transformation for Georeferencing?

texts:
There are several transformations that can be used when georeferencing a raster:

1st order Polynomial (affine)
2nd order Polynomial
3rd order Polynomial
Spline... etc

Are there any specific rules or rules of thumb, about which transformation should be used with a particular raster?

For example a  particular transformation should generally be used with aerial photographs & another should be used with satellite imagery while a third should be used while digitizing maps?

Is there any book which explain the mathematics behind this in detail?


texts: Using DE-9IM patterns in python plugin for QGIS?

texts:
I'd like to compare geometries using DE-9IM patterns in my QGIS plugin written in Python. I've tried to instantiate QgsGeos but it seems there's no such module in qgis.core, as reported by the following Python exception :

> from qgis.core import QgsGeos
Traceback (most recent call last):
  File "<input>", line 1, in <module>
ImportError: cannot import name QgsGeos

I could also use shapely but I don't know if this module is included in all QGIS distributions.

Any suggestion ?


texts: Relation Reference to simulate filtered drop-down lists

texts:
I have a layer of Land Uses that has a 3-tiered classification system (basically so that it can be broken down for different uses by other people).

I have managed to create Project Relations and use the Relation Reference widget to achieve the drop-downs in the forms -

Attribute Look-up Table:

Project Relations:

Relation Reference Widget:

Attribute Pop-Up:

(1) Is this the correct/intended methodology to achieve filtered drop-downs like this (where one attributes list of available values is dependant on a previous field)? Is just seems more logical to be using the filter in the Value Relation widget, except I can't get it to filter based on another field...!

(2) Is there a simple way in the pop-up to remove the table icon that brings up another form (assuming to add additional 'child' records to my parent feature)?
As that is not what I am doing in this case, it's misleading to users to have it there.


texts: Defining spatial operators (Equals, Disjoint, Touches, Within, Overlaps, Crosses, Intersects, Contains, Relate)?

texts:
Ok, I have to admit, I don't really delve too deep into the semantics of the spatial operators list below, I'm more of a user, I write software to drive the SQL most of the time and don't really think about it too much.

However, I have a situation where a spatial operation I'm doing is going dog slow, on a table that has great indexes on and has never caused a problem before now.

As a result, I'm trying to figure out which of the following is my best bet to use, that will give me the best performance, however, sometimes I find my searches don't return anything when I know they should, and sometimes return stuff when they shouldn't.

SO I'm reaching out, to ask if some one can give me the actual definition, such as criteria to match, fully within enclosing geometry, or within outer ring etc etc of each of the following operations:

Equals(Geom1, Geom2):int
Disjoint(Geom1, Geom2):int
Touches(Geom1, Geom2):int
Within(Geom1, Geom2):int
Overlaps(Geom1, Geom2):int
Crosses(Geom1, Geom2):int
Intersects(Geom1, Geom2):int
Contains(Geom1, Geom2):int
Relate(Geom1, Geom2):int

If you can give a concrete example of when they will and will not match too, that would be great.

For reference, I'm searching on a table of Linestrings, using a rectangular polygon, rotated to be pointing in the direction of travel (It's part of a traffic managment app)

For further reference, here is a sample of the SQL the app generates to perform a search:

SELECT recordID,AsBinary(geometry) AS geometry,Distance(GeomFromText('POINT(-1.84101 54.85078)',4326), geometry) AS distanceFromGps FROM linegrid WHERE Intersects(GeomFromText('POLYGON ((-1.8413149820810311 54.850782468607292, -1.8409507853094111 54.850952257034713, -1.8408279009723911 54.850864894077496, -1.8411920982612455 54.850695105650068, -1.8413149820810311 54.850782468607292))'), geometry)

You can see in this example I'm using intersects trying the others as mentioned doesn't quite give the results I expect, hence why I'm asking the question.

Please note, I'm not particularly bothered about actual faster/other ways of doing this at the moment, BUT if you do happen to spot a better way do feel free to shout up :-)

For now all I'm mostly trying to understand is the differences between each of the operations.

As for doing it faster/better I intend to open a different question for that at a later date.

For more reference, the Spatial engine this is working against is Spatilite 2.3.0, and I can't upgrade it due to the devices it's running on not having a more up to date build for it (Windows CE)

I'm going to choose mintix'es answer here as that provided a super easy to understand grid and description of all the various operations, which answered my question.

I'm still faced with getting performance on the routines up, but I do now have a few strategies to try on that one, I may still seek advice here from the group with regards to that.

@Vince deserves some credit too, as that white paper is definitely the definitive answer, but it's also quite challenging to read (I had to read through it several times) , for anyone wanting the real I am on these operations The Clementini paper is the way to go, if you need an easy to understand overview while developing solutions however, the FME link should be open in another browser tab.


texts: How to calculate the height of buildings from satellite/aerial images using MATLAB

texts:
I want to calculate the height of buildings from satellite/aerial images using MATLAB automatically with no or least possible human interaction.  I am not using LIDAR or any DEM information.  Can anyone help me out in this?

My project is simple but very difficult to be completed. There are some GIS tools available for determining height of buildings but they are proprietary software and not automatic. What I am trying to do is to read satellite images in MATLAB and calculate the heights of buildings present in the images. There is no compulsion of what type of images I should use i.e I can use aerial or satellite images, normal or stereo images or oblique images.
I hope this information is sufficient to answer my question, I would be really happy to provide further clarification if needed.

Thanks in advance.


texts:
I'm working with a lot of shapes and raster files in one map document (.mxd) in ArcMap and I would like to send the project to another person that works with me.

I would like to know if there is a way to add all files from a map document to a geodatabase, and if map packages can be opened in other GIS software packages, such as QGIS.


texts: Seeking solution to using PostgreSQL/PostGIS views in ArcGIS for Desktop?

texts:
We are currently trying to 'clean up' our geodatabase. For example: we have 30 datasets of districts, each of a different year (1986-2016). To make it more visible, we want to combine them all into 1 dataset where we use a view for the most recent one. There are some issues when we create a view. First, we cannot set a primary key (we use PostgreSQL together with PostGIS) in a view. A user must specify a primary key each time he drags the view into ArcMap. Second, for some reason, when creating a view, ArcGIS recognizes 2 datasets: a table and a query layer (see image). I have no idea why. In PgAdmin it looks like one tableview. Catalog also crashes when trying to open the table version.

Is there a way to set and save a primary key in ArcGIS, since it cannot be done using PostgreSQL? And what causes ArcGIS to show two datasets of the same view, a table and a query layer? If not, are there any alternatives to these problems?

My table definition:

CREATE TABLE gemeentegrenzen
(
  objectid serial NOT NULL,
  gem_nr integer,
  gem_naam character varying(30),
  jaar smallint,
  prov_code character varying(10),
  prov_naam character varying(28),
  cbs_code character varying(6),
  shape geometry,
  CONSTRAINT gemeentegrenzen_pkey PRIMARY KEY (objectid),
  CONSTRAINT enforce_srid_shape CHECK (st_srid(shape) = 28992)
)
WITH (
  OIDS=FALSE
);
ALTER TABLE gemeentegrenzen
  OWNER TO geo;

-- Index: sidx_54890_8

-- DROP INDEX sidx_54890_8;

CREATE INDEX sidx_54890_8
  ON gemeentegrenzen
  USING gist
  (shape);

My view definition:

CREATE OR REPLACE VIEW gemeenten_actueel AS
 SELECT gemeentegrenzen.objectid,
    gemeentegrenzen.gem_nr,
    gemeentegrenzen.gem_naam,
    gemeentegrenzen.jaar,
    gemeentegrenzen.prov_code,
    gemeentegrenzen.prov_naam,
    gemeentegrenzen.cbs_code,
    st_transform(gemeentegrenzen.shape, 28992)::geometry(Geometry,28992) AS shape
   FROM gemeentegrenzen
  WHERE gemeentegrenzen.jaar = 2016;

ALTER TABLE gemeenten_actueel
  OWNER TO geo;

Both are stored in a PostgreSQL 9.4 database with postgis 2.0. ArcGIS is version 10.3.


texts: Adding Correct and Incorrect as attributes automatically to a Line layer in Mapinfo

texts:
I have a Polyline layer that I am adding to. Two of the columns are Correct and Incorrect. When I add a new Polyline I want these two columns to be filled out automatically with the Correct and Incorrect values.

Does anyone know how to done this using MAPINFO?


texts: Why are WMS layers not loading on new computer?

texts:
I've created a QGIS project that makes use of WMS layers.  I'm sharing the project with many other computers. However when I load the QGIS project on those computers, the wms layers do not load.  even once I've added the wms server to the WMS layer manager it still doesn't load.  the only way that seems to work is to actaully load the WMS layers again into the project....obviously a pain since I have several layers from various sources.

Any ideas?  Is there a proper methodology for doing this?

thanks!


texts:
I am wondering how to place images on a line. For example, instead of a dotted or dashed line. I could include a symbol of a cloud a character (e.g. |) repeated along the line. it cannot show the image along with line.

  var style = new ol.style.Style({
      image: new ol.style.Icon(({
                     opacity: 1,
                     size:20,
                     src: 'https://i.imgsafe.org/73d1273.png'
                     })),
      stroke: new ol.style.Stroke({
      color: 'black',
      width: 4,
      lineDash: [5],
     })
     });

I did a lot of research, but I couldn't find any results.
Does it support OpenLayers? How can it be done if it supports it?

Code sample:
http://jsfiddle.net/jguxq4j0/8/


texts:
I have an image in projection Pulkovo 1942 / CS63 zone X2. When I try to open it in QGIS, it fails to detect the projection and asks to select the projection manually. Does QGIS handle this family of projections?

$ gdalsrsinfo '/home/user/Downloads/2DJI_0983.tif'

PROJ.4 : +proj=tmerc +lat_0=0.0833333333333333 +lon_0=26.5 +k=1 +x_0=2300000 +y_0=0 +ellps=krass +towgs84=23.92,-141.27,-80.9,0,0.35,0.82,-0.12 +units=m +no_defs

OGC WKT :
PROJCS["Pulkovo 1942 / CS63 zone X2",
    GEOGCS["Pulkovo 1942",
        DATUM["Pulkovo_1942",
            SPHEROID["Krassowsky 1940",6378245,298.2999999999998,
                AUTHORITY["EPSG","7024"]],
            AUTHORITY["EPSG","6284"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4284"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0.0833333333333333],
    PARAMETER["central_meridian",26.5],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",2300000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AUTHORITY["EPSG","7826"]]

Also I managed to reproject the image into EPSG:4326 using gdalwarp, but I'm not sure whether the reprojection went right. Do GDAL and PROJ work correctly with this projection? Maybe there is an open source software using which I could open the image to assure the reprojection accuracy?

QGIS 2.18 Las Palmas, GDAL 2.1.0, PROJ 4.8.0


texts: Is it possible to style a vector layer in Leaflet with multiple "style layers"? For example, could a polygon layer be styled with multiple borders, of different weights/colours/opacity? OpenLayers supports this, and it would be very useful for more complex cartographic styling.


texts: v.clean what input feature type for polygon?

texts:
I have a polygon layer where the Topology Checker detects hundreds of overlaps and tiny gaps. So I want to use the v.clean tool to clean the topology. For cleaning, I selected the tools bpol, rmarea and rmdupl as recommended in How to Clean Topology Error in QGIS
.
But what I have to select as input feature type? There's no option for polygon.


texts: Arcpy Python Snippet for SQL "LIKE" expression problem

texts:
I am trying to use arcpy.Select_analysis in python code with ArcGIS for straightforward geoprocessing. I have a field called "other tags" in a vector dataset which I want to select any record which contains the word "rail".
I can get it working interactively by typing the following code into the query builder:
"other_tags"LIKE ('%rail%')

When I copy it as a python snippet it gives me:
arcpy.Select_analysis(in_features=input, out_feature_class=output, where_clause='"other_tags"LIKE ('%rail%')\n')

and when I try and run this code directly rather than through the interface it doesn't work. The error is:
TypeError: not all arguments converted during string formatting.

I am not clear how to amend the where_clause to make it work.


texts: Distinguishing Python version when using IDLE with ArcGIS Pro and ArcGIS 10.x for Desktop architecture?

texts:
I have just started to work with Python 3.4.1 which comes with ArcGIS Pro 1.0.0, on my laptop, which also has Python 2.7.8 which comes with ArcGIS 10.3 for Desktop installed, and I would like to minimize the chances of running one version of Python when I intended and needed to use the other.

In Windows Explorer, I now have two Python folders: C:\Python27 and C:\Python34 which is to be expected.

On the Windows Start Menu under All Programs > ArcGIS there is a Python 2.7 folder with the IDLE (Python GUI) that I often use to start IDLE, so starting it from here should be fine.  I tested and it does start 2.7.8.

I had expected to see Python 3.4 folder next to Python 2.7 or maybe under All Programs > ArcGIS > ArcGIS Pro but it is in its own folder directly under All Programs.  It contains another IDLE (Python GUI) that starts 3.4.1 as expected.

However, things become a little confusing when I look at the list of programs on the Start Menu because there are two IDLE (Python GUI) entries to choose from:

By right-clicking on each and using the General tab of their Properties I have renamed them to be called IDLE (Python GUI) 2.7 and IDLE (Python GUI) 3.4 respectively so once again there is no confusion.

The other way that I usually start IDLE is by right-clicking on a *.py file and then choosing Edit with IDLE - but the Python version of IDLE I will get is now ambiguous.  By testing I find that it is 3.4.1.

Is there a way that I could turn my one choice of Edit with IDLE into these two choices?

Edit with IDLE 2.7; and
Edit with IDLE 3.4

I never edit with PythonWin so if I am limited to two programs that can be associated with *.py files for editing then I am happy to sacrifice Edit with PythonWin for Edit with Python 3.4.


texts:
I have used GeoRSS files successfully in openlayers v2.x but I am not able to find how to handle these GeoRSS files in Openlayers v3.x.

Please help if anybody has successfully tested the usage of GeoRSS in Openlayers version 3.


texts: how do I create a feature class from a table?

texts:
I have a database table that has the column and types defined. But it's a database table. I need a feature class so I can add spatial functionality, like add a SHAPE field. How do I create a feature class from a database table?

Update - I'm starting from a database table. The destination is an existing feature class in my enterprise SDE. I need to take some fields from the table and build a value for a SHAPE field. I also need to be able to use my fieldmappings to map fields from the source to the destination; hence the need to get this data into a feature class instead of a table.


texts:
I have a terrain, that looks like the below:

The geotiff file can be found here
I run the geotiff file with the taudem, with the threshold value =10
1. PitRemove.exe Site6.tif
2. D8FlowDir.exe -p Site6.tif -sd8 Site6sd8.tif -fel Site6fel.tif
3. aread8.exe -p Site6p.tif -ad8 Site6ad8.tif
4. Gridnet.exe -p Site6p.tif -plen Site6plen.tif -tlen Site6tlen.tif -gord Site6gord.tif
5. PeukerDouglas.exe -fel Site6fel.tif -ss Site6ss.tif
6. aread8.exe -p Site6p.tif -wg Site6ss.tif -ad8 Site6ssa.tif
7.  Threshold.exe -ssa Site6ssa.tif -src Site6src.tif -thresh 10
8.  Streamnet.exe -fel Site6fel.tif  -p Site6p.tif  -ad8 Site6ad8.tif -src Site6src.tif -ord Site6ord.tif -tree Site6tree.dat -coord Site6coord.dat -net Site6net.shp -w Site6w.tif

I would expect that all of the original cells in the original terrain should be contributing to the flow path, because there is no edge contamination, ie: no drainage is inwards from the boundaries or areas with no data values for elevation. For this terrain, the elevation levels peak at the center and then going downwards radiating from the center.
But when I examine the final output, ie: site6w.tif, it seems that there are quite a lot of missing cells ( the contributing cells are marked in dark, and since they are nowhere representing all of the possible contributing cells, we can thus infer that there are lots of cells that are missing), as shown below:

TauDEM seems unable to decide for a lot of cells, where the flow contribution is going to.
I also find that with or without the -nc parameter, the outcome is almost the same, there are still a lot of missing contributing cells.
I've no idea why lots of contributing cells are missing, any suggestion?


texts: Setting lower limit to interpolation using ArcGIS Desktop

texts:
I have rainfall data that I want to interpolate using spline. I have done Spline Interpolation, only that the values produce negative predicted values which is basically wrong. There should be no negative rainfall value.

How can I set lower limit to zero or set a true value limit?

I am using ArcGIS Desktop 10.0.


texts: Calculating value of string using array_foreach

texts:
I am relatively new to this, but I have an attribute table field (SPECIES_PR) where some fish species are delimited with a comma and I would like to evaluate each species in the string, then assign the highest value I've set among the group for each rated species (ARGR, BLTR, BRST, BURB) to a field vul.
I've gotten this far, but cannot think how I can link this to a max value expression.

EXPRESSION:
array_foreach(string_to_array("SPECIES_PR",','),
map_get(map('ARGR',1,'BLTR',1,'BRST',2,'BURB',2),@element))


texts: How do you toggle between layers that have been loaded using the FileLayerLoad.js?

texts:
I'm trying to toggle between geojson layers that were loaded using the filelayerload javascript. Right now I'm basically stuck on how to append the geojson layers to one variable once they have been loaded using filelayerload. My thought is by appending them to one variable ('overlayMaps' variable in the code below), I can just use L.controls.layer to set up the toggling. Here's my code so far:

        //Main map layer:
    var mainmap =L.tileLayer(
        'http://otile{s}.mqcdn.com/tiles/1.0.0/map/{z}/{x}/{y}.png', {
        attribution: '&copy; '+mapLink+'. Tiles courtesy of '+mapquestLink+mapquestPic,
        maxZoom: 19,
        subdomains: '1234',
        }).addTo(map);

    var baseMaps = {
        "Base Layer": mainmap
    };
    var overlayMaps = {
    };

    L.Control.FileLayerLoad.LABEL = '<i class="fa fa-folder-open"></i>';

    //Loads external geojson file:
     var fileload = new L.Control.fileLayerLoad({
        fileSizeLimit: 5000,
        fitBounds: true, //MOVE THE CENTER OF THE SCREEN
        layerOptions: {
        onEachFeature: function(feature, layer) {
            layer.bindPopup( "Name: " + feature.properties.Name);
        }
        }
    }).addTo(map);

    L.control.layers(baseMaps,overlayMaps).addTo(map);


texts: Generalizing (simplifying) polygon shapefile or feature class in ArcMap while preserving topology?

texts:
I came across Generalizing polygon file while maintaining topology in QGIS? which asks for a QGIS solution to my exact question. I'm trying to do this in ArcMap. I start with a polygon shapefile (WGS84).

I verified by importing it into a file geodatabase feature dataset topology that it does not have any gaps or overlaps. Then I run the Simplify Polygon ArcToolbox tool on it.

The output has visible gaps and overlaps between some polygons. I added that layer to the topology and it shows me the rule-violating areas.

I'm kinda an ArcGIS Topology noob though - one thing that has me confused is the Topology Toolbar remains entirely disabled. And even if that wasn't a problem I'm not keen to have to repair gaps and overlaps one by one interactively. Even if there aren't that many issues this time, there will be next time. Ideally I'd come up with a simplification with topology preservation workflow that is entirely automated. Is there a better way?


texts: Can I avoid decimal conversion calling CostDistance

texts:
Our client asked us to incorporate blockages and delays for movement of boats and ships on the local waterways.  In our one process we are creating polygons of the distance a boat can move at top speed for up to 90 minutes.

The change required us to look at several things and we decided to set the raster to a time-based seconds per meter.  In our test sample I have a point located on the waterway.  I have time zones established for how far a boat will move with no slowdowns.

ie.  Boat speed = 40 knots
     @15 minutes  10 nm
     @30 minutes  20 nm
     @45 minutes  30 nm
     @60 minutes  40 nm
     @90 minutes  60 nm
To define the seconds/meter we compute 40 * 1.94384449 rounded and converted to integer.

this yields a default raster value of 78 which is the number of seconds it takes for the boat to move one meter.  And when we convert nautical miles to meters the maximum distance a boat can travel at full speed is 8667360 meters.

All of these numbers are within a decent realm for storing as integer values in a cost raster, yet when I call CostDistance() it takes 33 hours to run and the result shows the 8667360 in scientific notation.  (which means it is being stored as floating point in the raster)

So my question is What can I do to force CostDistance to remain integer so it will not run like a dog and take 33 hours?


texts: Scale color ramp filtered on common layer columns, in QGIS

texts:
I'm creating an atlas in QGIS showing bus lines, tax fare zones, and travel flow lines.

I have a layer with bus lines - one feature per bus line.

Moreover, I have a layer with travel flow lines - there are many features per bus line.

Both layers have a common attribute column describing the bus line.

What I want is to show the travel flows lines that corresponds to a certain bus line.

Hence, I created an atlas via the atlas feature, where the coverage layer is featuring the bus lines.

In the travel flow line layer properties I made a rule based symbology, in order to only showing the travel flow lines that corresponds to each bus line.

In the filter I put:  "Rute" = attribute( @atlas_feature , 'Rute')

Now I want to style each travel flow line feature double. The line width in accordance to the amount of travellers per flow line (called 'rejsende_a') but also styled from a blue color ramp, the higher amount of travellers, the darker blue.

Here is my problem. I want QGIS to style the color of the lines from a expression, as you can see from the picture below.

I need to find the minimum and maximum values of the travel flow lines to the corresponding bus line that is the current atlas feature.

So far I've tried several expressions. The latest one is like:

ramp_color('Blues', scale_linear(Rejsende_a, minimum(CASE When "Rute" = attribute(@atlas_feature, "Rute") Then  "Rejsende_a"  end),maximum((CASE When "Rute" = attribute(@atlas_feature, "Rute") Then  "Rejsende_a"  end)), 0,1))

Be aware that:
Rute = bus line    and    Rejsende_a = travellers

But apparently the syntax is wrong.

Can anybody help out here? Can you use 'CASE WHEN' in the expression string? Or how can you tell the expression to choose only the feature from the travel flow lines that has the "Rute" column in common?


texts: Multi-modal public transport Routing using pgrouting or opentripplanner

texts:
I need multi-modal public transport routing for my iPhone app and also to solve a TSP problem.

I have been researching and I am unable to make a decision whether I should go with pgrouting or install Open Trip Planner.. Is there a chart comparison, what each of them offers and what is different about these softwares?