texts:
I want to intersect two layers (shapefiles) with QGIS 2.18.
For both layers, I tested validity with

QGIS Vector / Vector Tools / Check validity and
GRASS v.build.check (with default parameters)

both tests report no errors for both layers, but QGIS Intersection reports invalid geometry for one or more input features.

The two layers may have some boundaries, that should be equal, but may actually have slight differences (unprecise manual digitalisation of one layer) possibly resulting in several very small slices.

What do you recommend me to try?


texts: How to remove or disable a fundamental toolbar from ArcGIS Desktop?

texts:
I don't want hide a toolbar, rather I need to prevent a user from using some fundamental toolbars (not custom ones).

I have been reading anything about a file '.mxd' or '.mdx', but it's not clear to me yet.

I use ArcObjects SDK.


texts: How to set interval option for timeslider with same min and maxdate using Ol-ext:Timeline control?

texts:
I am trying to create a timeslider with Ol-extension timeline similiar to https://viglino.github.io/ol-ext/examples/storymap/map.control.timequake.html. All the sample codes have different min and maxDate with graduation option as day. The following code works

var tline = new ol.control.Timeline({
    graduation:'day',
    minDate: new Date('2020/5/1'),
    maxDate: new Date('2020/5/10')
  });
  map.addControl (tline);

However I need to create a timeslider with a 15 minute interval on a single day. I used the interval property  but nothing works. How to fix this?

 var tline = new ol.control.Timeline({
    interval:'15mn',
    minDate: new Date('2020/5/1'),
    maxDate: new Date('2020/5/1')
  });
  map.addControl (tline);


texts: How to handle multiple space separator in Grass?

texts:
I try to load a Lidar point cloud in .xyz format to Grass (in QGIS 2.0.1).

I tried to input the following options to the 'separator' field: tab, 'tab', '\t'. Just for the unlikely case that it worked I tried alse ' ', 'space', 'newline', 'comma' and ';'.

Two things have changed. First: empty space between two apostrophes didn't work yesterday, but it works now (I don't know why). Phrasing it as 'space' still doesn't work.

Second:  I found the data are actually separated by two or three spaces, not by a tab. When I set space as the separator, there is not the old wrong separator error, but Y coordinate seems wrong: it looks for it between the two spaces separating the fields. Writing '  ' (double space) to the separator field doesn't help. So I rephrased this question to fit my problem better: how to handle multiple space separator in Grass?

First few lines of the .xyz file look like this:

  3453939.450  5633900.030   361.670
  3453941.750  5633900.020   361.720

Other options are the default ones, just the "statistic" is "median" instead of "n".


texts: Configuring wps service in the lizmap web client

texts:
I can run lizmap web client after configuring the lizmap web client according to Updated - Installing Apache, QGIS Server, And Lizmap Web Client On Windows OS.  There is a way (Lizmap WPS Web Client) to configure it. But I don't know how to configure wps in the lizmap client. Can you give a specific example: including wps_url, wps_rootDirectories, redis_port and ows_url. I am a beginner about how to configure wms, wps service. And what is role of the qgis-wps4server?

[wps]
wps_url=http://wps:8080
wps_rootDirectories=/projects/wps
redis_host=localhost
redis_port=6379
redis_key_prefix=lzmwps
ows_url=http://map:8080


texts: I have a shapefile polygon that shows the extent of a state. I need to highlight only the boundary of the state. Can I do that in QGIS?


texts: Maintaining Clip Layer Attributes When Using Clip Tool in ArcGIS Desktop

texts:
I am looking to clip master map polygons to land owner boundaries. However, the normal clip tool does not allow me to identify which land owner the master map parcel belongs to once it has been clipped. The master map polygons only maintain their own attributes. I have tried a spatial join using 'completely within' operation, but this misses out some polygons.
The 'intersect' operation picks up several of the clipping boundaries where there are two adjacent to each other.
Is there another way to identify which landowner polygon each mastermap poly belongs to?


texts: Saving output layer using Python Toolbox?

texts:
Why when I try to save my layer in destination gdb (I pass path) does it save in default.gdb?

def getParameterInfo(self):
    """Define parameter definitions"""

    param0 = arcpy.Parameter(
        name="Input Features",
        displayName="Area",
        direction="Input",
        datatype="DEFeatureClass",
        parameterType="Required")

    param1 = arcpy.Parameter(
        name="Output Features",
        displayName="Output Features Class",
        direction="Output",
        datatype="DEFeatureClass",
        parameterType="Required")

    param2 = arcpy.Parameter(
        name="Input Number0",
        displayName="Cell Size Width",
        direction="Input",
        datatype="Double",
        parameterType="Required")

    param3 = arcpy.Parameter(
        name="Input Number1",
        displayName="Cell Size Height",
        direction="Input",
        datatype="Double",
        parameterType="Required")

    param4 = arcpy.Parameter(
        name="O",
        displayName="O",
        direction="Input",
        datatype="DEFeatureClass",
        parameterType="Required")

    parameters = [param0, param1, param2, param3, param4]

    return parameters

def isLicensed(self):
    """Set whether tool is licensed to execute."""
    return True

def updateParameters(self, parameters):
    """Modify the values and properties of parameters before internal
    validation is performed.  This method is called whenever a parameter
    has been changed."""
    return

def updateMessages(self, parameters):
    """Modify the messages created by internal validation for each tool
    parameter.  This method is called after internal validation."""
    return

def execute(self, parameters, messages):
    """The source code of the tool."""

    reload(sys)
    sys.setdefaultencoding('utf8')

    # teritory
    inputFeatureClass = parameters[0].valueAsText

    # my output
    outputFeatureClass = parameters[1].valueAsText
    outPath = os.path.basename(outputFeatureClass)

    arcpy.AddMessage(outputFeatureClass)

    cellWidth = parameters[2].value

    cellHeight = parameters[3].value

    o = parameters[4].valueAsText

    # my workspace require to fishnet
    # arcpy.env.workspace = outputFeatureClass

    # when I will be save file, each of them will be overwriting existing file about same name
    arcpy.env.overwriteOutput = True

    origin, yAxis, opposite = coordinates(inputFeatureClass)

    # create fishnet
    fishnet = arcpy.CreateFishnet_management("in_memory/Fishnet", origin, yAxis, cellWidth, cellHeight, '0', '0',
                                             opposite, 'LABELS', inputFeatureClass, 'POLYGON')
    # add field for id
    arcpy.AddField_management(fishnet, 'idField', 'LONG')

    # add id sequence
    with arcpy.da.UpdateCursor(fishnet, ['idField']) as cursor:
        i = 1
        for row in cursor:
            row[0] = i
            cursor.updateRow(row)
            i += 1
    del cursor

    # clip a final fishnet
    finalFishnet = arcpy.Clip_analysis(fishnet, inputFeatureClass,
                                       outPath)

    # create dict for results
    score = {}

    # adding every key value equal 0
    with arcpy.da.SearchCursor(finalFishnet, ['idField']) as cursor:
        for row in cursor:
            score[row[0]] = 0
    del cursor

    oFunction(o, finalFishnet, score)

    addResultFunction(finalFishnet, score)

    return

After execute


texts: Get counts in grid cell for multiple values

texts:
I have 500+ polygon shapefiles with nationwide coverage.  Each shapefile has possible record values of 0-4.  I want to get the value that intersects each grid centroid and then get a total count for all the shapefiles.

So my final output will be a 1 km grid with value fields 0,1,2,3,4 and counts for each value like below.

I'm limited to ArcGIS for software.

I know I can do this by leaving everything in vector form and creating a fishnet with label points, doing spatial joins and summary statistics, but its a bit unwieldy.  I'm not as familiar with the raster tools available with ArcGIS and can't find anything that fits my needs.  If I convert all the 500+ shapefiles to raster is there any easier way to do this?


texts: How can I avoid self-managing z-index on leaflet popups

texts:
UPDATED INLINE

Using Leaflet 0.7 (but I would like to move to 1.0 Beta 2) - NOTE, I do not have this problem in 1.0, beta 2.

I have a leaflet map with two tile layers (mapbox, and my own) one marker layer, and a GeoJSON layer. What I was finding was that since the GeoJSON layer is dynamically loaded (from a file contents) and that that happens after I had the marker layer, I the marker was getting hidden by z-index.

Ok, easy fix:

   theCar = L.marker(mCityCenter, {
         draggable: true,
         rotationOrigin: "32px 32px",
         rotationAngle: 90,
         zIndexOffset: 10000,
         icon: myIcon
     });
     carLayer = theCar;

Now, the reason I think I have to do this is that the GeoJSON that I load has lots of marker (GEOJSON points) which I have custom makers for. There can be 4K+ Markers. I think this is filling up the marker pane, and that why zIndexOffset works.

[Don't worry, the map is not that busy, the GeoJSON layer is actually the extension called Leaflet.Timeline, which only shows those markers based on the timline]. But if they they are not physically present in the GeoJSON layer, the timeline will not work]

But that has force me to also self manage the z-index of the pop-up associated with it. Since it was much lower than 10,000.

     var popup = L.popup();
     popup.setContent("CAR");
     popup.setLatLng(mCityCenter);
     carLayer.bindPopup(popup);
     popup.addTo(theMap);
     self.CarPopup = popup;
     $(self.CarPopup._container).css('z-index', 20000);
     theMap.closePopup(popup);

This is getting really ugly. Obviously, I am doing something wrong here. While it works, I really hate working against the paradigm.

Question: If hypothetically, I had 5 GeoJSON layers in Leaflet, would the vector data all sit in overlay pane and the markers sit in the marker pane. What happens if I have 4K markers in each pane, how would I keep the "top" level markers from interfering with the lower level ones?


texts: Finding largest visible area inside polygon using PyQGIS

texts:
I need to find the largest possible visible area inside many polygons. The polygons represent clearcuts inside forests.
For example if I stand here I can see 10 out of 14 hectares, nowhere else can I stand to see more.
The workflow I can think of now is:

Create a grid in each polygon

For each grid point:

Draw line of sight-lines from each point in the grid to all vertices
Intersect them with the polygon
Select by location the lines intersecting current grid point
Buffer - dissolve calculate area

This is going to take a long time.
Is there a tool for this?

To create sight lines:
polylyr = QgsProject.instance().mapLayersByName('polygon')[0]
g = [f.geometry() for f in polylyr.getFeatures()][0]
verts = [v for v in g.vertices()]

pointlyr = QgsProject.instance().mapLayersByName('Grid')[0]
point = QgsPoint([f.geometry() for f in pointlyr.getSelectedFeatures()][0].asPoint())

linelist = []
for vert in verts:
    Line = QgsGeometry.fromPolyline([vert, point])
    linelist.append(Line)

vl = QgsVectorLayer("LineString?crs={}&index=yes".format(pointlyr.crs().authid()), "myLayer", "memory")
provider = vl.dataProvider()

for L in linelist:
    f = QgsFeature()
    f.setGeometry(L)
    provider.addFeature(f)

QgsProject.instance().addMapLayer(vl)


texts:
I'm getting this error code with QGIS when I am using the tool Join attribute by location: I'm trying to join the attribute table of polygons into my vector points attribute table.

Processing algorithm…
    AlgorithmDialog
    Input parameters:
    { 'INPUT' : 'C:\Users\renau\Desktop\Geocoding_test_batch_Jakarta\3. spatial query on territory\Shape_FILE_UTM8_transformation\batch.shp', 'JOIN' : 'C:\Users\renau\Desktop\Geocoding_test_batch_Jakarta\3. spatial query on territory\Shape_FILE_UTM8_transformation\kelurahan_21148.shp', 'PREDICATE' : [0], 'JOIN_FIELDS' : [], 'METHOD' : 0, 'DISCARD_NONMATCHING' : False, 'OUTPUT' : 'memory:' }

Feature (61) has invalid geometry. Please fix the geometry or change the Processing setting to the "Ignore invalid input features" option.
Execution failed after 0.48 seconds

Loading resulting layers
Algorithm 'Join attributes by location' finished

How can i fix this error?

I can get also this error:

RuntimeError: wrapped C/C++ object of type AlgorithmDialog has been deleted
Traceback (most recent call last):
  File "C:/PROGRA~1/QGIS3~1.0/apps/qgis/./python/plugins\processing\gui\ProcessingToolbox.py", line 344, in executeAlgorithm
    if dlg.wasExecuted():
RuntimeError: wrapped C/C++ object of type AlgorithmDialog has been deleted

Python version: 3.6.0 (v3.6.0:41df79263a11, Dec 23 2016, 08:06:12) [MSC v.1900 64 bit (AMD64)]
QGIS version: 3.0.1-Girona Girona, a86bec25eb

Python Path:
C:/PROGRA~1/QGIS3~1.0/apps/qgis/./python
C:/Users/renau/AppData/Roaming/QGIS/QGIS3\profiles\default/python
C:/Users/renau/AppData/Roaming/QGIS/QGIS3\profiles\default/python/plugins
C:/PROGRA~1/QGIS3~1.0/apps/qgis/./python/plugins
C:\Program Files\QGIS 3.0\bin\python36.zip
C:\PROGRA~1\QGIS3~1.0\apps\Python36\DLLs
C:\PROGRA~1\QGIS3~1.0\apps\Python36\lib
C:\Program Files\QGIS 3.0\bin
C:\PROGRA~1\QGIS3~1.0\apps\Python36
C:\PROGRA~1\QGIS3~1.0\apps\Python36\lib\site-packages
C:/Users/renau/AppData/Roaming/QGIS/QGIS3\profiles\default/python
C:\Users\renau\AppData\Roaming\QGIS\QGIS3\profiles\default\python\plugins\mmqgis/forms


texts:
Its been years since I last had to code, so please excuse the following.

The idea here is that a user inputs their work and it checks for missing expected data. (I've taken out other code that works but is for other checks)

Simplified summary of what I'd like it to do:

Identify which data set (MASTER_ID) is being added.

Select all records from the masterTable that have that MASTER_ID.

Compare these selected records against userTable by Join based on KEY

If there is a record found in masterTable and not found in userTable, (NULL record found in joined userTable) output Key from masterTable for that row as an error.

-Key =unique string eg. 12345678-A123

-MASTER_ID =same for entire UserSubmitTable eg. 12345678

land = arcpy.GetParameterAsText(0)
area = arcpy.GetParameterAsText(1)

#Join Land & Area tables#
landList = "landList"
arcpy.JoinField_management(land, 'KEY', area, 'KEY', ['ID_land', 'land_NO','MASTER_ID'])
arcpy.Frequency_analysis(land, landList, ['KEY', 'land_NO','MASTER_ID'])

#get the active MASTER_ID as string #
for row in arcpy.da.SearchCursor(area, ['MASTER_ID']):
    if row[0] is not None:
        Active_ID = row[0]
    break

#query MASTER_ID #
queryList = "query_table"
arcpy.MakeQueryTable_management(landList, queryList, "NO_KEY_FIELD", "MASTER_ID = '{}'".format(Active_ID))

#search for errors #
for row in arcpy.da.SearchCursor(queryList, ['land_NO', 'MASTER_ID', 'KEY']):
    if row[0] is None and row[1] is not None:#
        arcpy.AddError(row[2] + " exist in reference but is missing")

I get:

  \management.py", line 7441, in MakeQueryTable raise e RuntimeError:
  Object: Error in executing tool


texts: Exporting KML from QGIS with relation (one-to-many) table as popup

texts:
Here is a challenging one for someone.

I have two tables in QGIS.

Holes (parent) - Which contains the borehole locations.
Strata (child) - which has many rows that contain the strata for each borehole

Although i can see this join in QGIS it would be really useful to be able to export a KML that has each point and a popup that has a table of the strata at each point...

Any idea?


texts: Changing Move Feature mechanism in QGIS?

texts:
I am currently using Qgis 3.2.1 and struggle with the move feature mechanism there. Its works like a select and drop function and it´s pretty akward to use. See the picture below: (moved features will be inserted on red dot position) Using the Snapping Tool doesn't work that well with that kind of movement either.

In QGIS 2.18.5 the move feature mechanism was much better. It was without this "dropping". While moving a feature it was following the mouse cursor.
At least for me serveral times faster to work with.

Is it possible to change the move feature mechanism in QGIS 3.2.1 so that it works again like in older versions ?


texts: Interactively selecting polygon in ArcMap and passing selected feature to ArcPy?

texts:
I currently have a Python 2.7 script that takes a hard-coded select by attribute from ArcMap 10.3.

arcpy.SelectLayerByAttribute_management("SoilLandscapePolygons_Layer", "NEW_SELECTION", "POLY_ID = 25984)

that then copies this polygon and its related tables into a new FGDB. What I want to do, is to have the user select a polygon and then have this polygon copied ... so that I can continue working on it in my script.  I'm not doing any editing what so ever but I need this one polygon to begin working through the various relationship tables to create a summary for the selected polygon.


texts: Calculating annual mean temperature from ERA5 data in Earth Engine

texts:
I am calculating bioclimatic variables using temperature data from ERA5 in Google Earth Engine. ERA5 has 488 images, one for mean monthly temperature for all months 1979-2019.

For calculating annual mean temperature (BIO1), do I take the mean across all 488 images as simple as:

var era5_meantemp = era5_complete.select(['mean_2m_air_temperature']);
var era5_meantemp_annual= era5_meantemp.mean(); //unit= kelvin

...or am I supposed to take the mean of the temp of all months in a year i.e., one annual mean temp for 1979, for 1980.... and have an image collection of annual mean temperature for all years i.e. 40 images (2019-1979)? If yes, how do I further calculate one value across all 40 years?


texts: mosaic dem with multiple coordinate systems

texts:
I'm struggling to mosaic multiple DEMs (~40) into one large raster. The DEMs are 30m NED rasters from the NRCS data gateway, for the New England region. They are in two different coordinate systems; NAD 1983 UTM Zone 18N and 19N.

I've read several posts debating whether or not DEMs should be projected or unprojected before running slope and hillshade tools. In the end, it sounds like running slope and hillshade on an unprojected mosaic of all DEMs would be most appropriate. That being said, I can't figure out how to create the mosaic dataset with DEMs from 2 different coordinate systems.

I have only a VERY basic understanding of ArcGIS (using ArcGIS 10.2).


texts:
I cannot connect to postgresql 9.0.5 with ArcSDE 10.1, using ArcCatalog 10.1 on Win 7 professional 64 bit.

I followed the steps described here, and installed the Pg-Agent from Stack builder. When I execute "Create Enterprise Geodatabase", I set "instance" to "localhost:5432", but I got this error:

Cannot connect to database because the database client software failed to load.  Be sure the database client software is installed and configured correctly.
Failed to execute (CreateEnterpriseGeodatabase).

I've made sure that the service is running and restarted the whole machine after these installations, but it still gives me the error.

Any help, please?


texts: Anyone have success with OpenLayers and ArcGISCache?

texts:
OPENLAYERS

Has anyone had success using OpenLayers and the ArcGISCache layer beyond the examples?  I'd like to mix tiles with different resolutions and although they seem to scale correctly they do not translate correctly.  Code to translate the tiles appears to be missing with the exception of a % translation necessary after the scale:

transformDiv: function(scale) {

    // scale the layer div

    this.div.style.width = 100 * scale + '%';
    this.div.style.height = 100 * scale + '%';

    // and translate the layer div as necessary

    var size = this.map.getSize();
    var lcX = parseInt(this.map.layerContainerDiv.style.left, 10);
    var lcY = parseInt(this.map.layerContainerDiv.style.top, 10);
    var x = (lcX - (size.w / 2.0)) * (scale - 1);
    var y = (lcY - (size.h / 2.0)) * (scale - 1);

    this.div.style.left = x + '%';
    this.div.style.top = y + '%';
},

ESRI

It appears that ArcGIS JSAPI has an even greater limitation...all the tiles must have exactly the same resolution.  Can someone from ESRI discuss/confirm this limitation?


texts: How does the offset values of LAS files relate to their coordinate system?

texts:
I have a conundrum concerning LiDAR data, the stored point offset, and the coordinate system.

I process large amounts of LAS data in TerraScan (TS), which sits on top of Microstation (MS). I receive LAS files that are created from raw data by CloudPro. Traditionally, when setting up a project, before importing points, I set the coordinate system up in TS, including any false easting/northing values.  I typically set these values to that of whatever coordinate system I happen to be working in (e.g. For NAD83 Virginia South State Plane US_Foot, false_easting: 11482916.667, false_northing: 3280833.333).

I just started working on a dataset in NAD83 UTM 15N meters (false easting: 500000). I set the TS false easting, and began setting up the project.  Unfortunately, after some investigation, the LAS files have an x,y,z offset (230247.299, 3550963.700, 0) that will not allow them to fit inside the bounds defined by the 500000m false easting. I understand this is due to MS's coordinate range of +/-2147484 of the defined coordinate origin x,y,z.  I also understand that an LAS file's coordinates take the scale and offset into consideration. From the ASPRS specification:

  X, Y, and Z offset:  The offset fields should be used to set the
  overall offset for the point records.  In general these numbers will
  be zero, but for certain cases the resolution of the point data may
  not be large enough for a given projection system.  However, it should
  always be assumed that these numbers are used.  So to scale a given X
  from the point record, take the point record X multiplied by the X
  scale factor, and then add the X offset.

This got me thinking, how (if any) does the offset of an LAS file relate to the coordinate system that it is projected to?


texts: Is there a specific disadvantage to storing coordinates in geographic coordinate systems?

texts:
If a coordinate is stored in a geographic coordinate system, it's quite easy to calculate the following things exactly:

Distances (geodesic)
Point-in-polygon
Areas
etc.

If you store the coordinates in a projected coordinate system, these things are calculated in a Euclidean way, which means there will be huge distortions if your map covers a larger area of the earth.

So, in my view, storing the coordinates in a geographic coordinate system has the

advantage that all calculations can be made exactly, and the
disadvantage that they must be projected before they can be rendered in a map.

Since the latter is something that can be done very fast on modern processors, I don't think it would ever be a huge disadvantage.

Are there any other disadvantages? If not, then why does one ever store geographic data in projected systems?

One answer to the second question could be that the projecting of coordinates was not always so cheap. On processors of the 70s, it could be better to store the coordinates in the way they are rendered, to avoid having to project them every time the map is rendered.


texts: QGIS Model Designer; save layer in geopackage with date

texts:
I have used the model designer to make a model in QGIS.
The model takes data different shapefiles, and runs a few algorithms.
The result is different layers loaded in an new project.
I manually store these layers in a geopackage.
Example of tablename:
week-37_high-detail-area
I would like to add a step to the model, so it stores the output as a new table in an existing geopackage. I would like the model to store the layer to the geopackage with the number of the week in the tablename of the geopackage.
Is it possible to add some sort of expression section in the tablename?


texts:
I am creating a QGIS plugin and want to list the layers loaded in a Combobox when writing self.dlg the error message object has no attribute appears. I am new in this domain. The code below is for the function.

from qgis.PyQt.QtCore import QSettings, QTranslator, QCoreApplication, Qt
from qgis.PyQt.QtGui import QIcon
from qgis.PyQt.QtWidgets import QAction
from .resources import *

from .mapping_utilities_dockwidget import MappingUtilitiesDockWidget
import os.path

  def loadLayers(self):
        self.dlg.cb_inputlayers.clear()
        layers = [layer for layer in  QgsProject.instance().mapLayers().values()]
        network_component=[]
        for layer in layers:
            if layer.type() ==QgsMapLayer.VectorLayer:
                network_component.append(layer.name())
        self.dlg.cb_inputlayers.addItems(network_component)


texts: Why does a LINESTRING change in length after union with a shorter but completely-overlapping LINESTRING?

texts:
I have two linestrings that I am trying to merge and determine their total length:

This is a contrived example: the green linestring B is fully-contained within the red linestring A (but this will not always be the case), and the length of the merged geometry should be the same as the length of linestring A.

In PostGIS, if I calculate the length of the red linestring A with

  select st_length(st_geomfromtext('(A geom)'), false)

I get a result of 0.6102 (arbitrary units).

If I then calculate the length of the union of the red linestring A and green linestring B, via:

  select st_length(

    st_union(

      st_geomfromtext('(A geom)'),

      st_geomfromtext('(B geom)')

    ), false)

I then get a result of 0.6328, even though the green linestring B completely overlaps with the red linestring A.

I've tried many variations of merging these two linestrings involving st_unaryunion, st_linemerge, st_collect, st_dump, st_removerepeatedpoints and no luck.

How can I reconcile this? There are situations where the linestrings to be merged will not overlap 100% and I'll need to measure the length of the total union geometry, but I have trouble trusting that number when it seems to fail in this scenario.

For what it's worth, linestring A is:

  MULTILINESTRING((-84.3381536299124 33.7648637741265,-84.337082 33.76486,-84.3370269 33.7648551,-84.3367653 33.7648604,-84.336237 33.764861,-84.3358725 33.7648538,-84.335818 33.764857,-84.335192 33.764848,-84.3350945 33.7648474,-84.335009 33.764844,-84.334164 33.764846,-84.333928 33.764846,-84.3337994 33.7648477,-84.332997 33.764858,-84.332176 33.76486,-84.331834 33.764852,-84.330605 33.764845,-84.330119 33.764848,-84.329495 33.76486,-84.328663 33.764858,-84.327896 33.764867,-84.3275311259472 33.764859190368))

linestring B is:

  MULTILINESTRING((-84.3372103956974 33.7648604521912,-84.337082 33.76486,-84.3370269 33.7648551,-84.3367653 33.7648604,-84.3365005784092 33.7648607006492))

EDIT: an added wrinkle:
The length of linestring B is 0.04079. If I measure the intersection between both linestrings:

  select st_length(st_intersection(st_geomfromtext('(geom A')),st_geomfromtext('(geom B)'), false)

I get a result of 0.01821, which implies that there is (0.04079 - 0.01821) = 0.02258 of its length that doesn't overlap, which is the difference measured between st_length(A) and st_length(union(A,B)).

EDIT 2: this is getting weirder! I ran ST_SymDifference to find the area where these lines didn't intersect/overlap. It looks like this:

So the point where the line doesn't cross the green dots is where both linestrings overlap, and that coincides to all of the points in linestring B that have identical matches in linestring A, which is every point except the first and last points of B.

So if I take the first point of linestring B (we'll call it P): POINT(-84.3372103956974 33.7648604521912) and see if it intersects with linestring A:

  select st_intersection(st_geomfromtext('(geom A)'), st_geomfromtext('POINT(-84.3372103956974 33.7648604521912)'))

I get "GEOMETRYCOLLECTION EMPTY" as a result.

If I then try to find the point on A that is closest to P, it is identical to P:

  select st_astext(st_closestpoint(st_geomfromtext('(geom A)'), st_geomfromtext('POINT(-84.3372103956974 33.7648604521912)')))

The result is "POINT(-84.3372103956974 33.7648604521912)"!

On top of that, this also fails:

  select st_intersection(st_geomfromtext('(geom A)'), st_closestpoint(st_geomfromtext('(geom A)'), st_geomfromtext('POINT(-84.3372103956974 33.7648604521912)')))

The result is "GEOMETRYCOLLECTION EMPTY", which seems quite broken since I'm requesting the intersection of a linestring with a point explicitly generated as being on that linestring. This query fails with the last point of B as well, but works with points 2, 3, and 4. The order of arguments (geom first, point second or vice versa) does not affect anything.


texts: GDAL Merge Py on ubuntu many tif tiles into one

texts:
Trying to batch mosaic 1000 Gtiff rasters (all in one folder) into 1 tif file via gdal merge on linux ubuntu.

Having trouble with the linux command to bring in all files in a folder into the output tiff.

Looking for sample code for use in gdal on ubuntu in the terminal.

For Ubuntu my best guess is

for f in *.tiff; do gdal_merge -of GTiff "$f"; done


texts: While Exporting from SDE, null geometry features attack?

texts:
I'm working on a layer in an SDE Database (ArcMap 10.1 SP1 for desktop (Build 3143)), I just export to a local working database on my workstation the layer I want to work on and load it back when done.

The problem is that when I export it seems that the records are exported fine (as I checked the tables) but on the map most of the points didn't appear.

The same layer appears normally on the ArcCatalog, but when exported the error happens.

The check geometry results "null geometry”, the Repair Geometry" didn't solve it.
Copying the Features to a new layer, result in the copying of the Features is that it copies the features with geometry and deletes the "null geometry".

Any suggestions?
Thanks


texts: Add postgis temporary table using pyqgis into canvas failed

texts:
I've tried to add my temporary table into QGIS canvas using pyqgis but it's failed. Here the code(view_tatanga is a temporary table that I've created using sql on pgadmin) :

uri = QgsDataSourceURI()
sql = "SELECT * FROM gis.view_tatanga"
uri.setConnection("localhost", "5433", "db_pbb", "postgres", "septin")
uri.setDataSource('', u'(%s\n)' % sql, 'geom', '', 'nop')
#uri.setDataSource("",sql,"geom","","nop")
tmadminkec = QgsVectorLayer(uri.uri(), 'view_tatanga', "postgres")
QgsMapLayerRegistry.instance().addMapLayer(tmadminkec)

I didn't get any error when I try to load a usual table (non temporary table). Any ideas ?


texts: How to keep ArcMap Toolbar in a nice order?

texts: I`m having an issue with my ArcMap. Every time I open it the tool bar is really messy, so I always have to start re-organising the tool bar first. Does anyone have an idea of how can I keep it in order without rearranging it every single time I use it?


texts: Where does QGIS store the list of server connections (WFS)?

texts: Over time I have collected a bit of a list of WFS server addresses. I would like to include the list of server connections in my overall backup. Therefore, I like to know where QGIS stores the list.


texts: How to cut/split multipolygon with blade (linestring) into two pieces using PostGIS?

texts:
I struggle to do split operation in postgis I tried to build a query but Im always fail.

Id like to cut a multipoligon with a blade (a linestring) into two pieces. These two objects have specific SRID of my country.

query:

CREATE TABLE Cut_result AS
SELECT ST_SPLIT(poligon,blade)
FROM (SELECT
ST_MakeLine(700264 126376, 780844 122481)::geometry(LINESTRING,23700) AS blade),
SELECT FROM county AS poligon;

ERROR

ERROR:  syntax error at or near "126376"
LINE 5: ST_MakeLine(700264 126376, 780844 122481)::geometry(LINESTRI...
                           ^
********** Error **********


texts: Intersect elevation interpolated from CAD file in ArcGIS

texts:
I have two datasets - contours and vector (water, roads,...)

Vector layers were used as breaklines to generate the contours. Then the contours were generalized a bit.

For quality check I need to find the elevation difference between those two datasets at their intersections.

If I use Intersect tool in ArcGIS and output it with both file attributes to a point file, I get the contour elevation nicely, since it is the same over all the element vertice Z values. But for vector data of roads or rivers every vertex has its own Z value and intersect only gives me the value of elements lowest Z value (if I have a line with two points and one is 3m and the other 4m and contour intersects exactly in the middle, the intersect elevation value should be 3.5m, but instead 3 is written to point file).

I thought of creating TIN from vector data set and adding value from there to the intersect points, but as soon as there is another vector object with different elevations close to intersect point that does not intersect with contour but crosses the other vector between its two vertices, it forces the TIN to follow it and with that the intersect value is not correct anymore.

Is there a tool that allows me to get the interpolated Z value between two vertices? It doesn't even have to be ESRI software.
I am not however very good at code.


texts: Problem in creating diagram with ArcGIS Diagrammer 10.2

texts: I want to create diagram using ArcGIS Diagrammer .I use esri transportation Geodatabase as sample. when i import XML file of the sample,the links of the relationships do not show in ArcGIS Diagrammer.How can i show the relationships?


texts: Use PostGIS functions from an other schema (than public)

texts:
Context

I have to set up a postgresql database with several schemas.

To separate things, I decided to install postgis elsewhere than the public schema, based on https://postgis.net/2017/11/07/tip-move-postgis-schema/ and Using schema other than public in PostGIS? for example.

I also installed other extensions, such as pgrouting or fuzzystrmatch or citext for example, which by default went into the main of my project schemas called main. And also postgis_topology which is set up into its own schema called "topology", if I'm not wrong.

For my own culture ; is it a good practice to centralise all extensions in one unique schema?

Error

Now, I want to use postgis functions in all my project schemas, but unfortunately, the postgis functions don't seem to be available as I went into this error while creating a database:

psql:create_db.sql:104: ERROR:  type "geometry" does not exist

What's really strange for me is that, OK, the script stops on this error, but when I manually create the table that precisely raises this error, for example with pgAdmin, it goes well!

Question

So why, during the execution of the script, it somehow doesn't find postgis when creating a table with a geometry column, but when creating the table in pgAdmin it does?

It sounds like an async process that is not finished...

Information

"PostgreSQL 10.12 (Ubuntu 10.12-0ubuntu0.18.04.1) on x86_64-pc-linux-gnu, compiled by gcc (Ubuntu 7.4.0-1ubuntu1~18.04.1) 7.4.0, 64-bit"

"2.4 USE_GEOS=1 USE_PROJ=1 USE_STATS=1"


texts: Leaflet popups with source's special characters

texts:
I've got a map where I integrate spatial point data from an ArcGIS Server. Now I want a popup triggered via a click on the point. I know how everything works and the popup is displayed.

But I've got the problem that some of the source's columns names are using German special characters (ü,ä,ö).

That leads to a strange behaviour when displaying the popup; instead of the data it shows the column's name. So how could I change the code's column tag to get the right information?

The code is published under this JSFiddle.

Any help is welcome.

BR


texts: Can Mapbox legends / overlay boxes be displayed differently when viewed on mobile devices?

texts: I created a webmap using Mapbox GL JS for my website. In a desktop browser, it looks exactly the way that I want it to look. When I load it on a mobile device, the information box and legend are too wide and take up too much of the map area. Is there a way to make them display smaller when loaded on a mobile device?


texts:
I recently discovered QGIS.  The software is very impressive, even if it is not always easy to know how to unlock the features.  This article provides good pointers, although I can't make it work.

I am trying to draw 90 minute isochrones across borders, but am stuck: ORS cannot handle more than 60 minutes, while TravelTime Platform does not cross borders - strange in Europe.

I also tried HQGIS, but for some reason nothing happens when I click the 'Calculate isochrones' button. QNEAT3 does not function either, even if I use a road shapefile in a non-WGS84 projection.

Any suggestions?


texts: QGIS print composer: Transparency in legend symbols

texts:
Is there a way to make the legend in the QGIS (v2.0.1) print composer reflect the same degree of transparency as set in the Style box in the main window? The legend symbol colours are shown with no transparency but these look very different to the transparent colours on the map...

Doing a search of this question on a different website I found the following answer:
"There are two transparency sliders in the style dialog. Altering the big
one at the top doesn't affect the transparency of the legend symbol or
the style preview in the dialog, but the small transparency slider below
the size unit selector does."

However, I only see one transparency slider in the style dialog. Or am I looking in the wrong place? Any ideas?

Thanks!


texts:
I've received these coordinates from a land surveyor and no matter which converter I try or where I insert these coordinates.

I really can't find them in Google Earth or any other map.

If anyone can please assist;


texts: Trimble Geo 7x Data Transfer from Pathfinder Office

texts: I have developed a data dictionary which i want to transfer to my company's newly acquired Geo 7X but GPS pathfinder office is having trouble connecting to the device even though i can see it in Window Mobile Center


texts: Can I use an existing landcover shapefile to classify a landsat image better than ROI approach?

texts: I am trying to teach myself how to develop a land surface temperature map for the first time using landsat 8 images in QGIS. I have been following an online tutorial but it seems that, when I get to the part where I need to classify landcover into different categories based on spectral signatures, it becomes very difficult to do it accurately. If I already have a shapefile of landcover that shows areas that are buildings, versus roads, versus vegetation etc. can I somehow use this?


texts: How to reuse code with a few inputs each time in mapbasic (sub or function)?

texts:
I  need to run a string of code multiple times in a mbx while just changing 2 parameters. I tried a sub and function but both failed. How can I do this.

The code is

Print "Starting DCDB_FeatureName"

strFieldName="Feat_Name"
strOutTable="DCDB_FeatureName"

Create Table "scratch_table" (Search Char(60),X Float,Y Float) file strWorkingFolder+"scratch_table.TAB" TYPE NATIVEX Charset "UTF-16" Locale "en-AU"
Create Map For scratch_table CoordSys Earth Projection 8, 116, "m", 147, 0, 0.9996, 500000, 10000000 Bounds (-7745844.29605, 2035.05684549) (8745844.29605, 19997964.9432)
Map From scratch_table
Set Map Layer 1 Editable On
'Change Group By
Create Object As Union From DCDB_Full_Copy Into Table scratch_table Group by strFieldName Data Search=FEAT_NAME,X=X,Y=Y
Commit Table scratch_table Interactive
'Check which rows to delete
select * from scratch_table where Search = "" into Selection
Delete from Selection
'Clean up values and add search type
Update scratch_table Set Search = SearchReplace(Search, ",", " ")+" [FeatureName]" DropIndex Auto
Commit Table scratch_table
Pack Table scratch_table Graphic Data Interactive
'make sure x and y are updated
Update scratch_table
    Set X = CentroidX(obj), Y = CentroidY(obj)
'set order
Select * From  scratch_table Order By Search Into Sorted
Commit Table scratch_table
'set name of table to save data to
Commit Table Sorted As strWorkingFolder+strOutTable+".TAB" TYPE NATIVEX Charset "UTF-16" Locale "en-AU" Interactive
Close Table scratch_table Interactive
Close Table Sorted Interactive

It is just

dim strFieldName as String
dim strOutTable as String

That needs to change on each call of this code.


texts:
I have an point location and the information given is:

Easting: 519122.06

Northing: 1290630.07

Original thought was it's UTM, but does not appear to be correct because I know the point needs to be in Florida. Is there another system similar to UTM?  I need to convert these to decimal degrees.


texts: QGIS Server on a Windows machine -where do the qgs files go?

texts:
I think I have installed QGIS Server correctly - I get

I have followed the process in https://docs.qgis.org/3.16/en/docs/server_manual/getting_started.html#installation-on-windows
I have entered the service capabilities under Project-->Properties-->QGIS Server and used Test Configuration with no problems but where do I put the project and data in so that I can see it as a WMS layer? It seems like section 2.3 is for Linux and not Windows systems.
UPDATE 1
When I try http://localhost:8080/cgi-bin/qgis_mapserv.fcgi.exe?map=D:\MGM\MapData\animals.qgs&
I get

In httpd.conf I have

Is there anything else to do?
UPDATE 2
On trying the edits as per the answer for Apache
My XAMPP

My httpd_qgis-ltr.conf
# Parse our QGIS Apache conf files
Include C:/OSGeo4W64/httpd.d/httpd_*.conf

LoadModule fcgid_module modules/mod_fcgid.so

DefaultInitEnv O4W_QT_PREFIX "C:\OSGeo4W64/apps/Qt5"
DefaultInitEnv O4W_QT_BINARIES "C:\OSGeo4W64/apps/Qt5/bin"
DefaultInitEnv O4W_QT_PLUGINS "C:\OSGeo4W64/apps/Qt5/plugins"
DefaultInitEnv O4W_QT_LIBRARIES "C:\OSGeo4W64/apps/Qt5/lib"
DefaultInitEnv O4W_QT_TRANSLATIONS "C:\OSGeo4W64/apps/Qt5/translations"
DefaultInitEnv O4W_QT_HEADERS "C:\OSGeo4W64/apps/Qt5/include"
DefaultInitEnv O4W_QT_DOC "C:\OSGeo4W64/apps/Qt5/doc"

DefaultInitEnv PATH "C:\OSGeo4W64\apps\qt5\bin;C:\OSGeo4W64\bin;C:\OSGeo4W64\apps\qgis-ltr\bin;C:\OSGeo4W64\apps\grass\@grasspath@\bin;C:\OSGeo4W64\apps\grass\@grasspath@\lib;C:\Windows\system32;C:\Windows;C:\Windows\System32\Wbem"
DefaultInitEnv QGIS_PREFIX_PATH "C:\OSGeo4W64\apps\qgis-ltr"
DefaultInitEnv QT_PLUGIN_PATH "C:\OSGeo4W64\apps\qgis-ltr\qtplugins;C:\OSGeo4W64\apps\qt5\plugins"
DefaultInitEnv TEMP "C:\Users\MANGOA~1\AppData\Local\Temp"
DefaultInitEnv PYTHONHOME "C:\OSGeo4W64\apps\Python37"
DefaultInitEnv PYTHONPATH "C:\OSGeo4W64\apps\Python37;C:\OSGeo4W64\apps\Python37\Scripts"

Alias /qgis-ltr/ C:\OSGeo4W64/apps/qgis-ltr/bin/

<Directory "C:\OSGeo4W64/apps/qgis-ltr/bin/">
    SetHandler fcgid-script
    Options ExecCGI
    # Order/Allow is for Apache 2.2
    #Order allow,deny
    #Allow from all
    # Require is for Apache 2.4
    Require all granted
</Directory>

The fcgi is correct


texts:
Is it possible to have a workable ArcGIS map within Excel?

Reason:  We are generating a request form for fire flows.  When a request comes in we need to identify exactly which hydrant(s) will be tested.  It would be simplest if within the Excel request form we could embed our ArcGIS map which has our cities fire hydrant locations, instead of having to open ArcGIS separately, snip the portion of the map, and copy it into Excel for the team responsible for the fire flow testing.

I know this isn't very technical, but I was looking for a tutorial, or for someone to point me in the right direction so I could set this up.


texts:
I have a button on my toolbar that when clicked, will run a particular model, I have copied the code from here
http://resources.esri.com/help/9.3/ArcGISEngine/com_cpp/samples/com/Geoprocessing/5eb001ba-1d6b-4596-b5f8-27c2ad7cb210.htm

But when i click on the button i get the following error.
-2147467259
Automation error
Unspecified error

My code is as follows

Private Sub LoadExcelToSDE_Click()

On Error GoTo Error:

' Create the geoprocessor
Dim pGp As IGeoProcessor
Set pGp = New GeoProcessor

'Set Overwriteoutput to True
pGp.OverwriteOutput = True

'Add the custom toolbox containing the model tool
'pGp.AddToolbox "H:\Corporate GIS and Web Mapping\Projects\CSY_City_Development   \ArcMapPrototype\City_Development.tbx"
pGp.AddToolbox "C:\Documents and Settings\HalilxS\Application Data\ESRI\ArcToolbox\My Toolboxes\City_Development.tbx"

'create the parameter array
Dim pParamArray As IVariantArray
Set pParamArray = New esriSystem.VarArray

'Populate array of paramaters
'pParamArray.Add "H:\Corporate GIS and Web Mapping\Projects\CSY_City_Development\Data\Developers\UK Developers & Contact Details.xlsx\developers$"
pParamArray.Add "C:\Temp\UK Developers & Contact Details.xlsx\developers"

Dim pResult As IGeoProcessorResult

'Execute the model
Set pResult = pGp.Execute("LoadDevelopers", pParamArray, Nothing)

'Get the returned tool messages
ReturnMessages pResult

Error:
MsgBox Err.Number & vbCrLf & Err.Description

End Sub
Public Sub ReturnMessages(ByVal messages As IGeoProcessorResult)

Dim i As Long
Dim message As String

For i = 0 To messages.MessageCount - 1
    message = messages.GetMessage(i)
    Debug.Print
Next

End Sub

Any ideas on how to fix this, I have tried to change the location of the toolbox and the excel spreadsheet as well but the error still occurs.

I'm using ArcGIS 9.3.1 and VBA


texts: I am trying to visualize a 3D .img file in QGIS 3.6 . I can add the .img file, as well as a .tiff overlay to my map. When I create a new 3D view, I am able to select the .img as the elevation (which appears to work fine), but the view outputs a large, single color (green or gray depending on my overlay) blob with no defined texture. Does anybody have a fix for this? Adjusting the 3D view resolution options does not work as a fix. I have inserted relevant images:


texts: Clipping a Raster into 5x5(25) equal tiles through the use of arcpy?

texts:
I am curious to find our if there is a way to create a python code that enables you to clip a raster into equal parts example: 5x5(25) tiles with the use of arcpy?

From what I can think of you would need to divide the lengths of the extents into 5 equal parts from the width and 5 equal parts from the length and then determine the coordinates from that.

Can anyone think of a way to perform this with arcpy? or possibly provide a better method than the one I am thinking of?

    import arcpy
arcpy.env.workspace = r"folder path"
outWorkspace = "folder path"
arcpy.env.OverwriteOutput = True

rasList = arcpy.ListRasters()
print "The rasters in the workspace are:"
for ras in rasList:
    print "-- " + ras

rasDescr = arcpy.Describe("raster name")
print "\nThe type of the data is : " + rasDescr.DataType
rasExtent = rasDescr.Extent

print "The extent of the data : "
print "\tX-Min: " + str(rasExtent.XMin) + " X-Max: " + \
str(rasExtent.XMax)
print "\tY-Min: " + str(rasExtent.YMin) + " Y-Max: " + \
str(rasExtent.YMax)

print arcpy.Usage("clip_management")
coords = ("X Min   Y Min   X Max   Y Max")
output = outWorkspace + "\\" + "output name"
arcpy.Clip_management(rasDescr.Name, coords, output)


texts: Adding Addtional Styles to GeoServer via Python

texts:
I am having difficulties in adding additional styles to GeoServer Layers using pythons Requests module.

The below code works perfectly up to the try block. It will successfully publish layer from an SQL Data Store as a WMS and I have tested this via QGIS.
But as soon as I run the code in the Try block the Layer is can no longer be found in the WMS getcapabilities_1.1.1.xml and I can’t access the layer via QGIS.

The Layer still exists if I go to the GeoServer GUI and it even has the correct Styles displaying in “Selected Styles” and if I go to Edit Layer>Publishing>Root Layer in Capabilities and select “Yes” on the radial buttons it will fix the problem. The layer will be available in QGIS and will display the correct symbology.

Code

for xmlCommand in xmlCommmands:
    xmlCommandSplit = xmlCommand.split('|')
    xmlLayer = xmlCommand.split('|')[0]
    featureTypeURL = 'http://localhost:8080/geoserver/rest/workspaces/OS_OpenData/datastores/OS_OpenData/featuretypes'
    headers = {'Content-type': 'text/xml', 'Accept': 'text/xml'}
    xmlInsert = '''
                <featureType>
                    <name>{}</name>
                </featureType>
            '''.format(xmlLayer)
    createRequest = requests.post(featureTypeURL, auth=auth, data=xmlInsert, headers=headers)
    if createRequest.status_code == 201:
        print("{} - Published".format(xmlLayer))
    else:
        print("\t{} - Not Published - CODE:{}".format(xmlLayer, createRequest.status_code))
        reason = createRequest.content.decode()
        parser.feed(reason)
    try:
        style = xmlCommand.split('|')[1]
        xmlStyle = '''
                <layer>
                    <styles>
                        {}
                        <workspace>OS_OpenData</workspace>
                    </styles>
                </layer>
               '''.format(str(style).lower())
        styleURL = 'http://localhost:8080/geoserver/rest/layers/OS_OpenData:{}'.format(xmlLayer)
        headers = {'Content-type': 'text/xml', 'Accept': 'text/xml'}
        styleRequest = requests.put(styleURL, auth=auth, data=xmlStyle, headers=headers)
        if styleRequest.status_code == 200:
            print("{} - Style Set ".format(xmlLayer))
            reason = styleRequest.content.decode()
            parser.feed(reason)
        else:
            print("\t{} - Style Not Set - CODE:{} ".format(xmlLayer, styleRequest.status_code))
            reason = styleRequest.content.decode()
            parser.feed(reason)
    except IndexError:
        pass

Example xmlCommand:

bl_high_water_polyline_201910|<name>bl_high_water_polyline</name>

Here is a example of the URL: http://localhost:8080/geoserver/rest/layers/OS_OpenData:bl_high_water_polyline_201910

Here is a example of the Payload:

<layer>
    <styles>
        <name>bl_high_water_polyline</name>
        <workspace>OS_OpenData</workspace>
    </styles>
</layer>


texts: Cannot properly import MODIS MOD07_L2 .hdf files into R using rgdal

texts:
I've been trying to accurately import and process some .hdf files from the MODIS Atmospheric Profile product (MOD07_L2) in R for several days now. Still, there's something going wrong during data import. The below code can be reproduced using one exemplary .hdf file (MOD07_L2.A2013001.0835.006.2013001192145.hdf) that can be downloaded from Dropbox.

library(rgdal)

# Extraction of metadata via `GDALinfo`
filename <- "MOD07_L2.A2013001.0835.006.2013001192145.hdf"
gdalinfo <- GDALinfo(filename, returnScaleOffset = FALSE)
metadata <- attr(gdalinfo, "subdsmdata")

# Extraction of SDS string for parameter 'Skin_Temperature' (formerly 'Surface_Temperature')
sds <- metadata[grep("Skin_Temperature", metadata)[1]]
sds <- sapply(strsplit(sds, "="), "[[", 2)

# Raster import via `readGDAL`
sds.rg <- readGDAL(sds)

So far, so good, but here comes the confusing part:

> summary(sds.rg$band1)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's
 -14870  -14850  -14850  -14840  -14840  -14820   53529

Considering the fact that Skin_Temperature has an official valid range from 150 to 350 K (see MOD07_L2:Format & Content), the mean would inherit a value of

> (-14840 - (-15000)) * 0.01
[1] 1.6

after considering the corresponding add_offset (-15000) and scale_factor (0.01). Note that we're still talking about Kelvin, not °C. Extracting SDS No. 8, i.e. Skin_Temperature, using

library(gdalUtils)
gdal_translate(filename, dst_dataset = "tmp.tif", sd_index = 8)

and opening the resulting file called "tmp.tif" in QGIS resulted in seemingly reliable values centered around 15000, i.e. roughly 27 °C. However, importing "tmp.tif" back into R using raster again resulted in values comparable to the ones shown above:

> summary(raster("tmp.tif"))
              tmp
Min.    -14867.31
1st Qu. -14848.13
Median  -14845.89
3rd Qu. -14840.53
Max.    -14819.93
NA's         0.00

I've been searching the internet and stumbled across similar problems related to rgdal. However, when I tried to cast toUnSigned on band 1 of my previously generated 'SpatialGridDataFrame', I received the following error message:

> toUnSigned(sds.rg$band1, 16)
Error in toUnSigned(sds.rg$band1, 16) : band not integer

Apparently, the data imported into R is not even of type integer (what it is supposed to be), but numeric:

> sds.rg$band1[1:5]
[1]        NA        NA -14839.40 -14840.25 -14839.26

Is there an apparent mistake in my code, or is there any point I miss when importing the .hdf and .tif files using rgdal? I would be extremely grateful for any kind of help.


texts: Error launching GRASS GUI : error message related to gcmd.py or gi_set.py

texts:
I am trying to run GRASS GIS for the first time. However, when I try to run GRASS (by clicking on the icon in the start-up menu), I get the following error message:
> Starting GRASS GIS... Traceback (most recent call last):   File
> "C:\Program Files\GRASS GIS 7.8\gui\wxpython\gis_set.py", line 43, in
> <module>
>     from core.gcmd import GMessage, GError, DecodeString, RunCommand   File "C:\Program Files\GRASS GIS 7.8\gui\wxpython\core\gcmd.py", line
> 43, in <module>
>     from win32file import ReadFile, WriteFile ImportError: DLL load failed: The specified module could not be found. ERROR: Error in GUI
> startup. See messages above (if any) and if necessary, please report
> this error to the GRASS developers. On systems with package manager,
> make sure you have the right GUI package, probably named grass-gui,
> installed. To run GRASS GIS in text mode use the --text flag. Use
> '--help' for further options
>      grass78 --help See also: https://grass.osgeo.org/grass78/manuals/helptext.html Exiting... Press
> any key to continue . . .

I am using Windows 10, 64-bit.
I no almost nothing about running code from a console so if that is the solution, please indicate the very first steps.
EDIT: This error is the same as reported by the bug report on github: https://github.com/OSGeo/grass/issues/1218.
Is the solution to install another version of GRASS or wait for an updated version without the bug?


texts:
I have a MXD with some layers that have parentheses in them. I'd like to rename them via script.

Is this a valid approach?

import re
map = arcpy.mapping.MapDocument("CURRENT")
layers = arcpy.mapping.ListLayers(map)
for layer in layers:
    if re.search("[()]", layer.name):
        newName = re.sub("[()]", "", layer.name)
        arcpy.Rename_management(layer, newName)

Rename_management is giving me an error:

  Runtime error : ERROR 000840:
  The value is not a Data Element.


texts:
I have been given a CSV which has a mixture of point and polygon data.  I am trying to open this up in ArcMap 10.5.1

I have no idea what I should be doing!  I have tried to look for answers online but pretty much everything is a python script that is just dealing with points, and looks like it isn't much data, compared to my mix of types and large set.

Any suggestions on how I can view this data in my map??


texts: How to create buttons with Javascript that set a new view (extent and zoom)?

texts:
I am creating a web map using CartoDB, and am having trouble creating buttons that will set a new view when clicked. I am a novice when it comes to coding but i do have some basic knowledge.

I'm thinking its going to involve creating a new visualization using CartoDB's .createvis command and a SQL statement but I'm not sure.

Any advice?


texts: How do I extract elevation data from a DEM using points I have in a CSV?

texts:
So I have a CSV with a lot of fields. The first 3 fields are the ID, longitude, and latitude. I have already converted to decimal degrees, and said to display XY data. (Is there something else I need to do before I can use these points?)

Now I need to know what Spatial Analyst tool to use to extract the elevation for each point, and how to use that tool. I tried extract by mask and Extract Multi Values to points, but neither worked.

Extract by mask left me with just one elevation for all 531 points. Extract by multi values failed to execute.

What do you think? I'm very inexperienced, so a step-by-step guide would help a lot.


texts: Is it possible to zoom out a plot of a map in R?

texts:
Sometimes a plot looks zoomed in how to arrange that?
i use the rgeos package mainly. When i do a buffer for example somewhere the plot tends to zoom in some part of the map which is not good since it cuts parts even from the buffer too in the plot.


texts:
I want to create a polygon with given coordinates from scratch. Doing this using the QGIS VectorLayer class gives no output. I used the example from the pyqgis cookbook, wich is creating point features and modified it to expect an polygon output . I also tried to create multiple points with this code by adding more QgsPointXY() elements but this also gives me no output.

from qgis.PyQt.QtCore import QVariant
vl = QgsVectorLayer("Polygon", "temp", "memory")

pr = vl.dataProvider()
#pr.addAttributes([QgsField("name", QVariant.String),
#                  QgsField("age",  QVariant.Int),
#                  QgsField("size", QVariant.Double)])
vl.updateFields()

f = QgsFeature()
f.setGeometry(QgsGeometry.fromPolygonXY(QgsPointXY(0,0),QgsPointXY(0,10),QgsPointXY(10,10),QgsPointXY(10,0)))
#f.setAttributes(["Ada L.", 2, 0.3])
pr.addFeature(f)
vl.updateExtents()
QgsProject.instance().addMapLayer(vl)

What should be done to make this code running?

My environment: Windows10, QGIS 3.10 Python 3.7


texts: Callout text not respecting space(s) between words in ArcGIS Pro Layout

texts:
I suspect that I may have inadvertently changed a setting somewhere in ArcGIS Pro 2.6.2 because when I use Callout (Sans Serif) as my Default Text Symbol any space(s) I use between words in the text I am trying to place are not appearing in the expected way.
To reproduce this try the following:

Start ArcGIS Pro and use Start without a template to create a new project
Insert a New Layout (I used A4,Portrait)
Use the dialog launcher of the Graphics and Text button group to open the Options dialog to the Text and Graphics tab

Click the Reset to Defaults button to try and clear any inadvertent settings
Click OK to close the Options dialog and accept those Default symbols
Use the Rectangle text button in the gallery of the Graphics and Text button group to drag out some text
Use the Element pane (which has Format Text as its title) to set the Text to be Project Site
Zoom the page so that the text element just added is easy to see, and there's room to add a few more below it.  The words Project Site have the expected space between them on the layout.

Use the dialog launcher of the Graphics and Text button group to open the Options dialog again, and change the Default Text Symbol to be Callout (Sans Serif)

Use the Rectangle text button in the gallery of the Graphics and Text button group to drag out some text below the first text element at about the same size
Use the Element pane (which has Format Text as its title) to set the Text to be Project Site

The expected space between the two words of Project Site is missing, and the "t" of Project even overlaps the "S" of Site.

Use the Rectangle text button in the gallery of the Graphics and Text button group to drag out a third text element below the first two at about the same size
Use the Element pane (which has Format Text as its title) to set the Text to be Project   Site (with three spaces between the two words)

The overlap between the two words of Project Site is now even more pronounced.
Why is the spacing between words not being respected?


texts: Printing result of custom model in Python using ModelBuilder and ArcPy?

texts:
I have a custom model (called TOOLNAME) in a custom toolbox developed in Model Builder on ArcGIS 10.1 Desktop. The model takes two input parameters (SingleFDI_D and dayindex) and outputs a Derived parameter called "MaxFDI_DateTime" which is data of String type.

The model runs well within the ModelBuilder on ArcGIS Desktop. It returns the result MaxFDI_DateTeime: 22/11/2017 12:00:00 AM as below.

I will need to run the model in Python.

How do I get the output of MaxFDI_DateTime out of the tool execution in Python?

arcpy.ImportToolbox(TOOL_DIR)
# The custom model should return a String, how to retrieve it here?
model_result = arcpy.{TOOLNAME}_{TOOLBOX_ALIAS}(SingleFDI_D, str(day))


texts:
I've run into an error when trying to open a Landsat-8 geotiff file in python, it seems like I am missing some tiff libraries but I'm not sure how to install them or where they need to go. Has anyone else encountered this issue when working with geotiff files?
I'm using the Anaconda 3.6 distribution of python on a windows 10 64 bit computer.
Python code:
from osgeo import gdal
path = "LC80040472015120LGN00_B4.TIF"
geotiff = gdal.Open(path)

Error Message:

Warning 1: LC80040472015120LGN00_B4.TIF:AdobeDeflate compression
support is not configured
ERROR 1: Cannot open TIFF file due to
missing codec.

Edit1:
Here is the index page on AWS for the image that I'm trying to open.
https://s3-us-west-2.amazonaws.com/landsat-pds/L8/004/047/LC80040472015120LGN00/index.html


texts: Issue with vertex level editing in ArcScene

texts:
I have a 3D data set of some buildings in ESRI Shapefile MultiPatch format. The buildings appear in 3D in ArcScene (10.2.1 on Windows 8.1 64bit) and I am able to move, rotate and scale buildings. However when I select a building in an editing session and I go to vertex editing, if I try to edit any of the vertices of my building I am presented with a 'Select a single editable feature to modify' prompt. I am unable to work out why I can't edit at vertex level for this set of buildings.

In addition clicking on the building appears to make it disappears until I select another building and the original building reappears. Anyone have any thoughts about this issue?


texts:
I have been working with ArcObjects Java for a short time and have not been able to find a solution to the following problems:

How to get the Raster object from a IRaster reference?

IRaster raster = rasterDataset.createDefaultRaster()

Raster rs = (Raster)raster

This gives classcastexception:  com.esri.arcgis.geodatabase.IRasterProxy cannot be cast to com.esri.arcgis.datasourcesraster.Raster

I have been unable to find the documentation of IRasterProxy. Apparently it is deprecated!!

How to over write pixel values of an existing raster and store it permanently?

I am following this example but to no avail. I can see the pixel array is altered but when I open the raster in ArcMap I see no change.

How to get Raster Info?

RasterDataset streamOrderRas = (RasterDataset)streamOrderRaster

Raster raster = new Raster()

raster.appendBands(streamOrderRas)

IRasterInfo rasterInfo = raster.getRasterInfo()

I get AutomationException: No such interface supported, stacktrace shows it goes inside com.esri.arcgis.datasourcesraster.Raster.getRasterInfo()

Please forgive if the questions appear stupid. If some one can show some working examples of the above problems, i'd be very grateful.


texts: Layer Order in Pop-up information - ArcGIS Online

texts:
Is there a way that I can choose the priority on how the pop-up window shows on ArcGIS Online?

For example, I have 2 layers named Active Country and Inactive Country. In my current data, I have 2 Nigeria polygons wherein each represent a program.

What I want to happen is to decide the pop-up display order of overlapped polygons.


texts: Tracklogs to single line OR several lines to one line?

texts:
Setting:

I have several trucks driving out of my factory delivering goods and coming back for more deliveries. These trucks are equipped with GPS so that they may return tracklogs to me. These tracklogs are in table format, so it is not hard for me to convert the logs to point shapes or line shapes. But after, say, a month I will have so many overlapping logs, it makes no sense to look at. Therefore I would like to "merge" all the roads together into one where ever they are close to each other. I would also like an attribute which tells me how many times each section has been driven. Finally I want to be able to add new tracklogs to the final result.

Link to image: http://postimage.org/image/1g1sfdnz8/

Tested:

I have tried using the "Merge Divided Roads (Cartography)" tool, but it does not give me the result i want. When using the "Merge Divided Roads (Cartography)" tool, it manages to merge the beginning of the roads, but as soon as they cross paths, it will not merge anymore. (I am guessing this has something to do with the tool being optimized for road networks and all I have are shape files.)

This is practically what "Merge Divided Roads (Cartography)" tool gives me:

Link to image: http://postimage.org/image/1g20p2eg4/

Questions:

Anybody know how to get two or more lines to become one line where they are close and almost parallel?
Would you also know how to get the count of how many times a section of the road has been driven?

Information:

I am using ArcGIS 10 with ArcInfo + some extensions


texts: How to create unique plane with PostGIS?

texts:
How to create a unique plane, with (x , y ,z) as center and a vector as normal, in PostGIS?

I already saw "Make3DBox" function, but it does not make a unique one.


texts: Making a buffer - degrees instead of meters

texts:
I just started working with GIS very recently.

I'm trying to create a buffer of 200m around the points of the layer "vlaamse venen naamloze laag". But in the distance box, the program will only allow me to create a buffer in degrees? I checked the project properties, but they are in meters.

I tried using this question, (gis.stackexchange.com/questions/97225/create-buffer-in-metres) but I can't figure it out.
Changing the CRS of the layers makes that the distance is put in meters but then i have a different problem. If I change them (e.g. to Lambert 72 or Lambert 2008), my layers dont overlap anymore, which they should.


texts: Computational time of well-known compression algorithm

texts:
I have googled for get information about computational time of some compression algorithms but information seems contradict each other. I need worst case performance of the following algorithm Sliding Window, Open Window, Threshold and STTrace.
All algorithms should have quadratic time. Do you confirm?


texts: Using Folium to visulize PCA image raise "Error expected type string, actual type long,actual value 1

texts:
This is how I get the PCs at the specified scale and in the specified region:

# Get the Principal Components (PCs)

pcImage = getPrincipalComponents(centered, scale, region) #this return PCA image
# Plot each PC as a new layer
for i in range(bandNames.length().getInfo()):
    band = pcImage.bandNames().get(i)
    map.addLayer(pcImage.select(['band']),{'min':-2,'max':2},'band')

Here is the pca functionin which am passing an ee image:

 # Principal Components (PC) in the region as a new image.
def getPrincipalComponents(centered, scale, region):
    #Collapse the bands of the image into a 1D array per pixel
    arrays = centered.toArray()
    #Compute the covariance of the bands within the region.
    covar=arrays.reduceRegion(
        reducer=ee.Reducer.centeredCovariance(),
        geometry=region,
        scale=scale,
        maxPixels=1e9)

    #Get the 'array' covariance result and cast to an array.
    # This represents the band-to-band covariance within the region.
    covarArray = ee.Array(covar.get('array'))
    #Perform an eigen analysis and slice apart the values and vectors.
    eigens = covarArray.eigen()
    # This is a P-length vector of Eigenvalues.
    eigenValues = eigens.slice(1, 0, 1)
    # This is a PxP matrix with eigenvectors in rows.
    eigenVectors = eigens.slice(1, 1)
    # Convert the array image to 2D arrays for matrix computations.
    arrayImage = arrays.toArray(1)
    #Left multiply the image array by the matrix of eigenvectors.
    principalComponents = ee.Image(eigenVectors).matrixMultiply(arrayImage)
    # Turn the square roots of the Eigenvalues into a P-band image.
    sdImage = ee.Image(eigenValues.sqrt())\
        .arrayProject([0]).arrayFlatten([getNewBandNames('sd')])
    # Turn the PCs into a P-band image, normalized by SD.
    # Throw out an an unneeded dimension, [[]] -> [].
    # Make the one band array image a multi-band image, [] -> image.
    # Normalize the PCs by their SDs.
    return principalComponents\
        .arrayProject([0])\
        .arrayFlatten([getNewBandNames('pc')])\
        .divide(sdImage)


texts: Conditional value maps in a comboBox using QGIS

texts:
I'm using a layer which has some attributes. The first attribute (lets say "Country", for example) has a comboBox and the user can freely choose which value he wants (For example, any country in North America). The second comboBox (lets say "Cities", for example), will display only the values which are related to the first attribute (only cities of the corresponded country). So, the values of the second comboBox will be somehow related to the first comboBox (not displaying cities which does not belong to the selected country). I followed the following website and I was able to do it as I wanted (using value relation): https://north-road.com/drill-down-cascading-forms/.
However, I had to create a whole new layer only for the sake of mapping "Countries" with its "Cities". Is there a way of doing it, without having to create a whole new layer only for mapping the values? Could I do it only using filters or something like this? (This will be useful for me, because I'm gonna have a special case, when I will have to map a very few values conditionally)
In the end, I would have something like this:


texts: Color blurring on GeoServer WMS layer loading using OpenLayers 5

texts:
I am using GeoServer 2.14.1 and OpenLayers 5 API to view the polygon data to my web map application using WMS service. I have used simple polygon style with black color.When zooming, the color of the data is fully occupied the screen (Chrome & Firefox) and then vanishing tile by tile as shown in the snapshot.

Before zooming

When zooming

Is this because of GeoServer style or do we need to configure anything at the layer publishing or do we need to tweak something in OpenLayers 5 CSS?


texts: Converting from dgn to shp or dwg files using QGIS?

texts:
I am a frequent recipient of microstation dgn files and don't like it because I don't have software to convert it to autocad. I do have ARCVIEW, but I can't edit the files. I have tried other open source programs to convert, but they all are lacking.

I'm wondering if anyone has used QGIS to convert from dgn to shp or dwg file types?

I'm skeptical because I noticed that QGIS supports dgn file types, but doesn't support dwg file formats.


texts: How to Access the Current Value of the Scale Factor of a Map as a Variable (eg $map_scale) in Print Layouts in QGIS

texts:
Does anybody know how to access the current value of the scale factor of a map as a variable (eg $map_scale) in print layouts in QGIS?
I'd love to put the scale factor in my maps...


texts:
I have QGIS with GRASS 7. I'm looking for the old style GRASS digitising tools.

Previously, when digitising boundaries, it was required to split any adjacent boundaries intended to be joined to. The split tool seems to be missing in GRASS 7.

The toolset was also included in earlier versions of QGIS and the two (GRASS and QGIS) worked well together. I can't digitize even simple areas using current versions of GRASS and QGIS.


texts:
I'm trying to generate some tiles from a very large dataset (5.8 million points). Each point represents a 20 x 20 grid square in the OSGB 1936 / National Grid CRS.

What I'd like is for each square to be coloured according to some property of the point. I've made very good progress by making my marker a square of size 20 map units, with a graduated colouring.

The issue I'm having (besides the time it takes to generate the tiles) is that the output tiles have very noticeable artefacts, that look like the kind of patterns you get when sampling an image at the wrong rate.

I think this is because the markers don't align to a perfect grid. Even in the QGIS viewer, I can see tiny light-coloured lines in between each marker.

There are two options I can think of:

Improve the sampling rate of QTiles
Fix the underlying issue of the markers not exactly lining up with the grid

Are there any ideas on how I can do either of these, or any other ways to improve the look of the tiles?


texts: Symbolize Rubberband Tracker as it tracks

texts:
I am trying to symbolize a rubberband tracker as it is drawing.  The default symbol is a simple black line.  When setting the symbol of the line in the ScreenDisplay method and and calling DrawPolyline, the symbol I have set does appear but this is after the draw has completed.

I would like to have the symbol I set be shown once a single point has been created (while its tracking the next point), not just when the construction of the line has completed.

This functionality seems to be possible as I can easily do it in ArcMap - I can select any symbol and then create a new line feature and it displays the selected symbol while I draw it.

Edit: I noticed that it doesn't refresh as cleanly as the built in ArcMap tools do either.  Strange, could there be some side effect happening?  I've added two clips of what I experience and what I would like to have happen and exactly how I would like the line to look.

Result from this code:

What I would like:

Relevant Code:

IScreenDisplay screenDisplay = ArcMap.Document.ActiveView.ScreenDisplay;
screenDisplay.StartDrawing(screenDisplay.hDC, (int)esriScreenCache.esriNoScreenCache);

// Custom symbol I create gets cast as ISymbol...
ISymbol symbol = (ISymbol)cartographicLineSymbol;

IRubberBand2 rubberBand = new RubberLineClass();
ICurve2 rubberLine = (ICurve2)rubberBand.TrackNew(screenDisplay, symbol);

screenDisplay.SetSymbol(symbol);
screenDisplay.DrawPolyline(rubberLine);
screenDisplay.FinishDrawing();


texts: Adding Z value to attribute table of shapefile from DEM?

texts: I have one shape file(polygon) downloaded from freely available sites. I want to add the Z value or height value into it. I have corresponding DEM file. The Z column has come after using the 'Add Surface Information' tool. But all the values are 0 in the attribute table and my DEM has elevation values within it. How can I include the Z column having values in it?


texts: QGIS 3.4 set transparency only for some layers below

texts:
I‘m working on a geomorphological map of a river and its sourroundings. As a base I have a topographic map and a DEM for the pre-mapping. The topographic map is not accurate enough in the river gulley, so I covered it with a polygon to give me the possibility to draw on the polygon instead of on the misleading topological map when I‘ll do my field work. But I set the opacity of my polygon to 100% it will also cover my DEM (but I‘d like to see through the DEM but not the topographic map). Is there a way set a layers opacity to different percentages for the different layers below?

Or do you have any other ideas how I can fix this problem, e.g. cutting the topographic map on the riverborders?


texts: Is there a hybrid enterprise GIS I can setup that will allow me to edit on ArcMap and have geoserver hosted in the cloud

texts:
I am looking for a hybrid enterprise GIS that will allow my company to continue using our ArcMap licenses but that doesn’t cost a fortune to do.
It seems like it’s either all open source or it’s all ESRI.
We need to be able to have a multi user/multi editor environment that will allow us to create web apps from.
Is it possible to have ArcMap for editing, GeoServer for serving the data to our editors, and use leaflet or ArcGIS online for web app development?


texts:
Sorry if this is a repeat, but I couldn't find an answer here.  And sorry if the answer is completely obvious, but I do appreciate your assistance with my first QGIS map.

When I try to use composer to create a printable PDF of a map I've made, clicking 'add new map', then drawing a rectangle, the map does not show up.  I can add and see a north arrow image or a scale bar, either on top of or to the side of the map rectangle.  I can add and adjust a frame around my map and I can change the background color.  I tried to show the grid and cannot.

I downloaded a map from www.naturalearthdata.com and composer seems to work fine for me, including moving the map, changing the scale, adding the north arrow image and the scale bar.

I'm using QGIS 1.8.0-Lisboa, with Windows 8.

Any suggestions?


texts: Plugin for fast export selected features from multiple layers

texts:
Is there any tool or plugin providing export of selected features from more than one layer on one click?

I need to create multiple exports from three layers. Feature from layer one, then feature from layer two and feature from layer three.

I need 3x click on selected layer, choose save as, click save only selected and select path and name. For 100 exports it is not very effective approach.


texts: Extracting data from one grid cells to another set of grid cells

texts:
I have my own point shapefile (each point has its lats n longs). I found a data which has lats and longs, and I imported it into ArcMap and created a shapefile which also has its own points (with lats and longs). Now I want to know whether there is a way for me to extract data from the new point shapefile I created to the older point shapefile. The thing is that the points of the two shapefiles are different (although some will accidentally match).

Is there a way I can export data from one of them to the other (especially from the new point shapefile to the old one)?


texts: ArcView 10.1 does not have this kind of topo map to "Add Data".  Where do I find this type of topo map, preferably free.  ArcGIS online does not seem to offer this type of topo, unless I am not seeing it on my side.


texts: How can I split a segment into two segments?

texts:
The problem is how to split segments or in other words how to insert new points into a polygon so that the segment where I created the new point is now two segments connected by said point. Every time I create a new point within a segment, the segment stays the same length as before and is not being split by the new point. (I don't know why the following snapshots have double columns on the left side but just ignore that. was not like this in the program). First I have the segment I want to split:

.

Then I go to edit mode for the segment layer, click on the knot tool and select the segment:

.

I double-click on the segment where i want a new point. The "new place for the point" appears:

.

I then go to pointlayer, activate the edit mode, click on "add new element" and create a new point near the segment:

.

After I activate the snap-function and define the radius, I place the point on the "x" with the knot-tool:

.

Now if I go on the info-pointer and click on the segment, the segment is still as one piece and is not split by the new point:

.

Why?


texts:
I have made a 1 / N union using the QGIS tool between an attribute table in CSV format and a shapefile using a common field(N).

The .csv attribute layer has the following field structure:

ID - N- T

1  - 1 -  T5

2 - 1  - T5

3 - 1 -  T5

and the shapefile layer the following:

N - F - G

1 - F1 - G1

2 - F2 -  G2

If I load the CSV layer from PostGIS I can edit the attributes but if I load from a folder it does not allow me to edit those attributes.

Is this edition only possible by loading the CSV layer from PostGIS?


texts: Layers not showing up in layers panel, but show in styling panel

texts:
I have a QGIS project that has stopped showing the layers in the layers panel and the map, but they show up in the layer styling panel:

I have had this issue before, and reported it on the QGIS Users mailing list - but did not get any responses.

Here is the QGIS file here, with the user/passwords removed:
https://www.dropbox.com/s/tfdbe02drj3uaua/Layers_Not_Showing.qgs?dl=0

Can anyone spot anything in the file that might be corrupted?


texts: Optimizing MySQL query to select all points with in a large polygon using spatial indexes

texts:
Re-posting it from stackoverflow as per their suggestion. Link to the thread in stackoverflow https://stackoverflow.com/questions/53182385/optimizing-mysql-query-to-select-all-points-with-in-polygon-using-spatial-indexe

Firstly, I admit that my experience with spatial functions is very minimal. I have a table in MySQL with 20 fields and 23549187 records that contain geographical data. One of the fields is 'point' which is of point data type and has spatial index on it. I have a query that selects all points within a polygon which looks like this,

select * from `table_name` where ST_CONTAINS(ST_GEOMFROMTEXT('POLYGON((151.186 -23.497,151.207 -23.505,151.178 -23.496,151.174 -23.49800000000001,151.176 -23.496,151.179 -23.49500000000002,151.186 -23.497))'), `point`)

This works well as the polygon is small. However, if the polygon gets massive, the execution times gets really slow and the slowest query until now ran for 15 mins. Adding the index had really helped to bring it down to 15 mins which otherwise would have taken close to an hour. Is there anything I can do here for further improvement.
This query will be run by a PHP script that runs as a daemon and I am worried if this slow queries will bring the MySQL server down.

show create table;

CREATE TABLE `table_name` (
  `id` int(10) unsigned NOT NULL AUTO_INCREMENT,
  `lat` float(12,6) DEFAULT NULL,
  `long` float(12,6) DEFAULT NULL,
  `point` point NOT NULL,
  PRIMARY KEY (`id`),
  KEY `lat` (`lat`,`long`),
  SPATIAL KEY `sp_index` (`point`)
) ENGINE=MyISAM AUTO_INCREMENT=47222773 DEFAULT CHARSET=utf8mb4

There are few more fields that I am not supposed to disclose it here however the filter won

Explain sql output for the slow query:

+----+-------------+--------+------+---------------+------+---------+------+----------+-------------+
| id | select_type | table  | type | possible_keys | key  | key_len | ref  | rows     | Extra       |
+----+-------------+--------+------+---------------+------+---------+------+----------+-------------+
|  1 | SIMPLE      | table_name | ALL  | NULL          | NULL | NULL    | NULL | 23549187 | Using where |
+----+-------------+--------+------+---------------+------+---------+------+----------+-------------+

Explain sql output for query with smaller polygons,

+----+-------------+--------+-------+---------------+----------+---------+------+------+-------------+
| id | select_type | table  | type  | possible_keys | key      | key_len | ref  | rows | Extra       |
+----+-------------+--------+-------+---------------+----------+---------+------+------+-------------+
|  1 | SIMPLE      | table_name | range | sp_index      | sp_index | 34      | NULL |    1 | Using where |
+----+-------------+--------+-------+---------------+----------+---------+------+------+-------------+

Looks like the biggest polygon does not use the index.


texts: How to pass custom data in GeoJSON point?

texts:
I have an application, which creates a GeoJSON stream on the server, which is then rendered on the client (using CesiumJS). GeoJSON is created using GeoTools.

I want to pass some custom data to the client about a particular point. This custom data will determine, how the client renders that point.

I add points using this code:

protected void addPoint(final double lon, final double lat,
    final GeometryBuilder builder, final SimpleFeatureBuilder fbuilder,
    final List<SimpleFeature> features) {
    final Point point = builder.point(lon, lat);
    fbuilder.add(point);
    final SimpleFeature feature = fbuilder.buildFeature(null);
    features.add(feature);
}

If I wanted to add a string property to GeoJSON, with, say, name id and value fooBar123.

How do I do it in Java (how can I modify the code above in order for id:fooBar123 to appear in the results GeoJSON string) ?


texts: ArcGIS Pro Dynamic Text Project Path abbreviation

texts:
I am using Arcgis Pro and have saved a template for a map which uses the default 'Project Path'  dynamic text option to display the folder the location of the project. The problem is some of the file paths being used in this template are exceedingly long and difficult to display.

Is there any option to abbreviate the file path in dynamic text to only show the last two or three components?

I.e. Instead of showing

Path C:\FolderA\ProjectB\OfficerC\Working\IndustryD\Analysis\Map.aprx

just show:
\IndustryD\Analysis\Map.aprx


texts: How to create a topologically correct dataset in PostGIS (with GRASS)?

texts:
I have shapefiles which I have loaded in postgis and now i am trying to use geometric function in postgis like ST_intersect and they fail because the data is not topological correct.

I have tried to use the function st_isvalid and others and could not correct the data. As a last resort i have loaded the data into grass which automatically builds topology and trying to export that layer into postgis fails with an error saying gid already exist.

Also, trying to export as a shapefile works but i encounter a problem when i load the shapefile back to postgis and try to do my clip as it complains that geometric column contains z column but geometry does not.

How best can i attend to these?


texts: ArcGIS Model - Iterate Row Selection - Multipal viewsheds

texts:
I've a CSV file with 106 rows with data suitable for a viewshed.

What I want to do is run an ArcGIS model where it works through the rows and produces 106 viewsheds for me (ideally the data could also be incorporated back into the original spreadsheet).

Initially I started playing with the 'Iterate Row Selection' but this does not seem to be working. It only produces one.

Does anyone have any ideas?

Here is an image as to how the model is set up. The Iterator is grouping by the ODI column.


texts: Query and display MySQL DB data from HTML

texts:
I want to link MySQL DB to the HTML and get the data according to the query below.

SELECT * FROM feature1 WHERE place = 'feature2'

php seems one of the way to do this, since I am one day old with php don't have much knowledge. I wrote the code below, but when run it gives me an error with no solid output. Is the code correct?

I want to query the DB using the above query and display data!!

<?php
$db = mysqli_connect('localhost','user','password','DB_name')
or die('Error connecting to MySQL server.');

> $query = "SELECT * FROM feature1";

mysqli_query($db, $query) or die('Error querying database..');
$result = mysqli_query($db, $query);
$row = mysqli_fetch_array($result);
while ($rows = mysql_fetch_array($query)) {
   $name = $rows['Name'];
   $address = $rows['Address'];
   $email = $rows['Email']

   echo "$name
$address
$email
$subject
$comment

";
?>


texts:
I have GPS tracks from 2 or more devices in multiple vehicles. At least one device reports the vehicle's position to a webui approximately every 5 seconds assuming that it has a cellular data connection. the second device logs its position to a file every second while it is functioning, however, it is not always functioning for the entire duration of the journey.

I would like to take these two tracks and merge them into a single track for analysis purposes. My naive first approach was to use the timestamps in the file and the webui to synchronize the tracks, however, there is no shared clock between the data stored from the two devices and I am not convinced that there is a constant delay in the webui log.

How can I merge the tracks, synchronize the timestamps, and interleave the points appropriately?


texts: Finding out if I've deviated from given route too much (computing similarity between routes)

texts:
I'm making an Android app and I'm trying to find out whether I've deviated from the given route too much (I'm getting location via GPS and the route is predefined).

I've tried calculating Fréchet distance but to no avail - it is only useful to compare two complete routes, but what if I need to compare my track (which is not yet complete) to a given track?


texts: Retain the rest of the Clip (Analysis) tool in ArcGIS

texts:
Which tool can I use to get or keep the rest of the output feature when I use the Clip (Analysis) tool in ArcMap?

For instance I need to calculate Blue minus Green.

Erase (Analysis) from Overlay will do.
That's what I was looking for.


texts: OpenLayers LineString: Not displaying the line (array of points)

texts:
I'm perplexed as to why tail_path is not drawing the line (can't see anything on the map) but the alert seems to be returning the correct data.

var path_travelled = new Array();
var path = json_coords.path_data;
for (i in path)
  {
    var lat = path[i].LINE_LAT;
    var long = path[i].LINE_LONG;
    path_travelled.push(new OpenLayers.Geometry.Point(long, lat));
    alert(path_travelled);
  }
var tail_path = new OpenLayers.Geometry.LineString(path_travelled);
vectors.addFeatures([tail_path]);

The alert for the path_travelled variable returns this (it adds a new POINT every time you click the OK button):

I'm using vectors.addFeatures (I have another variable I also use here but doesn't apply to this question so just removed it to avoid confusion).

Any input on the matter is appreciated.


texts: Can ArcGIS Transpose Fields tool output as numeric?

texts:
My table stores the same variable for each year in separate fields - ex. Pop2010, Pop2011, Pop2012. I used the Transpose Fields tool to transpose the years to long format. It transposes as expected except that the Value Field is stored as text.

Is there a way to set this as integer or double?


texts: Is it possible to download neighborhood boundaries from OpenStreetMap? I see on maps that they often indicate neighborhood names; is there a way to retrieve their coordinates and boundary polygon coordinates, much like the Zillow Neighborhood Boundaries data set?


texts: I'm using Mapnik in Python for batch creation of png files from my raster layers. <RasterSymbolizer> from XML style works perfect in my case. However I would like to put a small text annotation on each png file with the filename and maybe a small copyright note. Is it possible to do it within the Mapnik style?


texts: Aerial LiDAR: What is an acceptable resolution?

texts:
I am in the process of initiating a project involving laser scanning. The area to be scanned i s approximately 5000 km². We are initially going to use the data to create a solar potential application, but will also be obliged to use the raw data for DGM, Surface modelling, Flood risk assessments and even 3D City modelling.

We could order data with an accuracy of up to 10 points p/m². But usually for solar potential analysis not more than than 2 points p/m² is required. Having a more accurate data set sounds great, but I am also worried that my hardware won´t be able to handle these amounts of data.

I have recently created a DGM for the 5000km² with a point concentration of 1 point for 12m, and this was quite time consuming.

Has anyone had any experience working with similar data sets? what kind of hardware is able to handle these amounts of data? I have a new windows 7 machine with 12 cores and 8gb Ram...but is this enough?


texts: UI in QGIS for a UML class with multiple cardinality attribute

texts:
For example: a "way" class with a name[1:n] attribute. Mapping the UML class to SQL you can get something like:

create table way
(
   id serial not null primary key,
   geom geometry,
   other_attribute text
)

create table way_name
(
   way_id integer not null,
   name text not null,
   foreign key (way_id) references way (id) on update cascade on delete cascade,
   primary key way_name_pkey (way_id, name)
)

QGIS can automatic discovers the relationship between tables and you can set up a form for the "parent" entity (way) and include a subform-relation for the names. But it seems that creating a new way is a 2-step process in QGIS because when I edit a new way I can´t set its name. I have to close/accept the form,
save layers changes, open the form again and (now that the way exists in the database) set the way name. Is there any solution for 1-step new features editing?

Obviously, I can do it if no foreign key exists but I can lost the relational integrity.


texts: In QGIS 3.10, is there a way to click into selected features in the attribute table and make a selection from a selection, like in ArcGIS?


texts: Arcgis javascript - Set different symbol for each graphic in feature layer

texts:
I have a feature layer that contains 12 items with the same defined style. But I want to change each item to have a unique symbol when rendered.

var Layer = new FeatureLayer("http://...my map server.../0");

Layer.graphics.forEach(function (entry) {
    // loop through each item
    console.log(entry);
    console.log(entry.attributes.FID);

    var symbol = new PictureMarkerSymbol({ "url":"/images/icons/fid_" + entry.attributes.FID + ".png", "height":20, "width":20 });
    entry.setSymbol(symbol);
    var renderer = new esri.renderer.SimpleRenderer(symbol);
    Layer.setRenderer(renderer); // <-- this appears to override all previous items
});

map.addLayer(Layer);

I omitted some code to keep this simple; what I have now sets a symbol, but it's the symbol of the last entry that is looped for all items.

I need each entry to have a unique symbol defined by its FID (0-11)

I've search through ArcGIS's JavaScript API docs and can't find anything that addresses this.


texts: Problems registering a folder to ArcGIS for Server 10.2.1

texts:
I have ArcGIS for server 10.2.1 running and I have created an ArcGIS Server account.

This account has been given read/write access to the network at the office.

To verify this, I have logged on to a machine using the account which has been assigned to the server and connected some drive letters up to the locations of the folders I want to connect to.  I've also created folders and text files in there to make sure the account has write permissions.

I've also opened the "Configure ArcGIS Server account" utility to double check that the correct username has been assigned as the Server Account.

The problem that I'm having is that whenever I try to Register a folder to the ArcGIS server that isn't one of the local drives on the box, I get an error.

I get this regardles of wether I use the same path as the publisher path, or a different path.

But only when I try to map to a folder on the network.
Mapping to a folder on "C:" works fine.

Would anybody have any suggestions of what I could do next? Is there something really simple that I'm missing?


texts:
Do you know of a web site that has a json file I can download that contains the following data:

City Name or postcode
latitude in degrees
longitude in degrees
How many hours +- the city is from GMT (Ex: Boston is -5)

Nice to have but not required:
daylight savings time (Yes or No). Ex: Boston, Massachusetts uses daylight savings time and Phoenix, Arizona does not.

For now I'm looking to get the data for cities in the US. I will add other countries later if the web resource has multiple countries.


texts: Can't select using a custom function in QGIS 2.18.3

texts:
I'm following this tutorial on custom expressions and all goes well - my function shows up in the custom section and seems to work for the preview giving 21S for the first feature. When I extend the selection expression to GetUtmZone() = '21S' the output preview still says 1. But when I click on the select button no features are selected.

This works on 2.14.x but fails on 2.18.x.

The old style $GetUtmZone doesn't work in 2.14 or 2.18.


texts: Limit custom tiles to certain extent in OpenLayers

texts:
I want to fit my custom tiles at a certain region. For example let's say i have custom tiles for a certain country (example: USA) and I am using OSM or anything at all. How do i make sure that my tiles fit a certain extent. I tried putting in the extent but the layer does not load. Please help.

map = new Map({
      target: 'map',
      layers: [
        new TileLayer({
          source: new OSM({
            wrapX: false,
          }),
        }),
        new TileLayer({
          extent:[12679211,2560587,12703254,2544509],
          source: new XYZ({
            maxZoom: 16,
            projection: projection,
            tileSize: tileSize,
            tileUrlFunction: function(tileCoord) {
              console.log("--------------")
              console.log("before",tileCoord);
              console.log("after", " z = ",tileCoord[0] - 1," x = ",tileCoord[1]," y = ",(-tileCoord[2] - 1))
              return urlTemplate.replace('{z}', (tileCoord[0] - 1).toString())
                .replace('{x}', tileCoord[1].toString())
                .replace('{y}', (-tileCoord[2] - 1).toString());
            },
            wrapX: false
          })
        }),
      ],
      view: new View({
        center:  [12689317.435763245, 2552726.0298609314  ],
        zoom:13,
        minZoom: 0
      })
    });


texts:
How can i expand extent of a raster image with Python as the image shown in below. I have the coordinate information of a raster image such as Xmin Ymin Xmax Ymax,how can i expand its extent based on the specified size(Ex, expand 10px or 1px)? I know it can be easily implemented on AML through the command such as setwindow, con(isnull ...etc, but i have no idea on ArcPy.

Does anyone know about it?


texts: Multiple columns in TOC and Legend in ArcMap

texts:
I would like to give the user additional information next to the label in the TOC of a UniqueValueRenderer for each class. I can just add the numbers in the labeltext, but it does not align to proper columns.
Is there a way to achive an alignment or display an extra column for the classes of the UniqueValueRenderer (with arcobjects)?


texts: How to add label from header field name in ArcMap

texts:
I have some fields which contain numbers, when I added them as labels, it shows on the map as numbers. Is there any way to add header name of the field with the number?

see this http://rghost.net/8LftSZcc8/image.png