texts:
I have a map with some walking- and bike-routes and popups with a few details and a picture.

Now I want to set a marker on the first vertex of a GeoJSON polyline, but I can't find out how.

I'm new to Leaflet/Mapbox, and I put my map together from code snippets.

The map, as it is now, is available here.


texts: Installing third party python libraries in ArcGIS location

texts:
I'm running ArcMap 10.2.2 for desktop and have Python 2.7 installed as part of the software package. I'm trying to get the pyexcel and xlutils libraries installed. I ran a script to find the path of the arcpy module and figured this is where I would have to copy other libraries to in order to import them into my Python scripts. This is after reading the following article by ESRI:

https://www.esri.com/arcgis-blog/products/arcgis-desktop/analytics/a-simple-approach-for-including-3rd-party-python-libraries-with-your-scripts/

I now have both the pyexcel and xlutils folders in Python27/ArcGIS10.2/Lib/site-packages but am getting the generic import error stating the module does not exist when trying to import the modules.

Do I have to run pip in order to make these libraries available> I don't have admin rights to my machine and am dependent on IT staff to do anything requiring access to program files.


texts: Preventing symbol rotation of point character symbols combined with line features in ArcMap?

texts:
I have a line feature that needs to have a marker symbol at set distances along it. I want the markers to NOT rotate with the line orientation but can not find out how.
I am working with arcmap 10,4

I created a new marker symbol in style manager like this.

I then used symbol property editor to create a line symbol with a dashed line and marker line and marker line symbol.

This is the result I get. But I need the R labels to always be rotated to the same angle and correctly orientated to the data frame.


texts: GRASS 7.2.1 r.denoise/mdenoise error Windows 10

texts:
So I am trying to use the r.denoise add on in GRASS 7.2.1 to remove noise from a .dem elevation model.

Initially when running the add on I received the error:

ERROR: mdenoise required. Follow instructions in html manual page to install it (g.manual r.denoise).

So I downloaded mdenoise.exe from the r.denoise page and added it to PATH in my environmental variables after placing it in a OsGeo4W/Bin and also tried putting it into the Grass 7.2.1 bin folder. Yet the error still persists...

I noticed that with Linux systems they recommend compiling the mdenoise from source. Should I have to install something like gg++ to compile and then try again? Or am I missing something obvious?


texts:
I am curious to know on how to implement NNJoin into the python code. I have created a point shapefile using shapely. I have another shapefile which is a centerline and is used as a reference line. Now I want to find the shortest distance from each point to the centerline. The way which I have done so far is converted the centerline shapefile into points and used NNJoin to find the shortest distance and the nearest coordinates of the point on the centerline. I want to know as to how I can implement the NNJoin using python.

Is there a way to import NNJoin as a library and use it?

If not are there any alternatives to it which I can use?


texts: Customization required while downloading shapefile from Geoserver

texts:
By default, Geoserver outputs shapefile in ZIP file format. I need help with 2 points:

How can I avoid wfsrequest.txt file in shapefile, which have request query?
How can I choose attributes to download?


texts: Hosting DigitalGlobe stereo imagery on GeoServer?

texts:
I have a pretty fresh install of the latest Boundless OpenGeo stack (v4.5).  I have a requirement to host some stereo satellite imagery.  I'm using the sample imagery over China provided by DigitalGlobe (https://www.digitalglobe.com/product-samples).  The .zip file is quite large and contains a directory for panchromatic imagery, one for multispectral imagery, and a GIS_FILES directory that seems to just contain various footprint shapefiles.  Each of the imagery directories contain many supporting files: .ATT, .EPH, .GEO, .IMD, .RPB, .STE, .TIL...  DigitalGlobe has an Imagery Support Data Document (https://www.digitalglobe.com/sites/default/files/ISD_External.pdf) that explains the information in each file.  They all seem to provide information about the image pairs, overlaps, etc.

I'm confused how to properly host this data in GeoServer.

Should the two GeoTIFFs be added to independent stores in GeoServer?

Is there some other sort of mechanism that keeps the two images together or contains the metadata in the DG provided files?


texts:
I have a shapefile layer in QGIS containing Europes coastlines.

I want to extract only Norways coastline from this.

If I make a polygon layer that overlaps only Norway and then use the Vector selection "select by location" nothing happens.

My tries have resulted in that the entire layer gets selected.


texts: How to publish raster data in Arcgis Viewer for Flex?

texts:
I'm new to ArcGIS for Server and I would like to publish several raster layers (land surface temperature, NDVI etc.). I have uploaded raster to server, but I cannot add it to my ArcGIS viewer for flex. Is it possible?
If not, what is the best solution to publish raster layer that is located in ArcGIS server?


texts: How to find a point half way between two other points?

texts:
I'm working on implementing the OSM Mobile Binary format for map data, and because it uses a 16 bit delta for each point along a way, it cannot represent a way if there is a large distance between two points.

The documented solution is to inject new points into the way, to keep the gap smaller.

Here's the relevant section of my code:

// calculate the distance between this point and the previous point,
// multiplied by 1,000,000 so it can be represented as an integer.
int32_t realNodeLatDelta = ((point.latitude - lastPoint.latitude) * 1000000);
int32_t realNodeLonDelta = ((point.longitude - lastPoint.longitude) * 1000000);

// cast as a 16 bit int (reduces filesize by about half)
int16_t nodeLatDelta = realNodeLatDiff;
int16_t nodeLonDelta = realNodeLonDiff;

// check if the cast caused corruption
if (realNodeLatDelta != nodeLatDelta || realNodeLatDelta != nodeLatDelta) {
  // if this point is reached, we need to add new points in the way,
  // which can be represented in a 16 bit delta
}

I'm not sure how to calculate the position of new points along the way? I assume it needs to be a great circle?

Points far enough apart to cause this issue are rare, but when they do happen it's usually a political boundary or similar, which tend to be quite long, and need to be accurately represented.


texts: Changing (or just inverting) Basemap Colors in QGIS

texts:
I am building a vector (strava routes) map using QGIS and the Stamen Toner OSM Basemaps but want to edit it to change the colors of the baseman slightly. Is this possible? I see Arc seems to have a beta feature but not sure if there's a way to do so in QGIS. Is the basemap essentially just an image?
If easier, I've found the Stamen Toner map inverted is essentially exactly what I want (white streets on dark background). So I've been building the map in QGIS with my computer colors inverted, is there an easy way to just do that that wouldn't involve manually changing colors?


texts: Changing the background color of the Georeferencer in QGIS

texts:
Is there a way to change the background color of the Georeferencer in QGIS? What I mean is, can I somewhere change the color that does not belong to the image that is about to be georeferenced?

Let’s say I have a world file or I know the coordinates of the corners of my image, then, in the following example, I would have a hard time placing the points exactly as I cannot clearly differentiate between my image and the white background of the Georeferencer:

lower left corner

upper right corner

I am open to any suggestions, including Python.


texts:
I use the following nested loops to read the first point feature of a dataset and the rest points in that dataset. Then read the second point feature of a dataset and the rest points in that dataset and so on. The problem with these nested loops is that it reads the first point and the rest points only once. Why this happens?

Dim pFCursor As IFeatureCursor
Set pFCursor = pFeatureClass.Search(Nothing, True)

Dim pFCursor2 As IFeatureCursor
Set pFCursor2 = pFeatureClass.Search(Nothing, True)

Dim pFeature As IFeature
Dim pFeature2 As IFeature

Dim pPointA As IPoint
Dim pPointB As IPoint

Set pFeature = pFCursor.NextFeature
Do Until pFeature Is Nothing

Set pPointA = pFeature.Shape

MsgBox pFeature.Value(intPosFID) & "F1"

   Set pFeature2 = pFCursor2.NextFeature
   Do Until pFeature2 Is Nothing

         Set pPointB = pFeature2.Shape

         MsgBox pFeature2.Value(intPosFID) & "F2"

   Set pFeature2 = pFCursor2.NextFeature
   Loop

 Set pFeature = pFCursor.NextFeature
 Loop


texts: Error when using user-defined function on LAScatalog in parallel

texts:
Using lidR 3.1.3, I made a user-defined function for catalog_apply() to process a LAScatalog following this, but a it raises the following error :
 Erreur : pas de méthode pour 'filter_poly' applicable pour un objet de classe "c('LAScluster', 'Spatial')"

Here's the working code. Error appears when plan(multisession) is inserted.
library(lidR)
library(future)
library(sf)

filter_poly = function(las, ...)
{
  UseMethod("filter_poly", las)
}

filter_poly.LAS = function(las, polygons_shp, buffer, outClass)
{
  poly_buff = st_buffer(st_read(polygons_shp, quiet = TRUE), dist = buffer, endCapStyle = 'SQUARE', joinStyle = 'MITRE', mitreLimit = 3)
  las <- merge_spatial(las, poly_buff, "inpoly")
  las <- filter_poi(las, inpoly == FALSE)
  return(las)
}

filter_poly.LAScluster = function(las, polygons_shp, buffer, outClass)
{
  # The function is automatically fed with LAScluster objects
  # Here the input 'las' will a LAScluster

  las <- readLAS(las)                          # Read the LAScluster
  if (is.empty(las))
    return(NULL)              # Exit early (see documentation)

  las <- filter_poly(las, polygons_shp, buffer, outClass)      # Filter shp
  las <- filter_poi(las, buffer == 0)           # Don't forget to remove the buffer

  return(las)                                  # Return the filtered point cloud
}

filter_poly.LAScatalog  = function(las, polygons_shp, buffer, outClass)
{
  opt_select(las) <-  "*"         # Do not respect the select argument by overwriting it

  options <- list(
    need_output_file = TRUE,    # Throw an error if no output template is provided
    need_buffer = FALSE,         # Throw an error if buffer is 0
    automerge = TRUE)           # Automatically merge the output list (here into a LAScatalog)

  output <- catalog_apply(las, filter_poly, polygons_shp, buffer, outClass)
  return(output)
}

poly_file = "Z:/TRASH/polygons.shp"
buffer = 10
polyClass = 24

ctg <- readLAScatalog("Z:/TRASH/")
plot(ctg)
opt_output_files(ctg) <- paste0("Z:/TRASH/{ORIGINALFILENAME}_filtered")
opt_laz_compression(ctg) <- TRUE

output <- filter_poly(ctg, polygons_shp = poly_file, buffer = buffer, outClass = polyClass)


texts: Throwing E_FAIL errors so edit operation is aborted

texts:
I need to abort an edit operation somewhere, determined by a certain condition.

Heres the code I'm using:

    private const int E_FAIL = -2147467259;

    private void _Events_OnDeleteFeature(IObject obj)
    {
        _logger.Info("Iniciando handler de OnDeleteFeature...");
        try
        {
            // do we need to abort it?
            // this method returns a COMException with code = E_FAIL;
            if (DoINeedToAbort(obj, geaEditOperation.DELETE))
                return;
        }
        catch (COMException comEx)
        {
            throw comEx;
        }
    }

This is not showing me the appropriate message and it's not aborting the operation. BTW, this event is a IEditEvent, and it's being fired ok. The problem is to abort the operation, without calling _editor.AbortEditOperation(), which gives me an ugly message.

This is where I looked it in ESRIs manual (I know that the following applies to class extensions)

http://resources.esri.com/help/9.3/arcgisengine/ArcObjects/ao_start.htm#esriGeodatabase/IObjectClassEvents.htm

Any tips?


texts: Which version of C# should I learn for ArcObjects with .NET?

texts: I'm wanting to pick up ArcObjects and have decided on C# over VB as my language. I know that arcpy uses python version 2.7. If I had started learning python 3.2 to learn arcpy I'm sure I would have had issues. I don't want to have any extra issues learning C#. Which version of C# should I learn?


texts: Highlighting most recent point in time slider playback using ArcMap

texts:
I am trying to put together a video where spatial points from an animal are displayed.  However, many of these points are on top of other points and it is easy to lose track of what is happening.

Is there a way to highlight or otherwise change/draw attention to the most recent point being visualized on the map as the time slider plays?

I am using ArcGIS Desktop (ArcMap) 10.6.


texts: How to calculate correct angles for lat/lon points with the QGIS API?

texts:
I am trying to use qgis API to calculate the angle between two points by using function azimuth of class QgsPoint ( http://qgis.osgeo.org/api/classQgsPoint.html ), and my CRS is EPSG:4326 - WGS 84 , but I believe I am getting the angles wrong, as if QGIS were ignoring that the coordinates are Latitude and Longitude.
Could it be possibly true? If yes, how can I properly calculate the angle and why would one want to make such a confusing function?
Thanks in advance


texts:
I am looping through a csv, taking latlng values and trying to add there geometry to postgis enabled postgres table. xy values are int's

conn = psycopg2.connect("connected")
cur = conn.cursor()
cur.execute('''CREATE TABLE test pt_geom geometry);''')
conn.commit()
for row in csvfile:#sample code
    x = row[0]
    y = row[1]#etc...
    xy = (x, y)#tuple of the int xy's
    cur.execute('''insert into test
                        (pt_geom) values
                        (ST_SetSRID(ST_MakePoint(%s,%s),4326));''', (xy))
    conn.commit()

the script is running, connecting to the db, creating the table inserting the latlng but the values are coming up in a strange text format-I am not sure why that is..I have not found any info on this yet on the web. possible solutions? readsoning? the xy is an int value. is it spatially displayed? does the %s convert it to a string...anything would be helpful


texts: Calculating sequential values in ArcGIS Desktop?

texts:
Previous posts (such as Adding automated values to attribute field?) answer how to calculate sequential values in an attribute table.  The posts all seem to rely on the process described in this ESRI Tech Support page:

https://support.esri.com/en/technical-article/000011137

Here are the article's details:

  Create a new short integer field.
  Set the Parser to Python.
  Select Show Codeblock.
  Paste the following into the Pre-Logic Script Code:

rec=0
def autoIncrement():
 global rec
 pStart = 1
 pInterval = 1
 if (rec == 0):
  rec = pStart
 else:
  rec += pInterval
 return rec

  Paste the following code in the smaller box below the Pre-Logic Script Code:

autoIncrement()

  Click OK.

But at 10.5.1 I can achieve the same results with this single line in the field calculator:

!FID! + 1

The results are exactly the same as the more complex ESRI-supplied code, even if the the table is sorted on a non-FID field and/or there is a selected set of records.

What am I missing here?  Does the ESRI-supplied code do something that I'm not aware of?


texts: Connect two contour lines with calculated line length

texts:
Is it possible to make a plugin for QGIS like PEGGER for arcview 3.2?
Problem is, how to make path over contour lines with a gradient in percent?

Formula is: K= e*100/G*M ; K=10*100/0.08*10000=1.25m

K- distance between two contours
e- equidistance (vertical difference between the two contour lines)
G- Grade/Slope ( 0,08 = 8%)
M- map scale

When I set a grade by clicking with mouse line would be drawn.


texts: Count number of raster cells cover by shapefile

texts:
I have a raster image of Alberta province (1km * 1km resolution ), I put the other shapefile (i.e., grounds, hospital, schools etc) in QGIS . I want to calculate the number of raster cells (1km * 1km ) covered by ground, hospital etc.
Is there any way to find the how many raster cells covered by each required area?


texts: Line symbology based on a field for color AND size (linewidth)

texts:
It's possible in QGIS to add a symbology for lines which is based on a field.

This symbology may then use a method with either:

color
or
size

But is it possible to use both, i.e. colors AND size based on the same attribute field for classification?

It seems possible to define it with method "size" and then manually change the color (or the inverse; method "color" and specifying width manually) for each entry but it's rather tedious.


texts: Is there any method to choose the value of each pixel according to land cover

texts:
I have 2 images; one is land cover image (imgA) and monthly mean temperature image (imgB). I'd like to choose value of each pixel from imgB corresponding to imgA.
WhaT I expected is each pixel of certain land cover type should have different values in different area. I tried it using zonal statistic in ArcGIS and ERDAS. But zonal statistic give me only one value of each land cover (not in each pixel value). I'd like to keep all original pixel values of imgB according to each land cover. How can I do this using ArcGIS or ERDAS?


texts: Converting File Geodatabase Table to Access 2007 ACCDB table programmatically?

texts: Does anyone know of a method to convert a File Geodatabase Table into an Access 2007 Table programmatically.  I have been trying using the PYODBC library and a series of insert statements but am having difficultly with some of the data types and NULL.  Basically an unload to Access.


texts: I am currently doing a project where I want to use some form of mapping API (Google maps or an open source alternative) that can show the location of simulated wind turbines across Ireland. I am trying to understand what would be the most efficient protocol i could use to display information stored in a monitoring database onto the map as quickly as possible.


texts:
In previous versions of OpenLayers, we could filter WMS services with CQL. However, I can't filter my WMS services using CQL or other.

How do I do that in OpenLayers 3?

 var wmsLayer = new ol.layer.Tile({
    source: new ol.source.TileWMS({
        url: 'http://localhost:8282/geoserver/Isemri/wms',
        params: { 'LAYERS': 'Isemri:mahalleler' },
        serverType: 'geoserver',
        //query: 'OBJECTID_1=262'
        cql_filter: 'OBJECTID_1=237'
    })
});


texts: Converting MSSQL linked MapInfo Tab file to a shapefile using ogr2ogr

texts:
I have some tab files that are actually references to a view in a MS SQL Server table. It looks like

begin_metadata
"\IsReadOnly" = "FALSE"
"\MapInfo" = ""
"\MapInfo\TableID" = "055372be-523f-430a-b603-bbd462ad3efd"
"\DATALINK" = ""
"\DATALINK\ConnectionString" = "DSN=xxx_GIS;UID=sa;APP=MapInfo Professional®;WSID=AS00769;DATABASE=xxx_GIS"
"\DATALINK\Query" = "select ""assetID"", ""assetName"", ""ward"", ""locality"", ""GID"", ""location"", ""assetClass"", ""assetSubClass"", ""assetType"", ""hierarchy"", ""dimensionsLength"", ""MI_PRINX"", ""MI_STYLE"", ""OBJECT"" from ""NGSC_GIS"".""dbo"".""v"
"\DATALINK\Query\+1" = "w_ASSET_Bridge"""
end_metadata

In QGIS the view is like

<layer-tree-layer expanded="1" source="dbname='xxx_GIS' service='xxx_GIS' estimatedmetadata=true disableInvalidGeometryHandling='0' table=&quot;dbo&quot;.&quot;vw_ASSET_Bridge&quot; (SP_GEOMETRY) sql=" id="vw_ASSET_Bridge_1a8f7b4e_d11d_4a9e_b5ba_773aafe7014b" checked="Qt::Checked" providerKey="mssql" name="vw_ASSET_Bridge">
  <customproperties/>
</layer-tree-layer>
<layer-tree-layer expanded="1" source="../../GIS/Data/ASSET_Bridge.TAB" id="ASSET_Bridge_b0c1d436_6728_428f_a1dc_f22c15fe8c38" checked="Qt::Checked" providerKey="ogr" name="ASSET_Bridge">
  <customproperties/>
</layer-tree-layer>

Using the advice in https://alastaira.wordpress.com/ogr2ogr-patterns-for-sql-server/

I see the command should be something like

ogr2ogr -f “ESRI Shapefile” “C:\scratch\vw_ASSET_Bridge.shp” “MSSQL:server=xxx_GIS;database=xxx_GIS;trusted_connection=yes;” -sql “SELECT [ID],[GeoType],[GeoName],GEOMETRY::STGeomFromWKB([GeoShapes].STAsBinary(),28354) As [GeoShape] FROM dbo.vw_ASSET_Bridge WHERE [GeoShapes].STGeometryType() = ‘POINT'” -overwrite -a_srs “EPSG:28354”

But I am not sure what should be used for the tags in []


texts: Can .Net version of ArcGIS Server expose Oracle spatial data without ArcSDE, or is this limited to Java version?

texts:
I'm planning a new deployment of ArcGIS Server and I want to avoid using ArcSDE if possible. My vendor requires that we store our spatial data in Oracle and we prefer to use Microsoft server technologies whenever possible.

Our ideal environment is to host ArcGIS Server on Windows Server and expose the REST services over IIS.

Can the .Net version of ArcGIS Server expose Oracle spatial data without ArcSDE or is this limited to the Java version?


texts: QGIS - raster resampling having no effect=pixelization

texts:
In QGIS I cannot make the resampling setting (changing it from nearest neighbour to bilinear or cubic) to smooth the pixels of Sentinel-2 imagery, both as jp2 or tif. If I export the image, convert it to png, georeference it, then the resampling works and the image pixels are smoothed out. Irfanview has no problem resampling jp2 or tiff.
I cannot find anything online describing this problem. Am I missing something that allow the jp2/tif resampling on QGIS?


texts: creating a new layer/selection based on the conditions in Python script?

texts:
Is there any easy way how to write a python script, which would create a new layer based on location from a border line and "if" condition?

I have two layers A and B, both are made of several polygons. I need to select those polygons from both layers based on location, e.g. 1000 meters, they must touch each other or the border line, and based on comparing their attribute values of one selected column. Then i need to create a new selection or layer, if the difference between the values is big enough.

I mean something like "if column.A - column.B > x ^ polygon.A and polygon.B touch each other or the line/or the distance is 1000meters from the line, then create a new layer of those two polygons"


texts: Changing width of identify widget pop-up in Configurable Map Viewer?

texts:
How do I change the width of the identify tool pop-up in order to fit a table without having to scroll?

Image attached shows how it is at present. Note I want to change size of entire popup, not just the table.


texts: Measurement units not returned in project CRS units

texts:
I am attempting to calculate area of country polygons using the field calculator.

The layer in question is unprojected WGS84, but I have set the project CRS to World Equidistant Cylindrical (Sphere), which is listed with the projected CRSs. I've also enabled on-the-fly transformation.

From my limited understanding of the program (I am running QGIS 1.8.0), the measurements should be given in the units of the project CRS, and the proj4 string for the World Equidistant Cylindrical seems to indicate that its units are in meters.

Nonetheless, the measurements are off by orders of magnitude. For example, Afghanistan is significantly larger than 92 square meters, or even 92 sq km. I can only assume I am getting area in degrees. What must I do to change that?


texts: How to use locally served tiles for interactive maps instead of just a png image of tiles

texts:
I followed this link and it worked for me..

https://switch2osm.org/serving-tiles/manually-building-a-tile-server-14-04/

I can serve tiles locally but running it and by just typing the link like for instance - localhost/osm_tiles/0/0/0.png and it gives me an image with a particular zoom level, x and y axis.

Now I have this advanced localized openstreet map example I saw..
http://openlayers.org/en/v3.2.1/examples/localized-openstreetmap.html which is of course an advanced custom tile server. Taking the tutorial I got from switch2osm, is it possible for me to create that same kind of tile server?


texts: `folium.GeoJson` displays all polygons in GeoDataframe but not individual ones

texts:
folium.GeoJson displays all the counties of New York State:
import geopandas as gpd
gdf = gpd.read_file("Counties_Shoreline.shp")
m = folium.Map(location=[43.062776, -75.420884], tiles="cartodbpositron", zoom_start=7)
folium.GeoJson(data=gdf["geometry"]).add_to(m)
m

However, it fails to display a single county, outputting a map without borders when only the first row data, gdf["geometry"][0], is used:
m = folium.Map(location=[43.062776, -75.420884], tiles='cartodbpositron', zoom_start=7)
folium.GeoJson(data=gdf["geometry"][0]).add_to(m)
m

This despite the fact that the polygon appears fine, rendering inside Jupyter Notebook.
Data here: http://gis.ny.gov/gisdata/fileserver/?DSID=927&file=NYS_Civil_Boundaries.shp.zip
Shapefile: Counties_Shoreline.shp


texts: lat, lng to lng, lat in polygon WKT - axis flip

texts:
I have a sheet with the column with polygon information in WKT format. However, the data taken are from web mapping (OSM or google maps) and in the order of lat, lng. When I am using the solution from this thread (Visualising WKT geometry string in QGIS) which are what I need, I do not have another choice to switch x and y other than manually. But concerning the number of data, I need to find something more efficient.
Basically my polygon data is in the format

55.947686 -3.132133, 55.945932 -3.131232, 55.947878 -3.121362, 55.943168 -3.117714, 55.943841 -3.108273, 55.94771 -3.108745, 55.950161 -3.119216

but I need -

-3.132133 55.947686, -3.131232 55.945932, ....

There are so many items in the dataset so manual transfer is not an option.
Is there a tool somewhere (even outside QGIS environment) to do it for me?
Screenshot for an answer from BERA:


texts: Can anyone recommend any online courses for Python?  The MSc in GIS I hope to start later this year has Python programming modules, but you can substitute with an Environmental module if you so wish.  I would like to have a go at Python before I commit to that module within the MSc, as I have no previous programming experience


texts: Setting dashed stroke style when layers have features with sharing borders in QGIS

texts:
Layers like municipalities, countries,... have features sharing the border of their neighbours: they share the same vertices (except for the features at the edges off course).
I have a municipality layer, which I set with a dashed stroke style. However, for adjacent features, the borders are displayed as a solid line (see picture)

How can I display all borders with a dashed stroke style?
Now it seems to overlap.


texts: Converting xml -> dict -> xml using Python?

texts:
I want to convert xml (qml) file to dictionary, then convert this dictionary to new xml. (Also want to edit something, but first lets just convert it in both way without it.)

I'm using xmltodict module to xml->dict and then dicttoxml module,dom parse to string and file.write() to dict->xml.

    xml = open("path",  "r")
    org_xml = xml.read()
    dict_xml = xmltodict.parse(org_xml)

    new_xml = dicttoxml.dicttoxml(dict_xml, custom_root = 'qgis')
    dom = (parseString(new_xml)).toprettyxml()

    with open("path/new.qml", 'a') as file:
        file.write(dom.encode("utf-8"))

The problem is, that while converting I 'loose' something. First file is a categorized symbol style tempale and the second is single symbol.

Anyone tried conversion like this and can help?


texts: How do I install the las2ogr tool from liblas in Windows using the osgeo4w installer?

texts:
I was wondering if anyone has used the LAS2OGR tool?  I heard about this from a previous post.  The webpage gives the usage: http://www.liblas.org/utilities/las2ogr.html

las2ogr -i mydata.las -o points.shp -f "ESRI Shapefile"

But I am wondering what I run it in?  I am on a Windows & box with FWTools 2.4.7 and it is coming up as an unrecognized tool.  When I open OSG4W (gdal 1.9.2 released 2012/10/08) las2ogr is not listed in the tool set.  Where should I be accessing this from?  Do I need to download a specific .exe or .dll and save it to the GDAL folder?

Any help would be appreciated.

Thanks...


texts:
I'm very new to OpenLayers and all the GIS world.

I'm working with rasters (tiff that are converted to postgis sql)

I have this code:

        var feature;
        $.ajax({
            method: "POST",
            url: Routing.generate("getdpageometry"),
            dataType: 'json',
            data: { tipodpa:tipodpa, dpaid:object.id },
            async: false
        }).done(function( obj ) {
            feature = obj.feature;
        }.bind(this));

        var formatWKT = new ol.format.WKT();
        var extentFeature = formatWKT.readFeature(object.extent, {
            dataProjection: 'EPSG:3857',
            featureProjection: 'EPSG:3857'
        });
        var formatGeoJSON = new ol.format.GeoJSON();
        var feature = formatGeoJSON.readFeature(feature, {
            dataProjection: 'EPSG:32719',
            featureProjection: 'EPSG:32719'
        });
        feature.getGeometry().transform('EPSG:32719','EPSG:3857');
        Store.vectorSource.clear();
        Store.vectorSource.addFeature( feature );
        Store.map.getView().fit(extentFeature.getGeometry().getExtent(), Store.map.getSize());

This is part of a code that is showing a map but only a region (defined by an extent in the database) is colored, and adds some blue border to it, like this:

What I need to do is change that blue color, I need it transparent, I only need the border to be blue, and also I would like to change the border color too, so my question is:

Where should I be looking to change those colors?

I'm pretty confused, not sure if it's possible to change it like a parameter, or it's OL3 library default color (I've search into the code, I have changed some blue code colors but nothing happens), or it has to be related with the raster from the database, or where/how...

Using the chosen answer, now my code has:

    this.state.vectorLayer = new ol.layer.Vector({
        source: this.state.vectorSource,
        style: new ol.style.Style({
            fill: new ol.style.Fill({  color: 'rgba(0, 0, 0, 0.2)'  }),
            stroke: new ol.style.Stroke({ color: '#FF0000',  width: 2  }),
            image: new ol.style.Circle({
                radius: 7,
                fill: new ol.style.Fill({
                    color: '#ffcc33'
                })
            })
        })
    });

And the map:


texts:
I'm using this WMS service
    http://www.cartografiarl.regione.liguria.it/mapfiles/repertoriocartografico/LAVORIINCORSO/1633.asp?REQUEST=GetCapabilities&SERVICE=WMS&VERSION=1.3.0&

I can query it using GetFeatureInfo in a WMS request like this one:
    http://www.cartografiarl.regione.liguria.it/mapfiles/repertoriocartografico/LAVORIINCORSO/1633.asp?&LAYERS=L4313&QUERY_LAYERS=L4313&STYLES=&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetFeatureInfo&BBOX=1018299.271264%2C5479845.539429%2C1123018%2C5492228.338009&FEATURE_COUNT=10&HEIGHT=162&WIDTH=1370&FORMAT=image%2Fpng&INFO_FORMAT=application%2Fvnd.ogc.gml&SRS=EPSG%3A3857&X=891&Y=22&

The response is:

1086509.551130,5490410.633321 1086509.551130,5490410.633321 4 011023_002 Riccò Del Golfo Via Aurelia, Snc 03 immobile destinato autorimessa 1 21/10/05 altre associazioni a delinquere 03 destinato 19/07/2007 02 Comune istituzionali PRO_ciuffardi_paolo.pdf 011023_002_F.pdf 1086505.336252,5490406.481274 1086505.336252,5490406.481274 3 011023_001 Riccò Del Golfo Via Aurelia, Snc 03 immobile destinato appartamento 1 21/10/05 altre associazioni a delinquere 03 destinato 19/07/2007 02 Comune sociali PRO_ciuffardi_paolo.pdf 011023_001_F.pdf

Now I'd like to do the same request but for a bigger BBOX for example
    http://www.cartografiarl.regione.liguria.it/mapfiles/repertoriocartografico/LAVORIINCORSO/1633.asp?&LAYERS=L4313&QUERY_LAYERS=L4313&STYLES=&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetFeatureInfo&BBOX=832454%2C5423777%2C1118664%2C5564559&FEATURE_COUNT=100&HEIGHT=162&WIDTH=1370&FORMAT=image%2Fpng&INFO_FORMAT=application%2Fvnd.ogc.gml&SRS=EPSG%3A3857&X=891&Y=22&

I receive no response and no error ......Too big the bbox?


texts: Help with Geographic Weighted Regression: Condo Prices

texts:
I was working on a GWR Problem in ArcGIS to figure out the potential real estate values of Condos in a city. I thought I was following everything correctly but I keep receiving the same error over and over again.

I have attached a link to the dataset below.
The dependent variable is Finalprice. The independent (explanatory) variables are all the variables that are KM_* (distances to desirable features) and the DIST_TRECC (km to city center).

It keeps throwing the error that explains that there is a major model problem with multicollinearity. I also did the OLS regression analysis and a link to those results is below as well. Everything I read said that the numbers from the OLS analysis show that the GWR should be working. I am at a loss and not sure what to do next.

LINKS
DATASET: https://db.tt/9EUXfJy4
OLS Results: https://db.tt/LwS8pSIm


texts: Merge reversed linestring in correct order in postGIS

texts:
I want to merge linestring in the correct order in postGIS, but sometime they are reversed so something like this :

select astext(linemerge(ST_union(
ST_GeomFromText('LINESTRING(3 3, 2 2)'),
ST_GeomFromText('LINESTRING(1 1, 2 2)'))));

returns this :

"LINESTRING(1 1,2 2,3 3)"

But I want this :

"LINESTRING(3 3, 2 2, 1 1)"

I really don't know what to do,

Thanks.


texts:
I am doing this:

Access a few columns from a file geodatabase.

arcpy.env.workspace = path
datasetList = arcpy.ListTables("tableA") # Here putting the name of the table

 for dataset in datasetList:
      with arcpy.da.SearchCursor(dataset, ["tablename","query"]) as cur:
           for row in cur:
               print row[0]

RuntimeError: A column was specified that does not exist

This works fine if it reads from a certain table. If I put another table with the same mentioned columns as aliases but different names it gives the error. If aliases can't be used here then what is the usability because now the table has to be recreated with the correct names. It is a problem that it doesn't allow you to rename the columns.


texts: Simplifying land cover classification in QGIS?

texts:
I have recently completed a land use classification using over 100 ROIs for 5 different land cover classes and it has developed what looks like a highly specific spectral signature for each class. This has meant that the classification is almost like a series of disconnected coloured squares (see pic). It would be great if, for example, similar pixels surrounding or next to a different pixel could convert it to the same class so as to show more contiguous classes. Any ideas?
Green = forest
Purple = sub canopy cultivation
blue = degraded forest
cyan = monoculture
yellow = village
enter image description here


texts: FME - How to convert string represenation of date to a date format that SQL Server will accept?

texts:
I have dates in the form of 07-28-2010 I am extracting from DWG filenames. I need to push these to a DateTime formatted table in non-spatial SQL Server 2008. SQL Server doesn't like them as strings coming in from the ETL apparently. How can I set the datatype in the ETL from string to datetime? Is this even what I need to do? I'm pretty sure this is my issue based on past experience with pushing dates to SQL Server and this error from FME:

MS SQL Server Writer: 1 attribute value(s) failed conversion, and were written as NULL values


texts:
I was wondering if there is a way, or plugin, to save files (i.e. geojson, kml, shp, etc.) from QGIS directly into GitHub.

I know you can save it locally to a hard drive and copy it over, but I am looking for a way to save directly into a GitHub account.


texts: Is there a way to save and share topology rules in QGIS?

texts: I would like to create a set of topology rules and use them over and over again. I would also like to be able to share them with my team. Is there a way to do this in QGIS? I'm using QGIS 3.0 but I can switch to earlier versions of necessary.


texts: Where on-line can I find a GIS land use .lyr that shows standard colors for each land use code?

texts: I have been looking on-line for a GIS land use layer file that shows all the standard colors that represents each land use code.  I have not been successful in finding that.  Would you please direct me to where I can find this layer file on-line?


texts:
I'm an avid cyclist and want to reverse-engineer how Strava groups bike riders together.  Here is their method to determine if cyclists are riding together (they use time and lat/lon of a ride):
http://resources.esri.com/help/9.3/arcgisengine/java/gp_toolref/geoprocessing/proximity_analysis.htm

I want to create the algorithm and play with it in Tableau.  The only methods I can think of are very computation intensive.  The data has lat/lon every second.  I then put a proximity zone of 150m around each point and then compare each and every point with other riders' locations at the same time.  The result of each compare is 1 (yes) or 0 (no).  After completing the comparison of every point, determine if the rider was within the proximity zone 70% of the time.  With thousands of files to parse, this will take massive CPU power.

There must be a better way.  I really want to add a proximity zone around each route.  Then turn that route into some kind of "signature".  I would then overlay the signature of one route to another adjusted for time to determine if two people rode their bikes together.

Can anyone point me in the direction of solving this problem?


texts: Finding area covered by unmasked pixels in an image in Earth Engine

texts:
I want to calculate the area covered by forests in a region. Using Hansen image I have demarcated the forest regions in India. Now, how do I get the total area of forests?

 // FILTERING OUT FOREST COVER IN INDIA

var treeCanopyCover = hansen.clip(india).select('treecover2000');
var greater25 = treeCanopyCover.gte(25);
var treeCover_greater25 = treeCanopyCover.updateMask(greater25);

I know that I can use ee.Image.pixelArea(), and then subsequently use reduceRegion. But, I don't understand what scale I should give. In fact, I would like to know, what exactly are the roles of scale and maxPixels. Could someone enlighten me on this regard?

var areaImage = treeCover_greater25.multiply(ee.Image.pixelArea());

var stats = areaImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: india,
  scale: 100,
  maxPixels: 1e13
});


texts: ESRI.ArcGIS.Desktop.AddIns could not be found. Visual basic build error

texts:
I am trying to compile a project in visual basic into an Arc map add on.  I am receiving two warnings that are preventing the build from completing.  Both errors entail that the assembly ESRI.ArcGIS.Desktop.AddIns could not be found.  I received this project and source code by request from someone else.  From other responses on here to similar issues in compiling in Visual Basic I have made sure that the Register for COM interop is unchecked in the compile tab of the project properties and I have made sure the Copy Local and Specific Version in the References tab is set to false.

When I searched through my computer I was not able to find the ESRI.ArcGIS.Desktop.AddIns assembly in order browse and assign a reference path in Visual Basic.  It seems the it is not on my computer.  I have not had any experience with an assembly before, so my question is how do I get access to this ESRI.ArcGIS.Desktop.AddIns assembly?


texts: Watershed delineation using hydrological Python library

texts:
We would like to delineate watershed polygon for a given pour point using hydrological Python library.
Until now, we have tried these libraries:

pyshed
Grass (r.watershed)
Whitetoolbox

We have around 1200 pour points for which we would like to determine a watershed area.

Our generated flow direction for a given DEM seems valid but we do not know why we always get a strange watershed result.

Here is an example of what we have obtained for this pour point:

The calculated watershed area corresponds to this (the black pixel):

Did you already get this kind of strange results?

What could cause such results?

Could projection impact the calculation?

I do not know what we do wrong but we always obtains an invalid watershed result.

You can my drive folder URL where you find the DEM and the pour point used in my example.

https://drive.google.com/drive/folders/1LHkNDnDVNSLgx6Lw894mxL5XShHS1EWs?usp=sharing


texts: Export Google Earth Engine to Google Drive with service account (Python API)

texts:
I am  trying to build an application with Earth Engine's Python application. I first authenticate Earth Engine with a service account, then export the image to drive using the ee.batch.Export.image.toDrive function. However I cannot find a Google account that can access the drive folder in task.status(), after the task is completed.

{'state': 'COMPLETED', 'description': 'img5.9975525.997552', 'creation_timestamp_ms': 1605873230011, 'update_timestamp_ms': 1605873430624, 'start_timestamp_ms': 1605873339520, 'task_type': 'EXPORT_IMAGE', 'destination_uris': ['https://drive.google.com/#folders/1qFHOln6kdKPq6RJE_Vlc4FPSiBXVwu8W'], 'attempt': 1, 'id': 'H5BSMEH3SJG5ELP6FBLGXNK4', 'name': 'projects/earthengine-legacy/operations/H5BSMEH3SJG5ELP6FBLGXNK4'}

I have tried Enable G Suite domain-wide delegation for the service account, I also created a commercial account and added the service account to the Domain-wide delegation of the commercial account, but the folder is still unaccessible for both my personal account(which is the owner of the project) and the commercial account.
How can I access this folder?


texts: Calculating lat/long coordinates from DLS (Dominion Land Survey) Location

texts:
I am trying to calculate the legal subdivision center of a DLS location and even when I try the simplest calculation my numbers are always off from other calculators.

Take for example DLS location 1-1-1-1 W2. From what I understand being adjacent to the 102W meridian and adjacent to 49N as well it should be relatively easy to calculate This is also simple since I do not believe road allowances come into effect on this section as they begin on the top and left sides only AFAIK, and this section and LSD is the very bottom right side of the grid.

That said, even when I dumb down the calculations to go straight to the LSD level my numbers are off. An LSD is known to be 1/4 mile, and being in the bottom right section and bottom right legal division grid, it should be 1/4 mi west and 1/4 north north for the top left boundary from 49N/102W.

I used the rhumb line calculator here (http://www.movable-type.co.uk/scripts/latlong.html) to calculate .25 miles west (bearing 270), then .25 miles north (bearing 0) to get the the top left point of the polygon square.

I get the following top left corner coordinate: 49.0036182946, -102.0055151967

When I average that and 49, -102 I get: 49.0018091473, -102.0027575983

Yet the calculators I see are all extremely close to: 49.000976, -102.009228

So I must be missing something here because I cannot see how everyone else gets a number the same but not my calculation. Does anyone know what I am missing here?


texts: multiband raster: set no data values to 255,255,255,255

texts:
ArcMap 10.2:
I have a 4 band raster dataset (img-file) with no data values.
I want that the no data values appear in white (255,255,255) and therefore want to change it to these values.
I know how to do it with a 1 band raster as explained here:
https://desktop.arcgis.com/en/arcmap/10.3/tools/spatial-analyst-toolbox/h-changing-nodata-cells-to-a-value.htm
Identify the NoData locations with the IsNull tool:

    Input raster : InRas1

    Output raster : OutNullRas

Use Con on the results to give the original NoData areas a value:

    Input conditional raster : OutNullRas

    Expression : "Value = 1"

    Input true raster or constant value : 100

    Input false raster or constant value : InRas1

    Output raster : NullTo100_Ras

Is this working with a 4 band raster in a different way?


texts:
I am trying to divide a feature collection into a mosaic of different cells in order to aid with exporting over a large area. I found the following script which is mentioned in this stack overflow question to divide a geometry into bite-sized areas and I am trying to convert it to Python so that I can export data in an iterative fashion.
I have considered this question and this question, but neither seem to adequately capture the complexity of the above script. The issue pertains to how python treats scope and I have tried doing a (rather dumb) method of carrying variables as arguments through each  and every function, but this breaks down when we are forced to use .map()
Here is my naive approach to translating between the two (mostly converted using geemap:: js_snippet_to_py() function.
def tile(region, size_in_degrees):
    def tile_geometry(geometry):
        coords = ee.List(geometry.bounds().coordinates().get(0))
        min_point = ee.List(coords.get(0))
        max_point = ee.List(coords.get(2))

        def sequence(direction):
            start = min_point.get(direction)
            stop = ee.Number(max_point.get(direction))
            return ee.List.sequence(start, stop, size_in_degrees)

        def build_numbers(y):
            return ee.Feature(ee.Geometry.Rectangle([[ee.Number(x), ee.Number(y)], \
                [ee.Number(x).add(size_in_degrees).min(max_point.get(0)), \
                  ee.Number(y).add(size_in_degrees).min(max_point.get(1)) \
                ]]).intersection(geometry, 30))

        def func_nyp(x):
            return sequence(1).map(build_numbers)

        tiles = sequence(0).map(func_nyp).flatten()

    return ee.FeatureCollection(tiles).filterBounds(geometry)

    def tile_feature_collection(feature_collection):
        return ee.FeatureCollection(
          feature_collection.iterate(iterate_func))

    def iterate_func(feature, acc):
        return ee.FeatureCollection(acc).merge(tile_geometry(ee.Feature(feature).geometry())), ee.FeatureCollection([])

     if (isinstance(region, ee.FeatureCollection)):
        return tile_feature_collection(region)
    elif (isinstance(region, ee.Feature)):
        return tile_geometry(region.geometry())
    elif (isinstance(region,ee.Geometry)):
        return tile_geometry(region)
    else:
        print('region must be ee.FeatureCollection, ee.Feature, or ee.Geometry. Was ' + str(region))

When I run the above code, the sequence() function is not known to the global environment. However, when I try to extract each of the nested functions, I quickly run into issues surrounding variables that need to be carried through from others.


texts: Python Script to saves as Raster Img with new CRS, Map View Extent, and x,y Resolution

texts:
I am currently trying to write a python script that can do the same thing as when you right-click on a raster image and click save as. So I am still figuring out how to change the CRS for the raster image in python, choose the map view extent option and not the layer extent option, and change the layer resolution

.

I've tried:

        new_cord = iface.mapCanvas().extent()

    crsSrc =  layer_find.crs()
    crsDest =  QgsCoordinateReferenceSystem(102671, QgsCoordinateReferenceSystem.EpsgCrsId)  # WGS84 or UTM Zone

    renderer = layer_find.renderer()
    provider = layer_find.dataProvider()

    pipe = QgsRasterPipe()
    pipe.set(provider.clone())
    pipe.set(renderer.clone())

    file_writer = QgsRasterFileWriter(filename)
    file_writer.Mode(1)

    file_writer.writeRaster(pipe, provider.xSize(), provider.ySize(), new_cord, provider.crs())

But the file that is being written is around 7gb and just causing my QGIS to freeze.


texts: Eliminating decimals without approximation

texts:
Is there any way I can eliminate decimals without approximating the number?
Example:
 Value    | Desired output
 1.69     | 1
 2.12     | 2
 105.9999 | 105

I want to use the field calculator for this. Do you have any suggestions?


texts: Joining multiple records to single feature using QGIS

texts:
I'm trying to create a join/relate in QGIS where I have a shapefile of buildings and to that I would like to join a non-spatial table (.csv-file) containing people who work in each building. So I have multiple records in my table that I would like to be able to join to my single features in my shapefile.
I can only run a join in QGIS where the first record is joined to the shapefile feature (ie building) but subsequent records get removed.
Could someone let me know how to complete this join/relate in QGIS?


texts: Multipage atlas with different maps in every page

texts:
I have to print 2 pages with different maps on every page. Every map represents a city as polygone and the scale in every map is the same. The setting of the zoom and the scale is not a problem: I have set the scale and for the extent in the map extention I've used the fallowing expressions
xmin(geometry(getFeature('frame', 'gid', 1))) ymin(geometry(getFeature('frame', 'gid', 1))) xmax(geometry(getFeature('frame', 'gid', 1))) ymax(geometry(getFeature('frame', 'gid', 1)))
and set map extension.

Now I have to find a way to have 2 pages as showed in the fallowing picture

The part in red is the same in every page. The 2 pages have the same dimensions and the same legend. In every map there are the same layers but zoomed on different polygones.
In print-composer I can set the number of pages but in this way I have to duplicate the red part (legend and general info). I have tried to overlap the maps of the first page on the ones of the seconds and find a conditional function on the transparency of the map/page number (ex: if mapA is in page 1 transparency=0, if mapA is in page 2 transparency=100 and il map F and G are on  page 1 transparency=100, if they are on page 2 transparency=0) but trying to find a way for page 1 and 2 with atlas I didn't work it out. Is there a way to change the maps and its position and dimensions through different pages dinamically?


texts:
I have a .shp file, and I need to get the dual graph from this file (that is to say, each polygon centroid should become a vertex, and two vertices are linked if their associated polygons are adjacent).

I managed to get the centroids from the polygons of my file.

How can I add lines between the vertices that should indeed be connected?


texts: Take color from other layer with data-override/expression

texts:
QGIS 3.8, Win. 10. In addition to a question asked earlier, I want to set the color of text-buffers to adapt to the color from another layer in the background. I have a shapeburst-style background (blue) and on top of it the label with buffer. I would like to keep the buffer, but change the color (override) to exactly the same color as the background (which is on another layer): is that possible? Can I create an expression that takes the color dynamically form another layer?

In the screenshot, the label-buffer for name 3a should have a darker blue value on the left and a brighter one on the right.

I need the buffer to to mask the leader lines that continue below the label. So the buffer does not make sense as it has exactly the same color as the background, but it makes senst to cover the lines.


texts:
How can I setup a portable QGIS to run from a flash drive / USB-Stick?
In the userguide for QGIS 1.7 Wroclaw there is only mentioned the option of "--configpath". But to me it is not clear how to deal with this information.

Can someone post a step-by-step explanation!


texts:
I've been trying to add correct coordinates to my grid in GIS.

I have a map of the Iberian Peninsula using an Open source layer from OSM/Stamen (CRS: EPSG:3857. WGS 84 Pseudo Mercator) and while in the print composer, after adding the grid this is what iI get:

Now, these coordinates do not make any sense to me. I know that they are related to the Interval units and offset but I am quite unsure how to apply correct coordinates on my map. So please tell me how to correct the showing coordinates!

Might help to understand:

Let me know if you need further details!


texts: OGR writing MapInfo Tabfiles with column Indexes

texts:
I'm using ogr2ogr to convert MSSQL to MapInfo File and would like to add indexes to some of my MapInfo columns, is this possible?

This is what I have so far

set ogr2ogr="C:\Program Files\GDAL\ogr2ogr.exe"
set sqlcon="MSSQL:server=SQL;database=Spatial;User id=user;Password=password;MARS_Connection=yes;"

%ogr2ogr% -f "MapInfo File" test.tab %sqlcon% -sql "SELECT * FROM [dbo].[tableName]"


texts: Everything but 'property' greyed out in calculate geometry

texts:
In calculate geometry, I need to change the units to square metres, however all but 'property' is greyed out. (See image)

I have tried this with editing on and off, same result.

What should I be doing differently?


texts: Cloud Shadow interferring with NDWI Analysis

texts:
I'm working on surface water identification in Northern Kenya. I'm using Landsat 5TM and Landsat 8OLI data in ENVI 5. I've masked clouds using the QA band (although have also tried the fmask tool- with default settings), and then run NDWI (Green-SWIR)/(Green + SWIR). However where the fmask has not masked all of the clouds and particularly shadows, I now get a lot of areas classified as water when actually they're clouds as shown in screen shot below. The image of the left is the cloud masked NDWI image, the middle is the NDWI image with a threshold set of -0.08, and the image of the right is the true colour image. It shows how the cloud shadow has not been masked and results in incorrect identification of water.

I have tried using different thresholds but there are cloud shadow pixels in the NDWI image that have the same value as water pixels. I am aware that I may not be able to 100% mask all cloud and cloud shadow. How can I improve my identification as I have a lot of pixels falsy identified as water; Is there a better way to make a distinction between cloud shadow and water?

Many thanks in advance.


texts: Identifying certain polygon from overlapping set

texts:
To clean up a large shapefile, I am looking for a way to do the following:

I have the ID of many polygons which need getting rid of stored in a list. These polygons are villages and a second shapefile has the watersheds to each of these, and they overlap.If I delete a village polygon I need the data assigned to it to be appended to the next village downstream. The next village downstream will belong to the next larger overlapping watershed, or the 2nd smallest one to overlap with the original village. The villages and watersheds are linked by an ID number.
First step would be: identify the ID of the 2nd smallest entirely overlapping polygon from the watershed shapefile which overlaps with the village identified in the list.

The "Identify" tool looks like it could do the trick but I can't find the script to it, does that exist for Python?

If anyone with experience in python for ArcGIS (ArcMap 10.2 Advanced license) or have tips do let me know!


texts:
While plotting the gadm map of Russia in R, I have a problem at the 180 degree meridian, where part of the Russian region is plotted separately. You can see it on the plot using the following code (it's the edge left and part is on the right). Can anyone suggest how to avoid this problem?

Here is a short version of the code for my question:

require(sp)

rus<-url("http://www.gadm.org/data/rda/RUS_adm1.RData")
print(load(rus))
gadm$regions = as.factor(1:88)
spplot(gadm,"regions")


texts:
I have a very large map file (.TIF) and a large DTM file (.GRD).
When I add them to my layers in MapInfo I get a .TAB file.

When my Map / DTM are visible MapInfo works and responds really slow and sometimes crashes.

I wanted to know, what can be done so that I could work with large layers with better performance and in one piece?

I have Global Mapper if some kind of conversion is needed.


texts: Remove small interstices on the edges of a polygone

texts:
How can I remove the white interstices on the edges (and interior) of a polygone or polyline, as represented in the first figure below? I want to get smooth and continuous edges, capturing the global external shape, as represented by the red line in the second figure below. Also, I want to keep the big white interiors.

I've tried the following methods, none of them being sucessful:

convex hull (vector, geoprocessing tools)
v.generalize.simplify (GRASS)
Simplipy plugin

I used different generalisation algorithms, such as Douglas-Peuker, which I think is not suitable for this problem.

I work with QGIS 2.18.4 runing on Windows 10.


texts: Performing spatial queries with the SODA API in Python

texts:
I am trying to fetch data from data.cityofnewyork.us with the Python client sodapy. It generally works as it should, but I just can't figure out how to perform spatial queries. For example, I would like to get the simplified geometry data by leveraging simplify(), but it doesn't return a geometry at all. How can I get this to work?
from sodapy import Socrata

identifier_build = 'qjdb-bxqz'
client_build = Socrata(domain, None)

results_build = client_build.get(identifier_build,limit=5,
                                 select="CNSTRCT_YR,simplify(the_geom,0.01) AS the_geom",
                                 where ="FEAT_CODE==2100")


texts: How to start developing with Jeppesen C-MAP?

texts:
I've just downloaded all the contents from the FTP of Jeppesen for their product called C-MAP CM93.
I've been reading the documentation contained in the SDK folders, but I can't figure out how to start developing in C# or Java with it.
Does anybody know what can I do? Where I can find a DLL? Any information on the Web?

Thank you!


texts: Openlayers 2.13 not calling correct XYZ tiles from couchdb?

texts:
End goal is to be able to access a pouchdb of map tiles on phone off-line that has been synced from remote couchdb of map tiles.

Openlayers is not calling the correct document from the couchdb to test before using pouchdb.  The Z parameter is correct, but the X & Y parameters are not close to what they should be when examining the database or original mbtiles data.

Example:

In plain files converted with mb-util from mbtiles database for Revelstoke trails, tiles are: z=11,x=351,y=684,685
In contrast Openlayers calls: z=11,x=53-59,y=89-92 for about the same region

Done so far

Built a couchdb of map tiles from a mbtiles database with Mapproxy
Built with both SW origin & NW origin
Projections: 900913/3857
Tiles display directly in futon as tile attachment
Pouchdb syncs to couchdb
Openlayers requests tiles outside the current tiles set

My non-working Openlayers code using XYZ:

var bc_orm_trails_couchdbLayerXYZ = new OpenLayers.Layer.XYZ("BC ORM trails couchdb local",
    ["http://localhost:5984/tiles_sw/SWoriginMercator-${z}-${x}-${y}/tile"], {
    spericalMercator: true,
    displayProjection: geographic,
    maxScale: "1000",
    minScale: "220000",
    maxExtent: BC,
    isBaseLayer: false,
    displayInLayerSwitcher: true,
    visibility: true
});

trails_mobileMap.addLayer(bc_orm_trails_couchdbLayerXYZ);

Tried pulling from local plain file folders created with mb-util to simplify -> doesn't work

var bc_orm_trails_plaintilesLayerXYZ = new OpenLayers.Layer.XYZ("BC ORM trails plain tiles local",
        ["./data/tiles/trails/${z}/${x}/${y}.png"], {
    spericalMercator: true,
    displayProjection: geographic,
    maxScale: "1000",
    minScale: "220000",
    maxExtent: BC,
    isBaseLayer: false,
    displayInLayerSwitcher: true,
    visibility: true
});

trails_mobileMap.addLayer(bc_orm_trails_plaintilesLayerXYZ);

This works from a remote php TileServer with same mbtiles data:

var bc_orm_trails_mbtilesLayerTileServer = new OpenLayers.Layer.XYZ("BC ORM trails",
        ["http://tiles-alpha.xxxxxxxxxx.xx/trails.mbtiles/${z}/${x}/${y}.png"], {
    type: "png",
    projection: spherical_mercator,
    displayProjection: geographic,
    //maxScale: "52000",
    //minScale: "220000",
    isBaseLayer: false,
    displayInLayerSwitcher: true,
    visibility: false
});

trails_mobileMap.addLayer(bc_orm_trails_mbtilesLayerTileServer);

I think I either have a problem with origin (SW vs NW) or possible a projection issue but have run out of possible solutions.

Ideas or possible solutions?


texts: Creating Esri Grid and populating cell values using Python/GDAL or ArcObjects?


texts: GDAL 1.9.1 ogr2ogr KML to PostGIS not incuding ExtendedData

texts: I have a KML file utilizing a Schema in the document and ExtendedData in each placemark. The ExtendedData does not seem to pass in conversion to PostGIS or shapefile.


texts: Which one is the correct Plate Carrée projection in QGIS?

texts:
I need to create a Plate Carrée map and I want to do so in QGIS. However I encountered a projection glitch in the Print Composer when I use EPSG:32662:

QGIS says that EPSG:32662 is deprecated and this gisSE answer suggests that also its successor EPSG:32663 is outdated and that EPSG:4087 replaced it. Spatialreference.org says that ESRI:54001 is the "World Plate Carree". Neither of the last three can be found in QGIS' CRS settings. Simply using unprojected WGS84 / EPSG:4326 doesn't look correct for my task.

Which one is the correct Plate Carree for QGIS and how can I set it?


texts: Exporting geographic coordinates in QGis as projected (flat) coordinates

texts:
I've created some polygons with Google Earth which I'd like to convert to my flat map's coordinate system.

I've loaded a few shapefiles into QGis of a county's roads, parcels, etc. The coordinate system for these points is projected using the "NAD83(HARN) / South Carolina (ft)" or "EPSG:3361" standard.

Next I load a kml file into QGis created using Google Earth Pro which is composed of a Polygon of geographic (lon/lat) points.

QGis is clever and has perfectly overlaid the geographic-points polygon directly where it should be on the flat/projected-points polygon. I can hover my mouse around and see the projected coordinate points where my polygon's exterior ring resides.

I want to export the points of my geographic-polygon into the projected coordinate system that the rest of my qgis map is based on. This is because all of my county data shapefiles are projected, that and PostGIS doesn't seem to allow me to do spatial joins of two geographic polygons...

Anyway, I right click on the imported vector layer (Google Earth-gen'd polygon) and click 'Save As'. (I've already enabled 'On the Fly CRS Transformation')

The plan is to export the polygon into projected coordinates, so I choose 'Selected CRS' and locate my proj4 code which is for SRID:EPSG-3361.

Exporting this does not accomplish my goal. The exported kml still contains geographic coordinates. QGis seems very capable of doing what I want to do, its overlaid two layers from entirely different coordinate systems. How might I go about transforming this clearly-understood (by qgis) polygon into the projected coordinate system upon which it apparently lies?

UPDATE 1:
- I created a second layer whose chosen CRS is the projected one, redrew a similar polygon on top of the imported one, and tried to export it, using, again, the projected coordinate system. IT ALSO got exported into geographic coordinates!

UPDATE 2:
I exported a shapefile of the county this time and IT TOO got exported into geographic coordinates!!! This tells a lot about the nature of the problem. Everything is being exported into geographic coordinates. Could it just be the fact that I'm using kml?


texts: ArcMap 10.2 crashing when opening .mxd files

texts: One of my coworkers emailed me about some issues with opening .mxd files in ArcMap 10.2, and through further questioning she revealed it actually started when she updated to 10.1 and has continued after her update to 10.2. Apparently when certain mxd files are opened (I'm guessing ones that were created before her upgrade to 10.1) in ArcMap it crashes with an error (fig. 1). If MXD doctor is run on the files they then work fine, I have attached the Diagnostic Report(fig.2).I thought maybe if they are just older 9.xx documents this could be the culprit but I read that those mxds should still open but will just be saved as 10.x map docs.  Any idea looking at the errors what the problem may be?


texts:
I have a map in a bitmap image format (.bmp) where each region is assigned an RGB colour. What I want is a list of every adjacent region, for every region.

So if region1 is adjacent to region2 and region3, I want this output in a text/csv file:

region1, region2
region1, region3

I can also deal with the output where we have duplicate adjacency links:

region1, region2
region1, region3
region2, region1
region3, region1

And I can also deal with an adjacency matrix:

(whatever), region1, region2, region3
region1, 0, 1, 1
region2, 1, 0, 0
region3, 1, 0, 0

And I can deal with the regions being named after their colour number. (probably the simplest thing to do) So the actual output would probably look like:

hFF0000, h00FF00
hFF0000, h0000FF

(or)

hFF0000, h00FF00
hFF0000, h0000FF
h00FF00, hFF0000
h0000FF, hFF0000

(or)

(whatever), hFF0000, h00FF00, h0000FF
hFF0000, 0, 1, 1
h00FF00, 1, 0, 0
h0000FF, 1, 0, 0

The output will then be processed in Excel.

This is the map, but there are others like it that I'd also like to process once I got the automation down:

I feel like what I'm looking for is a very basic feature for a geography software, and I've tried to find it in QGIS, but I didn't find it and I am confused by all the vocabulary I don't understand. What I remember from my last geography class is that we talked about rocks and coloured maps.

Like, I know what's a raster and what's a vector, but the software is constantly asking me for a coordinate system, and I don't know what to answer because this isn't a real world map and I don't care about the GPS coordinate. Also, every raster operation seem to want to operate on a single band, the bands are named after the RGB colour, so are all of those processes only working with 1 colour component of my pixels?? I can't deal with that, because I need all 3 colour components to know which region is which.

This map is the regions map for the Total War Three Kingdoms video game, taking place in China circa 200 AD, and I want the adjacency dataset because I'm building a modelling Excel workbook, and some buildings give a critical effect to their adjacent regions.

So my precise questions are :

Is this doable in QGIS?
If not doable in QGIS, can another free software do it? (free as in money)
And finally, how can I do this in said software?


texts: Migrating Defense Solutions objects from ArcGIS 9.3.1 to 10.0 gives Unable to cast COM object error?

texts:
I've got an ArcGIS Engine application that uses the IMADtedLayer interface to access VPF DTED catalog layers and adds them to a map. The code ran fine under v9.3.1. I just updated to v10, and now I get:

Unable to cast COM object of type 'ESRI.ArcGIS.DefenseSolutions.MADtedLayerClass' to interface type 'ESRI.ArcGIS.DefenseSolutions.IMADtedLayer'. ... No such interface supported

Here's the relevant C# code:

// geodbFWS is the IFeatureWorkspace object that contains the catalog
// sCatName is the catalog name (string)
ITable catTable = geodbFWS.OpenTable(sCatName);

if (catTable != null)
{
    IMADtedLayer aMADted = new MADtedLayerClass();
    bool bSetup = aMADted.Setup(catTable);
}

The error occurs on the "bool bsetup" line.

I've tried substituting new MADtedLayer(); for new MADtedLayerClass(); - with the same result. I've also uninstalled all ArcGIS components, cleaned out the registry and re-installed. I get the same error.

I'm using MSVisualStudio 2010, and the intellisense recognizes the IMADtedLayer interface, the MADtedLayerClass, and the .Setup method - as well as the fact that it should return a boolean.

I'm running the Windows 7 Professional Operating System. Has anyone else encountered this problem? Any insights will be greatly appreciated!


texts: Writing QGIS Plugin or function to find/pan/zoom to coordinate generated from number that has geolocation basis for its creation?

texts:
In my area, PINs (Parcel Id Numbers) are 10-digit numbers built from the coordinate of the centroid of the parcel (polygon).  I have no coding experience, but I know what I need is fairly straightforward and hope someone can help me.

Note- the 10-digit PIN is NOT a coordinate, it is an encoded coordinate and the actual coordinate would need to parsed out.

What I wish to do is create an input field labeled 'PIN Search' where a user can input a 10-digit PIN and after hitting enter the program would determine the coordinate and, ideally, select the nearest polygon at that coordinate and zoom/pan to its extent.  I would be happy if it would just zoom/pan to the coordinate at some predetermined scale.

I suppose I need to learn to write some python code.  I have no idea how to create an input field, I know what to do with the data to generate a coordinate, but then I would like to select the nearest polygon to that location and pan/zoom to its extent.  I would be happy to just pan/zoom to the coordinate, though.

Is this something I would need to create a plugin for or can it be done in some sort of macro?

EDIT
--How a PIN becomes a coordinate
CRS = NAD83 / North Carolina (ftUS)  EPSG = 2264

extents: xmin = 400000, ymin = 400000, xmax = 1300000, ymax = 1300000

PIN is 10-digit number of the form x1 y1 x2 y2 x3 y3 x4 y4 x5 y5
The coordinate is derived as X=x1 x2 x3 x4 x5 x6, Y=y1 y2 y3 y4 y5 y6 where x6 and y6 are = 0 and if x1 = 0,1,2,or 3 then 10 is added and if y1 = 0,1,2,3 then 10 is added.  Examples:  PIN = 9876543210 generates coordinate x=975310, y=864200:  PIN = 0123456789 generates coordinate x= 1024680, y= 1135790


texts:
I have a shapefile of 200 lines (all connected) with the coordinates x and y and a DEM raster. I would like to give the z value to the lines by taking it from the raster.

I have read several discussion on how to do the same operation with a points shapefile but I was wondering if would be possible to have the same results with lines. Maybe applying the z-value to a line vertex, or calculating the average elevation of the line length?


texts: GRASS v.net.iso does not include all lines

texts:
I want to define catchment area for bus stops in small city. As a result of v.net.iso in QGIS with GRASS I don't get all highways in that city. Here are my steps:

importing above layers into GRASS mapset
running v.net from GRASS plugin in QGIS (input: highways layer, operation: "New point is placed on each node...", arc layer: 1, node layer: 2, threshold: blank, assign unique categories: not checked
running v.net from GRASS plugin in QGIS (vectorInput: resultOfPrevStep, pointVectorInput: bus stops layer, operation: “Connect still unconnected points to vector network by inserting new line(s)”, arc layer: 1, node layer: 2, threshold: 100, assign unique categories: checked
running v.net.iso from GRASS plugin in QGIS (inputArcs: resultOfPrevStep, pointVectorInput: resultOfPrevStep, categories: 100,200,300,400,500,600,700,800,900,1000)

As a result I get two layers, one with calculated lines (yellow to blue gradient: correct), a second with unconnected lines (pink). What I am doing wrong?


texts:
I am using a shell script that is called by a batch file. The shell script contains the GRASS scripts that I am using for my program. Every time I use the r.mask script, it displays the error r.mask: command not found.

Can r.mask be called from a batch file? How can I use r.mask using a batch file and shell script?


texts: Cursor Select by Location Append Value to List

texts:
I have a point layer of Manholes and a point layer of start points for StormSewer pipes. I need to determine where Manholes intersect a start point of a StormSewer pipe (I used the Vertices to Points tool), grab the FACILITYID of the Manhole and then update a field (UPSTREAMFID) in my StormSewer pipe layer (different from the StormSewer start point layer)
I believe my workflow should look something like this:

use Search Cursor to Select by Location where Manholes intersect Storm Sewer start points

use getValue to store the FACILITYID of the Manhole that intersects the Storm Sewer Start Point in a list

use UpdateCursor to update the StormSewer  UPSTREAMFID field with the Manhole facility ID
a) an if statement needs to be used to ensure that the correct StormSewer line is updated with the Manhole facility ID (the StormSewer start point layer has a FACILITYID field that matches the StormSewer line it originated from)
b) something like if StormSewer Start Point facility ID = StormSewer facility ID update UPSTREAMFID with Manhole facility ID

Here is the code I have so far (I am struggling at this point):
import arcpy
from arcpy import env
arcpy.env.overwriteOutput = True
env.workspace = r"\\cold\redirect\Profiles\fultor\desktop\GIS\Projects\Stormwater Python\Shapefiles"

Manhole = r"\\cold\redirect\Profiles\fultor\desktop\GIS\Projects\Stormwater Python\Shapefiles\Manhole.shp"
StormSewer = r"\\cold\redirect\Profiles\fultor\desktop\GIS\Projects\Stormwater Python\Shapefiles\StormSewer.shp"
SS_Start = r"\\cold\redirect\Profiles\fultor\desktop\GIS\Projects\Stormwater Python\Shapefiles\SS_StartPoints.shp"

Manhole2 = arcpy.MakeFeatureLayer_management(Manhole, "Manhole_lyr")
SS_Start2 = arcpy.MakeFeatureLayer_management(SS_Start, "SS_Start_lyr")
StormSewer2 = arcpy.MakeFeatureLayer_management(StormSewer, "StormSewer_lyr")

field1 = 'FACILITYID'
field2 = 'UPSTREAMFA'

# Create Feature Layers
print ("Layers created!")

#Get all the unique strings from the first layer
FIDlist = []
with arcpy.da.SearchCursor(SS_Start2, 'FACILITYID') as cursor:
    for row in cursor:
        FIDlist.append(row[0])
print (FIDlist)

#Select each feature in sequence and get intersecting feature, then update
for FID in FIDlist:
    arcpy.SelectLayerByAttribute_management(SS_Start2, 'NEW_SELECTION', "FACILITYID = '{}'".format(FID))
    arcpy.SelectLayerByLocation_management(Manhole2, 'INTERSECTS', SS_Start2)
    with arcpy.da.UpdateCursor(StormSewer2, field2) as cursor:
        for row in cursor:
            row[0] = FID
            cursor.updateRow(row)
    #Clear the selections before the next iteration
    arcpy.SelectLayerByAttribute_management(SS_Start2, 'CLEAR_SELECTION')
    arcpy.SelectLayerByAttribute_management(Manhole2, 'CLEAR_SELECTION')


texts:
I'm no GIS expert and I hope I can explain my problem in English.

I have these polygons – and a border. Sometimes I need this exact definition of the polygons but if you want to have an overview map, it will look much better, if the polygons are expanded and there is no empty space between.

My idea was to use a Voronoi diagram to fill the empty spaces. I cut to original the Voronoi diagram with polygons, then I cut every polygon that is outside the border. In this step, I had quite a lot problems but, after some tries I got a very good result by using the point on surface function to check if the polygon is outside the border.

But I'm stuck in the next step – join the original polygons with the new created.

I try to use several PostGIS functions, sometimes combined. Base function that I used is st_union. I combined it with st_buffer and ST_SnapToGrid or st_simplyfiy and st_makevalid because I read that this gives good results.

But I have some of these errors in it (see picture) and as a result I need exact one polygon. Most of the results are multi-polygons.

I even try a PostGIStopology example, but I don’t get it working.
I try to script something in PHP to merge/join the polygon.

Some of them even complex one are exactly what if want but a lot of them don’t find the correct next neighbor point and the direction of the line is wrong on some parts.

Then I try something else more drastic ^^ I use st_union to merge the polygons as good as possible. To get only one polygon filter all polygon by the area and only use the biggest one corresponding to the original polygons. That one often have some wired line on the surface. I export the polygon with st_astext and spilt then by the brackets inside and also use the biggest one. The result is okay in higher zoom level, but if you zoom in use the white not used space.

Do you have any other ideas? At the moment, I use QGIS to visualize the output and I see that QGIS has a lot of function that may be able to do what I want (but I don’t know if it’s scriptable even if it looks like it can do that).

Note: Sorry I have to delete some images, because I'm just allowed to link two images


texts: Create points based on location of known points in a layer

texts:
I'm using ArcGIS 3D Analyst 10.5.0

I have a layer of points, we will call them points1, and I want to create a layer of points2 based on a set distance and direction away from points1. The distance will be a constant, while the direction is data that is included in the attribute table of points1.

Is there is a way in ArcMap to do all of these points at once?

Rather than creating each point using create feature?

I don't know Python, but I do know R and could potentially script this in R if I had a general idea of what commands to use (not much experience with GIS in R).


texts: Can a mask layer filter out color, rendering the area covered by the mask into grayscale?

texts:
Here I have a transparent mask layer which highlights an area of interest by making the basemap outside the area of interest a lighter color.

What I'd like to do is have the masking layer change the basemap to grayscale outside the area of interest. Is this possible?

The basemap is Bing Satellite. I'd rather not get into downloading the basemap, converting to raster, clipping the raster, etc. etc. I'm hoping for a graphic solution.

My question is related to and inspired by this question.


texts:
I am using the NOTIFY signal of PostgreSQL intensively to update a PostGIS-layer in the map (see https://oslandia.com/en/2017/10/07/refresh-your-maps-from-postgresql/).

Additionally to the channel 'qgis', I am sending a payload with the latitude and longitude of the newly added point.

Now, additionally to the update of the layer, I would like to pan the map automatically to the location given in the payload of the NOTIFY signal.

Question:

How can I reuse the NOTIFY payload of the qgis-channel of PostgreSQL in QGIS (e.g. within a Python-based layer-action).

Looking for any hints.


texts: Failure when trying to write/update some features in FME

texts:
By means of FME 2021, I am trying to update one column in a Point feature class inside FileGDB (FME: writer). Updates are comming from another feature class (FME: reader). I dont have unique ID but I do have x and y coordinates which are unique and match between reader and writer. Even though coordinates are identical in reader and writer, ONLY for some features is not working and pops-up:

No rows matched UPDATE query to table/feature class 'points' where x =
95330918.98604 AND y = 9535458.351159

Let's say for 50 features it is working but for some reasons fails in 3 or 4 features!
settings I used:

Feature opeation: update
Match columns: x,y

Do you have any idea how to do this?


texts: Create Shapefile/Feature Class with CSV that has Shape values?

texts:
I have a CSV with no lat/long values, but it does have Shape.STArea and Shape.STLength values. Is there a way to make this CSV into a feature class or shapefile? This data is supposed to be a polygon.

I have access to ArcMap 10.4 or QGIS 2.18.3


texts: Is there an api for world terrain/biome data?

texts:
I'm looking for an api where I pass a lat/lng, and get "Desert" or "Mountains" or "Ocean" back. Just the general biome of the area.

Most of what I've seen involves getting tiles back, so you can display a map, but that's not very useful for my use-case. I'm trying to build a mobile app that changes ui based off where you currently are.


texts: Showing selected feature in Attribute Table (Show All Records)?

texts:
I have a bunch of records in my attribute table (some of which are duplicates) I can select them individually and click the 'show selected tab' to see them, however it would be much more productive for me if when I selected a data point it would automatically go to it under the 'show all records tab' in the attribute table. This way I could automatically see the duplicates around it (by ordering my attribute data numerically).

As it currently stands when I make the select I have to manually scroll up and down etc to find it in the show all records' tab of the attribute table.


texts: Configuring attachments for ArcGIS Local Government "Citizen Service Request" App

texts: Has anybody configured the Citizen Service Request app that ESRI's Local Government has developed? I have it ready to go, but can't get the attachments to work properly, either pics or documents. I can't find a space in the js code for the service to be specified either...


texts: How to avoid error in "SAGA Strahler order"-process when running it on a merged terrain model file

texts:
I am trying to create a catchment in QGIS 3.8 following this tutorial: https://www.youtube.com/watch?v=2fvy9xtAHc4
However, my terrain model is originally divided into 4 overlapping tif-files. I can follow the steps in the video if I only use one of the four tif-files, but I wish to merge the four terrain models into one and analyse the the whole area. I have therefor used Raster > Miscellaneous > Merge to merge the four tif-files into one tif-file. I set the input and output nodata-value to zero, and everything seems to be working fine. There are no error messages, and the output layer seems to be correct.
Then I do the "Fill sinks (wang & lui)" process on the merge terrain model, and still everything seems fine. There are no error messages, and the output layer seems to be correct.
Then, when I try to do the "SAGA Strahler order"-process, it seems like the execution finishes but there is no output from it and I get the following error:

Loading result layers
The following layers were not correctly generated

 L:/temp/MFP/Elvenett.sdat

 You can check the 'Log Messages Panel' in QGIS main window to find more information about the execution of the algorithm.

I have tried to run this operation many times now: in a local folder, in a path without any space characters, only with temporary files but I always get the same result. I haven't found any similar problem online, except maybe this unsolved one (could be the same problem as me, but maybe not): Fill sinks and Strahler Order error in QGIS
EDIT:
I think this can have something to do with the terrain file size. I tried clipping the terrain model into a rectangle, and that magically did the work. At first I thought it had something to do with the shape of the terrain model, maybe it had to be shaped as a rectangle or another "regular" shape.
However, when I applied this theory to another terrain model I am working with, it did not help. I tried to clip out a very small and very irregular shape of that terrain model, and that did the trick.
So, is there a size limit to terrain models in the "SAGA Strahler order"-process? My terrain model, before cutting it, is "only" 150 MB. If this is too big, then that is a major limitation to this process.


texts:
I want to know if we can post the versions based on events through coding.

How we do batch reconcile and post through c# code?

Is it possible for automatic replication also?


texts: Updating attribute table of shapefile on display automatically using ArcPy?

texts:
I have a shapefile attribute table on display. After some processing, I have to right click click Refresh to get updated version of table. I am wondering if there is a way using arcpy to update it automatically every time I make an edit to shapefile.

Note: I have a dialog box with values which I change based on my requirement and press recalculation which gives me updated attribute table.


texts: How to reproject a raster file in QGIS with datum transformation?

texts:
I want to reproject a tiff raster file from projection A (ETRS89/UTM Zone 32) to projection B (DHDN/Gauss-Kruger zone 2). The two projections have a different datum so I have to use a datum transformation (NTv2 BETa2007 in my case).
I want to create a new raster file with projection B because I want to use it in another software. So reprojection on-the-fly is not what I need.

I have tried Raster > Transform (or similar - I do not know the English menu) but I can not see how to define a datum transformation here.

How can I do a permanent reprojection of a raster file in Qgis with a datum transformation?


texts: How to configure QGIS Sextante to use SAGA?

texts:
I have been trying to set up SAGA to work using Sextante in QGIS.  I have followed the instructions to configure SAGA in Sextante, with the following input in the path for the SAGA folder C:\Program Files\SAGA-GIS\saga_gui.

When I then try to use a SAGA algorithm I get this message:Saga folder is not configured. Please configure it before running SAGA algorithms.

When I go back to the Sextante configuration options, the path for the Saga folder is empty.  There is no information under the Sextante history log.

Anyone got any ideas about how I can resolve this issue?


texts: R: How to get a filename list of raster layers in a geodatabase?

texts:
In R, how do we get a file list of rasters in a geodatabase folder?

We know the case of vector layers as follow.

library(rgdal)
fgdb = "C:/path/to/your/filegeodatabase.gdb"
subset(ogrDrivers(), grepl("GDB", name))
fc_list = ogrListLayers(fgdb)

However, for raster layers, if only we know the raster layer name, we can read it by specifying the layer name as follow.

library(arcgisbinding)
arc.check_product()
a <- as.raster(arc.raster(arc.open(paste(fgdb,"/rasterlayername", sep=""))))

For more effective reading, I want to get the raster layers name list like the vectors layers case. How to do it??

RasList <- arc.open(fgdb)@children$RasterDataset
a <- NULL
for(i in 1:length(RasList)){
a[[i]] <- as.raster(arc.raster(arc.open(paste(fgdb,"/",RasList[i], sep=""))))
}


texts: OSM2PO Routing - reroute when traffic jam

texts:
is possible to force osm2po to reroute when traffic jam or other problem is detected on the route ? Or is it possible to tell osm2po routing to avoid some node/way even when it is best route ?

I'm not Java programmer so I prefer non Java solution


texts:
I'm trying to merge some DEMs downloaded from USGS. When I used ArcGIS 9.2 there was an expression to be used in raster calculator to do that. I used it again but, unfortunately, it doesn't work in ArcGIS 10.2.

What can I do?

I do not have the Production Contouring toolbar.


texts:
I'm trying to automate a procedure which involves a Dissolve. I've got the layer stored as a variable (layer_shp) and am trying to dissolve it through Python, as such:

dissolve_process = processing.run("gdal:dissolve", layer_shp, "", "memory")
layer_dissolve_shp = dissolve_process["OUTPUT"]
layer_dissolve_shp.commitChanges()
QgsProject.instance().addMapLayer(layer_dissolve_shp)

This is giving us the following errors:

After the dissolve_process = processing.run("gdal:dissolve", layer_shp, "", "memory"):

Traceback (most recent call last):
File "C:\PROGRA~1\QGIS3~1.4\apps\Python37\lib\code.py", line 90, in runcode
exec(code, self.locals)
File "<input>", line 1, in <module>
File "C:/PROGRA~1/QGIS3~1.4/apps/qgis-ltr/./python/plugins\processing\tools\general.py", line 105, in run
return Processing.runAlgorithm(algOrName, parameters, onFinish, feedback, context)
File "C:/PROGRA~1/QGIS3~1.4/apps/qgis-ltr/./python/plugins\processing\core\Processing.py", line 130, in runAlgorithm
context = dataobjects.createContext(feedback)
File "C:/PROGRA~1/QGIS3~1.4/apps/qgis-ltr/./python/plugins\processing\tools\dataobjects.py", line 72, in createContext
context.setFeedback(feedback)
TypeError: QgsProcessingContext.setFeedback(): argument 1 has unexpected type 'str'

Then, of course, when I try and commit changes:

AttributeError: 'str' object has no attribute 'commitChanges'

Think this should be quite an easy one but I'm not finding much about doing a dissolve. This isn't dissolving on a field, and I'm wanting it dissolved to memory and then adding into the project. This is QGIS 3.4.5

EDIT - this is the code I ended up using:

layer_dissolve_shp = processing.run("native:dissolve", {'INPUT':'Spray_SHP','FIELD':[],'OUTPUT':'memory:'})["OUTPUT"]
QgsProject.instance().addMapLayer(layer_dissolve_shp)

The layer exists as a memory layer on the QGIS layers list as "Spray_SHP" (note - this is not a shapefile - we're doing everything in memory and not saving to disk). This will just get that layer and dissolve it to layer_dissolve_shp, then add it as 'output'.


texts: How to create contour lines from point data using QGIS?

texts:
I'm using QGIS 2.4.0 and I have a shapefile of points with elevation data and X and Y coordinates. I'm using windows 7 OS.

The idea is to create a contour line SHP using my point data.

My question is: How do I start working on that issue? Is there a function? Is there a step by step process?


texts:
I am not sure how this is done but I Have a shapefile with state borders for a country. What I would like is to sort of delete/dissolve the state borders and be left with only the country border.   In the data attribute table, there are only polygons for the states and none for the country.
How is this done in Arcinfo or Global Mapper 10?


texts: Disconnect layer.selectionChanged when project is being closed or layer removed

texts:
On QgsProject.readProject I add a slot to some layer's selectionChanged. I tried disconnecting these slots in QgsProject.layersWillBeRemoved but here the project seems already closed, so I can't retrieve the layers to disconnect the slots.

There is a secont signature layersWillBeRemoved(self, layers: Iterable[QgsMapLayer]) but I failed to connect to it.

@pyqtSlot(list)
def layersWillBeRemovedHandler(self, layers):
  #disconnect

...

project.layersWillBeRemoved.connect(self.layersWillBeRemovedList) #without [list]

throws decorated slot has no signature compatible with layersWillBeRemoved(QStringList).

project.layersWillBeRemoved[list].connect(self.layersWillBeRemovedHandler)

throws there is no matching overloaded signal whether the handler is decorated or not with @pygtSlot(list).

So how do I solve this? Do I need to disconnect them at all, when the project is closed (or the layer is removed from the project)?


texts: Transform polygon to equivalent circle, QGis

texts:
I have a polygon and I'd like to find the radius of an equivalent circle from a particular point within the polygon. The polygon would have a large (infinite) number of points evenly distributed accross it. I am working in QGis so if there's a solution for this it would be great.

Or, to put it another way: What is the average distance from a point in a polygon to every other point within the boundary of the polygon? A circle, I believe, would give the shortest distance, but a rectangle with one long edge and one short edge would give the longest distance.

Edit 1:

@Whuber - For your edit 2 this is what I picture. I can't see how this is not the area under b between x=1 and x=0.


texts: Cannot open .pkinfo or add saved maps from ArcGIS online in ArcMap

texts:
I am trying to use this basemap in ArcMap (v10.8): https://www.arcgis.com/home/item.html?id=33ea4550c8144e66847d902e4766c2f7
If I click on "open in ArcMap", it downloads an item.pkinfo file. But, opening this file just opens a blank map, with the appropriate layers in the table of contents, but with red "x"s next to their names. If I save this map on my online profile and try to open it from ArcMap via file -> ArcGIS online -> my maps and data, it just loads the hillshade but not the streets.
What's going on here? Any ideas?


texts: Removing WMS layers from mxds using ArcPy?

texts:
I'm trying to remove WMS layers from an mxd using arcpy. I can successfully remove "MapServer" type services using the below code, but not WMS?

Here is an example mxd (zipped) with two different WMS layers that I'm trying to remove: here

Also, I've posted on the ArcGIS Resources forum too: here

Code:

import arcpy, os

inFName = "wms.mxd"
outFName = "wms_removed.mxd"
pth = r"c:\temp"
inFULLPTH = os.path.join(pth,inFName)
outFULLPTH = os.path.join(pth,outFName)
mxd = arcpy.mapping.MapDocument(inFULLPTH)

for df in arcpy.mapping.ListDataFrames(mxd):
    LyrList = arcpy.mapping.ListLayers(mxd, "", df)

    for lyr in LyrList:
        if lyr.supports("SERVICEPROPERTIES"):
            if lyr.serviceProperties["ServiceType"] == "WMS":
                arcpy.mapping.RemoveLayer(df, lyr)

mxd.saveACopy(outFULLPTH)
del mxd

Any ideas? Is arcpy.mapping.RemoveLayer the best approach ? I have hundreds of mxds to crawl through and want to remove WMS layers from them.


texts: Why are conversion results from proj.4 inaccurate?

texts:
I want to convert a Gauss-Kruger coordinate (EPSG::31467) into a WGS84 coordinate (EPSG::4326) using proj.4. The problem is that the output from proj.4 is very, very inaccurate.

Consider the Gauss-Kruger coordinate (3509502, 5399989). I converted it with three different online converters. All results agreed to within four decimal places: (48.7373, 9.1281). Moreover, I verified with Google Maps that (48.7373, 9.1281) is indeed the correct WGS84 coordinate (I know where it should be).

However, proj.4 converts the Gauss-Kruger coordinate into (48.7384, 9.1292), more than 140 meters away from the correct solution. Why?

$ echo 3509502 5399989 | invproj +proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +ellps=bessel +datum=potsdam +units=m +no_defs -f "%.4f"
9.1292  48.7384

The results are the same for version 4.8.0 and the newest version 4.9.3.

The solution (from @JGH):

$ echo 3509502 5399989 | cs2cs -f %.12f \
> +proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +datum=potsdam +no_defs \
> +to +proj=longlat +datum=WGS84 +no_defs
9.128151599064  48.737368376923 51.846968705766


texts: Geopackage Error - 'is mounted and in WAL mode, this combination is not allowed'

texts:
I am trying to find a replacement for shapefiles when sharing geospatial data. I am testing the geopackage format and with QGIS 2.18.2, I can create a geopackage and add additional layers to the geopackage using the 'save as' function.

I have created geopackage files in previous versions of QGIS (2.14.9), which creates a single gpkg file. This is attractive over shapefiles because it is one file instead of many. However, when I create a geopackage using 2.18.2 I get additonal two files - gpkg-wal, gpkg-shm.

The issue is when I try to open the geopackage created in 2.18.2 using ArcCatelog 10.4 I get an error "underlying DBMS error"..."is mounted and in WAL mode, this combination is not allowed".

My question is how do I create a single geopackage file that is not in WAL mode using QGIS?


texts: How do I add TIFF files to one layer/raster catalog?

texts:
We have several TIFF images (already tagged with position information) we want to load as a single layer into the map. It's the first time I dig into this raster/topology thing, but at least I managed to create a layer with one image as an IRasterDataset:

Type factoryType = Type.GetTypeFromProgID("esriDataSourcesRaster.RasterWorkspaceFactory");
IWorkspaceFactory wsFactory = (IWorkspaceFactory)Activator.CreateInstance(factoryType);
IRasterWorkspace2 rasterWS = (IRasterWorkspace2)wsFactory.OpenFromFile(@"c:\sourcefiles\");
IRasterDataset rds = rasterWS.OpenRasterDataset(@"example.tif");
IRasterLayer rasterLayer = new RasterLayerClass();
rasterLayer.CreateFromDataset(rds);
_map.AddLayer(rasterLayer);

This works well but I need to add all TIFF files to this raster layer. So I assume that IRasterCatalog could suit my needs. But my attempt to use it fails already on creation of the catalog object (taken from sample):

Geoprocessor geoProcessor = new Geoprocessor();
CreateRasterCatalog createRasterCatalog = new ESRI.ArcGIS.DataManagementTools.CreateRasterCatalog();

//Set parameters
createRasterCatalog.out_path = outPath;
createRasterCatalog.out_name = catalogName;
createRasterCatalog.raster_spatial_reference = map.SpatialReference;
createRasterCatalog.spatial_reference = = map.SpatialReference;

//Execute the tool to create a raster catalog
geoProcessor.Execute(createRasterCatalog, null);

On Exceute() I get a COMException, which isn't very helpfull. I tried to set the SpatialReference to IGPCoordinateSystem objects. I made sure that the output path exists.
Still the COMException raises. Has anyone a clue for me?

Is this the right approach at all?


texts: Allocating house with 2 coordinates in Leaflet.GeoJSON?

texts:
I have 2 coordinates of the house.

Is it possible to allocate this house with other color in Leaflet?

Is it possible using GeoJSON?

Or need to draw a polygon?


texts: Polygon gradient fills / tint bands in QGIS

texts:
Is it possible to use a gradient fill on a polygon in QGIS?

I would like to have my polygons with a colour at the edges, fading to clear/white as it moves away from the edge. The majority of the centre of the polygon will be white or clear.

It is something people I know who use Illustrator can do. I wondered if there was a way to emulate it in QGIS.

Here is an example. I would like something like the purple shading.

Update:

This method suggested below does work, however it produces some annoying extra bits for my polygons. For example, see the image below:

I am unsure how to remove these. I cannot use a mask as I need to be able to see the base mapping and other data outside of this polygon on layers that are below this layer. Using a mask would obscure all of these layers.

I have included a link to Underdarks blog post about this, which gives more detail, and contains the code to set this up.

http://underdark.wordpress.com/2011/08/08/creating-a-gradient-fill-for-polygons-in-qgis/