texts: Why a shapefile of a geological data is not displaying the fault lines in QGIS?

texts:
I am a new user of QGIS 2.4. and can not get a shapefile of geological data to display correctly. This shapefile contains polygons, fault lines, dips, volcanic structures, shoreline types, geological boundaries between other layers. Characters set: utf8. Reference systems: UTM WGS 84, Zone 17 North.

This layer was digitized from the geological map of Panama, 1:250,000.

When using ArcGIS the fault lines are displayed, but not when using QGIS. This shape contains the following extensions: .shp, .shx, .dbf, .prj, .qix., shp.xml.

What might be the issue?


texts: How to access array inside of Leaflet LayerControl?

texts:
When you make a Leaflet LayerControl you first make an array of markers, then pass the array to the LayerControl constructor and add the LayerControl to the map, like this:

myLayerGroup = L.layerGroup(arMarker).addTo(map);

However, I want the markers in the Layer Group to change when the user pans the map. I have coded the appropriate event handler, like this:

map.on('dragend',this.onDragEnd);

I then pass the new extent of the map to the server side via Ajax, query the database and get back the markers that are within the new extent.

What I want to do next is add the new markers and delete the old markers that are not within the new extent. But in order to do that I have to access the array that I passed to the LayerControl. I don't want to just delete all the current markers and add the new markers because many of the new markers will also be the current markers.

Is there a way to access that array within the LayerControl, or do I just have to delete the LayerControl and start again?


texts: Joining QGIS extent table with shapefile

texts:
I used QGIS to find the extents of a list of features from a census boundaries shapefile.

I used Vector - Research tools - Polygon From Layer Extent, which creates a new shapefile with a table of extent values for each feature.

However, I need to join this new table with the old one so that I have the extent values as well as other necessary fields from the original source shapefile. And the new extent shapefile doesn't have a join field.

Is there a way to join a polygon extent layer with its original source layer? I can try adding a new column to both tables with ordered numbers, but are the extent tables produced with the same ordering as their source?


texts:
I am trying to convert gpkg to shapefile using ogr2ogr command with python.
After running this line of code I got a message saying:

ogr2ogr.main(["", "-f", "ESRI Shapefile", "out6.shp", "kucni_broj.gpkg"])

Warning 1: Value '20А' of field out6.kucni_broj parsed incompletely to integer 20.

Then I added ENCODING parameter and converted it again using:

ogr2ogr.main(["", "-f", "ESRI Shapefile", "out6.shp", "kucni_broj.gpkg", "-lco", "ENCODING=UTF8"])

After that I got message:

ERROR 1: Failed to create field name 'kucni_broj_id': cannot convert to UTF8

Is there any way I can fix this issue?


texts: Calculating annual mean precipitation from DAYMET_V3 data in Earth Engine

texts: I am calculating annual mean precipitation from DAYMET_V3 data in Google Earth Engine. Here is an example of mean annual temperature from ERA5. How do I have annual mean precipitation from DAYMET_V3 data like this example?


texts: GPS device taking pictures and saving heading + HDOP?

texts:
I'm looking for a device that can do all of the following:

Take picture
Save lat-long GPS position + heading (angle to North) of picture [so need tilt-compensated 3D compass], in exif or other format (save tilt as well if possible)
Save precision-related GPS (NMEA) fields, as many as possible (at least HDOP, other if possible)
If possible: display OSM maps (better if I can put my own tiles)

Ideally HDOP and the like should go into both exif of the pictures and gpx track, but post-processing (gpscorrelate or the like) is OK.

I would need to know this ASAP to buy the device and learn how to use it: it's for an 80-day alpine OSM mapping party, starting early June 2012! See here:
http://www.outdoormaps.org
http://wiki.openstreetmap.org/wiki/ViaAlpina

There are already some data on this long high-distance trail, but they lack accuracy estimates. Heading will make the pictures much more valuable.

Thanks!

Mayeul Kauffmann

(PS: I've seen this question: Device for Logging GPS and Accelerometer Data?  but I need heading and HDOP)


texts: Obtaining the data of one selected feature in QGIS using Python?

texts:
I have a problem/question regarding python and QGIS. Is it possible to get the data of one selected feature in QGIS using python, and if so how?

In my example I have a layer which shows some parcels for urban planning. If I now select one of the parcels I want to get the data for exactly this parcel. This data should be read out by an python script and then added to an automatically created PDF file. The last step is easy, but I didn't get the data from this selected feature.

It is possible to get data from the active layer, but is it also possible to get the data from one selected feature of this layer?

Hope you can help me.
Regards


texts: Display Shapefile with leaflet and layer control

texts:
I try to implement a leaflet pluging, for display a local hosted shapefile. The display of the shapefile work, but i want to add a layer control (for togle shapefile layer).
the plugin link : https://github.com/calvinmetcalf/shapefile-js

the demo link :http://leaflet.calvinmetcalf.com/#3/32.69/10.55

I want to implement the layer control on demo page

The shapefile "Fr_adm.zip" is displayed, but no layer control.


texts: Generating sequence of numbers based on existing attribute and spatial configuration of polygons using ArcGIS Desktop?

texts:
I have come across a very specific problem and I'm sure there must be a way to get around it with Python in ArcGIS, but can't figure it out on my own.

The situation is as follows:

The initial data consisted of a Fishnet set of equal sized 150x150m polygons stack side-by-side and covering a total area of 10.2x6.9km (that is 68x46 rectangular polygons = 3128 total polygons). Each polygon received a unique ID (labelled TID) based on a quasi-random pattern.
The next step required the creation of a Fishnet of equally sized 25x25m polygons for the given area, which resulted to the production of a fishnet of 408x276 polygons. Using the Spatial Join command, the TID was transferred to the new smaller rectangular polygons and the  resulting table has 112608 rows.
In the resulting table, there is a column titled TID, which is of type Long and carries the numbers from the underlying layer. More specifically the numbers run from 1 to 3128, and as each polygon was dissected to smaller ones there are now 36 rows associated with each number (therefore 3128*36= 112608 rows).
In the same table, there is also a column titled GID and is of type Text, as it is meant to carry the ID of each specific polygon, ideally using some form of concatenation per the following formula: TID + "-" + Sequential numbering of the specific rectangle. This column is currently empty.

What I want to do is the following (the result was created manually editing properties):

For this to work it needs the following:

Assign sequential number based on specific pattern (South to North for each column) for each TID.
Concatenate the TID with the sequential number.

Is this possible with Python and, if so, would it be possible for someone to guide me through?

Currently, I'm a novice learner of Python, hence most concepts are beyond my grasp.


texts: Programmatically populating shapefile field with polygon areas in OGR?

texts:
I'm trying to programmatically populate a shp field with polygon areas:

import sys
import ogr

ds = ogr.Open( 'tttttttttt.shp', update = 1 )
if ds is None:
    print "Open failed./n"
    sys.exit( 1 )

lyr = ds.GetLayerByName( "tttttttttt" )
lyr.ResetReading()

field_defn = ogr.FieldDefn( "Area", ogr.OFTReal )
lyr.CreateField(field_defn)

for i in lyr:
    # feat = lyr.GetFeature(i)
    geom = i.GetGeometryRef()
    area = geom.GetArea()
    print 'Area =', area
    lyr.SetFeature(i)
    i.SetField( "Area", area )

ds = None

But I get an empty values in the field...


texts:
Can anyone cast light on the accuracy assessment of the Mathew Hansen dataset?

What will be the proper way of using this data for research?

I am asking this since I know that any of the remote sensing classification work without ground truthing/accuracy assessment is useless - is my understanding correct?

I do not find anything on the accuracy assessment of this data on the internet. Even not the methodology for this data though a bit very uncomplete discussion at Methodology for Hansen classification of Global Forest Watch?


texts: Creating Definition Query tool in ModelBuilder that uses in-line variable

texts:
I'm having trouble using an in line variable derived from the Parse Path tool in ArcGIS ModelBuilder to create a definition query for a layer file.

I'm importing an xyz file to a feature class, and using Parse Path to obtain %Value%.
I would like to set a defintion query to something like the one shown below.  I have a short Python script that was used to manually enter the def query, but would much rather use %Value% to derive it.


texts: Filtering out irrelevant point features in QGIS2Web Leaflet export

texts:
I'll make some simple maps with data-points with Leaflet or QGIS2web but it seems like it is big drawing blue transparent dots at the points I have no data. I just want to see my data-points (like the orange ones in the example).

I even didn't find the blue-dots in the js, CSS or HTML code.

I want to make maps with different layers one can switch on/off. For that, I need "empty" data points that aren't shown.


texts: Preventing a new/duplicate column to bo created in a spatial join in QGIS?

texts:
I wanted to have an common (between two shapefile's) column filled in by certain criteria by creating a buffer first (based on that criteria) and the do some spatial join steps. However I noticed that the column's of the shapefile were added to my new buffer shapefile after doing the spatial join to the buffer.

I assumed that, because the columns did not share a name this was a problem. So I made sure to that the name of respective columns were to same. I redid the proces, but the same thing happened, so I deleted all fields that were not necessary (after making a copy of the shapefile). I redid it again and the same thing happened.

To illustrate

How can I prevent this, what I am doing wrong here?

Edit: Go here if you want to know what I originally wanted: Filling a field in a layer using a spatial relation with another layer/shapefile?


texts: Error using projectRaster in R: Error in if (maxy == miny) { : missing value where TRUE/FALSE needed

texts:
I am trying to reproject a raster from Cylindrical Equal-Area (Lambert). Central Meridian: -160. Datum: WGS 1984 to latitude and longitude But I'm getting multiple error messages.
The original raster looks like this when plotted with rasterVis:

My code is as follows:

library(raster)
#read raster
raster(r)

# Define new Proj.4 spatial reference
sr <- "+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0"

# Project Raster
projected_raster <- projectRaster(r, crs = sr)

#But I get the following errors:
Error in if (maxy == miny) { : missing value where TRUE/FALSE needed
In addition: Warning message:
In rgdal::rawTransform(projfrom, projto, nrow(xy), xy[, 1], xy[,  :
53 projected point(s) not finite

I am not sure whether this is a problem with the raster or my code?
Please help!


texts: Using some 'temporary variables' in QGIS Geometry generator symbol layer type expression?

texts:
Making use of the Geometry generator symbol layer type, I draw rectangles of dimensions @nv_bg_w (width) and @nv_bg_h (height) (project variables) to line features, either at coordinates text_x, text_y (attributes, if not NULL) or alternatively at the center of the line by the following expression:

geom_from_wkt(
    'POLYGON((' ||
    COALESCE("text_x", x(point_on_surface($geometry))) ||' '||  COALESCE("text_y", y(point_on_surface($geometry))) || ','||
    (to_real(COALESCE("text_x", x(point_on_surface($geometry)))+ @nv_bg_w )) ||' '||  COALESCE("text_y", y(point_on_surface($geometry))) || ','||
    (to_real(COALESCE("text_x", x(point_on_surface($geometry)))+ @nv_bg_w  )) ||' '||  (to_real(COALESCE("text_y", y(point_on_surface($geometry))))- @nv_bg_h ) || ','||
    COALESCE("text_x", x(point_on_surface($geometry))) ||' '||  (to_real(COALESCE("text_y", y(point_on_surface($geometry))))- @nv_bg_h ) || ','||
    COALESCE("text_x", x(point_on_surface($geometry))) ||' '||  COALESCE("text_y", y(point_on_surface($geometry)))|| '))'
)

As can be seen x(point_on_surface($geometry)) and y(point_on_surface($geometry)) occur very often. At least in this simple example, this makes the code harder to read than it would have to.

So my question is: Is there a way to store the latter expressions in some temporary variables, something like (pseudocode):

@mx=x(point_on_surface($geometry))
@my=y(point_on_surface($geometry))
geom_from_wkt(
'POLYGON((' ||
... #and so on


texts: Export RGB data (from a tiff image) in excel format

texts:
I have a tiff file having 3 bands R, G, B. I used the following program in R to obtain the pixel values for each band separately.

r<- raster("lena.tiff", band = 1)
writeRaster(r, filename = "lenadd.asc", format = "ascii", datatype= 'INT4U')

It gave me an ascii file containing the pixel values for R band. Similarly, by changing the band value to 2 and 3, I can get the other two files. However, I want to export the data in excel format. But I suppose, writeRaster command doesn't support xlsx format. So, can someone please suggest me any alternative of writeRaster command so that I can export my RGB data in excel format?

I am new to R so pardon me if the question is too trivial.

Edited :

class       : RasterLayer

band        : 1  (of  3  bands)

dimensions  : 512, 512, 262144  (nrow, ncol, ncell)

resolution  : 1, 1  (x, y)

extent      : 0, 512, 0, 512  (xmin, xmax, ymin, ymax)

coord. ref. : NA

data source : C:\Users\lenovo\Documents\lena.tiff

names       : lena

values      : 0, 255  (min, max)


texts: D3 can't render GeoJson file (Simple structure in Montreal, Canada) created in QGIS

texts:
I am attempting to render SVG and attach it to a div with id="map". I can do this successfully with GeoJson which I acquired from this tutorial: http://maptimeboston.github.io/d3-maptime/#/1

Excerpt:

var neighborhoods_json = {
"type": "FeatureCollection",

"features": [
{ "type": "Feature", "properties": { "Name": "Roslindale", "density": 5.5800 }, "geometry": { "type": "MultiPolygon", "coordinates": [ [ [ [ -71.1259, 42.2720 ], [ -71.1257, 42.2723 ], [ -71.1256, 42.2724 ] etc...

I am rendering this using the following code:

var width = 960, height = 720;
var colors = d3.scale.category20();

var svg = d3.select('#map').append('svg')
  .attr('width', width).attr('height', height);

var albersProjection = d3.geo.albers()
  .scale(190000).rotate([71.057,0]).center([0, 42.313])
  .translate([width/2, height/2]);

var geoPath = d3.geo.path().projection(albersProjection);

svg.selectAll("path")
  .data(neighborhoods_json.features)
  .enter()
  .append('path')
  .attr('fill', function(d, i) {
  return colors(i);})
  .attr('d', geoPath);

This works great, as Boston is in the USA and Albers is an appropriate choice of projection.

The data which I need to render, however, is a simple movie theater traced over a location in Montreal, Canada using QGIS ( http://imgur.com/Dddkci2 ). I saved the layer as GeoJson using a setting of EPSG:4326, WGS 84. This results in the following file being added to my project:

var theater_json = {
"type": "FeatureCollection",
"crs": { "type": "name", "properties": { "name":     "urn:ogc:def:crs:OGC:1.3:CRS84" } },

"features": [
{ "type": "Feature", "properties": { "id": 1, "d3Css": "l0r0" }, "geometry": { "type": "Polygon", "coordinates": [ [ [ -73.433681281940878, 45.445198056579976 ], [ -73.433681281940878, 45.445201593863992 ], [ -73.433916566065932, 45.445339547768 ], etc...

I am unable to render this in my project. The svg is appended, but regardless of whether I choose to append it as 'path' or 'polygon', the segments are all 0x0 and no visualization is created. I am confused whether or not a projection should be used, as I have seen suggestions of utilizing a projection of "null" in case the data was already outputted with a projection from QGIS. I have also seen suggestions that Canadian locations should use the mercator projection. In any case, here is how the code currently stands:

var svg2 = d3.select('#map').append('svg')
  .attr('width', width).attr('height', height);

var mercProjection = d3.geo.mercator();
var theaterPath = d3.geo.path().projection(mercProjection);

svg2.selectAll('path')
  .data(theater_json.features)
  .enter()
  .append('Polygon')
  .attr('fill', function(d, i) {
    return colors(i);})
  .attr('d', theaterPath);

Regardless of whether I use mercator, null, or any other variety of projection (stereographic, azimuthal equidistant, azimuthal equal area), I am never able to render any visualization. The paths/polygons, however, are still being created with segments of 0x0 in every case. Such an observation seems to suggest that the wrong projection is being chosen, as I observe the same behaviour when I choose the wrong projection for the Boston data.


texts:
While migrating to the newest stable Tomcat Server, it seems that latest versions of GeoServer (2.18, 2.19RC) is no more compatible with Tomcat 10, likely due to the fact that The Java package has changed from javax.servlet to jakarta.servlet as part of the migration to Servlet 5.0 API. See https://tomcat.apache.org/migration-10.html.
When deploying the Geoserver .war to Tomcat 10, we get
org.apache.catalina.core.StandardContext.listenerStart Error configuring application listener of class [org.geoserver.GeoserverInitStartupListener]
java.lang.NoClassDefFoundError: javax/servlet/ServletContextListener
Are there any plans that GeoServer will run with this setting in near future?
What alternatives do we have?


texts: QGIS strongly offsetted symbols for categorized layers not visible on layer list (layer panel)?

texts:
I am constructing a symbology for telecom analysis, which so far is mainly point and line features.

So far it is working very well except that I have to offset the symbols to be able to avoid superposition and cluttering, and that also is OK except that the offset has to be of such a magnitude that the symbols are not longer visible on the layers panel.

Since there are many symbols, it is very important that the user of the map can see quickly the meaning of the symbols but since the spacing is fixed on this panel all  that is visible in layer panel are the names of the categories and blanks where the symbol should be visible.

Is it possible to modify (increase) the spacing between categories in the layers panel?

I am using QGIS 2.2.


texts:
I already installed the version 1.9 of Qgis and when i runned e Script python on i could

this error   'NoneType' object has no attribute 'children'   however it works fine on qgis 1.8

could someone help me please !!


texts: How to Select polygons within a feature based on rules applied to FID

texts: I want to select all the polygons of a vector grid which have a FID multiple of N (e.g. N=5). Is that possible?


texts: I am looking for AIS Shipping Route Data. I have not had much luck finding it in either SHP format or even in CVS format does anyone have a real good source of finding this type of data? To add to my problems I need only a specific area of coverage? It is starting to look like I may need to purchase the Data.


texts: Reloading cache on open Attribute Table using ArcPy?

texts:
I used python to update a shapeFile through cursors.  The Shapefile is updated, but there could be an open Attribute Table.

Manually you can hit "reload cache" to update the Attribute Table and see the changes.  I would like the python script to reload the cache, so the user does not have to remember to do this.


texts: How can i get the center point for victoria

texts:
I am following this tutorial. But as per this tutorial the the center point is located in Netherlands. Now I need to change this for Victoria, Australia.
Could you please tell me how can i change it. Using BRutile Lib for sharpmap.

Transform transform = new Transform(new PointF(629816f, 6805085f), 1222.992452344f, this.Width, this.Height);


texts: Identifying Joshua trees with LiDAR data?

texts:
I'm working on a LiDAR project to determine where Joshua trees are located within a specified study area. Due to the vegetation cover being so sparse, there are really on 2 canopy species there, which are Joshua trees and cottonwoods. I believe this to be a relatively easy LiDAR analysis due to very limited species richness in the canopy. My approach has been to create a bare earth raster (DEM) and then a 1st return raster. I would then subtract the bare earth from the 1st return raster to create a vegetation raster. I would be able to easily remove any noise (e.g. power lines, buildings) by using a basemap for verification. Because the client wants to see all Joshua trees >=12ft, I would simply reclassify the vegetation raster. By doing this, I should be able to see all tree species, which should be Joshua trees, within my study area. With this methodology, I've only been able to create the highest bare earth point locations in the study area, which is not what I want and am confused as to why this is my output.

This is the methodology I've followed in ArcMap:

Create Bare Earth Layer

Create a las dataset of the selected study area with the Create LAS Dataset tool
Make a las dataset layer with this layer with the Make LAS Dataset Layer tool
  a. Select 2 (ground) from the Class Codes
Convert this layer to raster with the LAS Dataset to Raster tool.

Create Vegetation Layer

REPEAT STEPS 2 AND 3 AGAIN BUT SELECT 1ST RETURN UNDER Return Values (optional) WHEN USING THE MAKE LAS DATASET LAYER TOOL.
Subtract the Bare Earth raster from the 1st Return Raster with the Minus tool

 1st Return (raster) – Bare Earth (raster) = Vegetation Layer

Use the Reclassify tool to determine what is 12 ft and greater:

       Classification: Natural Breaks (Jenks)

       Classes: 2

       Break values: 3.66, 10.725098

Does anybody have any experience with this and might be able to provide some tips/pointers in where I might be going wrong?  If people know of better methodologies, I am open to ideas!


texts: ArcGIS Javascript BasemapToggle requires double click on iOS Cordova

texts:
I am developing a Cordova application using the ArcGIS Javascript SDK.  One minor annoyance I have run into that I can not quite figure out is why the BaseMapToggle & GeoLocate options require a double click on the iOS application.  The web demo shows that this works on a single click.  A minor issue, but not very initiative when most everything on mobile is a single touch event.  Here is my code:

//NOTE: gMapVaults is a multi-dimension array declared  above
var map,toggle,geoLocate;
require([
    "esri/map",
    "esri/dijit/BasemapToggle",
    "esri/dijit/LocateButton",
    "esri/geometry/Point",
     "esri/symbols/SimpleMarkerSymbol",
     "esri/symbols/Font",
     "esri/symbols/TextSymbol",
     "esri/symbols/SimpleLineSymbol",
     "esri/Color",
    "esri/graphic",
    "esri/layers/GraphicsLayer",
    "dojo/domReady!"
], function(
  Map, BasemapToggle, LocateButton, Point, SimpleMarkerSymbol, Font, TextSymbol, SimpleLineSymbol, Color, Graphic, GraphicsLayer
)  {
      map = new Map("map", {
        center: [gMapVaults[0][2],gMapVaults[0][1]],
        zoom: 18,
        basemap: "satellite"
      });

        toggle = new BasemapToggle({
            map: map,
            visible: true,
            basemap: "streets"
          }, "BasemapToggle");

        geoLocate = new LocateButton({
            map: map
          }, "LocateButton");

        toggle.startup();
        geoLocate.startup();

        map.on("load", function() {
            for (var i=0;i<gMapVaults.length;i++){
                    var p = new Point(gMapVaults[i][2], gMapVaults[i][1]);
                    var s = new SimpleMarkerSymbol(SimpleMarkerSymbol.STYLE_CIRCLE, 20,
                        new SimpleLineSymbol(SimpleLineSymbol.STYLE_SOLID,
                        new Color([255,0,0]), 1),
                        new Color([0,255,0,0.25]));
                    var g = new Graphic(p, s);
                    g.setAttributes({
                        name: gMapVaults[i][0]
                    });
                    map.graphics.add(g);

                    var t = new TextSymbol(gMapVaults[i][0]);
                    t.setOffset(0,12);
                    var p2 = new Point(gMapVaults[i][2], gMapVaults[i][1]);
                    var g2 = new Graphic(p2, t);
                    map.graphics.add(g2);
            }

        });

});


texts: Selling Story Maps via ArcGIS Marketplace?

texts:
I just went to the ArcGIS Marketplace to see if I could buy the latest and greatest Story Map from around the world.  This might be some ripping yarn supported by maps.

When I searched the ArcGIS Marketplace for 'story map' the only hit I got was an Esri Story Maps page telling me that:

  Creating story maps is included in ArcGIS Online with both ArcGIS
  Online Public Account and ArcGIS Online Subscription

Using the ArcGIS Marketplace, and only the capabilities included within an ArcGIS Online Subscription (which I have by virtue of my ArcGIS for Desktop maintenance), is it possible to create a Story Map, place a paywall on it, and potentially see revenue trickle in, or what might be a showstopper that would prevent it being that easy?

This is a simplified version of something that I would like to do.  Assume that:

a Story Map could be authored which is of enough interest to attract an audience who might be willing to pay to view it.
I am an authorized marketplace provider even though that is not currently the case.

This question is about trying to understand how to participate in the ArcGIS Marketplace by using a Story Map which is listed in My Contents of ArcGIS Online as a web map app, as the simplest example of an app that might be made available for sale.


texts: QGIS Advanced Python Field Calculator Syntax

texts:
I am attempting to calculate the difference between two azimuths using the following code in the QGIS Advanced Python Field Calculator :
def azimuthDifference():
Model_Azimuth = <Model_Azimuth>
RMSTraffic_Azimuth = <RMSTraffic_Azimuth>

diff = RMSTraffic_Azimuth - Model_Azimuth
absDiff = abs(diff)

if absDiff == 180:
    return absDiff
elif absDiff < 180:
    return diff
elif RMSTraffic_Azimuth > Model_Azimuth:
    return absDiff - 360
else:
    return 360 - absDiff
value = azimuthDifference()

When I run the model that it is a part of, I receive the following error message:
File "<string>", line 17, in <module>
File "<string>", line 5, in azimuthDifference
TypeError: unsupported operand type(s) for -: 'NoneType' and 'float'

When I set the Model_Azimuth and RMSTraffic_Azimuth variables to integers (e.g. Model_Azimuth = 45), the code executes without any error, so I believe my problem is related to pulling the attribute values for the feature using <Model_Azimuth>.


texts: Importing data from Excel and maintaining the time format

texts:
I am trying to import a file from Excel which is a GPS data set of consecutive points at 15 minute intervals. When I import the file and save it as a shapefile the time column is not read properly and GIS converts to 00:00:00 and then further down the column shifts to 01/01/1900.

I have tried to format the column in Excel as either hh:mm, hh:mm:ss, dd/mm/yyyy hh:mm or no number format and each has the same end result. The imported table shows the correct times but as soon as the data table is converted to a shapefile to be displayed the time format is lost.

I have never noticeably had this problem when importing data previously and I did not need to format the time. Now it seems that no formatting method will work properly.

How can I format the time column to retain the time values when imported and saved as a shapefile?


texts: Calculating Feature Overlap Area within Single Shapefile

texts:
I have looked throughout this site to find a method to calculate the overlap area of multiple polygons within a single shapefile but haven't found one that works for my case. I've tried the intersect tool but that doesn't seem to accomplish my task as it requires multiple input shapefiles.

Above I have included a image of my shapefile and the polygons. Each polygon represents an animal's territory and has several attributes assigned to it (age, sex, etc). In this example, there are some polygons (animal territories) that don't overlap and others that overlap multiple polygons. What I would like to do is to calculate the total overlap area (dark purple) within the shapefile. Eventually I would like to know the identities of all the animals that overlap each other. In the entire shapefile I have 300+ territories and multiple years of data so I need an approach that would be somewhat efficient.


texts: Valid values of mollweide projection in QGis

texts:
i just imported a vectorlayer (without any reference points) to qgis as mollweide projection to save it to longlat. What i am wondering about is a mathematical problem: Not every x,y coordinate on the screen is valid when i put it to the inverse formular for this projection. So how do i find the valid area so i can resize my vectorlayer properly?

Many thanks in advance an sorry for mistakes, im german
MX117


texts: Openlayers 3: Plotting WKT format polygons on single layer from array


texts: JFileDataStoreChooser.showOpenFile does not filter by extension

texts:
Running GeoTools 15.1 on Ubuntu

Using this code:

...
import org.geotools.swing.data.JFileDataStoreChooser;

public class GeoJson2Shape {

   public static void main(String[] args) throws Exception {

      File shpFile=null;

      System.out.println("Starting");

      File geojson = JFileDataStoreChooser.showOpenFile(".GeoJson", new  File("."), null);
...

I get this window:

No filtering but displaying all files. Is that OK?


texts: Importing shapefile problem, unknown geometry type in QGIS

texts:
I have used ArcGIS extensively, but am at a new job with only access to QGIS. I am trying to import  a few shapefiles and am having some problems. First, I opened the shapefiles in esris free browser, and eveything opened without a problem.

I open the files in QGIS and they are loaded browser panel showing that they are polygon and polyline, respectively, but no features are displayed on the screen. Further, when I click zoom to features, I am taken to the right place, which indicates that something is there.

I ran the check geometry validity and got an error of Unknown geometry Type 0. However, I looked in the metadata and it correctly displays its geometry type.

Data can be found at the link below:
Click "project documents" tab, and then download "other/2d benchmarking data"
I am looking at the shapefiles in "Test 5". I tried shapefiles in the other "test" folders and those seemed to be working fine:

http://evidence.environment-agency.gov.uk/FCERM/en/Default/FCRM/Project.aspx?ProjectID=669DA341-AE12-494F-9E53-29C6633FA0CD&PageId=a0fe6dfc-506a-452c-9bff-a7ec06b4e6b0


texts: What causes these ServiceArea artifacts?

texts:
I'm building a bus service coverage tool, using ESRI's javascript API with NetworkAnalyst Route and Service Area tools on ArcGIS Server 10.2.1 doing the heavy lifting.

Workflow:

Identify 2 or more 'stops'
Route between them, in order
Generate service areas for each stop

Most of the time this works well, but I'm getting artifacts in the service areas for some stop locations:

The green diamonds are stops.
The grey areas are service areas (returned from the ServiceArea task)
The pink lines are the 'output lines' (returned from the same ServiceArea task)

If I move stops 4 and 5 around, I get different artifacts:

(including one for stop 8, which I didn't move!)

If I set ServiceAreaParameters.overlapPolygons = true, I get no artifacts:

Does anyone know what causes these artifacts, and if there's some setting (other than overlapPolygons = false, because I explicitly want the result produced by the non-overlapping request) that stops it from happening? Is it just a bug in Server 10.2.1? Has it been fixed in a later version?

stRoutePlanner._serviceAreaTask = new ServiceAreaTask(stRoutePlanner._serviceAreaTaskURL);

var saParams = new ServiceAreaParameters();
saParams.facilities = new FeatureSet();

$.each(stRoutePlanner._stopLayer.graphics, function(idx, p) {
    saParams.facilities.features.push(p);
});

saParams.defaultBreaks = [500];
saParams.overlapLines = false;
saParams.overlapPolygons = false;
saParams.doNoLocateOnRestrictedElements = false;
saParams.outSpatialReference = new SpatialReference({ wkid:102100 });
saParams.outputLines = NATypes.OutputLine.TRUE_SHAPE_WITH_MEASURE;
saParams.outputPolygons = NATypes.OutputPolygon.DETAILED;
saParams.useHierarchy = false;
saParams.returnFacilities = false;
saParams.restrictionAttributes = [];
saParams.restrictUTurns = NATypes.UTurn.NO_BACKTRACK;

stRoutePlanner._serviceAreaTask.solve(saParams, stRoutePlanner.processCoverageResults, tris.errorHandler);


texts: Calculating Stream gradient in ModelBuilder - Shapefile iteration?

texts:
I'm working on some project (in ModelBuilder) and of its element is calculating stream gradient.

I have shapefile with rivers, each segment have specific length. I want to calculate stream gradient in my model using ModelBuilder, but I have some problem with that.

I've created points on both ends of each river segment using Feature Vertices to Points, then to those points added elevation from DEM using Extract Values To Points.

It looks like that:

Points referring to specific segment have the same Object_ID and Length (column form river shapefile):

Now I'd like to calculate stream gradient using formula:

(segmentStartElevation-segmentEndElevation)/segmentLength

In other words if Object_ID is the same for two points I have to take elevation (here RASTERVALU) of point with lower FID subtract from elevation of point with higher FID and divide by Length.

I've tried to figure it out, but don't know how to make those calculation.

Any ideas how to do that?


texts:
I know this is a very basic question, but I've been looking for an answer and haven't been able to find one, so here it goes. I am trying to build a long query to select features that have different values in the same field, in order to make it as agile as poossible. All in one selection.

Example: In a given attribute table, in the field "H1" select features that have values "31p", "34n", "61f" and "31q" in it.

"H1"  =  '31p'  OR  '34n'  OR  '61f'  OR  '31q'
"H1"  =  '31p'  OR  "H1"  =  '34n'  OR  "H1"  =  '61f'  OR  = '31q'

None of them work. Should I use "IS" logical operator instead of "="?


texts: Enabling LWGEOM for OGR tools (QGIS 2.8)?

texts:
I am doing some scripted operations on a Spatialite DB using OGRINFO, OGR2OGR etc. in the QGIS 2.8 OSGEO4W shell.

Unfortunately, LWGEOM functions such as ST_Split() are unavailable and lwgeom_version() returns NULL.  GEOS functions appear to work fine.

Is there an easy was to make LWGEOM functions available for the OGR tools?


texts: PyQgis Intersect line with polygons. Result is not Correct (some of polygons result not intersec the line)

texts:
I have a track (Line) and I want to know wich countrys/cities (Polygon) intersect this track.

I try in two ways:

a) With PyQgis - I use QgsSpatialIndex().intersects(track.geometry().boundingBox())

            for track in layer_track.getFeatures():
                 # Find Features
                 index = QgsSpatialIndex()
                 index = QgsSpatialIndex(layer_to_intersect.getFeatures())
                 find_result = index.intersects(track.geometry().boundingBox())

b) With Qgis console - Join Atributtes by Location

And the results was in both cases are:

As you can see in the image there are polygons that are not intersected by the line and appear in the search results.

The correct result should be as in the following image:

I understand what happens, and that is that to do the search, the line becomes a rectangle, and the result of the operation shows the polygons that intersect this rectangle.

And this is not the result of what I am looking for.

I am thinking of converting the line to points and seeing what points intersect with the polygons and thus have a satisfactory result.

When I work with Arcgis I use Spatial_Join to do this, and the results was satisfactory.

Any idea to help me?


texts: Converting Coord extend Lambert II (WGS84) to Geographical coordonate (lat lng) using ESRI API

texts:
I am using the esri API and I am facing an issue to convert Coordinate extend Lambert II to Geographical coordinate.

  For example :  x : 633951 y : 1959059.

I have try to use webMercatorUtils like that :

webMercatorUtils.xyToLngLat(X, Y);

But it returns the wrong coordinate.

I also tried to convert by the class Point like that :

let mypoint = new Point(this.bassinCoord[i].xLieu, this.bassinCoord[i].yLieu, new SpatialReference({ wkid: 27582}));

I use the wkid for France : cf : http://resources.esri.com/help/9.3/arcgisserver/apis/rest/pcs.html

Here is my full code :

        while (i < this.bassinCoord.length) {
          let mypoint = new Point(X, Y, new SpatialReference({ wkid: 27582}));
          let radius = radius;
          let normalizedVal = webMercatorUtils.xyToLngLat(X, Y);
          let circle = new Circle({
            center: mypoint,
            radius: radius
          });
          let graphic = new Graphic(circle, symbol);
          gl.add(graphic);
          i++;
        }

I do not want to use an external library other than Esri (if possible).


texts: Query astropy or skyfield to get planets positions into a QGIS layer

texts:
I want to have a layer in QGIS where I give a date/time input and it will query astropy to get the sun, moon and planets position at that specific moment.
https://docs.astropy.org/en/stable/coordinates/solarsystem.html
https://rhodesmill.org/skyfield/


texts: how to download World MODIS image from ESRI ImageServer?

texts: I need to download the following dataset that is available online in ESRI website: https://sampleserver3.arcgisonline.com/ArcGIS/rest/services/World/MODIS/ImageServer. I can view the map in arcmap as instructed but I need to download the original raster so that I can further process. This layer has 12 files of 12 different months.


texts:
I have a shapefile which contains 10.000 polygons and I want to know how many of them falls on a specific land use. Unfortunately, the land use map that I have to use is a raster map and not a vector map.

I created a mask from the raster map with 0-1 values for the land use I'm interested in, and I thought to convert the mask to a vector file: this way I can obtain the information that I want with simple vector processing tools.

However, I'm struggling to convert the raster to a polygon (maybe it's because it is a large raster, 10940 rows x 11404 columns with 20m pixel size). I've tried several methods, namely:

Raster to Vector tool in QGIS;
rasterToPolygons function in R;
gdal_polygonize using the OSGeo4W shell);

This far, 1 and 2 literally crashed my computer after a while; 3 is running since 20 hours and it's around 70% complete, but after all this time I don't even know what to expect and I certainly do not have the time to re-run the script.

Is there a faster way to convert this raster?


texts:
I would like to know how to remove clouds (create a cloud mask and its shadows) for Sentinel 2A data, using R (Windows).
I have seen Information about this on Python or Matlab, but I work so far only in R.


texts: QGIS plugin: Set the size of the plugin when loaded and change the name

texts:
I have two quick questions with the plugin I built:

(Red circle) How can I set the size of the plugin window when it loads into QGIS? The default size when loaded is a little bit too small and not all tabs are visible.
(Yellow circle) Is there any method to change the name of the plugin? I change the code of the Sava_attribute plugin (in the tutorial) to make my own but I would like the change the name of the plugin.


texts: draw centerline of road polygon in python

texts:
I'm trying to draw a center-line of this road element but I'm not able to make it with the PostGIS function ST_ApproximateMedialAxis()

Also tried the python centerline library but it giving me error

Code :-

for key in list(self.dct_polygons.keys()):
poly_geom = self.dct_polygons[key]
centerlineObj = Centerline.Centerline(poly_geom, self.dist)
self.dct_centerlines[key] = centerlineObj.createCenterline()

Error :-

vor = Voronoi(border)
File "qhull.pyx", line 2494, in scipy.spatial.qhull.Voronoi.__init__
File "qhull.pyx", line 354, in scipy.spatial.qhull._Qhull.__init__
scipy.spatial.qhull.QhullError: QH6154 Qhull precision error: Initial
simplex is flat (facet 1 is coplanar with the interior point)

For the below polygon this is not working

Even for few polygons I have been able to draw centerlines but which was not up to mark.


texts: Why is polygon used as range in GIS and maps?

texts:
I am enroute to create an app that finds nearest venues. I already have another SO question on this regard.

I don't have any experience with such features, and some of the concepts are not being very obvious to me. So, I have been thinking that I am not taking things on right direction.

I have seen in number of places with both GIS and in case of mobile towers, the range is calculated in polygons. I have also been reading a bit about Haversine Formula, and it appears it does not use Polygon, it is basically two-point distance (from my understanding).

Why use Polygon for reference ?
Why is circle not used for this ?
How does Haversine formula integrate with this polygons concept ?

#Copied : Original question in SO as refered by Nate


texts: What alternative to Voronoi polygons should I use for my QGIS problem? 2d Clustering?

texts:
I have a table with all the hotels in Paris. The columns are the name of the hotel, the geocode and the neighborhood (arrondissement) of Paris that hotel is in.

I've plotted all these points in GQIS and then selected all those hotels that have neighbourhood = 1st arrondissement (highlighted in yellow).

My problem is that it seems my data has some points that have the neighborhood value = 2nd arrondissement within the area where the 1st arrondissement is - as you can see there are orange points within the yellow points.

What would be the best alternative to Voronoi polygons in this case? Perhaps some kind of 2d clustering algorithm for which there is a QGIS plugin?

UPDATE: Here's where I am now:

So soon I'll have the whole of Paris covered with these overlapping polygons. My final question is, how do I take the 'averages' of these polygons so that I'm left with a series of non-overlapping, interlocking polygons (kind of a 'mosaic' with no spaces between the polygons) so I have complete coverage with no overlap or gaps?

I guess where two polygons are overlapping I want the border to me the average of the two.


texts: How to export layer to shapefile that includes all the data in ArcGIS Pro?

texts: I'm trying to export an entire layer in ArcGIS Pro into a shapefile that I can save. The steps I have taken are: right click on the layer under contents, navigate to data, click export features. Here I choose the location I want it to save, and the name of the shapefile. Then I click run. This all works fine. But when I open this shapefile there is only one row available to see under attribute table, which is wrong. Anyone know how to fix this?


texts: Converting X and Y to lat/long using UTM WGS-1984 48N projection system

texts: I am using street directory API and it is using the UTM WGS-1984 48N projection system. I only need to convert static values! I am mostly using JavaScript but at the moment I'm alright with any solutions!


texts: Obtain percentage of overlap with another multipolygon in R

texts:
I have two MULTYPOLYGONS, one where a row is a election district*year pair, and the other is a collection of flood events.
 district.df
      year ward_ons     cycle                       geometry
    1 2007       E1   NA-2007 POLYGON ((527370.8 183470.7...
    2 2008       E1 2007-2008 POLYGON ((528891.1 182192.6...
    3 2009       E2   NA-2009 POLYGON ((370294.2 414678.7...
    4 2010       E3   NA-2010 POLYGON ((375025.4 414992.1...
    5 2011       E3 2010-2011 POLYGON ((375150.8 410809.8...
    6 2018       E3 2011-2018 POLYGON ((373286.3 414364.5...
    7 2007       E4   NA-2007 POLYGON ((373168.6 411597.8...
    8 2010       E4 2007-2010 POLYGON ((374783.2 406209.4...

Flood data:
    flood.df
    Simple feature collection with 8 features and 2 fields
geometry type:  GEOMETRY
dimension:      XY
bbox:           xmin: 317656.2 ymin: 90783.2 xmax: 546460.6 ymax: 631125.7
projected CRS:  OSGB 1936 / British National Grid
  year            name                       geometry
1 2007      River 2007 POLYGON ((359637.7 268239.7...
2 2007       Tank 2007 POLYGON ((325444.1 92717.57...
3 2008  Yorkshire 2008 POLYGON ((318550.7 103058.8...
4 2009 Flood East 2009 POLYGON ((541472.6 112593, ...
5 2010  Occurence 2010 MULTIPOLYGON (((545863.4 11...
6 2012      Storm 2012 POLYGON ((473637.4 103927, ...
7 2011      Flood 2011 MULTIPOLYGON (((524617.6 42...
8 2017      River 2017 POLYGON ((393387.6 631125.7...

What I am trying to do is to get a column in the election district multypolygon which is the share of its area that overlaps with any of the polygons in the flood multipolygon with the same value for year.
This is what I tried (influenced by Extracting intersection areas in R)
# Create function
overlap.fraction <- function(election, flood, yeari){
  flood.year <- flood[flood$year == yeari,]
  election.year <- election[election$year == yeari, ]
for (i in 1:nrow(election.year)){
  int.i <- as_tibble(st_intersection(i, flood.year))
  int.i$affected <- st_area(int.i$geoms)
  affected.by.county(paste0(yeari)) <- int.i%>%
  dplyr::group_by(code) %>%
  dplyr::summarise(affected.area = sum(affected))
}
}

And then tried it
i.1990 <- overlap.fraction(district.df, flood.df, yeari = 1990)

But keep getting:
Error in UseMethod("st_intersection") : no applicable method for 'st_intersection' applied to an object of class "c('integer', 'numeric')"

I tried different alterations to the function and the loop inside it. Ideally, I would even get a function or loop that would do it at once instead of having to do it for every year, as is the case for what I am trying here. I am doing it this way because I am getting closer to getting it done this way, but I understand this is not the most effective way.


texts: Creating an Add-In to click on an ArcFM function

texts: I would like to create an add-in to press ArcFM save every five minutes. Is there any way to do this?


texts: Fitting known size polygons into irregular polygons using ArcGIS Desktop?

texts:
I have a task that sounds simple, but I just can't get through the last step. I need to find suitable areas where wind farms that have a size of 500x150 meters, could be placed. The restrictions are slope - best areas are form 0-2 degrees and 2-3 degrees are worth considered - and a 3 km buffer from the selected location. Data used here was SRTM as this is preliminary study. My workflow was to buffer the locations, clip the DEM, create slope, reclass into 3 categories and convert to shapefile.

My categories of reclass are as follows:

First category - 1 = 0-2 degrees - green;
Second category - 2 = 2-3 degrees - orange;
Third category - 3 = all other data - red.

The issue that I am facing is that I would like know many 500x150 meters polygons can I fit in the green area, without overlapping the orange and red areas and without overlapping the polygon itself. The end result should be a shp with the possible positions for the aforementioned 500x150 meters polygon.

I should also specify that the 500x150 meters polygon can have any position/orientation and that I have access to ArcGIS Desktop 10.6/ArcGIS Pro 2.2/Global Mapper 19 and QGIS and my coding skill relate to very simple tasks that I replicate from other users.

I had a look at Checking if polygon fits inside another polygon using ArcGIS Desktop? where a similar things is described, but with circles and using ETGeowizard, also an older version of Arcmap.

Also, Packing Polygons within polygon using ArcGIS Desktop? is almost what I need, but the answers are just to above my level.

EDIT: the purple area is an existing possible location drawn by hand. I am  trying to figure out the maximum number of 500x150m polygons that can be mosaic-ed into the green area. They don't necessarily have to touch each other.


texts: Feature Labels wrong with GeoJSON data source

texts:
EDIT:

I found out that the attribute table changes when zooming, my GeoJSON output does not change however.

I am working on a QGIS project where I've imported points as a GeoJSON data source from a web server. When I label the points ind QGIS and select the name field as the label field. I get the correct result initially.
When I zoom, the labels will change to some other field value from GeoJSON.
I've tried several QGIS versions and also changed the GeoJSON output for the field positions to be static.

Please help!!


texts: Use jsPDF to generate the PDF of a map on HTML using a GeoServer service

texts:
I've been bugging people here for two weeks now to no avail.

I have a project linking maps of a local GeoServer service, what I want is that when I press a button a PDF is downloaded of whatever is on the map. Quite simple, but I haven't been able to get it to work or get any help online.

This is my code

https://jsfiddle.net/dw8kb737/

This is what the HTML displays, and what the PDF should look like.

Looking at examples from all these sources

Export PDF from OpenLayers map with jsPDF Export PDF from OpenLayers map with jsPDF https://stackoverflow.com/questions/25698741/download-html-div-into-pdf-using-jspdf https://openlayers.org/en/latest/examples/export-pdf.html https://gist.github.com/ka7eh/88761650efd3425080035e8535230d15 https://github.com/MrRio/jsPDF/blob/master/examples/basic.html https://github.com/mpetroff/print-maps https://jsfiddle.net/crabbly/kL68ey5z/ http://jsfiddle.net/canvasjs/cm1qyk2L/

Did not yield any results whatsoever with me.

I've been working on this for days, why is it not working, what am I missing?

Obviously the code there doesn't work as it doesn't even have a button as I was trying a lot of different things.


texts: Reduce processing time in Erase function using Arcpy in IDLE

texts:
I am processing 18,000 feature classes (saved in a file GDB), each with ~500 records (polygon features) in ocean areas, and I am running a script to remove area overlapping with the coast, using the Erase function. The script is running successfully, but literally taking days to complete. Some files take seconds to process (I'm assuming this is due to the fact, which I've verified, that they are not overlapping with coastline) while some feature classes take over an hour to process.

I am working locally on my machine, running the script in IDLE, with ArcGIS 10.1 (all licenses) on a 32-bit Windows 7, 3.33GHz Intel DuoCore, 4GB RAM computer.

I'm learning the ropes with Python and would really appreciate tips on my script if there are steps I can take to speed up processing time. Thank you!

# Import modules
import arcpy
import sys
import os
from arcpy import env

arcpy.env.overwriteOutput = True
print arcpy.env.overwriteOutput

#Erase coastline from selected buffers.
workspace = "C:\\DIR\\points_buff\\"
arcpy.env.workspace = workspace
erase_output = workspace+"erased.gdb\\"
erasedirectorypath = "C:\\DIR\\points_buff\\selectedbuff.gdb"
fileList = []

for dirpath, dirnames, erases in arcpy.da.Walk(erasedirectorypath,
                                              datatype="FeatureClass",
                                              type="Polygon"):
    for erase in erases:
        fileList.append(os.path.join(dirpath, erase))
        filenameparts = erase.split('_')
        if int(filenameparts[1]) > int(1967):

            try:
        #Set variables for Erase tool
                in_features = erasedirectorypath+"\\"+erase
                coastline_feat = "C:\\DIR\\Basic_NE.gdb\\Coastline_Eastern_US"
                out_erase_feat = erase_output+erase[:-4]+"_er"
                arcpy.Erase_analysis(in_features, coastline_feat, out_erase_feat, '')
                print arcpy.GetMessages()
                print "Coastline has been erased from " + erase + " and saved to " + erase_output

            except:
                print arcpy.GetMessages()


texts: qgis2web Popup Hyperlink Based on Field with URL

texts:
I have a string field called PIN. And I have a string field called PIN_SOA which has URLs.

I'm creating a webmap using qgis2web 3.16. I have popups enabled, and of the two mentioned fields I only want to see PIN in the popup, not PIN_SOA. Here's how it looks now:

Is there a way to make PIN a hyperlink that opens the URL of the associated PIN_SOA field? For example, I want it to look like this:

I'm referencing these similar questions.
QGIS - Convert a text field to a clickable hyperlink, hotlink, or URL?,
Configuring URL link for popup window directly in qgis2web plugin?


texts: Carto DB - How to manipulate the drawing order of point data?

texts:
I’m doing a bubble/choropleth combo map with point data (cities) using CSS in CartoDB, and I want my biggest bubbles (that correspond with population size) to draw on the map first so that the map is easier to read and the infowindow feature is easier to use. I’ve assigned IDs to my data according the population size so that the city with the largest population has carto_id 1 and so forth (that’s how I would do it with ArcMap, because it draws ID 1 first). I’ve also tried to write my CSS so that the largest bubbles would be drawn first (I have little CSS experience). The CSS is below. I can’t figure out what order it’s drawing the bubbles in. How can I fix this?

#tm_db001_merge2{
  marker-fill: #F11810;
  marker-line-color: #FFF;
  marker-line-width: 0.5;
  marker-line-opacity: 0.6;
  marker-opacity: 0.6;
  marker-comp-op: screen;
  marker-placement: point;
  marker-type: ellipse;
  marker-allow-overlap: true;
  marker-clip: false;
  marker-multi-policy: largest;
}
#tm_db001_merge2 [ population <= 2880546] {
   marker-width: 62.0;
}
#tm_db001_merge2 [ population <= 2000000] {
   marker-width: 52.0;
}
#tm_db001_merge2 [ population <= 1000000] {
   marker-width: 43.0;
}
#tm_db001_merge2 [ population <= 750000] {
   marker-width: 35.0;
}
#tm_db001_merge2 [ population <= 500000] {
   marker-width: 28.0;
}
#tm_db001_merge2 [ population <= 250000] {
   marker-width: 22.0;
}
#tm_db001_merge2 [ population <= 100000] {
   marker-width: 17.0;
}
#tm_db001_merge2 [ population <= 75000] {
   marker-width: 13.0;
}
#tm_db001_merge2 [ population <= 50000] {
   marker-width: 10.0;
}
#tm_db001_merge2 [ population <= 25000] {
   marker-width: 8.0;
}

#tm_db001_merge2 [ percent_inner_city <= 100] {
   marker-fill: #e31a1c;
}
#tm_db001_merge2 [ percent_inner_city < 75] {
   marker-fill: #fd8d3c;
}
#tm_db001_merge2 [ percent_inner_city < 50] {
   marker-fill: #fecc5c;
}
#tm_db001_merge2 [ percent_inner_city < 25] {
   marker-fill: #ffffb2;
}


texts: Getting 64-bit background geoprocessing running in ArcGIS 10.5?

texts:
I've got fresh installation of ArcGIS 10.5 on my machine and was wondering from where and how can I get 64-bit background geoprocessing.

Help pages are not explaining anything regarding installation.

In previous versions it was (I think) located on the installation media (similarly as pointed to in this or this questions for 10.2):

My IT department informed me that this is not the case any longer and in some newer thread on geonet forum I found information that for version 10.3 this is downloadable from my.esri.com

(most likely with an institutional account..) but I cannot find any details regarding 10.5 version.

Has anyone managed to get and install it?


texts: Displaying multiple stacked layers in ArcMap

texts:
I am looking for a way using Layout View in ArcMap 10.3 to display multiple layers in this manner:

I would call this a "map stack", yet I couldn't find a way to do this or how the displaying method is called.


texts: Problem Importing shapefile into postgresql with simple geometries

texts:
I have a Line Type shapefiles to be imported to postgis database on my postgresql, I'm using the SHP loader, I checked the "Generate simple geometries instead of MULTI Geometries" options so that the data will be routable. But it gives me an error message like this:

Importing with configuration: jalan_tangsel_test1, public, geom, C:\ms4w\apps\coba1\data\test1\jalan_tangsel_test1.shp, mode=c, dump=1, simple=1, geography=0, index=1, shape=1, srid=0
Shapefile type: Arc
PostGIS type: LINESTRING[2]
We have a Multilinestring with 19 parts, can't use -S switch!
Shapefile import failed.

But if I unchecked the simple geometries options the shp import will be success, please help me

regards


texts: Definition query to extract last updated record

texts:
I am struggling to build a definition query using ArcGIS to select last updated record in relation to unique CADID (core data example below).

CADID is unique, updates recorded additional record in the FileGDB with same CADID with recent date and time stamp on it. I need to extract only unique record on (CADID) with last update.

Core dataset contains all records to keep track of historical updates but for current processing need only the latest updated record as follows.

I have heaps of datasets, and updates are happening to all datasets, manual selection is a real headache.

I applied different tricks without success. I can supply sample data to run the query.


texts: Get layer coordinates with mouse click ArcMAP 10.2

texts:
Im trying to get the layer coordinates with a mouse click in ArcMAP. I did this in QGIS and im not finding anything on how to get this working in ArcMAP. Here is the QGIS version:

def selectPoint(self):
        self.pointEmitter = QgsMapToolEmitPoint(self.iface.mapCanvas())
        QObject.connect(self.pointEmitter, SIGNAL("canvasClicked(const QgsPoint, Qt::MouseButton)"), self.selectNow)
        self.iface.mapCanvas().setMapTool(self.pointEmitter)

I got this working by creating a new plugin in plugin wizard and Choosing a "Tool" in the Toolbar.


texts: Show fields in infowindow from my custom layer

texts:
I want to display two fields from the POIData layer. What I get now is ${ToegangsNr} ${Baanvak} in text. Does anyone know how to do this in an DynamicMapServiceLayer?

POIData = new esri.layers.ArcGISDynamicMapServiceLayer("https://gis.strukton.com/ArcGIS/rest/services/POI_Data/MapServer", { id: "POI", "visible":false });
map.addLayer(POIData);
dojo.connect(map, "onClick", function(evt) {
    map.infoWindow.resize(175, 125);
    var symbol = new esri.symbol.SimpleMarkerSymbol().setStyle(esri.symbol.SimpleMarkerSymbol.STYLE_X).setSize(12);
    symbol.outline.setWidth(4).setColor("blue");
    dojo.connect(map, "onClick", function(evt) {
        map.graphics.clear();
        var graphic = new esri.Graphic(evt.mapPoint, symbol);
        map.graphics.add(graphic);
        map.infoWindow.setTitle("${ToegangsNr}");
        map.infoWindow.setContent("
Baanvak: ${Baanvak}
	");
        map.infoWindow.show(evt.screenPoint, map.getInfoWindowAnchor(evt.screenPoint));
    });
    dojo.connect(map.infoWindow, "onHide", function() {
        map.graphics.clear();
    });
});


texts: ArcGIS Calculate mean NDVI and assign it to polygons attribute table

texts:
I am new to ArcGIS.

I have a group of polygons with attribute tables and I have generated an NDVI map [Edit: by NDVI map i mean a File Geodatabase Raster Dataset with 10980x10980 resolution, 1 band, cell size 10x10, format FGDBR, 32bit, without compression with values from -1 to +1 ] from the whole area. I now need to calculate the mean NDVI value for each region of interest/polygon and then assign the value calculated to the tables accordingly.

How can I do that?


texts: Finding outermost border of set of geomertries (circles) using PostGIS?

texts:
I have a (large) set of geometries (circles) that overlap.!
I would like to find the nearest location outside of this group from my current location.
In this example my location is the black dot. In my scenario I would like to get out of this group of circles and find the edge / border of this group, do I can draw a point on (or just outside of it), this would be the red dot.

Is there any way to achieve this using postgis without first having to draw an infinite amount of points to test against?


texts: Making relationship class between feature class and Excel table with ArcGIS Desktop?

texts:
I have a feature class of parcels and an Excel table of other information related to those parcels. Goal: would like edits from an Excel table to get passed onto feature class. Problem: for a true ArcGIS "relationship class" both items have to reside in same geodatabase, and I cannot have the Excel table in the geodatabase.

What is the best workflow for relating an Excel table to an ArcMap Feature class where edits get passed onto one another, in either direction. I know "relates" are temporary until you export to new feature class. I am confused on how to have an Excel table linked to feature class where the relationship is permanent but the data is dynamic.


texts:
Sorry If you've read this before - I've tried to change the content so that it reflects what I need a little better. This question was marked as a duplicate although I read the the answer to the other question and it didn't answer my query, so was advised to re-post.

I've created thiessen polygons for point data I have in ArcGIS 10.1. Unfortunately the data locations are post-code centroids, and about half of the records are lost due to the fact that they fall in the same post-code centroid as other points.

What I'd like to do is offset same-location points, within the bounds of areas defined by polygons that are in a different layer. I don't want to generate points, as the sampling points already exist.

I would prefer the points to be set apart by area or by polygon length/width - I'd basically like them to be reasonably regularly spaced in relation to the polygons. This is so I can then create thiessen polygons from all of my data points, as it won't consider overlapping points. I would prefer that non-overlapping points be left where they are. I want each of the stacked points to be evenly dispersed within the bounds of each postcode, but the dispersed points to remain in a separate layer from the post code boundaries.

Otherwise if there are any other interpolation methods which would allow me to offset data points or cope with two points in the same location it would be great to know what they are. I expect this could be problematic as each of the points can have different data, which would then become quite difficult to display / visualize as they'd overlap.

I'm having huge difficulties with this and would really, really appreciate your help.

EDIT: I've found a new tool called Disperse Markers - this sounded promising until I discovered you can't limit dispersal by polygon bounds from another layer.


texts: QGIS 2.16.2 Windows 10 installation error: ImportError: DLL load failed: %1 is not a valid Win32 application

texts:
I have just installed QGIS 2.16.2 from binaries. But it can't load QGIS utilis because of this problem:

Couldn't load QGIS utils.
Python support will be disabled.

Traceback (most recent call last):
  File "", line 1, in
  File "C:/PROGRA~1/QGIS2~1.16/apps/qgis/./python\qgis\utils.py", line 21, in
    standard_library.install_aliases()
  File "C:\PROGRA~1\QGIS2~1.16\apps\Python27\lib\site-packages\future\standard_library\__init__.py", line 457, in install_aliases
    __import__(oldmodname)
  File "C:\PROGRA~1\QGIS2~1.16\apps\Python27\lib\site-packages\future\backports\__init__.py", line 17, in
    from .misc import (ceil,
  File "C:\PROGRA~1\QGIS2~1.16\apps\Python27\lib\site-packages\future\backports\misc.py", line 24, in
    from socket import getaddrinfo, SOCK_STREAM, error, socket
  File "C:\Python27\Lib\socket.py", line 47, in
    import _socket
ImportError: DLL load failed: %1 is not a valid Win32 application.

 Python Version:
2.7.5 (default, May 15 2013, 22:44:16) [MSC v.1500 64 bit (AMD64)]

QGIS Version:
2.16.2 'Nødebo', 979e1a2

Python path:
['C:/PROGRA~1/QGIS2~1.16/apps/qgis/./python', u'C:/Users/OPdaTe/.qgis2/python', u'C:/Users/OPdaTe/.qgis2/python/plugins', 'C:/PROGRA~1/QGIS2~1.16/apps/qgis/./python/plugins', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\matplotlib-1.3.1-py2.7-win-amd64.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\nose-1.3.3-py2.7.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\tornado-4.0.1-py2.7-win-amd64.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\backports.ssl_match_hostname-3.4.0.2-py2.7.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\certifi-14.05.14-py2.7.egg', 'C:\\Python27\\Lib', 'C:\\Python27\\DLLs', 'C:\\Python27\\Lib\\lib-tk', 'C:\\Users\\OPdaTe\\Documents\\Python Script', 'C:\\PROGRA~1\\QGIS2~1.16\\bin\\python27.zip', 'C:\\Program Files\\FreeCAD 0.16\\bin', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\DLLs', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\plat-win', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\lib-tk', 'C:\\PROGRA~1\\QGIS2~1.16\\bin', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\PIL', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\jinja2-2.7.2-py2.7.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\markupsafe-0.23-py2.7-win-amd64.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\pytz-2012j-py2.7.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\win32', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\win32\\lib', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\Pythonwin', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\Shapely-1.2.18-py2.7-win-amd64.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\wx-2.8-msw-unicode', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\xlrd-0.9.2-py2.7.egg', 'C:\\PROGRA~1\\QGIS2~1.16\\apps\\Python27\\lib\\site-packages\\xlwt-0.7.5-py2.7.egg']

Any idea how to fix it? I have already unistalled it and install it twice.


texts: Why are .las file sizes inflated when read using package 'lidR' via Rstudio?

texts:
I frequently use the package 'lidR' to do forest point cloud analysis and the package functions OK, just I noticed that every time I load in point clouds using Rstudio, the file sizes seem to be much larger than they should be. The code to load in the files is simple, such as:

library(lidR)

lashome <-('C:/Directoryname')
setwd(lashome)

#input las files
UASfile <- ('UAS_1cm_aligned.las')

#read las file
las = lidR::readLAS(UASfile, select = "*+")

The question is, why is the filesize reported as 1.8GB in Rstudio when the actual filesize is 695MB? Is this an Rstudio glitch that I should just ignore?


texts: using Python to calculate NDVI, output as a Geotiff

texts:
I need a script that can calculate NDVI from two separate input Geotiff files, then output the results as a Geotiff.  I wrote following script, most of which is borrowed from various on-line sources (I've never had any training with python and only minimal experience).  I know the calculation part works, because if I print out the ndvi array, I get reasonable values.  If I calculate statistics on the array, I also get reasonable results.  The script does create an output file called ndvi.tif, but all the values in the file are 0 ( and not 0.0).

I need to convert this into a batch script that can do this calculate on a large number of input files, but right now I can't figure out why it won't output calculated values to a file.  Any advice?

 +++++++++++++++++++++++++++++++++++++++++++++

 from osgeo import gdal
 import numpy as np
 from numpy import *
 g = gdal.Open("fhmod1.tif")
 red = g.ReadAsArray()
 g = gdal.Open("fhmod2.tif")
 nir = g.ReadAsArray()
 red = array(red, dtype = float)
 nir = array(nir, dtype = float)
 check = np.logical_and ( red > 1, nir > 1 )
 ndvi = np.where ( check,  (nir - red ) / ( nir + red ), -999 )
 geo = g.GetGeoTransform()
 proj = g.GetProjection()
 shape = red.shape
 driver = gdal.GetDriverByName("GTiff")
 dst_ds = driver.Create( "ndvi.tif", shape[1], shape[0], 1, gdal.GDT_Float32)
 dst_ds.SetGeoTransform( geo )
 dst_ds.SetProjection( proj )
 dst_ds.GetRasterBand(1).WriteArray(ndvi)


texts: Why is query layer rounding down coordinates to 6 digits?

texts:
We have a SQL table view with point locations and are creating a query layer from it. However, the resulting points do not plot to the same location as the latitudes & longitudes in the table. In the screenshot below, the green points are from the Go to XY toolbar with the table coordinates pasted in, and the red points are what resulted in the query layer.

After investigating a bit, I'm pretty sure it's not a GCS mismatch (they are both in WGS84), and the best I can tell is that the query layer is rounding the coordinates to 6 digits. Here is a comparison of the coordinates in the original table (on the left), and coordinates from calculating the geometry of the query layer after exporting it to a feature class. The calculated fields are formatted to have 6+ decimal places, so it's not just the display, and these being the actual plotted coordinates are supported by the Identify tool after zooming in an extreme amount.

Everything seems to be rounded to exactly 6 digits (except for point 203 which only had 5 in the original table). Why is this happening, and how can we create the query layer with the same precision as the coordinates in the original table?

I have no idea what information is relevant to this question, but field type for the coordinates in the original table is Double. We tried manually assigning the spatial reference when we created the query layer, but didn't see any options for the X, Y fields themselves.


texts: Using shp2pgsql instead of ogr2ogr to import shapefile to PostGIS?

texts:
I may very well be doing something wrong here, but:

If I import some shapefile to a PostGIS database using shp2pgsql, I have to first figure out the SRID/EPSG of that shapefile. I think this is, at minimum, a two step process. First I query the shapefile like this:

>ogrinfo -al -so someshapefile.shp

which returns the well-known text (wkt) projection information, but is a bit verbose and somewhat opaque [to me]. Something like:

GEOGCS["NAD83",
DATUM["North_American_Datum_1983",
    SPHEROID["GRS 1980",6378137,298.257222101,
        AUTHORITY["EPSG","7019"]],
    AUTHORITY["EPSG","6269"]],
PRIMEM["Greenwich",0,
    AUTHORITY["EPSG","8901"]],
UNIT["degree",0.01745329251994328,
    AUTHORITY["EPSG","9122"]],
AUTHORITY["EPSG","4269"]]

I then typically run the wkt info through a conversion tool like Prj2EPSG to find the EPSG/SRID.

At this point, I can import the shapefile using:

>shp2pgsql -I -s 4269 someshapefile.shp <schema>.<table> | psql -U <user> -d <dbname> -h <hostaddress> -p 5432

Note, I specify the SRID with the -s flag.

If, I run shp2pgsql without specifying the SRID, no projection is set and I think the geom column has to manually be updated to include a projection.

Alternatively, I can skip the lookup, and just use ogr2ogr:

>ogr2ogr -f "PostgreSQL" "PG:host=<hostaddress> user=<user> dbname=<dbname> password=<password>" "C:/shapefile.shp" -nln <schema>.<table>

which apparently sets the projection fine, presumably extracting it automatically from the source shapefile/prj.

Question

So what's the disadvantage to using ogr2ogr? Is there in fact a flag so shp2pgsql automatically extracts and set the right projection also? If not, why not?

Addendum

There's an interesting, perhaps slightly dated, comparative analysis of using ogr2ogr vs shp2pgsql available on naturalgis.pt. It demonstrates for their particular sample data, that ogr2ogr performs significantly better on small datasets, but that shp2pgsql performs slightly better on larger datasets.

I don't feel this provides a definitive answer. The codebases may have evolved, improving each one's performance. They've only tested a very small set of sample data. The representative "large" dataset isn't actually that large. Also, it is mainly discussing issues of performance, which certainly affects usability, but the original question is more interested in usability related user input requirements.


texts: How to apply a layer definition with a spatial query to a Dynamic layer?

texts:
Using the ArcGIS Server JavaScript API, it's possible to specify an attribute query in the form of a layerDefinition against a Dynamic layer.

For example, the following query will restrict the layer to show only features where the population is greater than 5,000,000:

var layerDefinitions = [];
layerDefinitions[0] = "POPULATION > 5000000";
dynamicMapServiceLayer.setLayerDefinitions(layerDefinitions);

Is it possible to specify a bounding box in the layer definition?

That is, apply a spatial query in addition to any attribute query. The query should be based directly on the layer's geometry field (not an attribute derived from the geometry). Assume the data is stored in WGS84 and is "geometryType": "esriGeometryPoint".


texts: Getting latitude / longitude in the QGIS Raster Calculator

texts:
I have raster files in the good and old EPSG:4326 (WGS-84) and I would like to create two other raster files with the exact same extent and resolution of it, but with the first raster having the the values of the latitude of each pixel while the second has the longitude.
It seems that in the past the QGIS Field Calculator accepted "$x" to represent the longitude, but today when I try to do it in the current QGIS Raster Calculator I get the message "Expression invalid" and cannot perform it.
So the final question is: How to represent latitude/longitude in the current (v 3.14) QGIS's Raster Calculator?


texts: Clipping map service by polygon feature in ArcGIS Pro

texts:
I have a map service published in a server. It actually is a raster, but published as map service.
Is it possible to clip the service and save the output to a disk?
If there is a way to save the map service to a disk, that would work, too.


texts: Export vector layer with labels to Google Earth as kmz

texts:
I have a vector file with contours which i am interested to export it as kmz to google earth with labels similar to the ones QGIS labeling function produces. However I read that Google earth does not offer such function and renaming each contour to the corresponding elevation through NameField is not enough for this case i work on. Is there any work around to this issue?

"This is a known limit. Google Earth only displays labels for point features (placemarks). To display labels for line or polygon features in Google Earth, midpoints or centroids for those feature classes must be created as a point layer in ArcMap and imported as a separate layer to Google Earth."

Reference:
https://support.esri.com/en/technical-article/000012154


texts: Using this site, we can see that the (approximate) centre of Southampton, UK has OSGB coordinates of 441847,112279.  However, looking at the UK grid, where each large square is (500 km)2, then these coordinates don't seem to follow: where does the leading number come from? The traditional 6 figure grid ref for the location is SU418122; given a 6 figure GR denotes a (100 m)2 area, this shows that the location is 41.8 km east and 12.2 km north of the SU grid square origin.  Does the extra figure not denote an order of magnitude larger in distance?


texts:
I got an Apache and MapServer bundle. On MapServer I made a WMS service and connect it to clients (Leaflet, QGIS etc.) with an URL like this:

http://localhost/cgi-bin/mapserv.exe?map=C:/ms4w/apps/example/polt.map

In MapServer documentation I found that Apache can convert this URL in something like this:

http://localhost/mywms

In httpd.conf I see that rewrite module is on

LoadModule rewrite_module modules/mod_rewrite.so

Also I had added this block

<Directory C:\ms4w\apps\htdocs>
  AllowOverride All
</Directory>

After I created .htaccess file in C:\ms4w\apps\htdocs

RewriteEngine on
RewriteRule wmsmap?(.*) /cgi-bin/mapserv?map=C:/ms4w/apps/example/polt.map&$1

After I restarted Apache. But when I try to open with URL

http://localhost/wmsmap?layer=Green&mode=map

I got

The requested URL /wmsmap was not found on this server.

Where is my mistake?

UPDATE

I had removed Directory block what I added in httpd.conf file. And edit default block so I got:

<Directory "/ms4w/Apache/htdocs">
  Options Indexes FollowSymLinks
  AllowOverride All
  #Require all granted
</Directory>

After I put .htaccess file in /ms4w/Apache/htdocs.
And now I got another error on same URL:

msSaveImage(): Unable to access file. Failed to create output file (/tmp/WMS15201604824704.png).

What permission do I need?


texts: Export numpy array as TIFF with rasterio fail with the folder name :No such file or directory

texts:
I experienced unexpected behaviour of rasterio and I couldn't understand yet where is my mistake.
I have opened rasterio TIFF, did some manipulations on it as numpy array and then export it as TIFF using rasterio.open and dst.write.
The problem is that the saving fails many times due to the name I give to the TIFF.
This is how I tried to do this:
        src=#original tiff with metadata
        numpy_res=#the array of the results after the analysis
        tiff_name='results/test7/2021-01-10/246/534514d2a593107sat31ecnvy4u0e1694k/results.tiff'

        with rasterio.open(tiff_name,
                           'w',
                           driver='GTiff',
                           height=numpy_res.shape[0],
                           width=numpy_res.shape[1],
                           count=1,
                           dtype=numpy_res.dtype,
                           crs=src.crs,
                           nodata=None, # change if data has nodata value
                           transform=src.transform) as dst:
            dst.write(numpy_res, 1)

this is fails with next error :

RasterioIOError: Attempt to create new tiff file
'results/test7/2021-01-10/246/534514d2a593107sat31ecnvy4u0e1694k/results.tiff'
failed: No such file or directory

now , I have to say that I do have folder called results, which is empty folder tha suppose to save the resulted array.
I have also try to change result to another name of folder that is not exists yet but that didn't work.
The weird thing is that if I  give it simple name without the directory, such as:
'test.tiff'

it does work, but then it does not save it to the folder I want.
Also, if I write the directory as shwn in the docomentrt, the:
/tmp/new.tif

it does not throw any error, though it also does not create new folder 'tmp' and the result cannot be found.
My end goal is to be able to save the results into the folder I specify.


texts: Error while exporting geometry to CSV with MMQGIS

texts:
I am importing a DIVA-GIS shapefile into QGIS. It has a lot of polygons on it. My aim is to group different groups of polygons together and dissolve their internal boundaries so that I end up with some larger polygons of amalgamated smaller ones. My understanding was that this is what the QGIS Vector > Geometry Tools > Multipart to Singlepart tool was for.

I have run this several times on the relevant polygons. It seems to run OK.

Then I go to MMQGIS > Export Geometry to CSV (the idea being to use the lat / long data in Tableau to create more map detail). This command has worked for me with the same DIVA-GIS data (just not a file that had had the Multipart to Single part command run on it).

I am getting the error below every time. Can anyone tell me what I am doing wrong / how / whether I can fix it?

An error has occurred while executing Python code:
AttributeError: 'QgsVectorFileWriter' object has no attribute 'writerow'
Traceback (most recent call last):
  File "C:/Users/Matt/AppData/Roaming/QGIS/QGIS3\profiles\default/python/plugins\mmqgis\mmqgis_dialogs.py", line 1594, in run
    attribute_file_name, field_delimiter, line_terminator, self.mmqgis_status_callback)
  File "C:/Users/Matt/AppData/Roaming/QGIS/QGIS3\profiles\default/python/plugins\mmqgis\mmqgis_library.py", line 3145, in mmqgis_geometry_export_to_csv
    node_file.writerow(row)
AttributeError: 'QgsVectorFileWriter' object has no attribute 'writerow'

Python version: 3.7.0 (v3.7.0:1bf9cc5093, Jun 27 2018, 04:59:51) [MSC v.1914 64 bit (AMD64)]
QGIS version: 3.10.0-A CoruÃ±a A CoruÃ±a, 6c816b4204
Python Path:
•   C:/OSGEO4~1/apps/qgis/./python
•   C:/Users/Matt/AppData/Roaming/QGIS/QGIS3\profiles\default/python
•   C:/Users/Matt/AppData/Roaming/QGIS/QGIS3\profiles\default/python/plugins
•   C:/OSGEO4~1/apps/qgis/./python/plugins
•   C:\OSGeo4W64\bin\python37.zip
•   C:\OSGEO4~1\apps\Python37\DLLs
•   C:\OSGEO4~1\apps\Python37\lib
•   C:\OSGeo4W64\bin
•   C:\OSGEO4~1\apps\Python37
•   C:\OSGEO4~1\apps\Python37\lib\site-packages
•   C:\OSGEO4~1\apps\Python37\lib\site-packages\win32
•   C:\OSGEO4~1\apps\Python37\lib\site-packages\win32\lib
•   C:\OSGEO4~1\apps\Python37\lib\site-packages\Pythonwin
•   C:/Users/Matt/AppData/Roaming/QGIS/QGIS3\profiles\default/python
C:\Users\Matt\AppData\Roaming\QGIS\QGIS3\profiles\default\python\plugins\mmqgis/forms


texts: Error when clipping raster with vector layer

texts:
I'm trying to clip a raster (GeoTIFF) with a vector layer in QGIS (version 2.18).

I get the error "ERROR 1: Attempt to create 0x0 dataset is illegal,sizes must be larger than zero." when I changed the resolution (x and y resolution = 250 meters, the same resolution of my input raster).

However, when I choose for "keep the resolution of the input raster", the program script indicates a resolution "-tr 0.002083333 0.002083333" and gives a clipped raster with only missing values ("nan").

My raster and vector polygon have the same coordinate reference system (WGS 84).


texts: Iterating and calculating with fields using PyQGIS

texts:
Inside my layer there are fields like this:
name    |    data       |    length
name_a  |    field_a_x  |    ...
name_b  |    field_b_x  |    ...
name_c  |    field_c_x  |    ...
name_d  |    field_d_x  |    ...

Is it possible to iterate only over the data-column, get the fields into a list and replace the data after calculating (or doing something else...)? The calculating-part will be easy but I'm having trouble understanding how to get and afterwards replace the data?
This is what I got so far...
from qgis.core import QgsProject

layer_name = "ExampleLayer"
layer = QgsProject.instance().mapLayersByName(layer_name)[0]
ftrs = layer.getFeatures()

list_old = []
for feat in layer.getFeatures():
    list_old.append(feat['data'])


texts: How to use QGIS Accuracy Assessment plugin?

texts:
I am new user and i want perform accuracy assessment in qgis. I have 1 classified map and another for validation. I have installed the plugin. what am i doing after?

Plugin: http://plugins.qgis.org/plugins/accassess/

Docu: https://bitbucket.org/jkibele/accassess/overview


texts: lat/lon of intersection between line of sight from a given lat/lon and altitude above ground

texts:
Given a starting latitude, longitude and altitude, and a line of sight defined by azimuth and elevation, I want to find the latitude and longitude (assuming spherical Earth, WGS84 or other) where the line of sight pierces a given altitude.

Background:

I have a set of data points from an all-sky camera (i.e. a wide-angle lens pointing at the sky). The data show auroral emissions. I have elevation and azimuth for each pixel (each array element), and I know the coordinates (latitude, longitude, altitude) of the camera. Based on the assumption that the auroral emissions are coming from a given altitude above ground (e.g. 150 km), I want to show the data on a map.

I aim to implement the calculation in Python (using numpy), and would like an efficient way of doing this (i.e. preferably not element-by-element), since the data arrays contain upwards of 200k pixels. The distance between pixels are on the order of a few km, which should give an idea of the precision needed.


texts: Combining two areas on a shapefile (two polygons) into one

texts: I am trying to combine two polygons or areas on a map into one.   I have added a field called region to try and dissolve the borders by region.  This did not work.  I am now trying to combine them using merge selected features in QGIS but I get the error:


texts: Assigning new values to raster attribute table by comparing two fields?

texts:
I have created a new raster from two raster datasets using combine operation in ArcGIS. The new attribute table has two fields from the input rasters. What I need is, to compare these two fields and assign meaningful values to a new field.

For example:

see the attached image (I don't know how to past a table in here)

... this list goes to few hundred.

“New Class” is the field that I need to assign the new values.

I have been using “code block” in field calculator using python. But it is too long that I have to write individual statements for each unique combination. I am looking for an efficient way other than my method.

Something like reading from a reclassification table?


texts: Creating expression using arcpy.SearchCursor?

texts:
Example:
On the basis of the table which is not defined by the number of columns and rows, you need to evaluate expressions as shown, using the arcpy.SearchCursor ....

First, create a list of fields, using the arcpy.ListFields
Then, using SearchCursor needs to calculate the defined an expression, for each column separately ...

I should be to existing code insert defined an expression....

 values = [] # store the sum values here.
 fields = arcpy.ListFields(fc, "*")

 # get the OID/FID field name to skip
 desc = arcpy.Describe(fc)
 if desc.hasOID:
     OIDname = desc.OIDFieldName.upper()
 else:
     OIDname = ""

 for field in fields:
     if field.name.upper() != OIDname: # skip the OID/FID field.
         if field.type in ("Double", "Integer", "Single"):
             # sum each suitable field, but not the NULL ones - they would be bad
             with arcpy.da.SearchCursor(fc,field.name,field.name + " is not NULL") as sCur:
                 thisValue = 0
                 for row in sCur:
                     thisValue += row[0]
                 values.append(thisValue) # this will be the inserted row
             fieldNameList.append(field.name)

 with arcpy.da.InsertCursor(fc, fieldNameList) as cur:
     cur.insertRow(values)

@Michael Miles-Stimson provided the preceding code.

Instead of the sum of the "thisValue" to create a new expression?


texts: Is there OpenLayers layer type that supports just labels?

texts:
Is there an OpenLayers layer type that supports just labels on a map without an underlying feature?
I think this is called an annotations layer in GIS.
If not, I'm not sure how to approach this and any advice would be appreciated. I thought styling a vector point layer so that the point itself would be invisible might be a good approach but then the user won't be able to select/edit/delete the feature, which is part of the requirement.
The application is a bit like the 'Create Map' functionality in Google Maps. The user can add markers, lines and polygons to the map using OpenLayers vector editing controls. A popup allows them to edit a label, change the styling and delete the feature. The map settings (zoom, centre, background layer) and features are saved to a SQL Server database and stored as geometry. This is sent to the database as WKT but the layers are loaded as GeoJSON from a handler on the website.
The client has asked for the ability for users to just place labels on the map (without a marker) and save these too. As I mentioned above, I thought about implementing this as a labelled point layer with the points styled to be invisible but this stops the user from selecting/editing existing features.
My question was whether there was an annotation layer type where the user can create annotations, select them and move them around on the map - similar to vector layers.
It looks like there isn't so how should I go about approaching this?


texts: Pyqgis: "a.geometry().intersects(b.geometry())" wouldn't find any intersections but should. Why?

texts:
I am teaching myself pyqgis and try to find those polygons in one layer that intersect with polygons in another layer. It seems very trivial, and on GIS.SE are several similar questions that are suggesting me (with my limited python/pyqgis-knowledge) that this should work:

#ausschluss and gebaeude are QgsVectorLayers with Polygons
#(shapefiles with the same crs and clearly intersecting)
gebaeude.setSelectedFeatures([]) #deselects everything
selections=[] #declares it is a list
for f in gebaeude.getFeatures():
    for a in ausschluss.getFeatures():
        if a.geometry().intersects(f.geometry()):
            #same result with .within() and even with a and f switched (the docs aren't that clear on that)
            selections.append( f.id() )
            break #only one or less intersection are possible
gebaeude.setSelectedFeatures(selections)
gebaeude.invertSelection()

#gives 0, 78137 an 78137 (no intersection found)
print len(selections)
print gebaeude.selectedFeatureCount()
print gebaeude.featureCount()

What is it, I am doing wrong?

Edit: Running QGis 2.8.3 on OS X (10.8.5) with the KyngChaos-Longterm-Built


texts: Creating line between two points with common field?

texts:
Does anyone know how I could create a line between two points that have a common code?

For example: I have the start and end point of a pipe as a lat/long for a network of pipes, each with a different ID. I want to create a line feature for each of those with the pipe ID which I can then calculate the length from.

Anyone know of any arcscript tools in ESRI software or QGIS tools I could use to accomplish this without coding?


texts: The source of and possible solutions for DEM anomalies?

texts:
I am working with a DSM (Digital Surface Model) created using stereo-paired spot 6/7 satellite imagery. This DSM, and metrics derived from it will be used to create a Predictive Ecosystem Map (PEM).

My problem is there appears to be two kinds of anomalies in the model. The first type appear as uniform "Bumps" throughout the DSM, the second appear as large "triangular" shaped polygons. These triangles mostly occur in relatively flat areas of the DSM, where-as the bumps appear almost everywhere else (See Photo B below).

Photo A shows the pan-sharpened true-color SPOT image of the study location, Photo B Shows the DSM viewed with Hillshade effect. Hillshade does exaggerate theses anomalies, however, the anomalies appear in other metrics created from the DSM.

My questions are:

1) Why are these anomalies not more uniform throughout the DSM if they are induced through sensor error, or what are the possible sources of these anomalies, and

2) What might be done to reduce the influence of these features when creating DSM metrics, and in the final analysis.

I am using ArcMap 10.3.1


texts:
I am trying to run the tool on this simple table using Point X, Point Y, Xt0 and Yt0 fields for calculating distances. But immediately as the tool runs, it stops with an error "an expected field was not found or could not be retrieved properly". This error is shown for every single feature in the results window. I have tried exporting to a new table, removing features with null XY coordinates, exporting to a geodatabase, restarting ArcMap...everything possible. Even google doesn't give me a solution.

What is wrong with my table?


texts: Using multiple address locators in ArcGIS

texts:
I am geocoding millions of addresses for my research (every registered voter in six states).

I am using two geocoders: one that is really fast and one that has fewer unmatched addresses when I geocode.

Is there a way to use ArcMap or ArcCatalog to use the second address locator on only the unmatched addresses from the first?

Right now I am importing the attribute table using the R programming language, exporting only the unmatched addresses as a CSV, re-importing them to ArcGIS, and geocoding, but I figure there has to be a faster way.


texts: Linking datasets for two polygon layers in ArcMap?

texts:
I have two polygon shp files that do not have a common ID.

I'm trying to link the datasets with a common ID based on spatial location.

Is there a tool/process to link features (common ID in the attribute table) from one polygon layer that are spatially located inside of features of the other polygon layer?

For example, imagine you have a parcel/property boundary layer, and a buildings layer.

Would there be a way to link all building features with the surrounding parcel/property boundary feature with a common ID.


texts:
I am experiencing a strange problem by using the GRASS toolbox on QGIS 2.14.3 Essen.

Even if in settings>options>CRS I select to use the project CRS when creating new layers, the results of GRASS toolbox tools, like r.watershed and r.fill.dir,  have an "user-defined CRS" with code 100000, even if the project CRS is EPSG 20137.

How can I solve this?


texts: Add and remove map layers using radio buttons

texts:
I've struggled to use radio buttons to toggle multiple map service layers. So far it's fine to add but I can't then add another map service and remove the previous at the same time - I just can't figure out the problem:

$("input[type='radio']").click(function () {

    var radioValue = $("input[name='exampleSwitch']:checked").val();
    if (radioValue == 'on') {

        map.addLayer(mapidurl);

    } else if (radioValue == 'off') {
        debugger
        map.removeLayer(mapidurl);
    };
}

The HTML:

      <input class="switch-input" id="switch3" type="radio" name="exampleSwitch">

Can anyone offer any guidance?


texts: Problems Connecting to WFS using ArcGIS Pro

texts:
I cannot connect to this Web Feature Service using ArcGIS Pro even though it is public:
https://maps.isric.org/mapserv?map=/map/wosis_latest.map
I also tried https://maps.isric.org/mapserv?map=/map/wosis_latest.map&service=WFS&request=GetCapabilities
error message:


texts: ArcGIS Online data fields "Create List" and "Delete" options not present

texts:
I have two fields present in multiple layers, and the "Create List" and "Delete" options are not present. Navigate to the item details page --> Data --> Fields
The fields in question are named Feature and Mapname.  The layers have been designated as Authoritative.
But when I made a copy of one of the layers that have this issue in ArcGIS pro and Shared as a Web Layer and the options are visible.   I do not want to reload all the layers that are affected by this as I do not know how or why the loss of these options occurred.


texts: Saved layout not showing in new Map of ArcGIS Pro

texts:
I am new to ArcGIS Pro and I am trying to use the same layout I used in Map 1 as in the figure below.

Now, I created another map (Map2), and I want to use the layout I saved ( I followed this link instructions to save layout3).

The next steps I took are the following

Insert Layout > Import Layout (The one I saved)

Upon right clicking on the map frame in the table of contents, I selected my new map (Map 2).

Then, this happens. The map does not appear.

Is it because it is not linked with any of the other objects?

Also, the spatial reference system is the same as the first map..

For reference the new map I made is this:

Is there a way to do this such that I do not repeat everything and just use the layout I saved so that both figures are consistent?


texts: How can I optimize layer filtering using setSubsetString or QgsFeatureRequest().setFilterFids(my_indexes)

texts:
i'm developping a plugin (under QGis 2.8.1) for traffic assignment where I want to show the results of my simulation at each time step. Right now I'm using Time Manager plugin but it gets very slow when my layer has hundreds of thousands of attributes. In my case I know exactly what feature IDs I want to show at each time step so I thought it would be easy to make it faster. Then I tried 2 solutions:

Here is what I tried: at each time step of my loop I set the ordered list of indexes of attributes to show (they are always ordered in my case)
indexes = list(j for j in range(index_start, index_end))
then I tried 2 solutions:

Use of

my_layer.setSubsetString("fid>=indexes[0] AND fid<=indexes[-1]")

Use of

request = QgsFeatureRequest().setFilterFids(indexes)
my_layer.getFeatures(request)

I use self.iface.mapCanvas().refresh() to refresh the view
Solution 1 works as it refresh the view with the right filtered features but it is even slower than using a SQL expression not based on the indexes.
Solution 2 is fast but the display of the layer doesn't change.
Any idea why?


texts: Using ArcGIS Online Base Maps without Esri Software?

texts:
Since the end of January 2011:

"ArcGIS Online map services are
available at no cost regardless of
use." - http://www.esri.com/software/arcgis/arcgisonline/standard-maps.html
However, the use of ArcGIS Web Mapping
APIs may require a subscription. - http://www.esri.com/software/arcgis/arcgisonline/common-questions.html

All the basemaps available from ESRI are accessible by REST services which do not require an API. There is already a new layer type in OpenLayers that can access these services.
The Terms of Service state:

In particular, you may not (or attempt
to) (a) access Services or resources
not made available to you through a
standard web browser or Esri products;

Update - these TOS may in fact refer to the esri.com site and "services" rather than the map services themselves.
So is it legal and within the terms of service to use any of these layers in any web application such as OpenLayers or through the Google or Bing mapping APIs as these are all used in a standard web browser?