texts: How can I query a geometry field in PostGIS for line segments within a given lat/long bounding box?

texts:
In google maps, I want to use map.getBounds() and take those values and query my PostGIS instance for all line segments that intersect that range. I am doing something wrong. I tried:

SELECT * FROM pipeline_denorm
WHERE pipeline_denorm.wkb_geometry &&
ST_Transform(
    ST_MakeEnvelope(
        -121.76986352563478, 37.85809027719253,
        -121.19308129907228, 38.07462286323802, 4326),
    2223);

based on a little research. Unfortunately no results are being returned and I am certain there is data that should be being returned. These are LineStrings. Basically want any record that contains a point within the bounding box.


texts: Map with OpenLayers 3 from GeoServer not showing but no error

texts:
References to localhost changed to mydomain.com otherwise Stack Exchange doesn't allow me to post.

I am using OpenLayers 3 to show a map. I can successfully show vectors on the map by using a GeoJSON or KML source from a URL. However, for larger maps with thousands of areas I'm trying to use GeoServer to serve an ImageWMS layer.

If I navigate to http://mydomain.com:8080/geoserver/wms?service=WMS&request=GetMap&layers=areas_for_type&viewparams=typecode:AdministrativeWard&width=660&height=330&format=application/openlayers&bbox=-180.0,-90.0,180.0,90.0&srs=EPSG:4326 then I can see the map in my browser. However, when I call it from code I can only see the underlying OSM tile. My snippet of JavaScript which renders the layer is:

var layers = [];
var tile = new ol.layer.Tile({ source: new ol.source.MapQuest({layer:'osm'}) });
layers[0] = tile;
var geoserverRoot = '@ConfigurationManager.AppSettings["GeoServerUrl"]';
var params = {
    LAYERS: 'areas_for_type',
    FORMAT: 'image/png',
    viewparams: 'typecode: ' + areaType
};
var wmsImage = new ol.layer.Image({
    opacity: 0.3,
    source: new ol.source.ImageWMS({
        url: geoserverRoot + 'wms',
        params: params
    })
});
layers[1] = wmsImage;
wmsImage.addEventListener("change", function() {
    map.getView().fitExtent(wmsImage.getExtent(), map.getSize());
});
var map = new ol.Map({
    target: 'map-container',
    layers: layers,
    view: new ol.View({
        center: [0, 0],
        zoom: 2
    })
});

there are no JavaScript errors and the GeoServer log shows:

2014-11-10 11:27:58,546 INFO [geoserver.wms] -
Request: getMap
    Filters = null
    Tiled = false
    Palette = null
    FeatureVersion = null
    SRS = urn:x-ogc:def:crs:EPSG:3857
    Styles = [StyleImpl[ name=polygon]]
    Layers = [org.geoserver.wms.MapLayerInfo@72adad35]
    MaxFeatures = null
    Crs = PROJCS["WGS 84 / Pseudo-Mercator",
  GEOGCS["WGS 84",
    DATUM["World Geodetic System 1984",
      SPHEROID["WGS 84", 6378137.0, 298.257223563, AUTHORITY["EPSG","7030"]],
      AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]],
    UNIT["degree", 0.017453292519943295],
    AXIS["Geodetic latitude", NORTH],
    AXIS["Geodetic longitude", EAST],
    AUTHORITY["EPSG","4326"]],
  PROJECTION["Popular Visualisation Pseudo Mercator", AUTHORITY["EPSG","1024"]],
  PARAMETER["semi_minor", 6378137.0],
  PARAMETER["latitude_of_origin", 0.0],
  PARAMETER["central_meridian", 0.0],
  PARAMETER["scale_factor", 1.0],
  PARAMETER["false_easting", 0.0],
  PARAMETER["false_northing", 0.0],
  UNIT["m", 1.0],
  AXIS["Easting", EAST],
  AXIS["Northing", NORTH],
  AUTHORITY["EPSG","3857"]]
    Bbox = SRSEnvelope[-3.2580518936273526E7 : 3.2580518936273526E7, -1.1740727544603072E7 : 1.1740727544603072E7]
    Env = {}
    Angle = 0.0
    CQLFilter = null
    Elevation = []
    FeatureId = null
    StartIndex = null
    ViewParams = [{TYPECODE= AdministrativeWard}]
    RemoteOwsType = null
    RemoteOwsURL = null
    FormatOptions = {}
    BgColor = java.awt.Color[r=255,g=255,b=255]
    Transparent = true
    SldBody = null
    Sld = null
    StyleUrl = null
    StyleFormat = sld
    ValidateSchema = false
    TilesOrigin = null
    Exceptions = SE_XML
    StyleBody = null
    SldVersion = null
    StyleVersion = null
    Height = 600
    Width = 1665
    Filter = null
    Time = []
    Format = image/png
    Buffer = 0
    Request = GetMap
    RawKvp = {CRS=EPSG:3857, BBOX=-32580518.936273526,-11740727.544603072,32580518.936273526,11740727.544603072, TRANSPARENT=true, VERSION=1.3.0, FORMAT=image/png, SERVICE=WMS, VIEWPARAMS=typecode: AdministrativeWard, HEIGHT=600, REQUEST=GetMap, LAYERS=areas_for_type, STYLES=, WIDTH=1665}
    Get = true
    BaseUrl = http://mydomain.com:8080/geoserver/
    RequestCharset = UTF-8
    Version = 1.3.0

Why is this not showing on the map if there are no errors?


texts:
I have been investigating using TopoJSON files to map floor plans for various buildings. As GeoJSON these files render perfectly, however after converting to TopoJSON they are missing a lot of the lines used to render them and the number of geometries is much less than the number of features in the GeoJSON.

I experimented by increasing the quantize parameter but after a certain point it doesn't add any more detail. Is it possible to convert from GeoJSON to TopoJSON without having my features quantized?


texts: How to specify a polygon's coordinate reference system in a CQL_FILTER geometric filter?

texts:
The CQL tutorial page contains an example showing how a polygon can be used in a DISJOINT query:

DISJOINT(the_geom, POLYGON((-90 40, -90 45, -60 45, -60 40, -90 40))

What is the syntax for specifying the coordinate reference system of the POLYGON object's coordinates?

The assumption is that the polygon's coordinates are in the same CRS as the layer being queried, but I'm hoping that if the polygon's CRS can be specified, GeoServer may be able to perform the calculation on-the-fly and return a meaningful result regardless of the CRS of the layer. (Please let me know if this isn't going to work, and I need to re-project the polygon before performing the query.)

The answer at https://gis.stackexchange.com/a/29630/3112 says

  call setSRS(POLYGON(...), 'EPSG:xywz')

but I can't get setSRS to work, and I can't see that listed anywhere else.


texts: Displaying Legend for Variable Width Line in QGIS

texts:
I am using QGIS v3.4 to display a simple line as a line of varying width using the variable line width attribute facility of QGIS (please see screenshot below).

The layer displays correctly within the map window however I am experiencing challenges in displaying the variable lines as a legend within the print layout window.


texts: QGIS is locking the clipboard from other application? (SitesPlus)

texts:
We are starting to use QGIS along with Sales Forecasting software Sites Plus by Gravitec. An ugly problem just popped up: Whenever QGIS is open, any attempts to calculate a data table in SitesPlus results in an error stating "Unable to open Clipboard" and the changes are not committed.

This happens even when just opening a blank QGIS project, without adding any data. If I close out QGIS, the problem goes away.

Is there anything that could be causing QGIS to lock up the clipboard? Anyone else having a similar issue?

I've tried clearing the clipboard before attempting edits to no avail. QGIS seems to lock it up once it's open.

Thanks!

Running QGIS 2.2 Valmiera & SitesPlus 2010 on a Win 7 machine
Have duplicated the error on multiple machines.


texts: How to append variables to shapefile in R?

texts:
I want to create a map of China that demonstrates the value of the ethnic heterogeneity index.
I want to generate the map as in this webpage, https://mgimond.github.io/Spatial/mapping-data-in-r.html
I gathered the .shp file of provinces from Hijmans, Robert J., University of California, Berkeley. Museum of Vertebrate Zoology. Second-level Administrative Divisions, China, 2015 [map]. No Scale Provided. Retrieved April 15, 2021, from https://geodata.lib.utexas.edu/catalog/stanford-nt024fn0432
But I do not know how to append "ethnic heterogeneity index" to this file?
The ethnic heterogeneity index is stored as a CSV. file and the common variable for both files is "province".
My current csv file does not contain the longitude and latitude information for provinces, is that the problem?
chinaprovincedf <- as(chinaprovince, "data.frame")
Error in as(chinaprovince, "data.frame") :
  internal problem in as(): “sf” is(object, "data.frame") is TRUE, but the metadata asserts that the 'is' relation is FALSE

I am not sure how to proceed.


texts: Interoperability Macromedia Freehand and Esri ArcGIS

texts:
I will receive a Adobe Macromedia FreeHand file, which contains multiple vector data that I have to export, so that they can be used in ArcGIS.

Does anyone know in which format FreeHand stores vectors by default and to which formats - that are interoperable with ArcGIS - they can be exported?

A list of ArcGIS interoperable foramts can be found here: http://www.esri.com/library/fliers/pdfs/data-interop-formats.pdf


texts: ArcGIS API for JavaScript Bookmark in separate Div

texts:
Using the ArcGIS API for JavaScript I need to create a number of different bookmarks and access them through the HTML in separate div tags.  I created this JS Fiddle in attempt.

I thought that if I defined separate bookmarks and bookmarkitems, it would work.  No success yet.  I think I am close, missing something simple/.

      var bookmark1 = new esri.dijit.Bookmarks({
      map: map,
      bookmarks: Bookmark_1

    }, dojo.byId('Bookmark1'));


texts: Converting AutoCAD drawings to QGIS shapefiles

texts: I am having problems converting an AutoCAD 2012 drawing (.dwg) file to a shapefile (.shp) in a QGIS project I am working on using the Coordinate Reference System NAD_1983_StatePlane_Indiana_East_FIPS_1301_Feet. My main problem is that the shapefiles are not located on the correct coordinates nor is object scale correct. What am I doing wrong here? I know there are questions regarding AutoCAD/QGIS conversions, but I am teaching myself and I need a little more help on the basics.


texts: Does the Chamberlin trimetric projection map (some) great circles to straight lines?

texts:
I know that, in general, the Chamberlin trimetric projection maps most great circles on the sphere to some sort of non-straight curve on the flat map.

But what about the 3 special great circles that go through each pair of the 3 control points of the Chamberlin trimetric projection?
Does the Chamberlin trimetric projection map those 3 great circles to straight lines on the flat map?


texts:
I'm using osmdroid and osmdroid bonus back and i can't cache my tiles for offline use.
I followed this url: https://github.com/osmdroid/osmdroid/wiki/Cache-Manager

This is the code that i use to cache tiles,

CacheManager cacheManager = new CacheManager(map);
BoundingBox bbox = new BoundingBox(33.298190, -8.572239, 33.148845, -8.393711);
cacheManager.downloadAreaAsync(getActivity(), bbox, 10, 15);

I don't know if i did it in the wrong way but when i run my project, my tiles are not downloaded, and the progressbar shows on for 2 secends and disappears with any error message.

I seek any help?


texts: Using filter functions in cql_filter throws "... cannot be cast to org.opengis.filter.Filter"

texts:
I'm trying to use CQL filters in WFS, like (with proper URL encoding):
&cql_filter=strEqualsIgnoreCase(name, 'Jane Doe')
&cql_filter=strMatches(name, '(?i).* doe')

But this tells me, e.g.:

java.lang.ClassCastException: org.geotools.filter.function.FilterFunction_strMatches cannot be cast to org.opengis.filter.Filter

The CQL and ECQL documentation claims:

CQL/ECQL can use any of the filter functions available in GeoServer.

The functions I'm trying to use are listed there, and I can indeed use some of the other functions, like the following works just fine:
&cql_filter=strToLowerCase(name) like '% doe'

However, many other functions throw an error like above.


texts: How can I get line features to appear as great circles on a map, instead of straight lines?


texts: GeoServer: Use oracle as raster store without spatial extension

texts: Is it possible to store raster data (e.g. GeoTIFFs) in an oracle table and access it from geoserver without using the oracle spatial extension/GeoRaster? I would store the images as blobs, but I don't find a way to access the rasta data from GeoServer and use it as a DataStore.


texts:
I am trying to extract all highway types for the world from OSM through Osmosis. My Osmosis command looks like below, however, after the extract, I find certain roads (especially, motorways and trunk types) are missing when I compare with the results against overpass-turbo. The output from overpass-turbo is complete for a region while what I get from the extract isn't.

osmosis --read-xml file=file/path --tf accept-ways highway=* --used-node --write-pgsql host="<host>" database="osm" user="<user>" password='<pw>'

Update:

This is how the overpass results (orange) and extracts I did (blue) for Boston seaport area look like. Results from overpass has everything in my extract and more as you can find in the image.

UPDATE

Actually, the problem was discovered to be a filter that is set to extract only street segments that has a name tag. Therefore, the issue is not related to osmosis or the source OSM data.


texts: How to create a specific svg file with data from a feature

texts:
I am working in c++/ QT4.8.5, with qgis api c++ 2.8.2

I have to display a dynamic icon for a feature.
some of the data change periodically and in order to display the change i need a svg that can get the data
I don't want to have a collection a svg file because there too possibilities
At this point i need only two parameters but the number of parameter will be more
My feature has two attributes speed and color each change
the speed determine the size of the arrow in different color
The color is set in the feature, it is the result of an algorythm.
 there is, for the moment , eigth colors and 3 type of speed ( slow, medium , fast )

do you have any idea on how can i proceed to create file an apply it on my symbol ?


texts:
What datums were used in North America before NAD27? I googled around on this without luck.

Also I found an old nautical chart of Lake Michigan, dated 1876, which states:

  Longitudes on this map depend on that of Fort Howard astronomical station adopted in 1874 as 88 00' 57.82". Latitudes are derived geodetically from Fort Howard.

Don't know what that means.


texts: Obtaining information from nearby points in shapefile using ArcGIS Desktop?

texts:
I am working with ArcGIS Desktop 10.2.

I have a point shapefile, something like cities, which contain city level information, population, area, average temperature in the city etc.

Is there a way to identify for each city the total population (or any other variable of interest) in the k nearest cities?

An example would be for US, I have a shapefile containing all the big cities, I would like to know what is the total population living in 5 (this parameter can change) nearest cities around New York. I would like to do it for all the cities in the shapefile.

I am in the process of coding it by first making a function for identifying the first k neighbors and then obtaining information from these neighbors. The idea to post this question is to find out that if there is a more straight forward way to do it via arc gis.


texts: Field calculator executing error 999999 field is not nullable?

texts:
I've succeeded running such code in the pycharm.
but I failed when I just copied the code to the code block of the field calculator.
I am quite puzzled.
What is wrong?

the code is:

def calVB(appeal,VB,ttmin):
    if(VB!=0):
        temp=appeal / (VB*ttmin)
    else:
        temp=0.1
    return temp

the expression is:

calVB( !facilities.appeal!, !ODMAX.VBFld!, !ODMAX.Total_Minu!)

arcpy code:

import arcpy
from arcpy import env
try:
    # set the environments
    env.workspace = "E:/Emily/code/result0610.gdb"
    env.qualifiedFieldNames = "UNQUALIFIED"

# Define script parameters
    inFeatures = "E:/Emily/code/point/test/ODMAX.shp"
    layerName = "veg_layer"
    newField = "influence"
    joinTable = "E:/Emily/code/point/test/facilities.dbf"
    joinField1="OriginID"
    joinField2="OBJECTID"
    #calcExpression = "(!population._F_1! / !ODMAX.Total_Minu!)+1"
    calcExpression = "calVB( !facilities.appeal!, !ODMAX.VBFld!, !ODMAX.Total_Minu!)"
    codeblock="""def calVB(appeal,VB,ttmin):
        if(VB!=0):
            temp=appeal / (VB*ttmin)
        else:
            temp=0.1
        return temp"""
    outFeature = "E:/Emily/code/point/test/"+layerName+".shp"

# Add the new field
    arcpy.AddField_management (inFeatures, newField, "double",field_is_nullable="NULLABLE")

# Create a feature layer from the vegtype featureclass
    arcpy.MakeFeatureLayer_management (inFeatures,  layerName)

# Join the feature layer to a table
    arcpy.AddJoin_management (layerName, joinField1, joinTable, joinField2)
    for fld in arcpy.ListFields(layerName):
        print fld.name
# Populate the newly created field with values from the joined table
    arcpy.CalculateField_management (layerName, newField, "0", "PYTHON")
    print "0"
    arcpy.CalculateField_management (layerName, newField, calcExpression, "PYTHON",codeblock)
    print "cal done!"
# Remove the join
    arcpy.RemoveJoin_management (layerName)
    print
    print
    for fld in arcpy.ListFields(layerName):
        print fld.name

# Copy the layer to a new permanent feature class
    arcpy.CopyFeatures_management (layerName, outFeature)

except Exception, e:

    import traceback, sys
    tb = sys.exc_info()[2]
    print "Line %i" % tb.tb_lineno
    print e.message

ERROR:

ERROR 999999: Error executing function.
The field is not nullable. [influence]
The field is not nullable. [influence]
Failed to execute (CalculateField).


texts: Using variable as fixed input value for Buffer in ArcGIS ModelBuilder?

texts:
I am using ArcGIS Desktop 10.0.

I want to show/create a buffer of a single point but the radius of that buffer depends on a simple formula (e.g. Radius = abc).

Instead of entering directly the radius, is there are way to just enter the constants on the ModelBuilder, use the result and connect it to the Buffer tool so that we can get that result?


texts: GeoServer GetFeatureInfo returns plain text result

texts:
I want to send GetFeatureInfo request to my GeoServer and get response as JSON.

My query string parameters are like this:

request:GetFeatureInfo
service:WMS
info_format:application/json
exceptions:text/javascript
version:1.1.0
...
...

But my response comes like this:

Results for FeatureType 'http://localhost:8080/geoserver/App:line':
--------------------------------------------
gid = 590545
geom = ...
name = ...
--------------------------------------------

I want to get response as JSON result.

I set exceptions:application/json but it did not work.

This was working before, but I moved GeoServer to another server. I think some settings are not set in a config file about JSON response.

My GetCapabilities response is like this:

  <GetFeatureInfo>
    <Format>text/plain</Format>
    <Format>application/vnd.ogc.gml</Format>
    <Format>application/vnd.ogc.gml/3.1.1</Format>
    <Format>text/html</Format>
    <Format>text/javascript</Format>
    <Format>application/json</Format>


texts: Installing python xlw and xlrd libraries for QGIS on Mac OSX with Python 2.6

texts:
I am trying to use the XYtools plugin for QGIS (because copy-paste is not working with just one of my files?!), and I'm concerned that I actually just need to update my Python but not sure if doing so would disrupt anything. I've read all the posts on

  install python xlw and xlrd libraries

errors, but most instructions seem to be for PC or assume Python version >= 2.7.

I installed QGIS about 4 months ago on a brand new machine running Yosimite, using the QGIS 2.12 built for Mac OS X Mountain Lion, Mavericks, Yosemite and El Capitan. Built with Qt Cocoa 4.8. The only Python I can find on my machine is version 2.6.

*update: In the meantime, I installed the Quick Export plugin and this gives me (almost) what I need. The strange thing is, even this appears to leave out several (seemingly) random entries in the table. I'm not entirely positive there isn't a problem with the project, as I have no problem copy-pasting from other projects. In any case, the question still stands (potentially for other plugins).


texts: Converting Points to Lines in ArcGIS Desktop and retaining attributes of both start and end points?

texts:
I have a table with animal location points. I need to calculate movement distances between points from the same animal. So I thought I would create lines from the points, which are linked by an animal ID code, then calculate the line distance of each movement.

I have a table structure similar to this:

FID; AnimalID; Date; UID; Weight

The AnimalID field links observations from the same animal, the date is the sequential ordering field, UID is a copy of the FID so that each record has a unique ID code and Weight is a measurement of the animal at the time.

Some animals were only located once (these get excluded by the creation of lines so is usually not an issue and I don't want to delete all of the single observations from the dataset every time I run this). The remaining animals are located multiple times.Of note is that some location points are identical for different individuals.

I am using ArcGIS Desktop 10.0 with spatial analyst extension.

At the end I need to have the following fields in the attribute table (or something like it):
AnimalID; Date 1; Date2; Weight1; Weight2; Distance

So I have looked at ETGeowizards Points to Line tool and this is a really good toolset. I specify the ID field as the AnimalID, the order field as Date and it allows me to attach the attributes of the UID field of the start and end point. It is almost what I need. HOWEVER, it creates a single polyline for each animal and I need to have different segments. While I can convert it to multipart, I lose the important attribute information (date, weight) that I can't link back to the original data.

I also tried the point to line tool in Geospatial Modelling Environment (GME) which is also almost suitable. This allows me to create multipart lines with the AnimalID transferred over. However, I can't include a UID field from each point, so I can't join the specific attributes of both points over. I tried to do a spatial join (one to one). However, as I have more than 2 points that overlap with some lines, this causes incorrect attributes to be transferred. A one-to-many spatial join causes replicates of the line (as there are two points to each line), plus the erroneous ones from overlapping points from other animals. It seems a nightmare to try and extract the correct ones from this.

As far as I can see, the Point to Line tool in GME is most appropriate if I could get it to copy a From_UID and a TO_UID fields over. Then I could do a join of the line layer to the original point layer -> export to shapefile -> do another join of the exported line layer with the original point layer (to get the attributes for the second point) -> export shapefile. Then I add a field -> calculate geometry for the distance calculation.

FYI this is what the GME Point to Line tool does:

  The tool allows a user to select a field name from the points layer to
  group the line segments, a field name to sort the groups (optional),
  and the point layer fields to be copied to the line layer. Line
  segments are created from two points with the same attribute in the
  group field. The attributes for the line segment are taken from the
  first point. The points DO NOT have to be sorted by group attribute
  and order.

I have no experience writing scripts but hopefully someone can think of a workaround.

I am used to an older version of ArcGIS, so I used ET Geowizard tools first. But yes, you can do this with the ArcGIS toolset Point to Lines tool. But still can't retain attributes or a to_ and from_ field.

Any more suggestions?


texts: Finding the intersection of an array of markers across multiple layer groups Leaflet

texts:
I have two challenges related to this. The first(1) is that I have created a bunch of variables that store the data needed for my markers and I cannot reliably add markers to layer groups when storing it this way. I can't figure out why. It is an array called 'parks'.

I wrote them out like previously but cannot use it this way any longer because I have hundreds of them:

var Lurparken = L.marker([59.3011142,17.9939973], {icon: greenIcon}).bindPopup("Lurparken"),
    Aspuddsparken    = L.marker([59.3075995,17.9944431], {icon: greenIcon}).bindPopup("Aspuddsparken"),
    Sparvagnsparken    = L.marker([59.3298868,18.0031605], {icon: greenIcon}).bindPopup("Spårvagnsparken"),
    Fredhallsparkens    = L.marker([59.3320485,18.0029481], {icon: greenIcon}).bindPopup("Fredhällsparkens Plaskdamm"),
    Uggleparken    = L.marker([59.3343715,18.0040208], {icon: greenIcon}).bindPopup("Uggleparken"),
    Observatorielundens    = L.marker([59.3413877,18.056007], {icon: greenIcon}).bindPopup("Observatorielundens Parklek"),
    Kronobergsparken    = L.marker([59.331715,18.0335563], {icon: greenIcon}).bindPopup("Parklek Kronobergsparken"),
    Drakenbergsparken    = L.marker([59.314187,18.038608], {icon: greenIcon}).bindPopup("Drakenbergsparken");

Instead I want to store them like this where the array after 'greenIcon' is a set of attributes. You can see they vary by each marker, sometimes there are 2 or 3 attributes but there can be up to 10 in my set of data. Example:

const parks = [
    ['Spårvagnsparken', 59.3298868, 18.0031605, 'http://leksplay.com/sparvagnsparken', greenIcon, ['wc', 'grill', 'fence']],
    ['Fredhällsparkens plaskdamm', 59.3320485, 18.0029481, 'http://leksplay.com/fredhallsparkensplaskdamm', greenIcon, ['wc', 'water']],
    ['Uggleparken', 59.3343715, 18.0040208, 'http://leksply.com/uggleparken', greenIcon, ['wc', 'pfence',]],
    ['Observatorielundens Parklek', 59.3413877, 18.056007, 'http://leksplay.com/observatorielundensparklek', greenIcon, ['wc', 'pfence', 'toddler']],
    ['Parklek Kronobergsparken', 59.331715, 18.0335563, 'http://leksplay.com/parklekkronobergsparken', greenIcon, ['wc', 'indoor', 'pfence', 'toddler']],
    ['Drakenbergsparken', 59.314187, 18.038608, 'http://leksplay.com/drakenbergsparken', greenIcon, ['wc', 'grill', 'indoor', 'pfence', 'toddler']]
];

Now the journey begins. I have managed (structuring the data in the first way above) to manually place each marker inside a layer group so that the standard Leaflet controller can display or remove. I also managed to use dummy groups in order to prevent the map from removing markers that belong to two groups when one of  the active groups is removed. There is a working example of this here (excuse my invalid SSL cert, launching a site is complicated).

I am new to JavaScript. I am looking for a way to (1) create markers from my data (array 'parks'), toggle those markers on with leaflet's regular controller i.e.:

var overlayMaps = {
       "string": obj
  };

Then (2) display the intersection of attributes (i.e. 'wc', 'fence', 'grill') on the map when a checkbox is checked.


texts:
I am a newer using GRASS. I have tried to open the data named "nc_spm_08". But it couldn't open because errors as fallows:

       __________  ___   __________    _______________

GRASS GIS homepage:                      http://grass.osgeo.org
This version running through:            Bash Shell (/bin/bash)
Help is available with the command:      g.manual -i
See the licence terms with:              g.version -c
See citation options with:               g.version -x
If required, restart the GUI with:       g.gui wxpython
When ready to quit enter:                exit

Launching <wxpython> GUI in the background, please wait...
GRASS 7.2.svn (nc_spm_08):~ > ERROR 1: dlopen(/Library/Application Support/GDAL/2.1/PlugIns/gdal_ECW_JP2ECW.dylib, 1): no suitable image found.  Did find:
    /Library/Application Support/GDAL/2.1/PlugIns/gdal_ECW_JP2ECW.dylib: mach-o, but wrong architecture
    /Library/Application Support/GDAL/2.1/PlugIns/gdal_ECW_JP2ECW.dylib: mach-o, but wrong architecture
ERROR 1: dlopen(/Library/Application Support/GDAL/2.1/PlugIns/gdal_ECW_JP2ECW.dylib, 1): no suitable image found.  Did find:
    /Library/Application Support/GDAL/2.1/PlugIns/gdal_ECW_JP2ECW.dylib: mach-o, but wrong architecture
    /Library/Application Support/GDAL/2.1/PlugIns/gdal_ECW_JP2ECW.dylib: mach-o, but wrong architecture
ERROR 1: dlopen(/Library/Application Support/GDAL/2.1/PlugIns/gdal_MG4Lidar.dylib, 1): Library not loaded: /Library/Application Support/GDAL/2.1/Libraries/liblti_lidar_dsdk.1.dylib
  Referenced from: /Library/Application Support/GDAL/2.1/PlugIns/gdal_MG4Lidar.dylib
  Reason: image not found
ERROR 1: dlopen(/Library/Application Support/GDAL/2.1/PlugIns/gdal_MG4Lidar.dylib, 1): Library not loaded: /Library/Application Support/GDAL/2.1/Libraries/liblti_lidar_dsdk.1.dylib
  Referenced from: /Library/Application Support/GDAL/2.1/PlugIns/gdal_MG4Lidar.dylib
  Reason: image not found
ERROR 1: dlopen(/Library/Application Support/GDAL/2.1/PlugIns/gdal_MrSID.dylib, 1): Library not loaded: /Library/Application Support/GDAL/2.1/Libraries/libltidsdk.9.dylib
  Referenced from: /Library/Application Support/GDAL/2.1/PlugIns/gdal_MrSID.dylib
  Reason: image not found
ERROR 1: dlopen(/Library/Application Support/GDAL/2.1/PlugIns/gdal_MrSID.dylib, 1): Library not loaded: /Library/Application Support/GDAL/2.1/Libraries/libltidsdk.9.dylib
  Referenced from: /Library/Application Support/GDAL/2.1/PlugIns/gdal_MrSID.dylib
  Reason: image not found
dyld: Library not loaded: /Library/Frameworks/GDAL.framework/Versions/1.11/GDAL
  Referenced from: /Applications/GRASS-7.2.app/Contents/MacOS/bin/g.version
  Reason: image not found
Unable to get GRASS version
dyld: Library not loaded: /Library/Frameworks/GDAL.framework/Versions/1.11/GDAL
  Referenced from: /Applications/GRASS-7.2.app/Contents/MacOS/bin/g.version
  Reason: image not found
Traceback (most recent call last):
  File "/Applications/GRASS-7.2.app/Contents/MacOS/scripts/g.extension", line 1726, in <module>
    version = grass_version['version'].split('.')
KeyError: 'version'
GRASS 7.2.svn (nc_spm_08):~ >

I have installed all frameworks so I don't know how to do.

Could anyone help me?


texts: Need Basic Viewer template for ArcGIS JSAPI 3.2/3.3

texts:
I am looking for basic viewer for ArcGIS JSAPI  but not  getting  proper sample/example.

I am going through ArcGIS resource also but looks like its not working.  Going through this question also.

Checking out layout example and getting confused that what is difference between Viewer template and layout?

So any help regarding this will be great!


texts: getExtent() in OpenLayers returns strange values

texts:
I'm trying to get the bounding box of the map area displayed by OpenLayers.

I tried to use map.getExtent() because it should be exactly what I need, but the returned values are "strange":

For example, if I center the map around Grenoble (France, 45.19 5.71) with zoom level 14, the getExtent returns : 634409.225396,5647894.800465,643008.391077,5652672.114732
I guess I'm doing something wrong, but I can't figure out what. The getExtent() is supposed to return a bound object with Lat/Lon values...

The map is created with :

var map = new OpenLayers.Map('map', {
  projection: new OpenLayers.Projection("EPSG:900913"),
  displayProjection: new OpenLayers.Projection("EPSG:4326")
});

I'm using OSM as my base layer.


texts: Extract Sentinel-1 "VH" values for each record of FeatureCollection in Google Earth Engine?

texts:
I am new in using GEE and for my research I have to extract from the Sentinel-1 ImageCollection (IC), the VH values for each record of my FeatureCollection (FC), which is basically a multipolygon (in which every row is a polygon). This has to be done for a 5-years period (2015-2020). I am using the 'reduceRegions' function over each image of IC (Code here: https://code.earthengine.google.com/e993dfbe08c29d0cbcc494f48799f5a6):

var sentinel1 = ee.ImageCollection('COPERNICUS/S1_GRD');

var S1_ic = sentinel1
                 //filtering by polarisation
                 .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VH'))
                 //filter by image collected in interferometric wide swath mode
                 .filter(ee.Filter.eq('instrumentMode', 'IW'))
                 //filter by the period
                 .filterDate("2020-12-25","2020-12-31")
                 .select('VH');

print(S1_ic);

tos = ee.FeatureCollection(S1_ic.map(function (img) {
  var rR = img.reduceRegions({
    collection: tos,
    reducer: ee.Reducer.mean(),
    scale: 10}
    );
    return rR;
})).flatten();

Export.table.toDrive({collection: tos, folder: 'Export GEE'});

What I want to obtain is a table where for all the days sensed by Sentinel-1 I get a VH value for each polygon of my FC, for example:
01/01/2015: AN_1 (it is the code of the polygon) -> '-1.3'(VH), AN_2 -> '-1.8'...
02/01/2015: AN_1 -> '-1.5', AN_2 -> '-0.8'... and so on.
I get serveral problems with the outcome table, since there are no VH values at all. Moreover, for each feature there are too many S1 data, not only for an entire day but also for each day of the period, which is quite strange since the temporal resolution should be of few days.
Instead of using reduceRegions, I tried to use sampleRegions, which gives me back the VH values but I get too many data per polygon anyway.
I hope somone can help me to resolve this.


texts: How to export attribute table with start and end coordinates of lines to Excel?

texts:
I am working in QGIS. I have a georeferenced line layer and I want to export the attribute table and the coordinates of start and end points of my lines to Excel.

Can I do this with QGIS?


texts: SVG for symbols in ArcGIS JavaScript API

texts:
I'm developing a web app based on ArcGIS for Javascript. My customer's problem requires using many different point symbols, some of them graphically quite complex. I use PictureMarkerSymbols to create my symbols.

First I've used PNG's for my PictureMarkerSymbols. I'm not quite happy with the drawing quality, they look kind of pixelated.

Now I'm thinking about changing to SVG. My first tests showed much clearer symbols. But drawback seems to be that export to PDF and PNG doesn't work anymore. Also, it seems to me that it takes more time to build the map, and I'm getting timeout errors now and then.

Is anyone experienced with SVG's for PictureMarkerSymbols and knows about the drawbacks?


texts: Calculate distance from a point to the border of polygon (US state)

texts:
I have a point, and I want to know if it is nearer than 100 miles from some US State border (California for example).

I assume that I would need to use this: US Regions State Boundaries. But if I compare my point with the points of a state polygon, I will only get a result if my point is near some of the polygon's points and not the border itself.

If this is possible, is possible to do something similar with ZIP Codes? For example: I want to get 33139, and every location on my array of locations that is near it.


texts: Label Symbology Pie Charts/Multiple Bar Graph

texts: I created bar graph and pie charts on every single feature on the basis of their values stored in attribute table but I can't label/plot the values on these graph/on the top of graphs. Currently this could be done by coverting the labels to annotation and then adjusting the location manually, This is really hectic and time taking way. Do we have any option to add/label values on the bar graph/ pie charts?


texts: Convert lat/lon to Mercator x/y, then measure distance

texts:
I'm attempting to convert latitude/longitude coordinates to X/Y coordinates on a 52,000 sq. ft. grid according to a simple Mercator projection.  I used the following Python code to accomplish this, adapted from this post.

from math import pi, log, tan

latitude    = 42.374994
longitude   = -71.075891

mapWidth    = 52000.000000
mapHeight   = 52000.000000

x = (longitude+180)*(mapWidth/360)
print "x", x

latRad = latitude*pi/180
print "latRad", latRad

mercN = log(tan((pi/4)+(latRad/2)))
print "mercN", mercN

y     = (mapHeight/2)-(mapWidth*mercN/(2*pi))
print "y",y

Using the examples of 30th Street Station and 69th Street Terminal in Philadelphia, I get the following results:

30th Street Station
#39.955641,-75.181178
x 15140.4965111
latRad 0.697357490195
mercN 0.761899319942
y 19694.4777688

69th Street Terminal
#39.962156,-75.258626
x 15129.3095778
latRad 0.697471198396
mercN 0.7620476662
y 19693.2500467

If I calculate the distances between these X/Y coordinates using the Pythagorean theorem, I get a value of

sqrt((15129.30958-15140.49651)^2 + (19693.25005-19694.47777)^2)
= sqrt(125.1474774)
= 11.25410054

I am assuming that in this case the result would take the same units as my mapWidth and mapHeight, in feet in this case.  But, a value of 11.25 feet is not reasonable.  When measured on Google Maps, these two sites are roughly 4 miles apart.  What am I missing?


texts: How can I symbolize a layer by multiple attributes (category and amount) in ArcGIS Pro 2.5?

texts:
I have a layer of properties that are categorized by type (text field - farm, house, business, etc.) and cost (numeric field). I want to assign each property type a different-shaped symbol (all of the same size) and to classify it to one of five cost categories by colour.
There are 120 properties, 5 cost classes to which they will be assigned, and 5 property types.
For example, all houses will be squares. Houses in cost class 1 will be green; cost class 2 will be blue; class 3 will be red; etc.
I am using ArcGIS Pro 2.5 and will be publishing the layer to Portal for ArcGIS.


texts: QGIS action that sets transparency in RGB raster

texts:
In the question: "GIS Action that loads a raster, creates pyramids and sets color" I have resolve the question for raster in B/W but now I need it to load RGB raster with white background that I would like transparent.

I think to need QgsRasterTransparency::alphaValue but I have some problem  how to use it in code. Any help is welcome

from PyQt4.QtCore import QFileInfo,QSettings
from qgis.core import QgsRasterLayer, QgsCoordinateReferenceSystem
import gdal
import os
from subprocess import call

#Define gdaladdo
gdaladdoFile = 'C:\\Program Files\\QGIS 2.18\\bin\\gdaladdo.exe'

s = QSettings()
oldValidation = s.value( "/Projections/defaultBehaviour" )
s.setValue( "/Projections/defaultBehaviour", "useProject" )
vl = QgsMapLayerRegistry.instance().mapLayersByName('catalogo_CTR')[0]
qgis.utils.iface.setActiveLayer(vl)
root = QgsProject.instance().layerTreeRoot()
group_name = "CTR_10000"
group = root.findGroup(group_name)
if group == None:
    group = root.addGroup("CTR_10000")
else:
    pass
fileName = 'D:\\GIS\\Dati\\IDT_CTRR_Raster\\[% A_CODICE %].tif'
fileInfo = QFileInfo(fileName)
baseName = '[% A_CODICE %]'
rlayer = QgsRasterLayer(fileName, baseName)

call([gdaladdoFile, '-ro', '--config', 'USE_RRD', 'YES', fileName, '2 4 8 16'])

raster_transparency  = rlayer.renderer().rasterTransparency()
ltr = QgsRasterTransparency.TransparentSingleValuePixel()
ltr2 =  QgsRasterTransparency.TransparentSingleValuePixel()
tr_list = []
ltr.min = 0  # Or another value
ltr.max = 0  # Or another value
ltr.percentTransparent = 100  # Or another value
ltr2.min = 8  # Or another value
ltr2.max = 8  # Or another value
ltr2.percentTransparent = 100
tr_list.append(ltr)
tr_list.append(ltr2)

rlayer.triggerRepaint()  # Tried with iface.mapCanvas().refresh(), but it didn't work

raster_transparency.setTransparentSingleValuePixelList(tr_list)
crs = QgsCoordinateReferenceSystem()
crs.createFromSrid(3003)
rlayer.setCrs(crs)
QgsMapLayerRegistry.instance().addMapLayer(rlayer, False)
group.insertChildNode(-1, QgsLayerTreeLayer(rlayer))

s.setValue( "/Projections/defaultBehaviour", oldValidation )


texts: Does ArcGIS online allow the use of ArcMap print templates?

texts:
The ArcGIS Server help page Publishing additional services for printing shows how to publish an MXD as a print template:

  Navigate to Toolboxes > System Toolboxes > Server Tools > Printing.
  Double-click the Export Web Map tool. Click Geoprocessing > Results to
  display the Results window. Right-click the ExportWebMap result and
  click Share As > Geoprocessing Service.

These instructions seem to apply to ArcGIS Server only. How can I use a print template, derived from an MXD, within ArcGIS Online?

(The question Setting up ArcGIS Print service for ArcGIS Online account? shows that ArcGIS Online supports a print service, but doesn't confirm that ArcGIS Online supports a custom template derived from ArcMap.)


texts: How to georeference a raster based on identical, already georeferenced raster

texts:
I have two identical rasters, one has been georeferenced, the other has not. I would like to georeference the second raster automatically using the first one. I do not want to create links by hand and the raster has been warped so simply shifting xy position does not work.

Let me know if you have any ideas... most comfortable working in model builder but would love any solution.

Thanks!


texts: ogr2ogr - function getsrid(geometry) missing error {PG to File}

texts:
I am issuing the following statement

ogr2ogr -a_srs EPSG:27700 -f "MapInfo file" "test.tab" PG:"host=xxx user=xxx dbname=xxx password=xxx active_schema=xxx" -sql "Select * from sale_region"

and I am getting the following message on completion:

  ERROR 1: ERROR: function getsrid(geometry) does not exist

  LINE 1: SELECT getsrid("wkb_geometry") FROM (Select * from sale_regi...

  HINT: No function matches the given name and argument types.  You
  might need to add explicit type casts."

All of the data appears to be in the generated file (I have also tried with exporting to SHP and Geojson with the same message appearing).

I think the error message has somethign to do with the SQL and/or that my geometry column is called wkb_geometry in the table.  So I'm obviously doing something wrong somewhere and am worried that some data may have been lost in translation.  If anyone can put my mind at rest or tell me whats wrong with the code, it'd be very much appreciated.

I did try:

ogr2ogr -a_srs EPSG:27700 -f "MapInfo file" "test.tab" PG:"host=xxx user=xxx dbname=xxx password=xxx active_schema=xxx" "sale_region"

and although it completed with no errors at the command line, it produced an empty MapInfo table.


texts: Auto-Update Layers in ArcCatalog from modified layers from .DWG

texts:
let's say that we have 2 layers in autocad .dwg :
1 layer specify the roads and 1 layer specify buldings.

Now if we create in ArcCatalog (Import-> Feature class single..) and load the .dwg and creating a layer (polyline) with name "Buildings" is any option after creating the polyline to be modified (auto-updated) if we will add more buildings in the .dwg file?

P.S actually i delete the layer (polyline) in arccatalog and recreate them from scratch.

Thanks in Advance


texts: I’m creating layers with qgis 3.4 raster calculator. How can I go back and modify the  calcs or do I have to start over every time?


texts: Reading values from dynamically created Qlineedit boxes using PyQGIS?

texts:
I am creating dynamic labels and line edit boxes based on the number of fields in the shapefile and populating the lineedit boxes using the feature attributes using the below code. Each lineedit box has the field name as object name.

I need to know how to read the current value of each dynamically created lineedit boxes.

    def addLayerNames(self):
     for i in reversed(range(self.dlg.verticalLayout.count())):
        self.dlg.verticalLayout.itemAt(i).widget().setParent(None)
     layers = QgsMapLayerRegistry.instance().mapLayers().values()
     layer_list = []
     for layer in layers:
        if layer.type() == QgsMapLayer.VectorLayer:
            layer_list.append(layer.name())
     self.dlg.comboBox.addItems(layer_list)
     selectedLayerIndex = self.dlg.comboBox.currentIndex()
     selectedLayer = layers[selectedLayerIndex]
     fields = selectedLayer.pendingFields()
     fieldnames = [field.name() for field in fields]
     for i,field in enumerate(fieldnames):
        self.labels = QLabel()
        self.linedit = QLineEdit()
        self.labels.setText(field)
        self.linedit.setFixedWidth(100)
        self.dlg.verticalLayout.addWidget(self.labels)
        self.dlg.verticalLayout.addWidget(self.linedit)
        request = QgsFeatureRequest().setFilterFid(0)
        layer.setSelectedFeatures([0])
        self.iface.actionZoomToSelected().trigger()
        feat = layer.getFeatures(request).next()
        self.linedit.setText(str(feat[i]))
        self.linedit.setObjectName(field)
        self.dlg.pushButton.setEnabled(False)


texts: GDAL/OGR ogrinfo for shapefiles in python

texts:
This might be a dumb question, but I wasn't able to find the answer anywhere.

Reading the metadata of a raster file in Python using GDAL is as easy as:

result = str(gdal.Info(raster.tif))
print result

In this case I want to do the same for a shapefile. I know I can use ogrinfo -ro -so -al foo.shp in the shell, but I have not found a way to do the same only in Python within a single line.

I'm looking for something like:

result = str(ogr.Info(shapefile.shp))
print result

Am I missing something obvious or is there really no function like this?


texts: Point Symbolizer() in rundemo.cpp file Giving error

texts:
The code is given Below:

roads2_rule_11.append(point_symbolizer("abc.png"));

The error is given below:

rundemo.cpp: In function âint main(int, char**)â:

rundemo.cpp:142:57: error: no matching function for call to âmapnik::point_symbolizer::point_symbolizer(const char [8])â
rundemo.cpp:142:57: note: candidates are:
/usr/local/include/mapnik/point_symbolizer.hpp:47:5: note: mapnik::point_symbolizer::point_symbolizer(const mapnik::point_symbolizer&)
/usr/local/include/mapnik/point_symbolizer.hpp:47:5: note:   no known conversion for argument 1 from âconst char [8]â to âconst mapnik::point_symbolizer&â
/usr/local/include/mapnik/point_symbolizer.hpp:46:5: note: mapnik::point_symbolizer::point_symbolizer(mapnik::path_expression_ptr)
/usr/local/include/mapnik/point_symbolizer.hpp:46:5: note:   no known conversion for argument 1 from âconst char [8]â to âmapnik::path_expression_ptr {aka boost::shared_ptr, mapnik::attribute> > >}â
/usr/local/include/mapnik/point_symbolizer.hpp:45:5: note: mapnik::point_symbolizer::point_symbolizer()

I have an image that i want to show at each point .How to input point_symbolizer()?


texts: Aligning and overlaying maps from different providers

texts:
I am quite new to this whole GIS domain, so I was hoping someone could shed some light into a problem I have been trying to solve for a while.
I consume data from 2 different providers:
http://www.bom.gov.au/australia/meteye/
https://absstats.maps.arcgis.com/apps/MapSeries/index.html?appid=d48afe3794c24eb1a9cb531ee2408bee
My goal is to align / overlay the maps between these two providers so I can feed it to some other analytical tools.
I had a look in the server configurations for each one:
https://tiles.arcgis.com/tiles/v8Kimc579yljmjSP/arcgis/rest/services/Population_Grid/MapServer?f=json
http://www.bom.gov.au/australia/meteye/meteye-config.json?_rls=598b5b1,
and noticed that absstats has a set number of zoom levels, while meteye seems more flexible by manipulating the BBOX, WIDTH and HEIGHT.
I started ingesting the meteye data with the following single map for the whole country: http://wvs1.bom.gov.au/mapcache/meteye?TRANSPARENT=true&FORMAT=image%2Fpng&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap&STYLES=&layers=IDZ73031&TIMESTEP=0&SRS=EPSG%3A4326&BBOX=108,-48,162,0&WIDTH=512&HEIGHT=512&mapScale=0.15625, but I might need to zoom in to get a higher resolution in both.
If I select the third zoom level in abstats, for example, how can I manipulate the meteye REST call to match the same area and zoom?


texts: Is it possible to rearrange layers on the fly?

texts:
I'm involved in a project where one of the main features is that the users will have the ability to select wich layers will be available for him to view on the map, the layers then should be loaded into the map and the toc components.

The project is web based and it will be done using Visual Studio .NET 2008 with Web ADF and Argcis Server 9.3.

Does anyone knows how this can be done? Thanks.


texts: How to edit a featureService from ArcGIS for Android?

texts:
I've been trying to edit data in a FeatureService, hosted on my own ArcGIS Server (I'm not using ArcGIS Online).

I have created a map in my ArcGIS Server's Mobile Content Directory, using the following JSON:

{
    "version" : "1.1",
    "operationalLayers" : [{
            "url" : "http://example.com/ArcGIS/rest/services/Bridges/FeatureServer/0",
            "visibility" : true,
            "opacity" : 0.75,
            "title" : "Bridges"
        }
    ],
    "baseMap" : {
        "baseMapLayers" : [{
                "url" : "http://services.arcgisonline.com/ArcGIS/rest/services/World_Topo_Map/MapServer"
            }
        ],
        "title" : "World Topographic Map"
    }
}

I've reached upto this point, by using Dave Bouwman's Blog and this ESRI Documentation.

I can add the map to the ArcGIS App on my Android Device. The data shows up, but the Edit tool (the little pencil Icon) does not show up.

Using this video from ESRI, I can get an editable map via ArcGIS online, but I want an editable map without ArcGIS Online.

Is it possible to edit a feature service directly, without using ArcGIS Online?


texts: Iterative Mean of Means Infinite Loop for calculating H-t Index

texts:
I'm trying to create a function that will calculate an Ht-Index for Quantifying the Fractal or Scaling Structure of Geographic Features (but am new to R)

which basically just finds nested means of continuous data (raster in this case) and separates them into 2 groups: heads and tails (values above and below the means). The Ht-Index is the number of times (-1) the 'heads' group can be subdivided without being greater than some threshold (e.g., 45%) of the distribution of both groups. I'm not sure how to efficiently do this in R but here is my failed attempt which creates an infinite loop iterating over the same input.

require(raster)

foo <- raster('/foobar.img')

hindex <- function(x){
iCounter <-0
repeat{
    head <- x[x>mean(x)]  # need to replace x each iteration for head and tail
    tail <- x[x<mean(x)]
    headpct <- length(head)/sum(length(head),length(tail))
    if (headpct >.45) break # If head is > 45% of values stop
        else
            if (length(head)== 1) break # if head is only one value stop
    x <- head
    iCounter <- iCounter+1
    }
return(iCounter-1)
}


texts: Converting Canadian census geocode to postal code?

texts:
I am looking to convert the geocode on the 2011 Canadian Census data file to postal code for free. The geocode looks like the following - 0010002.00.

I understand that Canada Post sells a PCCF file that will accomplish this.

geocoder.ca has file that will convert from long and lat to postal code but I can't see how to convert using the geocode in the census file.

Apologies if this is answered elsewhere in this forum, but I didn't see it on a search.


texts: How to upload a group of Shapefiles into QGIS?

texts:
I have a folder with 211 Shapefiles that I want to upload as a group into my project.

I have been looking for a way to upload all the files simultaneously, but I have failed to find one.

Is there a way to upload all these files in one go?


texts:
I am doing a work which I need to sample random points for many layers (classification of land use). Then, I'm trying to do a loop with my several layers and applied a functions of random pointos. I would like to use Random Points... functions (tool that I find in Menu>>Vector>>Research Tools>> Random Poits) because there I would choose 20 points pet class. This tool allows to make this.

It function (Random points...) belongs to fTools plugin, however I didn't find it functions in python console. All information that I saw in internet wasn't help me (or I couldn't understand).

How can I put this tool (Random points...) in my script? Is it possible?


texts: Artefacts when merging different resolution DEM

texts:
I am trying to merge a 10m DEM with a 30m DEM (10m DEM is a sub area of the 30m DEM). Both DEMs have the same projection and nodata value. I create a .vrt file:
gdalbuildvrt test.vrt *.tif -resolution highest -r cubicspline

That works fine, but at the border from high to low resolution, I get this artefacts:

Offset is not consistent:

10m data only doesn't produce any artefacts at the image border:

The first image shows hillshades I generated from the test.vrt. I couldn't find any solution to avoid getting this line at the transition between 10m and 30m data. The 10m data comes from an airborne laser, the 30m from satellite radar. So I guess only because of this different sources, the values are different and produce a gap.
Is there a possibility to tell GDAL to blend the datasets into each other (only at the outline of the 10m dataset)?


texts: QGIS 2.16 Composer - Georeference tif/pdf to import to Avenza PDF map

texts:
I want to create a geo-referenced tif or pdf to import to Avenza PDF map.

The Feature list for QGIS 2.16 states that composer creates a geo-referenced  file by default.

  "The existing option to create a world file has been separated from
  the map selection for georeferencing. The new behaviour is to always
  georeference outputs, and only create the separate world file if that
  option is checked." - reference
  https://www.qgis.org/en/site/forusers/visualchangelog216/index.html#feature-georeference-outputs-eg-pdf-from-composer

I have been able to create a georeferenced image but it is of no use and will not load into my map because it is in pixel coords.

C:\output>gdalinfo r4_composer_ix.tif
Driver: GTiff/GeoTIFF
Files: r4_composer_ix.tif
Size is 7015, 4960
Coordinate System is `'
Metadata:
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  TIFFTAG_XRESOLUTION=299
  TIFFTAG_YRESOLUTION=299
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 4960.0)
Upper Right ( 7015.0,    0.0)
Lower Right ( 7015.0, 4960.0)
Center      ( 3507.5, 2480.0)
Band 1 Block=7015x1 Type=Byte, ColorInterp=Red
Band 2 Block=7015x1 Type=Byte, ColorInterp=Green
Band 3 Block=7015x1 Type=Byte, ColorInterp=Blue
Band 4 Block=7015x1 Type=Byte, ColorInterp=Undefined

Is there a trick to getting the output to be georeferenced in the EPSG the project is using?

I also cannot get a '.tfw' file created when I toggle 'Save World File' On.


texts: Symbol Rotation for all features along line in ArcGIS Pro

texts:
I have a waterline with accesories set along the pipes and I would like to rotate the symbols towards the pipe direction (Paralel, perpendicular, etc Depending on the fitting or device type). For that, I created a Simbol Rotation field which does it by writting the angle, and it works fine but I have to do, object by object individually and this will take ages for all the network. It is any way to do it automatically?
Example: Take all Tees on the network and rotate them all in one step depending on the direction of the pipes


texts: How to align (georeference) images using QGIS?

texts:
I'm completely new to all this.

I'm trying to make a map of Han dynasty-era China for a mod for a game called Crusader Kings 2.

I've downloaded QGIS and Natural Earth 2, but I need this image to align with the rest of the map so I can use it as a guide for province boundaries.

I tried georeferencing the image, and I managed to get it added as a layer, but even with on-the-fly projection conversion, QGIS refused to show it at the same time as any other layer.

This was the result.

I'm also open to suggestions on any simpler methods I could use if QGIS is overkill for what I'm trying to do.


texts:
I want to change the shape extent of some shape files, I would like to know how to do it in QGIS, since I rasterized some shapes and when I used the calculator to add them, the resulting raster had deleted zones with null data...
What i'm trying to do is to add the boolean rasters of the distribution of species into one, to obtain a single one that shows the # of species in a region....

The first image shows the resulting raster (after I use raster calc to add every raster), the second one is the raster I made of the original shape file I used to clip all of my shapes for the region. When I inspect the resulting raster, I get null values (the blue colored section, so the image seems incomplete...


texts: Using 'Split Features' tool to split a polyline at a specific node with QGIS

texts:
In the example shown in the picture below I would like to split the light blue polyline into two polylines on a specific node (circled with number 1) with QGIS (version 2.12.3-Lyon).

Given that the layer is already editable, and given that I have set the snapping options at 1 map unit tolerance and snap it to vertex, which tool should I use, number 2 or number 3? (I get the names in Spanish and the picture seems easier to identify for anyone).

How do I select the node (number one) once the splitting tool is active? So far, if I just click on the node, the entire polyline disappears.

I understand the explanations already given by @Knightshound and @Joseph in answers to question 'Polyline split at node in QGIS' but need more detail into how to correctly proceed.


texts: What apps offer the best ROI for city government?

texts:
I work in a government agency. At a town hall.

I'm looking for ideas and suggestions involving technology and GIS maps. What are the trends that are up in the public sector?

What is the relative return on investment (ROI) for applications in these areas: health, transit, education, social welfare, sport and recreation, culture, youth, housing, public safety, urban mobility, sewage, environment, disaster management, environmental management, urban economy, rural economy, management PMI , community governance.


texts: Google maps layer as the base layer and WMS as overlays shows empty map

texts:
I have a problem of showing my base layer which is Google Maps. The browser console reports

  "Uncaught TypeError: undefined is not a function" in     var
  google_streets = new OpenLayers.Layer.Google(     "Google Streets",
  {numZoomLevels: 20}    );.

This is my code:

        var map;
        var untiled;
        var tiled;
        var pureCoverage = false;
        var format;

        var googleMercator = new OpenLayers.Projection("EPSG:4269");
        var baseProjection = new OpenLayers.Projection("EPSG:4326");

        OpenLayers.IMAGE_RELOAD_ATTEMPTS = 5;

        OpenLayers.DOTS_PER_INCH = 25.4 / 0.28;

        function init(){

            format = 'image/png';
            if(pureCoverage) {
                document.getElementById('filterType').disabled = true;
                document.getElementById('filter').disabled = true;
                document.getElementById('antialiasSelector').disabled = true;
                document.getElementById('updateFilterButton').disabled = true;
                document.getElementById('resetFilterButton').disabled = true;
                document.getElementById('jpeg').selected = true;
                format = "image/jpeg";
            }

            var bounds = new OpenLayers.Bounds(
                -179.147340000341, -14.37879999982,
                179.778480000052, 83.1480859363171
            );
            var options = {
                controls: [],
                maxExtent: bounds,
                maxResolution: 1.4020539843765352,
                projection: googleMercator,
        displayProjection: baseProjection,
                units: 'degrees'
            };
            map = new OpenLayers.Map('map', options);

            // setup tiled layer
            tiled = new OpenLayers.Layer.WMS(
                "US_Canada:northus - Tiled", "http://localhost:8080/geoserver/US_Canada/wms",
                {
                    LAYERS: 'US_Canada:northus',
                    STYLES: '',
                    format: format,
                    tiled: true,
                    transparent: true,
                    tilesOrigin : map.maxExtent.left + ',' + map.maxExtent.bottom
                },
                {
                    buffer: 0,
                    displayOutsideMaxExtent: true,
                    isBaseLayer: false,
                    yx : {'EPSG:4269' : true}
                }
            );

            // setup single tiled layer
            untiled = new OpenLayers.Layer.WMS(
                "US_Canada:northus - Untiled", "http://localhost:8080/geoserver/US_Canada/wms",
                {
                    LAYERS: 'US_Canada:northus',
                    STYLES: '',
                    format: format,
                    transparent: true
                },
                {
                   singleTile: true,
                   ratio: 1,
                   isBaseLayer: false,
                   yx : {'EPSG:4269' : true}
                }
            );
            var google_streets = new OpenLayers.Layer.Google(
                "Google Streets",  {numZoomLevels: 20}
            );

            map.addLayers([ google_streets, untiled, tiled]);

            map.addControl(new OpenLayers.Control.PanZoomBar({
                position: new OpenLayers.Pixel(2, 15)
            }));
            map.addControl(new OpenLayers.Control.Navigation());
            map.addControl(new OpenLayers.Control.Scale($('scale')));
            map.addControl(new OpenLayers.Control.MousePosition({element: $('location')}));
}


texts: Seeking Open Source Alternatives for FME

texts:
Do you know any good open source alternatives for FME?

Searching for some free but powerfull alternatives for geodata management, workflows (Python APIs), etc.


texts: Symbolize External ESRI GeoJSON by Local GeoJSON

texts:
I'm working in Jupyter Notebook and attempting to symbolize an external ESRI GeoJSON service (now an output available from ArcServer 10.5) by an attribute of a local GeoJSON, and return a choropleth map.  The premise is one of combining customer data with a publicly available service.

If the two GeoJSON files are both local, I have no problem loading them with geopandas and achieving the desired result on a folium map.  However the utilization of an ESRI GeoJSON service url is proving difficult.

Anyone know of a way to do this with an external GeoJSON?

Am I missing a transformation of the ESRI GeoJSON output?

Code:

    import geopandas as gpd
    import folium
    import requests
    import json

    mymap = folium.map([42.65, -73.75],
        zoom_start=7,
        tiles = 'cartodbpositron')

    geo = gpd.read_file('local.geojson')

    url = json.loads(requests.get('https://gisservices:domain/arcgis/rest/services/...f=geojson').text)

    mymap = choropleth(geo_path = url,
            fill_color='YlGn', fill_opacity=0.5, line_opacity=0.5,
            threshold_scale = [50,100,250,500,750,900],
            data = geo,  #data source I want to pull the attribute from
            key_on='feature.properties.NAME', #key field from external GeoJSON
            columns = ['COUNTY', 'Value']  #key field from local GeoJSON followed by the field I want to symbolize by
            )

Examining the data and the request I get a 200 response and what looks like a properly formatted geojson.

    {'crs': {'properties': {'name': 'EPSG:4326'}, 'type': 'name'},
        'features': [{'geometry': {'coordinates': [[[-75.88982572053291,
             42.723834549294345],
             [-75.88998066913965, 42.725786369710484],
             [-75.89004324847727, 42.72657462151952],
             [-75.8900956058853, 42.72723409146341],
             ...

However the final portion of the code returns an error that the geo_path for the folium choropleth contains an invalid file:
"TypeError: invalid file: {'type': 'FeatureCollection', 'crs':... "


texts: ArcGIS Javascript API - Update Renderer when a Feature Layer record is selected

texts:
I'm using ArcGIS JavaScript API to display a map with a feature layer with multiple records.

I wrote code to display a pop up when a record is selected on the map. Is it possible to change the color of the selected record to some other color by updating the renderer when its clicked on and persist the change?

I know the Javascript Api highlights the boundary when a record is clicked on but it doesnt persist, when another record is clicked the previous clicked record go backs to its original color.


texts:
I have list of geometry type in a database or:

List<com.vividsolutions.jts.geom.Polygon> polygons

Is there any JTS Topology Suite method or PostGIS function which would output intersecting groups as unions of Polygon?

For example, if I provide list of polygons(boxes) like below, I would expect to get list of unions: one made of 1, 2, 10 boxes, another of 5,7, etc.


texts: Distance and time fields are null in the output of Matrix From Layer algorithm of ORS Tools plugin in QGIS

texts:
I tried to calculate the travel distance and the travel time using Matrix From Layer algorithm of the ORS Tools plugin in QGIS, but in the attribute table of the output, the duration and distance columns are null.
How can I solve this problem?


texts: OpenLayers doesn't display US Census-based GeoJSON Layer

texts:
I'm getting started with OpenLayers and all was going well until I tried to use my own GeoJSON file.

I started with the Cartographic Boundary Shapefiles - Counties from the Census website, loaded it into http://www.mapshaper.org/ and exported it as a GeoJSON file.  However when I use that as my source (either a straight Vector layer or as an Image layer) it doesn't display anything. No errors- just a empty window with zoom controls. I suspect it may have something to do with projections so I made it the only layer, but it didn't help.

   <script type="text/javascript">
  var map = new ol.Map({
    target: 'map',
    layers: [
      //new ol.layer.Tile({
      //  source: new ol.source.Stamen({layer: 'toner-lite'})
      //}),
      new ol.layer.Vector({
        title: 'US Counties',
        source: new ol.source.Vector({url: 'data/counties.geojson', format: new ol.format.GeoJSON() }),
        style: new ol.style.Style({
          fill: new ol.style.Fill({
            color: 'blue'
          }),
          stroke: new ol.style.Stroke({
            color: 'olive',
            width: 10
          })
        })
      })
    ],
    view: new ol.View({
      center: "center",
      zoom: 1
    })
  });
</script>


texts: Are Date columns in feature classes stored as python datetime objects?

texts:
I have a Date column called readings_measurement_tstamp stored in my feature class. I would like to add a String column to this feature class and Field Calculate it to be be just the %H:%M:%S component of readings_measurement_tstamp. After adding the new column I calculated it according to the following expression:
    datetime.datetime.strftime( !readings_measurement_tstamp! ,'%H:%M:%S')

I was under the impression that Date columns can be referenced like python datetime objects. Therefore, I was very confused when running the field calculator with this expression gave the following error: Descriptor 'strftime' requires a 'datetime.date' object but received an 'str'

The error implies to me that readings_measurement_tstamp is being stored as text, not as a Date. Any ideas what I need to do in order to be able to access readings_measurement_tstamp as a datetime object? Thanks!


texts: Finding lowest elevation point among points common in two vectors for each polygon

texts:
I have two shapefiles

Administrative boundary
Stream extract flow (using r.stream.extract)

and Raster DEM file

I want to find lowest elevation point in each administrative polygon and this point should also be present on stream extract flow

Image 1 shows when stream_flow vector file is taken with polygon shapefile along with DEM

Image 2 shows when stream_flow raster file is taken with polygon shapefile along with DEM


texts: Why is ArcPy script to apply symbology not working?

texts:
This is the first python script I have tried to write and I can't find out why it is not working. All I want to do is apply symbology from one layer to another.

Working Code:

# Name: ApplySym.py
# Purpose: apply the symbology from one layer to another

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

env.workspace = "C:\\Users\\jrutledge\\Desktop\\ArcGIS"

inputLayer = "TEST.lyr"

symbologyLayer = "LTSig_Template.lyr"

arcpy.ApplySymbologyFromLayer_management(inputLayer, symbologyLayer)

Error Message:

    Traceback (most recent call last):
  File "C:\Users\jrutledge\Desktop\Well DB Run 2_2015 (plus GAMA)\ArcMaps2015_2\Python\ImportSymbology.py", line 11, in <module>
    inputLayer = arcpy.Layer(r"C:\Users\jrutledge\Desktop\Well DB Run 2_2015 (plus GAMA)\ArcMaps2015_2\MapOutput\ByWell\ShapeFiles\T00010_SOURCE_TEMPERATURE_C_C_ Events.lyr")
AttributeError: 'module' object has no attribute 'Layer'

Code:

# Purpose: apply the symbology from one layer to another

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

env.workspace = r"C:\Users\jrutledge\Desktop\Well DB Run 2_2015 (plus GAMA)\ArcMaps2015_2\MapOutput\ByWell\ShapeFiles"

inputLayer = arcpy.mapping.Layer(r"C:\Users\jrutledge\Desktop\Well DB Run 2_2015 (plus GAMA)\ArcMaps2015_2\MapOutput\ByWell\ShapeFiles\T00010_SOURCE_TEMPERATURE_C_C_ Events.lyr")

symbologyLayer = arcpy.mapping.Layer(r"C:\Users\jrutledge\Desktop\Well DB Run 2_2015 (plus GAMA)\ArcMaps2015_2\MapOutput\ByWell\ShapeFiles\T71850_NITRATE_LT_DIR.lyr")

arcpy.ApplySymbologyFromLayer_management(inputLayer, symbologyLayer)

Old Code(original attempt): RuntimeError: Object: Error in executing tool

# Name: ApplySym.py
# Purpose: apply the symbology from one layer to another

# Import system modules
import arcpy
from arcpy import env

# Set the current workspace
env.workspace = "C:\Users\jrutledge\Desktop\Well DB Run 2_2015 (plus GAMA)\ArcMaps2015_2\MapOutput\ByWell\ShapeFiles"

# Set layer to apply symbology to
inputLayerFile = arcpy.mapping.Layer(r"C:\Users\jrutledge\Desktop\Well DB Run 2_2015 (plus GAMA)\ArcMaps2015_2\MapOutput\ByWell\ShapeFiles\T00010_SOURCE_TEMPERATURE_C_C_ Events.lyr")
inputLayer = arcpy.mapping.ListLayers(inputLayerFile, "T00010_SOURCE_TEMPERATURE_C_C_ Events")

# Set layer that output symbology will be based on
symbologyLayerFile = arcpy.mapping.Layer(r"C:\Users\jrutledge\Desktop\Well DB Run 2_2015 (plus GAMA)\ArcMaps2015_2\MapOutput\ByWell\ShapeFiles\T71850_NITRATE_LT_DIR.lyr")
symbologyLayer = arcpy.mapping.ListLayers(symbologyLayerFile, "T71850_NITRATE_LT_DIR")

# Apply the symbology from the symbology layer to the input layer
arcpy.ApplySymbologyFromLayer_management (inputLayer, symbologyLayer)


texts: Qgis select distinct by default on Print Composer

texts:
I have a filter on Print Composer that selects "id_mov" and shows figures one by one into the map. This "id_mov" is not unique so different id can have the same "id_mov".
Qgis select distinct by default and for the same "id_mov" I can't find the way to show them all together in the same map composition. It could be like saying in SQL "no select distinct"


texts: I have an ArcGIS StoryMap (Cascade variant) that I would like to load an image into that someone can pan around. The image is of an old Sanborn map so zooming in and out would be a good idea. The current service that hosts the image allows me to embed it using a link (see here: http://dc.lib.unc.edu/cdm/ref/collection/ncmaps/id/9084), but it also includes all the other page information.


texts: Does a tool exist to visualize custom OSM changesets?

texts: I have custom changesets, not yet applied to the OSM map. I want to be able to visualize them before uploading them. A tool like achavi (https://wiki.openstreetmap.org/wiki/Achavi#Changeset_Viewer) is what I am looking for. Achavi however only works with history from OSM and not with custom changesets, not yet applied to OSM.


texts: A nice display of hillshade and my geologic map

texts:
i want to create a nice hillshade map covered by a grid or RGB map, with similar colour results  in ArcMap like this:

Is there a way to achieve this in ArcGis without transparency settings which fades colours, and not QGIS, or GlobalMapper?

More like a blend in, described in various threads but with ArcGIS Desktop 10.


texts: Deploying OpenGeo Suite App on Apache Tomcat?

texts:
I am a beginner trying to develop a web app using Opengeo Suite. I've manage to test it localy on my personnal computer and I would like to deploy it to the web. It seems like Tomcat Apache is what I need to work with but I can't seem to figure how to get started.

I know it is a "wide" question but I just need to know where to start, I'm usually good at learning by myself and figuring things out but I've been reading and searching on the web for a couple of days and can't really find how to do this. I have used ngrok (https://ngrok.com/) to share my app quickly but only the layout and base layer would show up (the overlay layers are missing).

Can someone give me an idea of where to start, or the major steps I should take in order to make my app available on the web?


texts: Calculating Kernel Density output values based on parameter value in ArcMap?

texts:
I am generating a map where I have sentiment value for geo-tagged data. I want to show a heat map based on the sentiment of location. I am using kernel density to generate this map.

I am having trouble in understanding how does the 'kernel density' function calculates output density values when I am entering the sentiment as the parameter?

What role does the radius and area parameter play in the result?

Is kernel density the right option to generate such a map?

I am using ArcMap 10.3.


texts: Intersections where streets have different names

texts:
In Overpass Turbo, I want to find all intersections where the intersecting streets have different names.
I imagine the query would be based on this---
[bbox:{{bbox}}];
way[highway][highway!~"^(footway|path)$"]->.a;
way[highway][highway!~"^(footway|path)$"]->.b;
node(w.a)(w.b);
out;

But I want it to exclude any 'junction' that is merely a single way that has been split up.
Essentially I can't figure out how to query node(w.a)(w.b) where [name of .a]!=[name of .b] (sorry for the pseudocode).


texts: GIS confidentiality, how to convert address to block level

texts:
I have a list of addresses in a spreadsheet to use with ArcGIS and geocode using the provided ESRI street level geolocator.

But I intend to make point information viewable for the public, thus I want to change sensitive address information to the block level.
EG: 2686 Workman Street, Los Angeles CA 90031 ->
recoded to the block level as 2600-2700 Workman Street, Los Angeles, CA 90031

What can I do to accomplish this other than manually rewriting hundreds of entries?


texts: Cloud and cloud shadow identification with MERIS

texts:
I am looking at MERIS images to use for marine applications (Chl-a. SPM) and would like to confidently remove cloud and cloud shadow pixels from the data.

I have a few ideas:
Using the blue band reflectance. Setting a threshold of 0.215 for Band 3 (490nm) reflectance seems to give a good approximation. See image.

Calculating a general pixel brightness and whiteness is another idea.
It looks like MERIS doesn't have a SWIR band which would have be useful for cloud identification.

Using the Oxygen absorption band. I don't have experience with this and don't have Band 11 in my output that I am working with. However it looks like irradiance of ~15 mW/(m^2*sr*nm) is a good threshold.

As for cloud shadow I'm not quite sure? Cloud top height from oxygen band and then geometric projection of cloud shadow sounds like the best bet.

Has anybody has any experience doing cloud and cloud shadow identification with MERIS?

P.S. There isn't a MERIS tag so I tagged with MODIS, perhaps there are some similar approaches.


texts: Adding analysis privileges to ArcGIS Online application?

texts:
I'd like to use ArcGIS Online analysis services in my web app application.

When I try accessing the service: http://analysis.arcgis.com/arcgis/rest/services/tasks/GPServer?f=pjson&token=

I get this error:

{
  "error" : {
    "code" : 401,
    "message" : "You do not have permissions to use this resource.",
    "details" : [ ]
  }
}

When I check the privileges of my ArcGIS Online app using: https://.maps.arcgis.com/sharing/rest/portals/self?f=pjson&token=

I don't see analysis privilege in appInfo section:

"appInfo": {
    "appId": "-------",
    "itemId": "----------",
    "appOwner": "----",
    "orgId": "------",
    "appTitle": "myApp",
    "privileges": [
      "premium:user:demographics",
      "premium:user:elevation",
      "premium:user:geocode",
      "premium:user:geoenrichment",
      "premium:user:networkanalysis"
    ]
  }

How can I add analysis privileges to my app?


texts:
I have shapefile (the 48MB just below the table here) showing borders of different levels in administrative hierarchies (e.g. borders within communities, between communities, and between federal states). I want to display data of communities now with QGIS. If I do that in a straightforward way by importing an Excel file containing the values of the communities and joining that, it nicely depicts the different levels, but also the borders within the communities. Is there a way to hide these inside-borders?

I could work around that problem in this special case by simply using a different shapefile, but I ran into that problem several times now and want to know whether there is a more elegant solution than using different shape files for different administrative levels.


texts: How to specify a label "weight" in GeoServer SLD to control label priority?

texts:
I'm labelling a cities layer using GeoServer SLD. I need the larger city labels to draw first, then the smaller city labels to draw if there is enough room. There is a pop_band field which I'm using to distinguish between large and small cities.

In ArcMap this would be achieved by assigning a higher weight (or importance) to the larger city label class. Does this concept exist in GeoServer using SLD?

For example, the large cities like Sydney and Brisbane (in grey) should be labelled here before labelling the smaller cities like Toowoomba and Wollongong:

Here is a snippet of the relevant SLD:

<Rule>
  <Name>LargePop</Name>
  <Title>Large cities</Title>
  <ogc:Filter>
    <ogc:PropertyIsLessThan>
      <ogc:PropertyName>pop_band</ogc:PropertyName>
      <ogc:Literal>2</ogc:Literal>
    </ogc:PropertyIsLessThan>
  </ogc:Filter>
  <PointSymbolizer>
    <Graphic>

        <WellKnownName>circle</WellKnownName>
        <Fill>
          <CssParameter name="fill">#727272</CssParameter>
        </Fill>

      <Size>18</Size>
    </Graphic>
  </PointSymbolizer>
  <TextSymbolizer>

   <ogc:PropertyName>town_name</ogc:PropertyName>

   <CssParameter name="font-family">Arial</CssParameter>
   <CssParameter name="font-size">12</CssParameter>
   <CssParameter name="font-style">normal</CssParameter>
   <CssParameter name="font-weight">bold</CssParameter>

 <LabelPlacement>
   <PointPlacement>
     <AnchorPoint>
       <AnchorPointX>0.5</AnchorPointX>
       <AnchorPointY>0.0</AnchorPointY>
     </AnchorPoint>
   </PointPlacement>
 </LabelPlacement>
 <Fill>
   <CssParameter name="fill">#000000</CssParameter>
 </Fill>
</TextSymbolizer>
</Rule>
<Rule>
  <Name>MediumPop</Name>
  <Title>Medium cities</Title>
  <ogc:Filter>
    <ogc:And>
      <ogc:PropertyIsGreaterThanOrEqualTo>
        <ogc:PropertyName>pop_band</ogc:PropertyName>
        <ogc:Literal>3</ogc:Literal>
      </ogc:PropertyIsGreaterThanOrEqualTo>
      <ogc:PropertyIsLessThanOrEqualTo>
        <ogc:PropertyName>pop_band</ogc:PropertyName>
        <ogc:Literal>4</ogc:Literal>
      </ogc:PropertyIsLessThanOrEqualTo>
    </ogc:And>
  </ogc:Filter>
  <PointSymbolizer>
    <Graphic>

        <WellKnownName>circle</WellKnownName>
        <Fill>
          <CssParameter name="fill">#0033CC</CssParameter>
        </Fill>

      <Size>12</Size>
    </Graphic>
  </PointSymbolizer>
  <TextSymbolizer>

   <ogc:PropertyName>town_name</ogc:PropertyName>

   <CssParameter name="font-family">Arial</CssParameter>
   <CssParameter name="font-size">10</CssParameter>
   <CssParameter name="font-style">normal</CssParameter>

 <LabelPlacement>
   <PointPlacement>
     <AnchorPoint>
       <AnchorPointX>0.5</AnchorPointX>
       <AnchorPointY>0.0</AnchorPointY>
     </AnchorPoint>
   </PointPlacement>
 </LabelPlacement>
 <Fill>
   <CssParameter name="fill">#000000</CssParameter>
 </Fill>
</TextSymbolizer>
</Rule>

I'm aware that I can use <MaxScaleDenominator> to control whether the smaller cities draw at all, but this doesn't address the issue here since it doesn't assign a greater importance to the larger city label layer.


texts:
I am trying to generate a map with terrain and elevation in QGIS and as a beginner I have followed tutorials
on how to generate contours and hill shades succesfully. However, I don't know what are the next steps in order for me to
apply style/colors on based on the elevation. As for my project I am using Global Elevation GeoTIFF from https://www.iscgm.org/gmd/.
I am not sure though if this is the right base map/file that I should use to produce what I want, or if this produces accurate elevation
data. The map I need to achieve is a terrain + elevation map just like what SkyVector uses. Also, if I am able to do this, will there be a way to stitch or collage
multiple layers from different sections of a map since processing the entire map takes forever, whilst for a certain country only takes a minute.


texts: Creating Drop Down List of Layers for ModelBuilder Tool?

texts:
In ModelBuilder I require the user select 1 of 5 feature class layers as the initial feature class to run the rest of the model on. The selected layer is to be then fed into a clipping tool as the Input Feature

How do I create a drop-down list of layers in the tool selected from ones I have set the file path for?

I want to have 5 feature class layers selectable as the input for the clipping tool.

ArcMap 10.2


texts: Python OGR: Features of layer inside another layer

texts:
I'm trying to find features from a layer that are within another layer. To do so, I first created a geometry with the geometries from the filtering layer (created by adding geometry or union), and then iterated over the filtered layers.
outLayer = source.CreateLayer(layer.GetName())
for _ in range(layer.GetFeatureCount()):
    feature = layer.GetNextFeature()
    geometry = feature.GetGeometryRef()
    if geometry.Within(geomFilter):
         outLayer.CreateFeature(feature)

The resulting layer has what I want but is very slow : it took 10s per feature with the added geometry, and 4s for the union.
While https://gis.stackexchange.com/a/59264/181819 proposes a solution for points in polygons, it doesn't work for other geometries (multipoints, linestring...).
I was wondering if there was a better way to do it.


texts: Programmatically styling MapInfo Table without MapBasic?

texts:
I'm doing some research into how I can convert a GML file to MapInfo TAB and then style it using the same script.

I have read about a column called MI_STYLE, however my understanding of this is that it is only used when working with spatial databases (I could be wrong here though!).

I cannot use MapBasic for this task, but do know how to use python and OGR. I also understand the 'language' MapInfo uses to style tables (e.g. Pen, brush etc.)

Does anyone happen to have an understanding on how MapInfo records this information and how it can be modified?


texts: How to get coordinates of all polygon drawn on mapbox

texts:
I am new to GIS world and looking to get coordinates of the polygon drawn on mapbox map , so that i can pass the coordinates to my python program. In below program I am trying to capture all co-ordinate in the array "polyCoord" and display first value when calculate button is hit, just to see if data is getting captured

<!DOCTYPE html>

    <meta charset='utf-8' />
    <title></title>
    <meta name='viewport' content='initial-scale=1,maximum-scale=1,user-scalable=no' />

    <link href='https://api.tiles.mapbox.com/mapbox-gl-js/v0.30.0/mapbox-gl.css' rel='stylesheet' />

<link rel='stylesheet' href='https://api.mapbox.com/mapbox-gl-js/plugins/mapbox-gl-draw/v0.16.0/mapbox-gl-draw.css' type='text/css'/>

Draw a polygon using the draw tools.

Calculate area

At the moment, on click of calculate nothing shows up, probably statement " var polyCoord = turf.coordAll(data);" is not valid.


texts: Loop through a list and concatenate a string in Arcpy

texts:
I'm making a program where the user enters some codes (GetParameterAsText) that are going to be used in a SQL expression to select by atributes. I tried to find similar topics but didn't find a solution. This is part of a bigger program, but I'm going to show only the part with the issue.

ingresados = "10703,208,407"
splitIngresados = ingresados.split(",")
try:
    for i in splitIngresados:
        SQLexpresion = "PRM_IDENTIFICA LIKE '" + i + "%' OR"
        SQLexpresion += i
    print SQLexpresion
except:
    print "Se ha detectado un error que impide ejecutar el script"
    print (arcpy.GetMessage(2))

I need to get an expression that looks like this:

"PRM_IDENTIFICA LIKE '10703%' OR PRM_IDENTIFICA LIKE '208%' OR PRM_IDENTIFICA LIKE '407%'"

There is another thing. The user can enter multiple codes, I need a way to count them so I know when to end or close the SQL (So I don't write another "OR" in the end).


texts: How does Pix4D align overlaping images based on one GPS point

texts: Images that can be imported into pix3d have only one GPS position labeled in center of the picuture. I am wondering how does this tool make mosaic (create one big picture with stiching together many overlaping images). Does it do with information gathered from other metadata such as yaw, pitch, roll?


texts:
I am currently recieving this error:
TypeError: unsupported operand type(s) for -: 'float' and 'NoneType' arcgis

for line 22:

invert_out = updateToRow[2]

Could anyone explain why this is happening?

import arcpy as ARCPY
from arcpy import env as ENV
import SSUtilities as UTILS
"""
===============================================================================
SLOPE EQUATION"""

def PipeSlope(input_table, ID_Number , invertIN , invertOUT , ShapeL , slope):
    slopeResults = []
    negativeSlopeList = []
#---------------------------------------
    slopeResults.append(ID_Number)
    slopeResults.append(invertIN)
    slopeResults.append(invertOUT)
    slopeResults.append(ShapeL)
    slopeResults.append(slope)
    negativeSlopeList = []
    print "slopeResults" , slopeResults
    with ARCPY.da.UpdateCursor(input_table, slopeResults) as slopeUpdateCursor:
        for updateToRow in slopeUpdateCursor:
            ID_report  = updateToRow[0]
            invert_in  = float(updateToRow[1])
"""TypeError: unsupported operand type(s) for -: 'float' and 'NoneType' arcgis"""
            invert_out = updateToRow[2]
            shape_length = float(updateToRow[3])
            newSlope = (invert_in - invert_out)/shape_length
            updateToRow[3] = newSlope
            slopeUpdateCursor.updateRow(updateToRow)
            if newSlope < 0:
                negativeSlopeList.append("ID NUMBER:")
                negativeSlopeList.append(ID_report)
                negativeSlopeList.append("SLOPE:")
                negativeSlopeList.append(newSlope)
    return negativeSlopeList

"""
===============================================================================
===============================================================================
INPUT: Script tool arguments (Global Variables)                                """
Update_Table = ARCPY.GetParameterAsText(0)
ID_Number = ARCPY.GetParameterAsText(1)                                         #    1st Parameter: Select Update Table
invertIN = ARCPY.GetParameterAsText(2)                                          #    1st Parameter: Select Update Table
invertOUT = ARCPY.GetParameterAsText(3)                                         #    1st Parameter: Select Update Table
ShapeL = ARCPY.GetParameterAsText(4)                                            #    1st Parameter: Select Update Table
slope = ARCPY.GetParameterAsText(5)                                             #    1st Parameter: Select Update Table

"""
===============================================================================
Environment(s)
==============================================================================="""
desc = ARCPY.Describe(Update_Table)
workSpace = r"" + desc.path
ENV.workspace = workSpace
"""
===============================================================================
FUNCTIONS:                                                                     """
negativeSlopeList = PipeSlope(Update_Table , ID_Number , invertIN , invertOUT , ShapeL , slope)
"""
===============================================================================
OUTPUT TABLE CREATION: Updated .
===============================================================================
# Create Output Text Table
# So the input and output can be reported in the tool report window
"""
header = "Update Structure Attributes in Pipe Table:"

row1  =  "| SUCCESSFUL AUTOMATION OF SLOPE CALCULATION!"
row2  =  "|               "
row3  =  "|List of Incomplete Equations:   "
total =  [row1 , row2 , row3]
for negRow in negativeSlopeList:
    total.append(negRow)

tableOut = UTILS.outputTextTable(total,header=header,pad=1)
ARCPY.AddMessage(tableOut)


texts: QGIS 3.0 - Random Points w/in Polygon Expression

texts:
I have a polygon layer with 3 non-intersecting polygons(i.e., id=1, id=2, id=3). I'd like to randomly distribute 10 points within each of these polygon for 30 total, to be used as survey starting-points.

Qgis prompts to enter an expression that would define this; however, I'm pretty unfamiliar with the language required. If I remember correctly, previously versions allowed you to just specify the number of points per polygon?

The documentation on qgis is pretty scant here, so I may be missing something very obvious.


texts:
I am an amateur trying to produce a DEM that matches a map from the 1960s. I have succeeded in georeferencing the map and aligning it with a DEM, but I do not know how to make the DEM emulate the "cut out" portion of the map image. I call it the "cut out" portion because, being an amateur, I do not know the proper term.
Here is the relevant part of the map:

The actual location of these islands is some ways south of the main island (Kyushu, Japan). Is there a best practice for making DEMs that handle such a non-contiguous scenario?  The only solution I can think of is to manually make the edits in Photoshop, but this feels too ham-fisted. I am using QGIS.


texts:
I would like to create a button using Python Add-ins that allows the user to select a point on the map then immediately bring up a user prompt checkbox that was created in a separate toolbox. I thought what I wrote would work, but when I click on any point on the map, nothing happens. Is there something wrong with the way my script was written?

This is it:

import arcpy
import pythonaddins

arcpy.env.workspace = "E:/GIS/File_Organization/Data/Points/Points.gdb"
arcpy.env.overwriteOutput = True

class SelectWithinARadius(object):
    """Implementation for Button_from_Point_addin.SelectWithinARadius (Button)"""
    def __init__(self):
        self.enabled = True
        self.shape = "NONE"
    def onMouseDownMap(self, x, y, button, shift):
        mxd = arcpy.mapping.MapDocument("CURRENT") #designate map document
        pointGeom = arcpy.PointGeometry(arcpy.Point(x, y), mxd.activeDataFrame.spatialReference) #create point from lat/long
        arcpy.MakeFeatureLayer_management(pointGeom, "point_lyr") #make feature layer from point
        pythonaddins.GPToolDialog(r'E:\GIS\File_Organization\Data\Checkbox_Toolbox\Checkbox.tbx', 'Feature_Selection') #add selection toolbox for user input

    arcpy.RefreshActiveView()

If not, was I not installing it correctly? I can also provide the code for the toolbox that I am trying to link through the GPToolDialog if that would help as well.


texts: Converting text file with coordinates to database table in ArcGIS Desktop with Python?

texts:
How do I convert a text file with point coordinates to database tables in the dBase (*.dbf) format in Python?

I am using ArcGIS Desktop.


texts:
Is torque supposed to be part of Cartodb visualizations by default?  It doesn't show up for me. All I get is simple, choropeth, category, bubble.

This post and video

(http://blog.cartodb.com/post/66687861735/torque-is-live-try-it-on-your-cartodb-maps-today)

make it seem like there's nothing to install or setup so I'm wondering why it's not working for me.

I'm just using a free account, but the website above says Cartodb has released it to all users so I don't think that should matter.


texts: Labels dropping off PDF when exported in QGIS with EasyCustomLabeling plugin?

texts:
I'm having an issue with two particular labels on two separate layers in QGIS. The label layers having being created using the EasyCustomLabeling plugin.

I'm creating an A0 printout with an OS base map and c.100 points with labels(All created with the EasyCustomLabeling plugin) The only solution I have at the moment is to drop the DPI on the export resolution in the composer. However at 1:7500 @ A0 the DPI is now so low that it is causing the base map to appear blurry when printed. Even when the DPI is dropped to 120 one of the two problem labels will still not appear.

I have ensured that Show all Labels for this layer (including colliding labels) is checked as well as setting the weighting on the Discourage labels from covering features to High Weight for all the layers with labels.

All of my labels are horizontal and none are curved.


texts: rtree intersection returns only zeros when filled with via stream from generator function

texts:
I have trouble getting anything useful out of a rtree (Rtree==0.8.2) intersection. Instead of returning the indices/IDs that I passed into the index on creation, I get zeros.

I have some data:

for element in data:
    print(element)

POINT (497582.829135131 5419474.690534886)
POINT (565856.4271351419 5988338.580648975)
POINT (539734.212266123 5614400.468566706)
POINT (683523.3673841552 5398102.521964455)
POINT (633258.6526692112 5405182.75073655)
POINT (353675.3945317736 5747638.443753068)
POINT (480313.0993177936 5324783.150275591)
POINT (476802.8443284882 5390481.26880122)
POINT (392880.6204936574 5829710.142526797)
POINT (483370.4752331344 5792650.901288186)

I have an area of interest that spans over some of those points:

aoi = data[4].buffer(100000)
# GeometryCollection((aoi, box(*aoi.bounds), MultiPoint(data)))  # image below

I create an index based on rtree's stream load example and the documentation on the __init__ method:

def generator_function(geometries):
    for i, geometry in enumerate(geometries):
    # None instead of objects as I don't want to store the objects in the index
        yield (i, geometry.bounds, None)
idx = index.Index(generator_function(data))

When I query the index with the aoi's bounds, I get two hits but their IDs are both "0".

candidates = list(idx.intersection(aoi.bounds))
print("We got {} results".format(len(candidates)))
print("Their indexes are: {}".format(candidates))

We got 2 results
Their indexes are: [0, 0]

My problem: I expected to get the IDs as I passed them from the generator function. Why do I not get those?

If I fill the index step by step instead, I get the result I expected. I see no real logical difference in what should end up in the index though.

idx = index.Index()
for i, element in enumerate(data):
    idx.insert(i, element.bounds)

candidates = list(idx.intersection(aoi.bounds))
print("We got {} results".format(len(candidates)))
print("Their indexes are: {}".format(candidates))

We got 2 results
Their indexes are: [3, 4]

# GeometryCollection((box(*aoi.bounds), MultiPoint([data[i] for i in candidates])))  # Image below


texts: Visualizing 3D focal mechanism (beach ball) in QGIS

texts:
Currently, the qgis2threejs plug-in in QGIS allows users to use a 3d model as symbology for point vector data. I made a 3d focal mechanism template (beach ball) and used it in qgis2threejs:

You can rotate the beach ball around 3 intrinsic rotational axes (meaning the axes move after each rotation; they are fixed with the beach ball), with X (EW axis) as Pitch, Y (NS axis) as Roll, and Z (Up-Down axis) as Yaw. These rotations are applied in X-Y-Z order

Do anybody here have an idea on how to convert focal mechanism data (e.g. dip, strike of fault plane and slip vector; or orientations of the P,T and N axes) into the the correct X, Y, and Z angles (I've done some reading and these angles are referred to as Tait-Bryan rotation angles)? If we can do this, then we can finally plot focal mechanisms in QGIS.


texts: Adding custom processing provider into QGIS

texts:
I am trying to add a new processing algorithm provider into QGIS 2.18. In the picture below there are GDAL/OGR, GRASS GIS, etc. I want to add a provider called My Provider.

So I read the doc here: https://docs.qgis.org/2.18/en/docs/pyqgis_developer_cookbook/processing.html
and saw:

So I followed:

Then I reopened QGIS and installed the plugin, but no new provider appears. I read this thread: How to see the 'example provider' in processing toolbox of QGIS?. But this looks like this is for version 2.8, my exampleprovider is located in QGIS_2_18\apps\qgis\python\plugins\processing\algs\exampleprovider and I dont see which import is pointing to the wrong file.

What should I modify to make my provider appear in the processing toolbox?


texts: Run model from python .. module object not callable

texts:
Having issues to run a model from python, below the actual code , always return
module object not callable.
ArcMap 10.1 using arcPy esri python 64-bit

Any ideas why?

    import arcpy

        try:
            arcpy.ImportToolbox("c:/test/geocode.tbx")
            arcpy.Model_41(<param>,<param>)
        except Exception as e:
            print e.message

error >> 'module' object is not callable

Removing try/except as suggested.

    import arcpy

    arcpy.ImportToolbox("c:/test/geocode.tbx")
    arcpy.Model_41("<param>,<param>")

error >>TypeError: 'module' object is not callable