texts: Calculating sinuosity of a path with GME in ArcGIS 10?

texts:
Does anyone have any experience with pathmetrics in GME?

I would ideally like to use the old Hawths's tool "line metrics" for calculating sinuosity of a daily path but I have ArcGIS 10. I downloaded GME & R and have run some things successfully (e.g. min convex polygon), but have had lots of problems trying to run "path metrics" in GME. It basically won't run and I don't know why. I get an output that says 0 files created.

Also, the format of my 'time' column in the original xls file doesn't seem to get copied when I create a polyline - it ends up with no time codes (shows 0:00:00 in every row). this means I can't easily calculate the time intervals (and so the rate of path movement).
Is there a specific time format that is needed?

Lastly, imagining that I do get it to run, how do I calculate the path sinuosity from the segment lengths, turn angles etc..?

Thanks!!


texts: Splitting overlapping points with same symbology in ArcGIS for Desktop?

texts: I'm using a basic license of ArcGIS 10.2.2 for Desktop.  I have a point feature class layer and I uploaded it twice for representing two features. First is a type of person and second is the performances level of this person. Sometimes they are overlapped.  How can I split a bit between points of each layer so I can see all of them?  I can't use xy offset because they have not necessarily a different symbology &and I have only basic license.


texts:
I have a web ADF application which connects to a map service (either via WMS or Internet Service using http). I'd like to use DCOM connections (faster as the documentation states), so I would like to use a 'Local Service' for the map service..

but only WMS or Internet Service types work, when I choose the Local Service type the map disappears and my app is no longer able to retrieve the map

Has anyone ever experienced such an issue?


texts: Getting features attributes on feature selection in QGIS python?

texts:
I want to know how to get features from a layer on click event.. using python(m trying to create a plugin).. the below code always show 0 output

layer = iface.activeLayer()
selection = layer.selectedFeatures()
print len(selection)
for feature in selection:
    print feature['NAME']


texts: GMAP.net markers are not shown in the correct position

texts:
I have to show on a gmap control in c # more markers to make use of the code below.

The strange thing is that if I display only one marker is displayed correctly, while if I have more than one marker they are not displayed correctly ie the position is wrong ! How can I solve?

C: Code:

 private void buttonCarica_Click(object sender, EventArgs e)
        {
            if (comboBoxUtenti.Text != "")
            {
                MapControlRisorse.MapProvider = GMapProviders.GoogleMap;

                Utente utemp = new Utente(comboBoxUtenti.Text);
                MarcaTempo tm = new MarcaTempo(utemp);
                DataTable dt = tm.CaricaRisorsa(dateTimePickerDataRisorse.Value.Date);
                if (dt.Rows.Count == 0)
                {
                    MessageBox.Show("\n ATTENZIONE:  " + dateTimePickerDataRisorse.Value.Date);
                }
                else
                {
                    //int conta = 0;
                    GMapOverlay markersOverlay = new GMapOverlay("markers");
                    foreach (DataRow dr in dt.Rows)
                    {
                        Double lat = Convert.ToDouble(dr["Latitudine"].ToString());
                        Double lon = Convert.ToDouble(dr["Longitudine"].ToString());
                        MessageBox.Show("Lat: " + lat + " \n Lon: " + lon);
                        MapControlRisorse.Position = new GMap.NET.PointLatLng(lat, lon);
                        GMarkerGoogle marker = null;
                        GMarkerGoogleType markerType = new GMarkerGoogleType();
                        var point = new GMap.NET.PointLatLng(lat, lon);

                        //Verifico se l'utente ha segnato un ingresso o un uscita
                        //Se è entrato il marker sarà verde se è uscito sarà rosso
                        if (dr["Stato"].ToString() == "Ingresso")
                        {
                            markerType = GMarkerGoogleType.green;
                            // MapControlRisorse.Overlays.Add(markersOverlay);
                        }
                        else if (dr["Stato"].ToString() == "Uscita")
                        {
                            markerType = GMarkerGoogleType.red;

                            // MapControlRisorse.Overlays.Add(markersOverlay);
                        }
                        marker = new GMarkerGoogle(point, markerType);
                        marker.ToolTipText = "" + dr["Utente"].ToString();
                        markersOverlay.Markers.Add(marker);
                    }

                    MapControlRisorse.Overlays.Add(markersOverlay);
                    MapControlRisorse.MinZoom = 5;
                    MapControlRisorse.MaxZoom = 100;
                    MapControlRisorse.Zoom = 7;
                    MapControlRisorse.MouseWheelZoomEnabled = true;
                }
            }
        }

Image of Location:


texts: Mapbox-gl offline openmaptile font not shown

texts:
I'be downloaded a mbtile file related to rtl countries and used it as an offline source for Mapbox-gl server.

All stuff are shown correctly except labels that have rtl characters. Even ltr characters are shown in correct format.

Can anyone give some solution to this problem or even some technical clues?


texts:
I am relatively new to QGIS. I use conditional formatting for various problems. In the qgis query editor, I would like to define a string variable for further usage, so I wouldn't have to retype the column name countles times.

Example:
A = 'B' -- doesn't work

CASE
    WHEN A LIKE   '<1,25'  THEN  color_rgb(1,2,3)
END
-- A should be a column name


texts:
I've got a table called "MargamParkMap" - which is made up of polygons. One of the attributes in the table is "authors".

I've created a view with the aim of merging the polygons by author.
What I've got is this:

SELECT
row_number() OVER (ORDER BY "MargamParkMap".author) AS qgisid,
st_union(ST_MakeValid(ST_Multi("MargamParkMap".geom))) AS geom,
"MargamParkMap".author
FROM "MargamPark"."MargamParkMap"

GROUP BY "MargamParkMap".author;

So - row_number() is generating a new column with a unique ID, and the geom column is made up of the union of the original geom column.

If I run the query on it's own, it returns a lovely little table with all the right stuff in.

However, in QGIS - when I go to the view there are now 2 of them:

that is a really small picture...

Anyways - one has the correct geometry type, multipolygon, and SRID, the other doesn't...Is there any way of only showing the correct view?

(geometry_columns view shows the srid as 0, and type as geometry)

Also, is there a way of automatically selecting the Feature ID as the qgisid?

(PostGIS version 2.3.0, QGIS version 2.16.3)


texts:
In a GIS application I'm developing there's an option to show two types of grids. First is a standard lat/lon grid. Second is a crosshair with marks and distances in four directions away from the center.

.

Since English isn't my native language I wonder how to better call these two types in menus, etc. to distinguish between them?


texts: Error 000358 with SelectLayerByAttribute_managment

texts:
I get this error when 'arcpy.SelectLayerByAttribute_management' statement is executed.
I tried many different ways to format my whereclause. I had printed it and it looks fine.
But the error still there.

Anyone had similar problem?

fc_layer is a feature layer, that created in the code from a feature class in file geodatabase.

  Runtime error  Traceback (most recent call last):   File "",
  line 6, in    File "c:\program files
  (x86)\arcgis\desktop10.3\arcpy\arcpy\management.py", line 7221, in
  SelectLayerByAttribute     raise e ExecuteError: ERROR 000358: Invalid
  expression Failed to execute (SelectLayerByAttribute).

My code:

field1 = ['OBJECTID']
with arcpy.da.UpdateCursor('fc_layer', field1) as cursor:
     for row in cursor:
         currentID = row[0]
         whereClause1 = ''' "OBJECTID" = {0} '''.format(currentID)
         arcpy.SelectLayerByAttribute_management("fc_layer", "NEW_SELECTION", whereClause1)
         arcpy.SelectLayerByLocation_management("fc_layer", "BOUNDARY_TOUCHES", "fc_layer", "","ADD_TO_SELECTION", "NOT_INVERT")


texts: Rename images of ImageCollection to add band name at the end on Google Earth Engine

texts:
Helllo everyone, I am trying to batch download an image collection using Rodrigo E. Principe's tool but I want to download each band seperately by image.select function.
When I download a single band, I get the name of the image, but I want to add name of the band to the image id.
Example of the id I get while downloading band 'B4':
20190731T083609_20190731T084610_T36TVK
Want I want to achieve while downloading band 'B4':
20190731T083609_20190731T084610_T36TVK_B4
This is the code I use for download. Thanks in advance!

var tile = ee.Geometry.MultiPolygon((([31.8168273934012,40.644794803776],[31.8338683044504, 39.655752572233],[33.1137799145259 ,39.6615512269476],[33.1154430676193, 40.6507988109166])))

var geom= ee.Geometry.Rectangle([31.5, 40.4, 33.1, 39.4])

// Load the Sentinel-2 ImageCollection.
var sentinel2 = ee.ImageCollection("COPERNICUS/S2");

//ImageCollection filtered by date and region
    var filtered = ee.ImageCollection(sentinel2
   .filterDate("2019-07-27", "2019-08-03").filterBounds(geometry));

Map.addLayer(filtered, {bands: ['B4', 'B3', 'B2'], max: 2000}, 'filtered')
Map.centerObject(filtered)

//Clip layer
Map.addLayer(table,{},'140bin')

//clipping Function
var clipimage = function(image) {
return image.clip(table);
};

//Function to Select band
var redband =function(image){
  return image.select('B4')
}

//Select band to download
var red = filtered.map(redband)
Map.addLayer(red,{ max: 2000}, 'red')
Map.centerObject(filtered)

//how many images
var count = filtered.size()
print("Count", count)

//Single Band Batch Download
batch.Download.ImageCollection.toDrive(filtered, '140binclip',   {scale: 10, region: geom.getInfo()})```


texts:
Situation:

We have installed Google Earth on the whole computer network.
We have one KML with all our projects.

Need:

Distribute the KML on all those Google Earth instances.

Requirements:

-Distribution of the KML through a Networklink object for easy updating.

-Batch update all those computers with that networklink.

-The networklink should be placed in the my places section.

Possible solutions:

-Commandline -> as far as I know there are no commandline parameters?

-Relocate the directory where myplaces.kml is stored to a writeable central location? -> Users are users, this would quickly degrade myplaces.kml to a junkyard.

-Relocate the directory where myplaces.kml is stored to a read only central location? -> This would degrade Google Earth as a viewer? Or only as a editor in temporary places?

Are there any other options or remarks on the options researched?

Tnx


texts:
My app returns standard GeoJSON, always as a FeatureCollection.
When I display the GeoJSON in my front end, I want to zoom to show just the features in the file. I'd like to include an enclosing boundary box to make that easier:
{
    "boundaries": [12.3456, 78.9012, 34.5678, 90.1234],
    "type": "FeatureCollection",
    "features": [{
        "type": "Feature",
        "properties": {
            "name": "value"
        },
        "geometry": {
            "type": "Point",
            "coordinates": [...]
        }
    }]
}

Is it acceptable to non-standard fields like this? Is this valid GeoJSON or will third-party tools choke on it?


texts: Defining the CRS of r.watershed and other GRASS functions when used in QGIS

texts:
I am delineating some catchments from a DEM terrain file using the GRASS functions r.watershed and r.water.outlet, however, I prefer to do it in QGIS as it saves me time by not having to create a GRASS project and as I prefer the QGIS interface. I do it both through the GUI and through the Python console, however, very often the outputs are not in the same CRS of my QGIS project.
Is there any way of defining the output CRS of GRASS functions within QGIS without having to create a new GRASS Mapset?
Maybe in the advanced parameters?
Or the only way is by defining a Grass Mapset in QGIS?


texts: Editing WFS with max feature limit in QGis

texts:
I'm serving a large layer (500 thousand polygons) using Geoserver WFS and what to be able to edit it using QGis. I've set a feature limit on the layer so I wont overload GeoServer with the request.

I can successully load and edit the WFS layer in QGIS but as soon as I reach the WFS feature limit, QGis will not bring any new polygons as I pan through the map canvas, even when I hit the refresh button in QGis. I was under the impression that QGIS would keep on updating the WFS layer whenever I moved around and that the feature limit in GeoServer would be a limit to a single request and not a total request. Am I getting things wrong?


texts:
I am new to Python and am having trouble figuring out how to reclassify over 500 rasters so that any values are equal to 1 and no data is equal to 0.

Here is my code with the error code at the bottom:

import arcpy
from arcpy import env
arcpy.CheckOutExtension("Spatial")
from arcpy.sa import *

#Set environment settings
env.workspace = "E:/GIS_2017/Shapefiles"

#Set local variables
inputDir = "E:/GIS_2017/Shapefiles/TestRasters/projected"
outputDir= "E:/GIS_2017/Shapefiles/TestRasters/projected/re1"
outputsuffix = "re"

rasList = arcpy.ListRasters()
for raster in rasList:
    inRaster = raster
    remap = RemapValue([[11002,1],[11004,1],[11005,1],[11006,1],[11007,1],[11009,1],[11015,1],[11013,1],[11014,1],[11016,1],[11019,1],[11029,1],[11030,1],[11036,1],[11043,1],[11044,1],[11046,1],[11048,1],[11049,1],[11050,1],[12039,1],[12040,1],[12042,1],[12123,1],[80002,1],[80003,1],[80004,1]])
    outReclassify = Reclassify(inRaster, "Value", remap, "NODATA")
    outRasterName = '{}{}'.format(raster, outputsuffix)
    outReclassify.save(env.path.join(outputDir, outRasterName)
    print "outRasterName Reclassified Successfully!"

Parsing error SyntaxError: invalid syntax (line 21)

is there a better way to do this?
Each raster layer only has one value which is one of the codes specified in the remap.


texts: GRASS v.generalize method=douglas not working as expected

texts:
I'm frustrated with v.generalize trying to simplify polygons. Those polygons stem from ArcGIS's raster to vector function (imported into GRASS from Shapefile), so they are pretty zigzag like this:

Trying to simplify them with

v.generalize input=polygons -c method=douglas threshold=0.1 layer=-1 output=polygons_dp

Leads to a weird result like this:

Some isles and even partially borders remain untouched. Tried v.clean, v.build.polylines before simplifying, but same result.


texts: How to install wps geoserver extension to Centos6.5

texts:
I installed OpenGeo Suit on CentOS 6.5 and the version of geoserver is 2.5-SNAPSHOT.
I made a yum update and the version of OpenGeo Suit is 4.1, geoserver version remain  2.5-SNAPSHOT.

I would like to install wps extension but I am a very new centos user and don't know the path to put the extensions files.
I also want to ask about the path to put GED SDI Portal v.3 to work in my local-host.


texts: Checking feature overlap using ArcPy cursors

texts:
I have a set of geometries that I want to check whether they touch another geometry. I have created a tuple "tiffs" where I want position 0 when touches is true. tif[0] and tif1 are strings. The following code processes with no errors but no output in the print statements at the end. Confirming in the image below, there is overlap with many features and both features are in the same projection. Why am I getting no output?

list_tif = []
with arcpy.da.SearchCursor(boundary, ['SHAPE@']) as search:
    for r in search:
        for tif in tiffs:
            with arcpy.da.SearchCursor(tif[1], ['SHAPE@']) as cursor:
                for row in cursor:
                    if row[0].touches(r[0]) == True:
                        list_tif.append(tif[0])
for p in list_tif:
    print p


texts:
I am working on QGIS 3.4.9 with Mac OS Mojave 10.14.5
I have two DEM raster files and I would like to merge them. Unfortunately, an error occurs and I don't know what to do.


texts: Creating Shapely MultiPolygons from shapefile MultiPolygons

texts:
I'm reading a shapefile using Fiona:

with fiona.open('data/boroughs/boroughs_n.shp') as source:
   mpolys = [p for p in source]
candidate = polys[0]['geometry']

this gives me a dict containing a list of lists of coordinates, and type 'MultiPolygon' However, I'm not sure how to use the coordinates data to create a Shapely MultiPolygon, since it requires a sequence of exterior ring and hole list tuples: [((a1, ..., aM), [(b1, ..., bN), ...]), ...].

Is there a 'correct' way of doing this, using Fiona and Shapely?


texts: Linux alternatives to visualize and analyze LiDAR datasets?

texts:
For a long time I have been using softwares based in Windows to visualize and analyze LiDAR data sets for forest application. Recently, I have started to move all my work to Ubuntu platform, but I'm still looking for good softwares to work with LiDAR data sets.

Related to Windows the best free software experienced by me were:

Just for visualization:

a) Quick Terrain Reader is capable of opening pre-built digital elevation models (DEMs) and point clouds and allows users to freely move through the terrain in a fast and intuitive way.

b) FugroViewer is a robust, easy-to-use freeware designed to help users make the most of their geospatial data.

c) PointVue LE is a FREE 3-D LIDAR visualization tool which can be used to visualize LIDAR data in ASPRS LAS Version 1.1 format.
For analyze and visualize:

a) FUSION/LDV provides a fast, efficient and flexible access to LiDAR and terrain data sets.

Since I moved to Ubuntu, it was quite difficult to substitute the software I was used to use. Up to now, I found some alternatives but not so good as the Windows ones:

For visualization:

a) LiDAR visualization was implemented as an out-of-core multi resolution point cloud renderer. The renderer is able to visualize the largest LiDAR scans we currently have, containing up to 11.3 billion (11.3 * 109) sample points, at interactive frame rates using a fixed-size memory cache.

b) LAG is a tool for visualisation, inspection and classification of LiDAR point clouds. It currently supports LAS and ASCII file formats.

Unfortunately, I couldn't find anything for processing like FUSION.

Has somebody experienced good free software to Linux related to LiDAR files analysis?

Do not just list other software, but explain why you use and recommend it!


texts: Converting ECI coordinate system to ECEF coordinate system?

texts:
I have start and end date/time of TLE from satellite and so a matrix of ECI coordinate (122*4).

Generally, i want convert from ECI to LLA, and i know that first step is convert from ECI to ECEF, and then convert from ECEF to LLA. I implement convert from ECEF to LLA in C++.

In Matlab, there is a function to convert from ECI to LLA: eci2lla().

How to convert from ECI to ECEF with this date/time in C++?


texts: Get coordinate of second to last vertex for python ArcGIS

texts:
I need to obtain the x and y attributes of the second to last vertex for each segment and add to the attribute table using Arcpy (or some other python module).
To get the Y co-orodinate of the first vertex I can use !Shape.FirstPoint.Y!

arcpy.CalculateField_management(Flow_Paths, "Y_Cord_A", "!Shape.FirstPoint.Y!", "PYTHON_9.3", "")

In QGiS I can simply call Xat(-2) and Yat(-2) to get the second to last co-ordinate, is there an equivalent command for ArcPy?


texts: Calculate Topographic wetness index in SAGA

texts:
i am currently working on a project which involves the quantification of the wetness in certain areas. However it produces quite different results from what i expected. I'll describe the steps i have taken.
To do so i used the ASTER global DEMv2 and clipped it to the parts to be examined. I have to clip the raster, because using the original DEM results in various errors saying "memory allocation failed".

Then i calculated the slope and Catchement area via the SAGA function Catchment Area (Mass-Flux-Method, parameter: flow-spit method area cell). I also tried the Parallel and recursive method, however the results are mostly the same.
Then i calculated the Topographic Wetness Index, but the results look very mixed (somewhere between a rainbow and random fuzzy bits) and values are always negative.

I expected that the resulting raster should look like and display index values like the following image in the pdf PDF, page 14.

Has someone here experience with the quantification of topographic wetness in an area and could provide a tutorial how to do it right?
I would appreciate using SAGA or QGIS, but i also have ArcGis 9 at Hand.

By the way:
If i just use the calculate the SAGA wetness index with the clipped elevation data above, SAGA always keeps crashing (i believe because of the large no-data areas, because smaller smarts work)


texts: Linear extrapolation of points using QGIS

texts:
I am currently working with hurricane data from NOAA’s IBTrACS. Each gathered data point is separated by three hours. The points and lines can be found in two different data sets on the IBTrACS website. I am trying to add points on the path for every hour. Since storms are nonlinear and change speeds greatly over three hours, I cannot use Densify as that function adds points based on distance.
I will be assuming between each point there is a linear trend. This linear trend also applies to the change in wind speeds (which is an attribute in the data). Shown here is the buffer of the radius of the max wind speed, the data collection points, and the hurricane track, all are  of Hurricane Maria. I want to be able to show a constant change in all the features if possible (i.e. the radius of max wind speed for three-hour segment changes from t0=50kts to t3=80kts, but with the hour intervals it would have attributes of t0=50, t1=60, t2=70 and then t3=80).
I have access to QGIS 3.6. I have referenced Extrapolating a line in PostGIS but do not know how to apply it to QGIS.


texts: Calling my Python file from the QGIS Python console?

texts:
Situation: I have a Python file called 'Hello.py' and another Python file called 'File1.py'. Hello.py imports File1.py.

Hello.py code is:

print 'hello people'
import File1
File1

File1.py code is:

print 'hello world'

The files are located:

My files is located C:\Program Files\QGIS 2.18\PyPack\

Using 'OSGeo4W.bat' provided in the QGIS directory, I execute the following commands:

C:\Program Files\QGIS 2.18>cd C:\Program Files\QGIS 2.18\PyPack

python Hello.py

And the output is:

hello people
hello world

What do I type in this window to run Hello.py, which gives the output?

I have tried:

1 - dragging the 'Hello.py' file into the QGIS Python console

2 - the following command:

execfile(u'C:\Program Files\QGIS 2.18\PyPack\Hello.py'.encode('mbcs'))

3 - I tried to open 'Hello.py' in the QGIS python editor and run it with the same output.

And the outputs were:

hello people
Traceback (most recent call last):
File "<input>", line 1, in <module>
File "C:\Program Files\QGIS 2.18\PyPack\Hello.py", line 2, in <module>
import File1
File "C:/PROGRA~1/QGIS2~1.18/apps/qgis/./python\qgis\utils.py", line 607, in _import
mod = _builtin_import(name, globals, locals, fromlist, level)
ImportError: No module named File1


texts:
I'm trying to compile CGAL libraries on fedora 27
but when I compile them using cmake ,it show's the following message:

CMake Error: Unable to open check cache file for write. /home/xx/Documents/CGAL-4.11/build/CMakeFiles/cmake.check_cache
CMake Deprecation Warning at CMakeLists.txt:22 (cmake_policy):
  The OLD behavior for policy CMP0042 will be removed from a future version
  of CMake.

  The cmake-policies(7) manual explains that the OLD behaviors of all
  policies are deprecated and that a policy should be set to OLD only under
  specific short-term circumstances.  Projects should be ported to the NEW
  behavior and not rely on setting a policy to OLD.

== Setting paths ==
Build CGAL from release in directory CGAL-4.11
Packagenames: CGAL-4.11
== Setting paths (DONE) ==

== Generate version files ==
CMake Error at CMakeLists.txt:275 (message):
  File VERSION or include/CGAL/version.h are missing, required by
  release-build.

any help ,please ?


texts: Add new value combination to unique values, many fields symbology in Arc

texts:
We have a zoning layer symbolized with unique values, many fields. We use three layers and get a combination of 50 - 60 symbols, at a guess. How can I add a new symbol for a new combination of values in the three fields without undoing all the existing symbology and starting from scratch. Thanks.

ArcGIS 10.0 sp1


texts: Create bitmap from OSM data and GeoTIFF in Python

texts:
I want to create a bitmap (as a 2D numpy array) for a given satellite image in GeoTIFF format. The bitmap should have the same dimensions as the GeoTIFF and should get the value 1 at every pixel where the GeoTIFF has water.

My current approach would be:

Get the bounds of GeoTIFF.
Query the Overpass API for all features with "natural=water".
Initialise a numpy array with same dimensions as GeoTIFF.
Translate Lat/Lon from OSM features to indexes in numpy array.
Use indexes from step 4 to set entries in numpy array to 1.

My problem lies at step 5. Water areas are often closed ways which is simply a list of nodes. This means that in my numpy array I would have to identify all indices which lie in the polygon defined by the given nodes. This sounds like a problem suitable for Shapely but I looked at the documentation for Polygons and couldn't find anything.

I could loop through every possible index in my numpy array and check if that point lies in my polygon but that seems to be a bit cumbersome for me.

So my question is: Is there an easy solution to the problem in step 5?

Does this whole process seem reasonable or is there a more straightforward way to create a bitmap for the occurrence of water in a GeoTIFF?


texts: How can I display my Excel file in Arcmap 10.4?

texts:
I am currently working on a project. I have created an excel spreadsheet:

So, I need to attach the spreadsheet into ArcMap 10.4 and show the states population percentage. I have already created a floating layer. I know that I need to join the two layers, but I am not sure what to base the join on. At first I did it based on state name and nothing happen. So, I need help joining the layers and I also need help projecting what is on my spreadsheet onto arcmap once the join is done.

Do I join or relate?


texts: I have a set of points (one to many) in polygon. The middle polygon has a lot of concentrated and correlated points. Looking to see if there is way or an algorithm within QGIS to reduce the number of points, for lack of a better word reverse interpolations.


texts:
I use GeoServer war 2.4.5 with Tomcat. When I try to add a WMS data source  from add store, it gives:

"Connection test failed: sun.security.validator.ValidatorException: PKIX path building failed: sun.security.provider.certpath.SunCertPathBuilderException: unable to find valid certification path to requested target" and wms Capabilities URL is "https://"

How can I fix this problem?


texts: Importing .bag (bathymetry) files into QGIS 3.x

texts:
I downloaded bathymetric data from the NOAA Bathymetric Data Viewer website: https://maps.ngdc.noaa.gov/viewers/bathymetry/
I've tried to import and view several .bag files from this website within QGIS with no luck. I have tried versions 3.10, 3.12, and 3.14. The files I tried to view were downloaded from: https://www.ngdc.noaa.gov/nos/H12001-H14000/H12053.html
I was able to import and view these files just fine in QGIS on a Mac. I, however, wasn't able to on a PC with a Windows 10 OS. I read some github forums where this seems to be a Windows issue: https://github.com/OSGeo/gdal/issues/1428
I was wondering if anyone knew a workaround to import the .bag data into QGIS to either fix the gdal issue if it is one or import the data using the python console, etc.
Here is what I see in QGIS:

It loads up the data in the correct geographic loaction and lets me change the symbology but nothing is viewable.


texts: Determining accessibility of U.S. street address

texts:
I'd like to determine whether a given U.S. street address is publicly accessible to vehicle/foot traffic. This would include addresses on both publicly and privately maintained roads that are not within a gated community or another access-restricted area.

For instance, 8541 Lake Windham Ave, Orlando, FL 32829 is in a gated community. The openstreetmap tag for this portion of the road is tagged as access=private. However, since OSM data does not include address ranges for each segment (AFAIK), I'm not sure how to associate the address with the private portion of the road.

The U.S. Census Tiger data does include address ranges. However, as far as I can tell from the documentation, there's nothing to indicate accessibility of the segment.

Is there a freely-available way to determine whether a U.S. street address is on a publicly accessible road?


texts:
The vector data displays, but its down near null island, I have been trying this for a couple hours now and nothing I try has worked.
Can you see why this code refuses to place the data on the right location?

Source is data EPSG:27700 var proj is EPSG:3857
I have swapped dataProjection & featureProjection which has zero effect.

var vectorSource = new ol.source.Vector({
    format: new ol.format.GeoJSON(),
    loader: function(extent, resolution, projection) {
        var proj = projection.getCode();
        console.log(proj);
        var url = "site.com?typeName=my:layer=WFS&version=1.0.0&request=GetFeature&outputFormat=application%2Fjson&bbox=" + extent.join() + ",EPSG:3857"
        var xhr = new XMLHttpRequest();
        xhr.open('GET', url);
        var onError = function() {
            vectorSource.removeLoadedExtent(extent);
        }
        xhr.onerror = onError;
        xhr.onload = function() {
            if (xhr.status == 200) {
                vectorSource.addFeatures(
                    vectorSource.getFormat().readFeatures(xhr.responseText,{
                        dataProjection: 'EPSG:27700',
                        featureProjection: 'EPSG:3857'
                    }));
            } else {
                onError();
            }
        }
        xhr.send();
    },
    strategy: ol.loadingstrategy.bbox
});


texts: Shortest 3D distance (depth, z) between an underground point and the surface (dtm)

texts: I have a point P1 (see sketch) which is located underneath the surface (DTM). My aim is 1) to find the shortest 3D distance between P1 and the surface (d1 in sketch) and 2) the surface location (P2 in sketch) where the shortest 3D distance leads to. How can I accomplish this using QGIS or Grass or a similiar spatial open source GIS tool?


texts: Fitting circles to crescent shaped clusters in raster data automatically

texts:
I have mobile lidar scanner (MLS) that was collected in a forest with big, mature trees.

My research involves fitting circles to tree stems to estimate their diameters at different heights.

Once I normalized the cloud and segmented out the trees, I sliced the point cloud into horizontal sections (like floors in a building), then converted them into geotiffs. The results look something like this in a good case:

As it stands, I'm using the 'draw circle from 3 points tool' in qgis shape digitizing toolbar to manually fit circles to each cluster. The result looks something like this:

Is there any intuitive way to automate this process? I find that doing it manually it is difficult to be consistent in addition to being extremely tedious. It would save me tons of time and effort if this could be avoided!

I think there may be a way to accomplish this via something like 'contour detection' in opencv, but I dont know how to go about integrating it with qgis.


texts: How to read MapInfo .tab files in ArcPy?

texts:
I want to use arcpy (the ESRI install of python 2.65) to analyse MapInfo TAB files...can this be done directly or do I need to convert them to shp files...if so would it be possible to do something like

tabs=glob.glob('*.tab') and then
for tab in tabs
 some function to convert tab to shp tab_shp1.shp

  (maybe ogr2ogr -f "ESRI Shapefile" mydata.shp mydata.tab) not sure how
  to call an external command line.

 shpObj=arcpy.Describe(tab_shp1)
 etc

I know there is functionality in FWTOOLS via OGR but can I access this via the ESRI python 2.65 install?

We have both tab and shp data and I am trying to create a master spatial index of all the datasets we have. I have written the code for ESRI formats but want to include the tab's if possible.


texts: Layer Style - save location in mssql-server

texts:
Where is the save location of a QGIS-layer-style in a mssql-server?

I am user in a larger mssql-server-database with read and write permissions in a specific schema. On saving the style for a layer in the mssql-database, I get this message:

"Unable to save layer style. It's not possible to create the destination table on the database. Maybe this is due to table permissions. Please contact your database admin"

Will the destination table for the layer style be saved in my user-schema or in a different place (e.g. dbo)?


texts: how to increase geowebcache performance with openlayers

texts:
I have installed geowebcache and cached polyline data. And I am using my geowebcache data on OpenLayers map.

        var polyline = new OpenLayers.Layer.WMS(
            "Lines",
            "http://localhost:8080/geowebcache/service/wms",
            { layers: "GIS:pipes", transparent: true, TILED: 'true', gridSet: "EPSG:94000" },
            { isBaseLayer: false}
        );

When I start the application my map tiles(using TMS) are coming from remote server fast. But my cache tiles(local area) are coming after my map tiles.

When I zoom in map, my wms polyline is appearing blur about 4-5 seconds then appearing clear.

I think it should be very fast because of local area and geowebcache. My server upload speed is 4MBit. And my servers are linuc OS on virtual box.

Is there any point about performance that I forget? Do you think?


texts: Understanding CadNSDI versus GCDB PLSS datasets difference and accuracy?

texts:
I recently submitted some coordinates I acquired via ArcGIS to the BLM to go out into the field and pinpoint some locations.

When I did the work at the computer I used ArcGIS and leveraged the BLM's CadNSDI dataset to pinpoint PLSS section lines.  When I heard back from the BLM, they claimed that my points were not accurate to footage calls in the field.

I started doing some digging and realized that there are two standards utilized by the BLM for PLSS data.  These two datasets are the CadNSDI(which I used) and the GCDB, which based on the information at the following link, seems to be an older standard:

https://sites.google.com/site/1corpsofdiscovery/bureau-of-land-management/blm-gcdb

Why do these two datasets not line up, and which one is more accurate for comparing footage calls from section lines in the field?

Both of the datasets are in NAD83 but do not line up with errors in the range of 20 feet(I am projecting into a State Plane projection). I can't seem to wrap my head around what is happening here.  From my perspective it seems like the BLM serves out two datasets that don't line up and then when you submit coordinates to them they say they are inaccurate when using the newer standard.


texts: How to calculate most efficient route passing through all the houses in the world?

texts:
I'm new to GIS.

I need some help in determining the best or most efficient route, using a flying sleigh, through all the houses in the world. One of my coworkers told me that this site would be the best place to ask, because I would find lots of helpful GIS experts.

I will need some guidance about what software to use, where to get the data, and how to process it. Since I had some extra expenses this month, I would prefer some Open Source solutions.

Thank you all!

PS: I'm in a bit of a hurry, as I need this for tomorrow!


texts:
I want to list all cities in lubelskie vovoidenship, which is relation 130919.
I also want to get latLong data for each city. To perform this I do:

SELECT
    name, ST_AsLatLonText(way), place
FROM
    planet_osm_point
WHERE
    place IN ('city', 'town', 'village', 'hamlet')
    AND ST_Contains((SELECT way FROM planet_osm_polygon WHERE osm_id = -130919), way)

Query runs ok, but problem lays within the results.
For example:

...
"Lublin";"15°1'48.000"S 45°23'24.000"E";"city"
...

The coordinates for every city are wrong. What should I do?

UPDATE

Data is imported using osm2pgsql:

osm2pgsql -U postgres -W -d gis2 --create ~/Pobrane/poland-latest.osm

Data downloaded from: http://download.geofabrik.de/europe/poland.html

Result using ST_AsEWKT():

"Lublin";"SRID=900913;POINT(2512485.39 6665744.97)";"city"


texts: Filter geojson file to Get last 24 Hours Data

texts:
i am using mapbox platform to showing the map in my site.

geojson is having time attribute
Ex: "time":1510301955145
I want to show the map of last 24 hours.


texts:
By mistake, I've dragged the digitising toolbar partly offscreen in QGIS 16 (Hanover), and can no longer grab the handle to fetch it back. There seems to be no way in the menus of dragging it back into position using the customisation/menu options. Can anyone help/advise please. I've already tried uninstalling and reloading QGIS, but the toolbar is still overlaps the edge of the screen!
Many thanks!!


texts: Problem with publishing GIS resources in ArcGIS Server

texts: I have ArcGIS 9.3 on Win 7 sp1 integrated with VB.Net 2008. I have connected my server in ArcCatalog and publish mxd to ArcGIS Server without any problem, put after I open ArcGIS Server Manager and publish GIS services to GIS resources it gives me this error:


texts: Merging more than 9 input files with lasmerge?

texts:
Is there a possibility to feed more than 9 LAS files to lasmerge tool? The user interface only allows 9 LAS files from LAStools toolboxes in QGIS or ArcGIS, but at liblas.org it’s described as there can also be a text file holding the file names to be merged.

Unfortunately it’s not said how to format the text file and my attempts in Python with the absolute path of each file per line fails.


texts: How does the sDNA calculate betweenness?

texts:
I want to calculate the betweenness value by my hand, and then compare with sDNA calculated. But I found that I can't get the value same as what sdna calculated.

Here's my case:

Assume all the links length of this network is 500 meters. While calcualte the betweenness value for the link A with raidus 600 meters, totally 7 links should be considered(please refer the bold links in the picture).

I calcuate the betweenness as the following:

betweeness value of A = 1/2*6+9+1/3 = 12.3333

But the sDNA calculate the betweeness 6.3333.   Please refer the following sDNA picture.

How does the sDNA get that result and why?

My test environment:
- ArcGIS Desktop 10.2.2
- sDNA 3.4.5


texts: Plotting State Border without Water using R and Raster?

texts:
Suppose I want to plot New York or Ohio using R:

library('ggplot2')
library('map')
library('raster')

map <- raster::getData("GADM",country='USA',level=1)
map <- map[map$NAME_1 %in% c('Ohio'),]
map <- map_data(map)

ggplot() +
  theme_void() +
  geom_polygon(data = map,
           aes(x=long, y = lat, group = group),
           fill="grey", alpha=0.3)

Here are the results:

How can I not display extra land in the Ohio plot, and only display land for New York?

It looks like it plots the state borders, see google maps for comparison.


texts: How to preserve original direction of NHD segments when creating route using gp.CreateRoutes_lr tool?

texts:
I am trying to figure out how to preserve the original line direction of NHD High Resolution data that I am selecting out and creating a route with. I'm using ArcGIS 9.3.1 and this is a python script.

I am selecting out some of the NHD streams by the ReachCode field into a new feature class (I have also tried this with a feature layer in the context of the issue I'm posting about but creating a separate feature class allowed me to debug further), and then am taking that feature class and dissolving it. I then take that dissolved NHD section and am creating a route with the gp.CreateRoutes_lr tool.

I've narrowed down the issue to the step when the script creates the route and it seems that the direction of the line feature gets flipped when I create the route (I have definitively verified this is happening). Is this what is supposed to happen? I've never worked with the linear referencing tools before this so am unsure if this is what should be happening.

Each of the line segments have a defined FROM and TO node direction which is important for streams in terms of correct flow direction. The reason the direction matters to me for my script is that I am creating the vertices via the gp.FeatureVerticesToPoints_management tool get the elevations of the FROM and TO node to determine the stream slope for points located along the NHD section. This impacts the slope calculation obviously.

I just compared a site that had a negative slope (with assumption that a negative slope implies an instance where the route direction got flipped as I describe) with a site that had a positive slope (with assumption that a positive slope indicates an instance where the route direction did not get flipped). Both assumptions were true, and it seems for the former the direction got flipped, and the latter it did not. So it seems that this tool does not flip the directions consistently...

Is there any way to preserve the original direction of the NHD segments when I create a route? Or, what would influence the direction of the output route line feature?


texts: Allocation error by deleting symbol layer

texts:
In my plugin I have a loop where a Symbol is created in first run through.
During next steps the symbollayerlist should be filled.

In creating with standard constructor always a symbollayer is added. So I tried to delete the symbollayer in first time. But this causes a memory allocation error.

Then I tried to give the constructor an empty layerlist created by myself, so should I can delete the layers by myself, but also the constructor addes automatically an layer and I get the error by deleting it.

Here is a snippet of my code:

QgsSymbolV2 * p_Symbol(0);
bool oFirst=true;

...

if (oFirst)
{
        QgsSymbolLayerV2List * sLyrs= new QgsSymbolLayerV2List();
        qDebug()<<QString::number(sLyrs->count());
        QgsFillSymbolV2* p_FSymbol = new QgsFillSymbolV2(*sLyrs);
        qDebug()<<QString::number(p_FSymbol->symbolLayerCount());
        p_Symbol = dynamic_cast <QgsSymbolV2 *>(p_FSymbol);
        qDebug()<<QString::number((int)p_Symbol,16);
        qDebug()<<QString::number(p_Symbol->symbolLayerCount());
        if (p_Symbol->symbolLayerCount()==1)
                p_Symbol->deleteSymbolLayer(0);
        oFirst=false;
}
...

Output:

"0"
"1"
"1847bca0"
"1"
First-chance exception at 0x77e4e43e in qgis-bin.exe: 0xC0000005: Access violation reading location 0x4d46be13.
Unhandled exception at 0x77e4e43e in qgis-bin.exe: 0xC0000005: Access violation reading location 0x4d46be13.

What could be the reason?
I use QGIS 2.18.11 and Visual Studio 2010 on Win7.
In prior Version QGIS 2.14.3 and Visual Studio 2008 this problem doesn't occur.


texts: Error 000358 or Error 999999 invalid SQL

texts:
i've been searching arround all the answers in the questions of Error 358 and beyond, but i doesn't fix my problems, i tried every post and so on the problems stills there.

I make two options in my code, one make the error 358 appears, and the other option make the error 999999 with invalid SQL. If i put the code "by hand" in the SelectLayerByAttribute_management it works, but i need to be an input by the user. The query its the same by hand and by the script and i dont know how to make it run.

Sorry for the english, not my maternal lenguaje.

The code is:

Zona_trab=raw_input("Introduzca la ruta de la GDB/carpeta sobre el que se trabajara")
arcpy.env.workspace = Zona_trab
CapaObj=raw_input("Introduzca el nombre de la capa de tuberias")
CapaArq=raw_input("Introduzca el nombre de la capa de arquetas")
puntos = raw_input("Introduzca las arquetas con datos de permalogs, separadas por espacios: ")
query = ''
first = True
for punto in puntos.split():
    if not first:
        query = query + 'OR '
    query = query + '\\"FID\\" = ' + punto + ' '
    first = False
query = '" '+query+' "'
TempCarp="C:\\Trabajo Carto\\TempPrograma"
arcpy.MakeFeatureLayer_management(CapaObj,"capa tuberias")
arcpy.MakeFeatureLayer_management(CapaArq, "capa puntos")
arcpy.SelectLayerByAttribute_management("capa puntos","NEW_SELECTION", query)
arcpy.SelectLayerByLocation_management("capa tuberias", "WITHIN_A_DISTANCE","capa puntos","1 Meters", "NEW_SELECTION")
arcpy.FeatureClassToFeatureClass_conversion("capa tuberias",TempCarp,"Final.shp","","","")

The second option, wich give me a error 999999 its deleting the arcpy.SelectLayerByAttribute_management("capa puntos","NEW_SELECTION", query) and making the selection in the arcpy.MakeFeatureLayer_management(CapaArq, "capa puntos", query) this makes the error goes in arcpy.SelectLayerByLocation_management("capa tuberias", "WITHIN_A_DISTANCE","capa puntos","1 Meters", "NEW_SELECTION").

I hope get some response because i'm stuck in there :(


texts: Generate custom WKT Coordinate System definition from local / WGS84 pairs of coordinates

texts:
I am wondering if it possible to create a custom WKT for coordinate transformations based on pairs of WGS84 lat/lon and the corresponding local N/E values. So like a minimum 2 pairs?

If so is it possible to do with a tool (like pyproj/PROJ4/GeoPandas) and not a graphical tool.


texts: Convert an imported land surface temperature feature collection to an image


texts: Selecting nearest 300 buildings (polygons/centroids) from list of points using ArcGIS for Desktop?

texts:
I have been trying to look for an easy/quick way to select a certain amount of polygons (layer 1) from a 90 points shapefile (layer 2).

The idea is to identify the nearest 300 buildings from each of the 90 vaccination sites. Selecting these buildings based on a distance is easy but based on the number of features around it seems much more complicated.

I could not find any good way (I tried "Near", "Make Closest Facility", "Neighborhood location"...but they don't seem to do what I'm looking for).

I can't find any tool doing this. Doing this for only 1 point is ok but for a larger number of points is not obvious.

Any suggestion?

I'm using ArcGIS 10.3 for Desktop.


texts: Copy and paste feature using ArcGIS Online editing?

texts: Wondering if there is a way to copy and paste a feature (geometry and attributes) using ArcGIS online editing.  I have searched the internet a bit and have been unable to find any way to do it.


texts: PostGIS point table with dynamic content needs to have schema change: keep the change in 1 table OR separate them?

texts:
I have a web app with PostGIS point table as its main back-end data. Through this app, points can be :

inserted (by clicking map AND filling form)
updated (by moving the point location AND/OR by filling form)
deleted
shown by LeafLet.JS

Currently there's a need to alter the point table schema from the existing ('old') schema to the 'new' schema. The main motive is that I am gonna introduce several new mandatory fields. Accordingly the web app needs to be changed to cover this.

Upon the new release of the app:

New data goes in through a new form that complies to the new schema
Updating old data also through the new form (that complies to the new schema)
Hence the old form is no more
Both 'old' and 'new' points got rendered to the same Leaflet.js map (with ST_Intersection operator)

There's 2 option on how the change is implemented in DB:

A) Keep both 'schema' in the same table

B) Store them in different table

In details:

A) Keep both 'schema' in the same table, only distinguished by a flag (field) with 2 possible values. Rows existed before the change have all new fields NULL. When old rows updated they gonna get their flag switched from 'old' to 'new'. New rows will always get flag that represent 'new'.

B) Store them in different table. Here there's no need to have flagging field. The 'old' table will have zero-growth. When updated, rows from 'old' table got moved to the 'new' table. New data goes to 'new' table.

Plan B is preferable for its clear distinction between the 2 schema.

Plan A is -i think- more friendly for the ST_Intersection (as i think that the Intersecting Leaflet BBOX against 2 table costs alot).

Now:

Do I put all consideration correct?
Which one is better: A vs B ?


texts:
I'm looking for a query that returns overlapping areas of two tables with polygon geometries. Using row_number() over() the entries of the gid column are not unique in the case of dumped multipolygons. So I get some troubles when loading the layer into QGIS.

CREATE TABLE poly_intersection AS
SELECT
row_number() over() AS gid,
(ST_Dump(ST_Intersection(poly1.geom, poly2.geom))).geom::geometry(polygon, SRID) AS geom
 FROM poly1, poly2
 WHERE ST_Overlaps(poly1.geom, poly2.geom) = 't' OR ST_Contains(poly1.geom, poly2.geom) = 't' OR ST_Contains(poly2.geom, poly1.geom) = 't';

Meanwhile I've modified my query. But I wonder if there is a simpler solution.

CREATE TABLE poly_intersection AS
SELECT
row_number() over() AS gid,
g.geom::geometry(polygon, SRID) AS geom
FROM
(SELECT
(ST_Dump(ST_Intersection(poly1.geom, poly2.geom))).geom
FROM
poly1, poly2
WHERE
ST_Overlaps(poly1.geom, poly2.geom) = 't' OR ST_Contains(poly1.geom, poly2.geom) = 't' OR ST_Contains(poly1.geom, poly2.geom) = 't')
AS g;


texts: Calculate the area in sqkm across UTM zones using arcpy cursors

texts:
I am trying to update the area field in sqkm of my table using Search and Update cursor from the arcpy library. The area needs to be updated using the different utm zone. I intersected my shp (that I called "outIntersectPath") with the utm zone shp. Now my intersected shp has a column with the UTM ZONE that goes from 1 to 60 (integer) and another column NTH_STH that represents the Northern ("N") and the Southern ("S") hemisphere.

This is the code:

def utmConversion(outIntersectPath):
 # calculate Area per UTM Zones
 count = 1
 utmNorthZone = 32600
 utmSouthZone = 32700
 field = ["SHAPE@AREA", "BASIN_AREA_SQKM", "UTM_ZONES", "NTH_STH"]
 arcpy.AddField_management(outIntersectPath, "BASIN_AREA_SQKM", "DOUBLE")
 print "Field {0} has been added to the {1}".format(field[1], outIntersectPath)

 while count < 61:

    with arcpy.da.SearchCursor(outIntersectPath, field) as cursor:
        for row in cursor:
            utmZone = row[2]
            hemisphere = row[3]
            if utmZone == count and hemisphere == "N":
                utmNorthZone += 1
                count += 1
                # print utmZone
                # print utmNorthZone
                out_coordinate_system = arcpy.SpatialReference(utmNorthZone)
                print out_coordinate_system.name

                with arcpy.da.UpdateCursor(outIntersectPath, field, "", out_coordinate_system) as ncursor:
                    for nrow in ncursor:
                        nrow[1] = nrow[0] * 10 ** 6
                        ncursor.updateRow(nrow)
                print "Field {0} has been updated to {1}.".format(field[1], outIntersectPath)
            elif utmZone == count and hemisphere == "S":
                utmSouthZone += 1
                count += 1
                out_coordinate_system = arcpy.SpatialReference(utmSouthZone)
                print out_coordinate_system.name
                with arcpy.da.UpdateCursor(outIntersectPath, field, "", out_coordinate_system) as scursor:
                    for srow in scursor:
                        srow[1] = srow[0] * 10 ** 6
                        scursor.updateRow(srow)
                print "Field {0} has been added to {1}.".format(field[1], outIntersectPath)
            else:
                print "The Process failed..."

What I am doing wrong?.. it is not calculating the area per each UTM Zone.


texts:
I have a polygon file, taken from an Oracle database, that looks something like this:

   ID FMC_AREA_TYPE_NAME DESCRIPTION
1 162            HARBOUR     Esbjerg
2 164            HARBOUR    Thyborøn
                                                                                POLYGON
1 MDSYS.SDO_GEOMETRY(2003, 8307, NULL, MDSYS.SDO_ELEM_INFO_ARRAY(1, 1003, 1), MDSYS.SDO_ORDINATE_ARRAY(8.38, 55.46, 8.45, 55.46, 8.45, 55.49, 8.38, 55.49, 8.38, 55.46))
2    MDSYS.SDO_GEOMETRY(2003, 8307, NULL, MDSYS.SDO_ELEM_INFO_ARRAY(1, 1003, 1), MDSYS.SDO_ORDINATE_ARRAY(8.2, 56.68, 8.26, 56.68, 8.26, 56.72, 8.2, 56.72, 8.2, 56.68))
                                                                                POLYGON_AS_CLOB     NAME
1 8.380000 55.460000 \n8.450000 55.460000 \n8.450000 55.490000 \n8.380000 55.490000 \n8.380000 55.460000 \n  Esbjerg
2 8.200000 56.680000 \n8.260000 56.680000 \n8.260000 56.720000 \n8.200000 56.720000 \n8.200000 56.680000 \n Thyborøn

How can i make that into polygons, preferably using R. I have never worked with Oracle before. I tried using the sf package.

st_as_sfc(h$POLYGON_AS_CLOB)

But it seems that sf it doesnt support that geometry type:

OGR: Unsupported geometry type

#

Used the answer from Spaceman, and got this. Not pretty, but it works.

clob <- as.character(h$POLYGON_AS_CLOB)

mp <- list()
for(y in 1:nrow(h)){
x <- st_polygon(list(matrix(as.numeric(unlist(strsplit(strsplit(clob,"\n")[[y]]," "))),ncol=2, byrow=TRUE)))

mp = list(x, unlist(mp,recursive=FALSE))
}

mp1 = st_multipolygon(mp)

xx <- as(mp1, "Spatial")

proj4string(xx) = "+proj=longlat +ellps=WGS84"

mapview(xx)

havne <- st_as_sf(xx)

havne$SI_HARB <- 1

st_write(havne, "Q:/Havne.shp", delete_layer = T)


texts:
I created my own marker symbols (arrows) for point features. Now I want to rotate them with the Rotation-Tool from the Editor-Toolbar, but the Arrows always snap back to their original horizontal position.

What am I missing?

To rotate the arrow i proceeded as follows:

start editing session and select the feature layer containing the arrow symbols
select the arrow-symbol i want to rotate
select the rotation-tool
select the arrow and drag it in the wanted position

--> and then it just snaps back to its original position...

News: When I select rotation for onw of the attribute table columns in the layer properties, rotating seems to work but nuw the symbol just disappears!


texts: ogr2ogr update failing to update shapefile int PostGIS

texts:
I have a bunch of shapefiles that I am looping in with a bash script to load and update files but the ogr2ogr fails. I am running gdal 1.92 on ubuntu 14.04.

for file in `ls *.shp`;do ogr2ogr -progress -append -update -skipfailures -a_srs "EPSG:4326" -nlt MULTIPOLYGON -f "PostgreSQL" PG:"dbname=data_test" $file -sql "SELECT gid,name,id from ${file%.*}" -lco SCHEMA=test;done

  ERROR 1: Layer test.boundary already exists, CreateLayer failed.Use the layer creation option OVERWRITE=YES to replace it.


texts: What's an analogy for a Network Link (KML)?

texts:
I'm trying to communicate how a product could/should/does behave with Network Links but am often getting tripped up by explaining the difference between the "link" and the "fetched content" -- which is a very important distinction when it comes to saving KML/KMZ, for example.

I'm aware of the shortcut to target file analogy, but I'm looking for something for more of the non-technical folks.

Any ideas of a real-world (preferrably non-technical) analogy to Network Links?


texts: Calculate sum of one attribute of a point shapefile that distributed in polygon shapefile

texts:
I have two shapefiles, one of them is polygon and other is the point. Each of points have attributes that one of them is population. I want to calculate population size in all of polygons.

How can I do that?


texts: QGIS connection to MSSQL Server issue - no geometry columns table found

texts:
I know that this question has been asked before but I believe my case is a little different... We do have a geometry columns table as well as a spatial references table in our spatial database, however I still have following error message.

Like I said, we have a geometry table and a spatial references table for QGIS to look at.

When I add a new connection, ONLY these newly created feature classes in our database are recognised as spatial data. All the old entries are recognised as non-spatial tables.

Following is the details of above two tables in geometry columns table

I am puzzled. Anyone could shed some light please?


texts: Polygon not matching world street map in ArcMap

texts: I have a polygon layer that is not matching the World Street map by a few hundred meters in ArcMap. How can I correct this? I have tried projecting the polygon shapefile to WGS 1984 and WGS_1984_Web_Mercator Auxiliary Sphere but with both, I get: ERROR 000151: The geotransformation is not valid.


texts: How to iterate all input folders (containing coverages) using Python with ArcGIS 9.3?

texts:
I'm trying to write python scripts to batch convert thousands of coverages to shapefiles in arcgis 9.3. The coverages files are in different folders:

d:/area_A/building
d:/area_B/building
d:/area_C/building
...

and the output file will be merged or unioned into a single shapefile:

d:/output/building_merged.shp

I'm stuck in iterating all the input folders..here is what I've got so far:

import arcgisscripting
gp = arcgisscripting.create()

try:
    gp.Workspace = "d:/area_A/building"
    gp.FeatureClassToShapefile("polygon","d:/output")
except:
    print gp.GetMessages()


texts: Different number of observations when using SampleRegions

texts:
I have a code in GEE which creates two images from sentinel-1 (SAR) and from sentinel-2 (NDVI) , creates from them paired image and then creates scatter plot from the two bands in the image. In order to create the scatterplot, I know I need to run kind of training and i'm a bit confused about it and about my results.

This is the code I have that seemed to be working:

// make an image for the two variables
var pairedImage =  ee.ImageCollection.fromImages([imageNDVIcor,SARreproject]).toBands().rename(["NDVI","SAR"]);
print("pairedImage",pairedImage);
// Generate a sample of points within the region
var sample = pairedImage.sampleRegions(geometry, null, 150);
print("sample",sample);
// Generate chart from sample
var chart = ui.Chart.feature.byFeature(sample, 'NDVI', 'SAR')
    .setChartType('ScatterChart');
print("chart",chart);

The point is that I get results that look nice but I can't understand, for example:
this is suppose to be my point? what is this? what are those values?

var sample = pairedImage.sampleRegions(geometry, null, 500);

and this is with:
var sample = pairedImage.sampleRegions(geometry, null, 150);:

and this is when it's 1000:

How come using less points lead to more ponts in the correlation? I feel confused and i'm afraid that I don't understand something very basic with the sample points.

My end goal is to understand why I have gotten those two different charts when 500 has lesspoints than 100  and 1000 had the least. what'shapenning here?


texts: Creating temporary buffer layer in QGIS?

texts:
I have been able to select features using an attribute from a Point vector file:

expr = QgsExpression("\"NAME\" = 'PESCOE'")
collegePoint = collegeLayer.getFeatures( QgsFeatureRequest (expr))

#get ids
ids = [i.id() for i in collegePoint]

#select the features
collegeLayer.setSelectedFeatures(ids)

Now for each of the selected features I wish to create a temporary buffer layer (not as presented in Buffering in pyQGIS?, as it is not quite clear).

How do I achieve this?


texts: Extract road defects from terrestrial LiDAR, in Erdas 2012 or ArcGIS10.1

texts:
I'd like to extract road defects from a lidar point cloud. So, wherever a road surface drops or raises by a few mm relative to 100 cm around it.

I know it's possible but I just can't find or figure out a workflow for it.

Erdas Webinars has the workflow for version 2011, but this is very different in 2013 as it uses a raster .img file instead. I also looked at the workflow files from Imagine Help.

In ESRI it's a similar issue.

Any advise on workflows for this? We own both software so I would like to manually extract defects from a terrestrial LiDAR scan using either of these.


texts: How to make change to layer in ArcMap with python/ArcObjects?

texts:
In ArcMap 10.1 from the python window I am attempting to rotate a symbol for a point.  After I run my code I can look at the symbol properties and the rotation expression is set on the symbol, however the symbol in the data frame is not rotated.

What am I missing to get the points in the data frame to show the point rotated?

(also worth noting, if I click on "OK" in the expression builder and rotate windows, then click "apply" in the layer properties window, then "OK" the symbols are then displayed in the data frame properly)

import snippets101 as snippets
snippets.GetStandaloneModules()
snippets.GetDesktopModules()
import comtypes.gen.esriCarto as esriCarto
import comtypes.gen.esriFramework as esriFramework
import comtypes.gen.esriArcMapUI as esriArcMapUI
import comtypes.gen.stdole as stdole
import comtypes.gen.esriDisplay as esriDisplay
from snippets101 import CType
from snippets101 import NewObj

pApp = snippets.GetCurrentApp()
pMxDoc = CType(pApp.Document, esriArcMapUI.IMxDocument)
pMap = pMxDoc.FocusMap
pAV = CType(pMap, esriCarto.IActiveView)

pFeatureLayer = CType(pMxDoc.FocusMap.Layer(0), esriCarto.IFeatureLayer)
pGeoFeatureLayer = CType(pFeatureLayer, esriCarto.IGeoFeatureLayer)
pRotationRenderer = CType(pGeoFeatureLayer.Renderer,esriCarto.IRotationRenderer2)
pRotationRenderer.SetSymbolRotation3DFlags = (4)
pRotationRenderer.SetSymbolRotation3DExpressions('','','[SDO_ORIENTATION] * -1')
pRotationRenderer.SymbolRotation3DRotationTypeZ = 0
pAV.refresh()


texts: Identifying fewest locations from which to survey patch of habitat using QGIS viewshed?

texts:
We have a species which we know lives in a specific habitat. We have mapped that habitat across the landscape. We plan to put a VHF collar (radio-tracking collar) on some individuals. These collars work best if they have line of sight to a handheld receiver (up to a maximum of c.5km). We plan to go to certain locations in the landscape to scan for the collared individuals.

These scanning locations are usually picked intuitively, by going to the highest points that are easily accessible (e.g. on a road). It is a very time consuming process, and this seems pretty inefficient.

My method went like this: I am using QGIS Viewshed Analysis. I turned a map of roads into nodes. I then performed a viewshed analysis of those few thousand nodes, and overlaid it with the habitat map. This gives me the total quantity of habitat which is viewable from a road.

Is it possible to identify the amount of habitat which is viewable from each location/node individually?
Can we then identify the locations which see the largest amount of habitat with the smallest overlap?


texts: How to combine neighboring parcels to find maximum stream length per parcel combination - ArcMap

texts: I am typically assigned with the task of finding parcels in a watershed which have the highest linear footage of stream running through them. For this I can simply intersect my stream layer with my parcels layer and dissolve by parcel ID to generate a table with stream length by individual parcel. What I'd like to do is create a model or script that will go through every parcel in the dataset, and connect the (let's say 5 for example) most adjacent parcels, and compute linear feet stream for each combination of 5 adjacent parcels. I am stumped on how to go through each parcel and create a feature for the closet 5 parcels for each parcel in the dataset. FYI these parcel datasets typically consist of tens of thousands of parcels. Any ideas on how to go about this?


texts:
I have a bunch of sensors at fixed locations. I want to push this data into QGIS through an API. I am not able to find a way to push data through API.

Looking for some help or examples preferably in Python.


texts: Extrapolating surface in ArcGIS Desktop?

texts:
I'm trying to calculate the volume of a landform under water in ArcGIS and here's what I decided to do:

Convert the raster of the bathymetry that I had to a gridded point dataset using Raster to Point.
Clip out the points contained within the landform (leaving me with a somewhat large hole). I actually have multiple landforms, some creating smaller holes, but the biggest with a surface area of over 300,000 sq m.
Re-interpolate the surface sans landform.
Use the Cut and Fill tool to compare the original surface to the interpolated surface (without the landform) in order to get the volume.

Makes sense to me --- the trouble is, I've been trying to find the right interpolation method to capture a surface sans landform. I think my hole was way too big for IDW, because as soon as the contours hit the re-interpolated hole they shot off in funky directions. I did a TIN that did fairly well, but I'm not sure how accurate or believable a TIN can be. Lastly I tried kriging, but processing was taking in excess of 24 hours, so I cancelled the process.

Can Anybody give me some tips on the best interpolation method for this task in Arc?


texts: arcpy custom script tool to find and replace floats and ints?

texts:
I created a custom script tool that can find and replace text. See the image link for what my tool interface looks like too, and I'll point out the find and replace fields are both set to string in their data types. I'm really not experienced in custom script tools, and would like to know how to make some "buttons" or checkboxes of sorts to choose between finding and replacing strings or floats or integers. Would that be a tool validation-side issue or a script-side issue?

import arcpy

FC = arcpy.GetParameterAsText(0)
oldText = arcpy.GetParameterAsText(1)
replaceText = arcpy.GetParameterAsText(2)
targetField = arcpy.GetParameterAsText(3)

queryString = targetField + ' = ' + "'" + oldText + "'"

try:

    with arcpy.da.UpdateCursor(FC, (targetField), queryString) as cursor:
        for row in cursor:
            if row[0] == oldText:
            row[0] = replaceText
            cursor.updateRow(row)

            arcpy.AddMessage("String where " + queryString + " found and replaced with " + replaceText)

except:
    arcpy.AddMessage(arcpy.GetMessages())
    arcpy.AddError("Failed to find and replace where " + queryString)


texts: python-gdal write a GeoTIFF in colour from a data array

texts:
I have a 14400 x 7200 (360 x 180 at .025 deg resolution) numpy array containing data which I want to make into a GeoTIFF. I have done the following which creates a GeoTIFF fine:

# Set file vars
output_file = "out.tif"

# Create gtif
driver = gdal.GetDriverByName("GTiff")
dst_ds = driver.Create(output_file, WIDTH, HEIGHT, 1, gdal.GDT_Byte )

# top left x, w-e pixel resolution, rotation, top left y, rotation, n-s pixel resolution
dst_ds.SetGeoTransform( [ -180, 0.025, 0, 90, 0, -0.025 ] )

# set the reference info
srs = osr.SpatialReference()
srs.SetWellKnownGeogCS("WGS84")
dst_ds.SetProjection( srs.ExportToWkt() )

# write the band
dst_ds.GetRasterBand(1).SetNoDataValue(255) #I set my nodata values in array to be 255
dst_ds.GetRasterBand(1).WriteArray(data_array)

But the resulting tiff when opened looks to be binary in colour (black & white). How can I assign a colourmap to this image? Even just in grayscale? My values have a small range (1.7 - 1.9). Should I just remap this to 0-255?

Edit:

I ended up scaling my data array from ~ 1.7-1.9 to 0-255 and that worked nicely.

import numpy as np
data_array_scaled = np.interp(data_array, (data_array.min(), data_array.max()), (0, 255))

This gives me grayscale only though -- how do I generate in colour?

Edit 2: Updated code

import numpy as np
import os
from osgeo import osr, gdal

data_array = np.load('fini.tif.npy')
data_array = data_array[::-1]
data_array_scaled = np.interp(data_array, (data_array.min(), data_array.max()), (0, 255))
data_array_scaled[data_array_scaled == 156.91450597804007] = -9999

r = data_array_scaled

RES    = 0.025
WIDTH  = data_array_scaled.shape[1]
HEIGHT = data_array_scaled.shape[0]

output_file = "out.tif"

# Create GeoTIFF
driver = gdal.GetDriverByName("GTiff")
dst_ds = driver.Create(output_file, WIDTH, HEIGHT, 1, gdal.GDT_Byte)

# Upper Left x, Eeast-West px resolution, rotation, Upper Left y, rotation, North-South px resolution
dst_ds.SetGeoTransform( [ -180, 0.025, 0, 90, 0, -0.025 ] )

# Set CRS
srs = osr.SpatialReference()
srs.SetWellKnownGeogCS("WGS84")
dst_ds.SetProjection( srs.ExportToWkt() )

# Write the band
dst_ds.GetRasterBand(1).SetNoDataValue(-9999) #optional if no-data transparent
dst_ds.GetRasterBand(1).WriteArray(r)

band = dst_ds.GetRasterBand(1)
colors = gdal.ColorTable()

colors.SetColorEntry(1, (112, 153, 89))
colors.SetColorEntry(2, (242, 238, 162))
colors.SetColorEntry(3, (242, 206, 133))
colors.SetColorEntry(4, (194, 140, 124))
colors.SetColorEntry(5, (214, 193, 156))

band.SetRasterColorTable(colors)
band.SetRasterColorInterpretation(gdal.GCI_PaletteIndex)
del band, dst_ds


texts: QGIS Georeferencing - removing co-ord bubbles from .tiff file layer

texts:
I'm new to QGIS ...

I've georeferenced a raster image using the georeferencing tool. That went well. However, when the new TIFF layer opened, it still had the yellow co-ord bubbles. How do I remove these as I can create a map and they are not a good look?


texts:
I'm trying to open an WFS from my GeoServer in QGIS. But with very little success.  I'm getting an error:
"Error when parsing GetFeature response : Error: not well-formed (invalid token) on line 1"

Searching this error indicates it could be an encoding issue...
My data is on a PostGIS database with encoding = Latin1, and my GeoServer is set to UTF-8. Could be the problem, but changing GeoServer global settings to "ISO-8859-1" which should be equal to latin1, does not help.

Any ideas that might fix it?

postgis version: 9.4.10
geroserver version : 2.11.1
qgis version 2.18.7


texts: How to use thematic layer for supervised classification in eCognition

texts: I'm classifying vegetation in eCognition. From ArcMap I already have a shapefile with training polygons that I would like to use for this but so far I was only able to insert this as a thematic layer and use for the segmentation process. Also, is eCognition able to handle the attribute table information where the classes for the polygons are stored or will I have to split the shapefile with all the information into individual shapefiles?


texts: Ability in QGIS modeler to just run a step or subsequent steps in a model

texts:
In the current QGIS LTR is it possible to just run one step or all subsequent steps in a model if the previous steps have all been run and are not changed?

Sometimes steps take hours to run so if they haven't changed, I'd like to just run a step for testing rather than the whole thing without having to split the model into multiple different models.

If you deactivate an upstream process then the downstream ones also get deactivated and if you disconnect it then you have add an input and delete this once tested -all of which are time consuming.


texts: Why is SelectLayerByAttribute_management claiming "Column Does Not Exist"?

texts:
The solution to this question is not working for me.

My error message:

  File "SumAreas.py", line 24, in
  arcpy.SelectLayerByAttribute_management(fc, "NEW_SELECTION, '"CLASS_CODE" = "A"')
  File "c:\program files (x86)\arcgis\desktop10.2\arcpy\arcpy\management.py", line 6688, in
  SelectLayerByAttribute raise e
  ExecuteError: ERROR 000358:
  Invalid expression Failed to execute (SelectLayerByAttribute).

Further message in Python Window's side window:

ERROR 000358: Invalid expression
A column was specified that does not exist.

I am trying to select based on a column ("CLASS_CODE") which is of "String" type.

My line 24 (referenced in Error Message above):

arcpy.SelectLayerByAttribute_management(fc, "NEW_SELECTION", '"CLASS_CODE" = "A"')

I can confirm that feature class "fc" definitely has a column called "CLASS_CODE".

So if the column exists, and is spelled correctly, why am I getting this error?


texts: How to load(preload) Geojson with JS before start another script?

texts:
I need to load my Geojson file before start my script .
I need to load file before start this script :

var geojsonPrecitany = false;
var geojson = $.getJSON("./tematicke_vrstvy/tematicka_vrstva.geojson", function() {
    geojsonPrecitany = true;
});

Which checked that geojson is loaded ( if loaded geojsonPrecitany = true ).

OR

I need the solution how to load geojson file before start script with map (in Openlayers v3)

Because the work with geojson in js is asynchronous so i need to have the file peloaded before i start work with them in map script.


texts: I have pair of images: can I generate digital elevation model from that pair of images using QGIS? if it's possible, can anyone provide the way to do that?


texts: How to specify a password with shp2pgsql using CentOS?

texts:
I'm developing a python script that uses shp2sql via subprocess.call().
The string I'm using works fine on Linux Mint:
shp2pgsql -s 4326 ../data/download.shp ccd_hoods_new | psql -h hostname -d databasename -U username -P password -q

However, when I run the script on CentOS, I'm prompted for a password, which presents a new problem when deploying and automating via CRON.
The only documentation I can find is this cheatsheet, which confirms that I'm doing things correctly.
But why on the CentOS machine would it not be recognizing the -P flag?
Both servers are running the same script that points to the same host running PostgreSQL 9.6 and PostGIS 2.4.


texts: Adding today's date into attribute table using QGIS Processing modeler?

texts:
I am looking for a way to use the Processing modeler to add a date field to an output dataset in QGIS 2.6.1

I have tried using the field calculator algorithm without success. (Also the advanced python field calculator without success).

I have a vector layer and I add a field calculator algorithm and set the field type to Date, field length to 30, precision to 0, create new field and set the formula to $now
I have also tried datetime.date.today()

The error I get is;

  Prepare algorithm 0: Field calculator
  Running Field calculator [1/1]
  Parameters: INPUT_LAYER =C:/Users/lee/Desktop/TEMP/test.shp, FIELD_NAME =loadDate, FIELD_TYPE =3, FIELD_LENGTH =30, FIELD_PRECISION =0, NEW_FIELD =True, FORMULA =$now
  Failed
  Error executing algorithm 0 underlying C/C++ object has been deleted

(I am looking to automate the process of adding a 'loadDate' field to features so that users know when a feature was loaded into a postgis table. I want to repeat it over a number of layers, and do it on a regular basis hence the use of the modeler).


texts:
I'm trying to understand how data is stored in the GeoTIFF format. I understand that the file begins with header, where are stored information
which says that the file is GeoTIFF file. I also know that GeoTIFF is only container, so which types of data it can store? When I have multispectral image,
how are the values of the pixels stored? Something like this?
data_array=[[value1,value2,value3],[value1,value2,value3],[value1,value2,value3]...]
It would be great to describe it on these pictures or similar pictures. I can't find these informations anywhere. Can you explain it to me?


texts: Raster:area() in R before and after projecting data

texts:
I have the following snippet of code:

tracts<-readOGR('data/.')
proj4string(tracts)
"+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs"

tracts_projected <- spTransform(tracts, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))
proj4string(tracts_projected)
"+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0"

If I take the min area of both, the values are different, even though both spatialpolygondataframes are the same except for the difference in projection.

min(area(tracts)/1000000)
 1.221094e-05

min(area(tracts_projected)/1000000)
9.505325e-06

Why is this so? And more importantly which is the better way to get area


texts: What is wrong with how I am using ST_Within?

texts:
I am learning postgis 2.1. I created a postgis table from a shapefile and imported it into the DB. It made a table like this:

db=# \d districts
     ...
     geog       | geography(MultiPolygon,4326) |

I also created a table "callsforservice" importing tabular data from a CSV, which I then converted into points with a SRID = 4326.

db=# update callsforservice set latlongpointdetails = ST_SetSRID(latlongpoint, 4326);

I want to know which callsforservice.latlongpointdetails are included in which districts.geogs

I tried with this query, which takes very, very, very long time. The limit 1 is only to shorten the time the query is running. After 2 min or so I kill it:

select * from callsforservice, districts where ST_Within(geometry(districts.geog), callsforservice.latlongpointdetails) limit 1;

If I add an additional condition to the where clause (where districts.district='1-1') -- this query returns immediately with 0 rows, even though districts does have a row with this ID.

Do I have the syntax right for finding out which calls came in what districts?

The docs for ST_within say "the source geometries must both be of the same coordinate projection, having the same SRID." If I examine the two tables with the \dt command then the geog column in districts shows type geography(MultiPolygon,4326) but the latlongpointdetails column in callsforservices only shows type geometry -- without parenthesis like in the districts table. Maybe this is the reason for the failures of ST_Within? Where am I going wrong?


texts: How to read vector datasets using GDAL library

texts:
Up to now, I have used GDAL in order to read raster information and geotiffs and show them with the help of OpenGL. Now I have a vector dataset, and by that, I mean I have the files with types:

DBF File , PRJ File , SBN File , SBX File , SHP File , LOCK File , XML Document , SHX File

I want to read the vector information (e.g vertexes of the geometries and geometrical and topological relationships between feature) and then show them using OpenGL.

What part of GDAL library should I study?
Can you introduce me a tutorial link that has fully covered what I
want to read from the vector dataset?


texts: How do I tell GRASS where to find Python?

texts:
For some reason, when I start GRASS on my Mac it says that it has found Python 2.6.7. However, /usr/bin/python version is 2.7.2 and I can't find a lower version anywhere on my system. Potentially relevant system info follows:

Mac OS X 10.8.3

GRASS version 6.4.2 (2012) and all required frameworks installed with binaries from http://www.kyngchaos.com/software/grass in the past few days, so everything should be as up-to-date as I know how to make it.

I suppose this is really a two-part question.

How can I tell where GRASS is finding Python?
How can I tell GRASS where it should be finding Python?


texts:
I need to reorder the fields in a tab file from within QGIS modeler. Table Manager plugin allows this but I can't seem to find a way to do this in modeler. It needs to be run on thousands of files.

Can we import/use plugins like table manager in the graphical-modeler or in the python window?

@Olle - I don't see the same options. How do I add them to my toolbox? They are active in processing-->options


texts:
I am trying to interface NMEA output of multiple GPS units transmitting wirelessly to a single PC, each GPS position can then be shown on georefance map of that area.

My question

1.) How to interface multiple GPS data with one PC, since one GPS is interfaced automatically in QGIS?

2.) How to display each GPS position simultanously on one map?

Thanks


texts: Invert the legend scale representation [in R environment]

texts:
Through the following code,

cols <- (heat.colors(255))
plot(raster_sum3, col=cols[round(seq(1,255,length.out=10))], legend=TRUE) # raster of the study area
plot(site_3003, add=TRUE) # site boundaries

I plot a raster obtaining the attached image

I want to reverse the legend scale, in order to place lower values above and higher values in the bottom of the scale. Despite several researches, I did not find the appropriate syntax to reach my aim.


texts:
For my research project at university I would like to geocode about 100-200 addresses in mostly rural areas in Australia.

I'm having problems finding a suitable geocoding online service. Services that are based on OSM are probably not detailed enough. Services that use Google for their reference data cannot be used as it would violate Google's Terms of Service.

Can anyone recommend a service, or does anyone have experience in Australia?


texts: Selecting Features and Averaging Attributes

texts:
I currently have two point shapefiles, SHP A with approximately 75 points and SHP B with approximately 2000 points.  What I would like to do is:

1) For each point in SHP A I would like to find the closest 10-15 points in SHP B

2) Find the average of Field1 for the selected SHP B points

3) Place the average Field1 for the closest 10-15 points into attribute table for each point in SHP A

I can easily comprehend and complete the selection process of this task, and I'm sure I could complete it in a very long way, but I am wondering if there is any tools or processes to do this in an easier manner.


texts: Getting data of area covered by polygon in ArcGIS Desktop?

texts: We have a layer that contains land parcel information such as area for each parcel, and a unique parcel ID. We also have multiple polygons that cover some parcels. We're trying to figure out how to determine which parcels are covered by the polygon and which parcels are not covered by the polygon.


texts: Use GDAL 2.4 in Geoserver ImageI/O-Ext extension on Windows

texts:
Is it possible to use any version of GDAL newer than 1.9.2 in Geoserver on Windows (64bit) to support the ImageI/O-Ext for serving JP2 images?

In 2 days of test I could have the 1.9.2 version running.. ..but any other attempt to use other pakages failed.


texts: Postgis pgAdmin ST_area returning different to other GIS's

texts:
I have a shapefile which is returning a different area value in pgAdmin than in FME, Cadcorp SIS and QGIS and I don't know why it is different in pgAdmin.

The setting used in pgAdmin are:

SELECT ST_Area(GEOMETRY) FROM greenbelt_shp

the shapefile is available at:

https://drive.google.com/open?id=1AuXPKUoLermKBiEbDrxp9YHIG-58AHQr

Can someone can tell me why it is returning a different value?


texts: How to create an outline of a ship given the geom of the GPS aerial, ship dimensions, true heading

texts:
I have a table of AIS positions for a vessel and its true heading. I also know the distance of the GPS aerial to the bow, and the starboard side, and the overall length and width of the ship.

I want to know how to take the known point and create a second point 200m in the direction of the true heading, create a third point at +90 degrees to the true heading at the length of the distance from the aerial to the starboard side, create 4th point at 180 degrees to the true heading at the distance of total length of the vessel, etc until I have the points for each corner of the vessel and then create a line string to join the points and create an outline of the vessel pointing in the direction of its true heading.

I want to use QGIS / PostGIS to solve this.


texts: Make Raster Layer with Overlapping Images

texts: I am trying to set up a template which will display different maps based on layer definitions. One of the layers is a raster layer. Is it possible to combine a number of overlapping images into a single layer so specific images can be selected with a definition query?


texts: GPS Track Editor to crop tracks with select by clicking start and end points

texts:
I'm trying to take a bunch of gps tracks of mountain bike rides and turn this in to a trail map. I want a single track for each trail.

I'm looking for an editor that will allow me to select a segment of a track by clicking the start of the segment, clicking the end of the segment and the editor will select the segment between these points. I'm not a big fan of dragging a box or polygon around my segments since they tend to wind around and sometimes overlap.

Ideally, the editor would be able to display multiple tracks at a same time. I've already got qgis and arcgis available, so if I can select in the method I prefer using that software, I'm fine with that.

Currently, I've just been using the website ridewithgps to cut these segments out, but that requires a lot of saving as separate routes and exporting chunks. However, all I have to do to select the segment I want is click on the elevation profile at the start point, and then at the end point and copy that segment to a new route which I then export as a gpx file.

Is there any desktop software for mac or pc that will allow me to cut out segments in this fashion or a way to do this in qgis or arcgis?

Here's an example of how I select a segment in ridewithgps.

I hover over the trail or the elevation profile and click the start point, then I click the end point on either the trail or the elevation profile and that segment (not including any overlapping portions) is selected. I can tweak the extents by sliding the selection on the elevation profile and then I can then copy that to a separate track. Most of the time it takes me 4 clicks or less to break out one track. I'm looking for similar behavior or another process that works quickly in a desktop app.

Ultimately, I'm going to upload the segments I copy as individual trails on mtbproject and use for making printed maps with labels for each segment.


texts: Converting double values into time in QGIS's Field calculator

texts:
How can I achieve the following result (in column "output", string) in QGIS attribute table?

The values in "duration" were calculated with some other formulas and represent a certain time duration given in double data type.
So 2.584 actually means 2 hours and 0.584 hours (=35 minutes).
Is there a way/expression to get the "output" column directly from the "duration" column in QGIS Field Calculator?


texts: PyQGIS/QGIS2.18 : Get currently active map tool

texts:
I would like to get the name of the map tool currently in use in QGIS, so that my plugin can react accordingly. The idea would be to do something a bit like akbargumbira's answer to that other question, but for map tools instead of map layers. It should look like that :

if QGisMapTool.instance().value()=="Identify":
    do something
elif QGisMapTool.instance().value()=="Measure":
    do something else
(...)

I've searched on Stack and in the API documentation, but I have been unable to find the appropriate method so far.


texts: Summing values of grid moving upstream along flow paths in GRASS

texts:
I would like to sum the values of a grid, moving upstream along flow paths. In other words, I have a grid A with a value for each cell. I need to build a new grid, B, in which the value of each cell should be the sum of the A values of all the cells located downstream along the flow path.
How can this be done in GRASS?


texts: How can I remove "red lines" that become permanent to my map canvas while trying to digitize different shapes. The program just leaves the outline of what I'm trying to draw. It still remains even if I close all of the layers and empty everything. I've encountered this error numerous of times when creating circles, lines etc.


texts: Iterating over layers to export data as a CSV

texts:
I'm trying to export particular columns of data from the attribute table to a CSV, I've added in a functionality to remove the excess commas to make this easier so that I end up with two neat columns of all the data I need. I can tell the first for loop is working as its exporting the data, but it is no longer removing the commas, in the second half of the code. I am very, very new to all this.
This is what I'm working with:
layers = iface.mapCanvas().layers()

with open (r'C:\Users\James\Desktop\Computer science\Python\code\HER_results.csv', 'w') as file:

    for layer in layers:
        layerType = layer.type()
        if layerType == QgsMapLayer.VectorLayer:

            for getDat in layer.getFeatures():
                line = '%s, %s\n'% (getDat['Name'], getDat['PrefRef'])
                file.write(line)


texts: How to aggregate values of a shapefile according to another shapefile in QGIS?

texts:
I have a shapefile that has fine resolution, say 1km by 1km. There are several attributes that I am interested. Let's call this shapefile S1.

I have another shapefile, which is the administrative shapefile, say, all the towns in a state. Let's call this shapefile ST.

ST and S1 have the same spatial span.

So, for a town in ST, there are several pixels in S1 that fall into the town.

I would like to find out the average values of an attribute in all the pixels in S1 that fall into every town in ST.

I tried to Google solutions but I do not know the correct terminology for this type of operation and did not find useful tutorials.

Could you provide some hints or help?


texts: Using Stacked Labels with $length and single quote in QGIS?

texts:
I am using QGIS and making stacked labels.

Using $length, how do I get a ' (single quote) to show up after the length ie 386'?

My current expression is:

run  +  '\n'  +  "sequence"   ||  '\n'  ||  round($length,0)


texts: Writing Make Table Query script in ArcGIS Desktop

texts:
I am trying to get the make table query running in a script. IF I simply export the script from a model, it works without a problem. However, I do have to change the variable names, as this will be part of a bigger script. As soon as I change any variable names, it won't work anymore. See example:

# Import arcpy module
import arcpy

# set outputs to be overwritten
arcpy.env.overwriteOutput = True

#WORKS
# Local variables:
AA_allDomains = "D:\\Student\\DomainsToTables\\featureDomains.gdb\\AA_allDomains"
AB_BusBar_Domains = "D:\\Student\\DomainsToTables\\featureDomains.gdb\\AB_BusBar_Domains"

QueryTable1 = "QueryTable1"
# Process: Make Query Table
inList = [AA_allDomains, AB_BusBar_Domains]
##fieldList = "AA_allDomains.Code #;AA_allDomains.Desc #;AA_allDomains.DomainName #;AB_BusBar_Domains.dataset #;AB_BusBar_Domains.featureclass #;AB_BusBar_Domains.field #;AB_BusBar_Domains.domainName #"
fieldList = [["AA_allDomains.Code", 'CODE'],["AB_BusBar_Domains.dataset", 'DATASET']]
expression = "AB_BusBar_Domains.domainName = AA_allDomains.DomainName"
arcpy.MakeQueryTable_management(inList, QueryTable1, "USE_KEY_FIELDS", "", fieldList, expression)

#DOESN'T WORK
# Local variables:
in1 = "D:\\Student\\DomainsToTables\\featureDomains.gdb\\AA_allDomains"
in2 = "D:\\Student\\DomainsToTables\\featureDomains.gdb\\AB_BusBar_Domains"

QueryTable1 = "QueryTable1"
# Process: Make Query Table
inList = [in1, in2]
##fieldList = "in1.Code #; in1.Desc #; in1.DomainName #; in2.dataset #; in2.featureclass #; in2.field #; in2.domainName #"
fieldList = [["in1.Code", 'CODE'],["in2.dataset", 'DATASET']]
expression = "in2.domainName = in1.DomainName"
arcpy.MakeQueryTable_management(inList, QueryTable1, "USE_KEY_FIELDS", "", "fieldList", expression)

I already tried to play around with different DB types, different ways of accessing the fields but nothing worked out so far.


texts: QGraphicsTextItem's position not updated?

texts:
I'm using
QgisInterface->mapCanvas()->scene()->addItem(QGraphicsTextItem)
to add a QGraphicsTextItem to the map canvas. The problem is that when I use the pan tool to move the display the added QGraphicsTextItem is not moving !
How can I make it update its position whenever the pan tool is used on the map ?


texts: Problem with PostGIS ogc_fdw: unknown SRS

texts:
I am trying to connect a WFS from within PostGIS using the ogr_fdw foreign data wrapper
[https://github.com/pramsey/pgsql-ogr-fdw].

With most servers this works great!
With one specific server I get the error:

  EPSG PCS/GCS code 404000 not found in EPSG support files.  Is this a valid
  EPSG coordinate system?

The SQL code is like this:

CREATE SERVER kvb_wfs FOREIGN DATA WRAPPER ogr_fdw
OPTIONS (
   datasource 'WFS:online-service.kvb-koeln.de/geoserver/ows',
   format 'WFS' );

I tracked it down and figured out, that the problem comes from the GeoServer which actually has this code defined as being a wildcard crs.

I tried to put the code in the gdal file (psc.csv) and the proj4 file (epsg). No success.

What am I doing wrong?


texts: Can't see multiple layers when adding a WMS

texts:
I am using QGIS 3.8.2.

I have a WMS (https://openwms.statkart.no/skwms1/wms.topo4.graatone?) that I want to add, and in Arcgis it lists several layers (33 total) in the table of contents. When I add it in QGIS, only one is showing, it doesnt matter if I select all or only the parent layer in the WMS. Is there some boxes that needs to be ticked ahead of the operation?


texts: I want to add search widget to my cmv application, but nothing is displayed when I click on search, I have copied and paste all the files from github in my viewer folder and I started by creating a folder called widgets then I pasted the search folder, I also pasted searchwidget.js in config folder.


texts:
I am testing CARTO Airship. I have a query regarding airship widget default text. There are some default text values in all widgets like "There is no data to display" or "NO DATA AVAILABLE". I also need to change those text to my local language. I could not find any attribute keys like "heading", "description" etc. for this in the docs.

Is it possible to change those text values?


texts: Subset netCDF based on last valid variable by level

texts: I have a netCDF file (in this LINK) which contains monthly average values for 4 parameters at 11 different depth levels (0-200m). I would like to extract the monthly value for temperature "votemper" at the maximum valid depth (i.e. when the water is 20 m deep there are only NA vlaues for greater depths) however, the last valid depth varies spatially based on the water bathymetry. I have got as far as turning the netCDF into an array in R but am having trouble with how to extract just the values from the maximum valid depth value.


texts: Backwards compatible python scripts in ArcGIS

texts:
As I move through my small collection of python scripts written for arcgis 9.3 I'm discovering that many things work with little or no changes in arcgis 10. For example this pattern "just works":

import arcgisscripting
gp = arcgisscripting.create(9.3)
gp.workspace = 'd:/scratch'

even though the new way is:

import arcpy
from arcpy import env
env.workspace = 'd:/scratch'

So my question is, what it is the best structure to follow so that scripts are only written in the v10 way when they need functionality not available in earlier versions? Or perhaps better phrased as, what is the best structure for backwards compatiblity? I don't want to make my stuff automatically unavailable to people on earlier versions, just when it wouldn't work anyway.


texts: Selecting elements by grouping with one attribute and specific rule using QGIS

texts:
I have a question about a specific query for a Shapefile-Layer in QGIS 3.4.

I have a Point-Shapefile with the name point.shp. The Shapefile has a lot of different points with a large amount of attributes. But two attributes are interesting for us. The first attribute has the name 'Zone', and the second is called 'Direction'. Both attributes have the data type string.

This picture visualize the structure of the Shape:

Now, we have a rule: in every zone, the direction number must be unique. In the example, all points in zone 9 an 22 are legit. But we have a problem in zone 10. Direction 1 occurs twice. Because of that, all points in zone 10 are invalid.

Can I select all the points, where one zone has multiple directions with identical numbers?

Maybe, when I selected the points with the incorrect rule, I can reverse the selection to export the valid values.


texts: MRT (MODIS Reprojection Tool) error: 'Input file does not exist'

texts:
EDIT: Solved, see my answer.

I've downloaded 8 MOD44B granules from Reverb (using the FTP script method), and am trying to use MRT to mosaic them. I'm falling at the first hurdle though - when I add the 4 files (4 from 2000, 4 from 2007) I get an error saying firstly One of the input files does not exist or the inputs are not the same data product type then Input file does not exist!. If I add files individually, I get an error saying Error in Module: ReadHrFile, Opening Input Header File.

Has any one had this error before, or know how to fix it? I've updated Java and reinstalled MRT.

Versions:

MRT: 4.1 (March 2011)
Windows: 7
Java: 6 (Updated to 7, still doesn't work)

I would use the online interface, but I need MOD44b which isn't available with that for some reason.


texts: Join spatial point data with multiple polygon data using R


texts: Can the ArcGIS solution for uneven buffering be achieved in QGIS? version 2 - uneven tree canopies

texts:
I want to create tree canopies with uneven canopy extents.
A previous solution is below:
Representing uneven tree canopy spread around point in QGIS
A solution for uneven buffering has been provided for ARCGIS - link provided below. The method could also potentially be used for uneven tree canopies. However, I use QGIS.
The solution for uneven buffers (as in ArcGIS) may generate a more attractive representative tree canopy. Can this method also be undertaken in QGIS?
A representation of the ArcGIS solution

Buffering point with different distance toward N, E, S, W using ArcGIS for Desktop?