texts: Easting and northing conversion into lat and long gives wrong information python

texts:
This is my easting and northing
315846.356,5808849.713

p2 = Proj(proj="utm",zone='55H')
x,y=315846.356 ,5808849.713
lat,long=p2(x,y,inverse=True)

I got Longitude  144.29314340365048 and Lattitude 52.39886078782927
The lat is correct but long information is not
The coordinates are from MELBOURNE, AUSTRALIA
zone information is obtained from google earth


texts: Speed up spatial query on large table on psql

texts:
How to speed up spatial query on large table of bytea.

I tried the suggestion from this website to create an index using b-tree GIST but the index is exceed the maximum size of 8191.
Any suggestions?


texts:
I have the problem with low quality raster. Now it looks like this:

a busy cat http://iv.pl/images/53039057684486451032.png

I would like to change this to looks like this:

a busy cat http://iv.pl/images/87768106012364755171.png

I'm using ArcGIS, ArcScan. After some work with Photoshop and reclassify with ArcGIS I get something like that.

a busy cat http://iv.pl/images/35342487136630589409.png

And after ArcScan I get:

a busy cat http://iv.pl/images/98057626549525851657.png

I don't really know if I'm using ArcScan correctly or the problem is with low quality raster. Can you help me?


texts: Is there a big (if any) performance difference between Geoserver SQL Layer and OpenLayers geojson Layer?

texts:
I have a table with multiline geography items in SQL Server 2008 (Routes).
I want to show them on OL Google Map but I'm having trouble creating SQL Data source in geoserver (something related to TCP/IP although its opened and enabled in SQL Config).
As a workaround , I've create a web service which brings GeoJSON from the SQL and add it as a layer.

My question is , How "bad" performance wise is my solution if I want to show couple of routes that each is made from around 7k points ?


texts: Why is using Bing Maps with QGIS free but not with ArcGIS?

texts:
Since the end of 2013 a free use of Bing Maps in ArcGIS is no longer possible. You have to obtain a Bing Maps key directly from Microsoft.

In QGIS (with OpenLayers plugin)  it is possible to use Bing Maps (and also Google Maps) without a Bing Maps key and without a fee.

Why it is possible to use Bing Maps with QGIS for free while you have to pay for it with ArcGIS (Desktop)?

ESRI source: http://downloads2.esri.com/ArcGISOnline/docs/bing_faq1.pdf


texts: OD Cost Matrix snapping distance using ArcGIS Network Analyst?

texts:
I am using the ArcGIS Network Analysis - OD Cost Matrix layer tool and I was wondering if anyone knew of a way to incorporate the snapping distance of the origin and destination to my total distance value.  Right now, the Origin snaps to the closest line feature in my network and the same for the destination point, and the tool only gives me a distance between these 2 snapping points, but I would like to include the 2 snap distances as well.

Does anyone know if this is possible?

I know I could use the Near command ontop of this tool to get the total distance, but I was hoping there was a setting in Network Analyst that could solve my problem.


texts: How to extract ALL bands from multiband raster in ArcGIS?

texts:
I have seen many questions about how to extract a singe band from a multiband raster in ArcGIS, but how do you extract ALL of the bands and save them as separate tifs?  I want to automate this process using Arcpy because I have many bands in one multiband raster.  I also want to be able to do a batch import to separate bands from multiple multiband rasters.

I'm guessing I can use CopyRaster, but I'm not sure how to use it when I need to make copies of everything!

import arcpy
import os
from arcpy.sa import Raster

#setup workspace
arcpy.env.workspace = "Y:\Commons\Kelly_Carly_Aaron_Project\lvl 3\OLCI-Sentinel3\Florida"
arcpy.env.overwriteOutput= True

#setup output directory
out_dir= "Y://Personal//crobbins//Final_Project//Separate_bands_results//"

#list all rasters in workspace
RasterList = arcpy.ListRasters()

#get list of bands from a multiband raster
def list_bands(raster_workspace):

    raster_workspace = arcpy.env.workspace
    RasterList = arcpy.ListRasters()
    bands = [Raster(os.path.join(raster_workspace, b)) for b in RasterList]

    print bands
list_bands("Y:\\Commons\\Kelly_Carly_Aaron_Project\\lvl 3\\OLCI-Sentinel3//Florida\\Copy_of_2017.0527.tif\\")

This doesn't work.  It just lists the names of all the multiband tiffs I have in Y:\Commons\Kelly_Carly_Aaron_Project\lvl 3\OLCI-Sentinel3\Florida.


texts:
I am running ArcMap 10.2 and have been trying to complete a fairly large viewshed analysis and have chosen to use the visibility tool in order to be able to produce an optional above-ground-level (AGL) layer which will be the focus of further analysis.

The visibility is from a layer consisting of 2 points and the raster is a DSM. I have run several repetitions of this tool (both 3D Analyst & Spatial Analyst) and every time the process completes successfully producing a viable viewshed, but there is never any AGL to be found. There are no errors accompanying the completion of the tool either.

Why would this be? How do I get the tool to produce an AGL? I have tried using only one viewer point as opposed to two without any difference in the outcome.


texts:
I've just read the new post at Medium about Sentinel-2 color correction with Sentinel Playground. It says that

  The default gain factor is set to 2.5, which corresponds to white-point reflectance of 40% (a pixel with 40% or more reflectance in all three bands will appear white in the resulting image).

I get very similar visual result by stretching contrast in QGIS [0, 4095] --> [0, 255] for each band. I've chosen these ranges based on histogram and radiometric resolution of 12 bits. Please, help me understand what operations are performed with the image when setting gain factor at Sentinel Playground?


texts: Connecting river segments with R riverdist package

texts:
I have imported river courses from OpenStreetMap into QGIS in order to recreate the Rhine river system for analysis in R with the riverdist package.

In order to calculate distances between sampling spots the river segments should be connected, of course. I can check this with

topologydots(rivers=Rhine)

And see that many segments, well almost all, remain unconnected (red dots)

Using

cleanup(Rhine)

as suggested in the riverdist tutorial asks me several questions but does not connect the river segments.

Do I have to connect these segments one by one in QGIS or is there any solution within the riverdist package?


texts: How to add 'top' and 'bottom' directions with leaflet marker label

texts:
I need to add a label to marker from bottom. But there are no any value to add for direction option in bindLabel method.

 var marker = L.marker([ latitude, longitude ], {
         icon : iconInfo,
         iconAnchor: [10, 10],
         labelAnchor: [6, 0]
     }).bindLabel(displayingText, {
         direction: 'right',
         noHide : true
     }).addTo(leaf_map);

It only have right|left|auto options only.
How can I achieve this with leaflet marker ?


texts: Actual slope between two points from raster slope value

texts:
We have our own raster map and I have raster pixel slope value on map which is calculated using formula as shown in the link below
Pixel Slope Calculation, which is used to show different color as per the pixel slope value. So for a given point I can have these slope value anywhere on the map.

Earlier for calculating straight line slope between two points, I used to read elevation difference and then calculate straight line distance between those two points and finally dividing these two values. But I think with these pixel slope value my calculation can be minimize.

Is it possible to calculate the same value using pixel slope values of these points? If yes then how?


texts: Compile QGIS / Can't find Qt Cryptographic Architecture

texts:
I'm trying to compile QGIS-master 2.12 for the first time following the INSTALL guide (https://github.com/qgis/QGIS/blob/master/INSTALL).
I'm stucked while compiling with the "CMake Error at cmake/FindQCA.cmake:64".

Here is my environment:

Windows XP SP3
Qt 4.8.7
Qt Creator 3.6.0
Microsoft Visual Studio 10.0
Windows SDK for Windows 7 and .NET Framework 4
CMake 3.0.2
OpenSSL 1.0.2f
OSGeo4W (expat, fcgi, gdal, grass, gsl-devel, iconv, pyqt4, qt4-devel, qwt5-devel-qt4, sip, spatialite, libspatialindex-devel, python-qscintilla)
Cygwin (bison, flex, git)
I tried to install QCA from the source but I can't get it work.

Any advice on how to go ahead?

Here is the log when I try to compile QGIS in QtCreator:

-- The C compiler identification is MSVC 16.0.30319.1

-- The CXX compiler identification is MSVC 16.0.30319.1

-- Check for working C compiler: C:/Program Files/Microsoft Visual Studio 10.0/VC/bin/cl.exe

-- Check for working C compiler: C:/Program Files/Microsoft Visual Studio 10.0/VC/bin/cl.exe -- works

-- Detecting C compiler ABI info

-- Detecting C compiler ABI info - done

-- Check for working CXX compiler: C:/Program Files/Microsoft Visual Studio 10.0/VC/bin/cl.exe

-- Check for working CXX compiler: C:/Program Files/Microsoft Visual Studio 10.0/VC/bin/cl.exe -- works

-- Detecting CXX compiler ABI info

-- Detecting CXX compiler ABI info - done

-- QGIS version: 2.13.0 Master (21300)

-- Could not find GRASS 6

-- Could not find GRASS 7

-- Found Proj: C:/OSGeo4W/lib/proj.lib

-- Found GEOS: C:/OSGeo4W/lib/geos_c.lib (3.5.0)

-- Found GDAL: C:/OSGeo4W/lib/gdal_i.lib (2.0.2)

-- Found Expat: C:/OSGeo4W/lib/libexpat.lib

-- Found Spatialindex: C:/OSGeo4W/lib/spatialindex_i.lib

-- Found Qwt: C:/OSGeo4W/lib/qwt5.lib (5.2.1)

-- Found Sqlite3: C:/OSGeo4W/lib/sqlite3_i.lib

-- Found PostgreSQL: C:/OSGeo4W/lib/libpq.lib

-- Found SpatiaLite: C:/OSGeo4W/lib/spatialite_i.lib

-- Could NOT find ZLIB (missing:  ZLIB_LIBRARY ZLIB_INCLUDE_DIR)

-- Could NOT find PNG (missing:  PNG_LIBRARY PNG_PNG_INCLUDE_DIR)

-- Could NOT find JPEG (missing:  JPEG_LIBRARY)

-- Could NOT find TIFF (missing:  TIFF_LIBRARY) (found version "4.0.2")

-- Could NOT find ZLIB (missing:  ZLIB_LIBRARY ZLIB_INCLUDE_DIR)

-- Found Qt version: 4.8.6

-- Touch support disabled

-- Found QScintilla2: C:/OSGeo4W/lib/qscintilla2.lib (2.8.4)

CMake Error at cmake/FindQCA.cmake:64 (message):   Could not find QCA Call Stack (most recent call first):   CMakeLists.txt:312 (FIND_PACKAGE)

-- Configuring incomplete, errors occurred! See also "//mac/Home/Desktop/build-QGIS-master-Desktop-Default/CMakeFiles/CMakeOutput.log". See also "//mac/Home/Desktop/build-QGIS-master-Desktop-Default/CMakeFiles/CMakeError.log".


texts: CMV app - files not accessible within network; cloudflare & @font-face

texts:
I am having a connectivity issue. My network will not allow me to connect to:
1. //cdnjs.cloudflare.com/ajax/libs/proj4js/2.3.12/proj4.js
2. CSS3114: @font-face failed OpenType embedding permission check. Permission must be Installable

Where can I reference these files within the application? The WAB widgets will not work now.


texts:
In QGIS 3.0, none of my GDAL raster tools work (toold in Raster Analysis, Raster Conversion, etc). I've been able to find workarounds because many of these tools are available in other toolboxes.

Example: The Processing Toolbox has a "Slope" tool. Here is a screenshot of my parameters:

Adding to the confusion, my log is not helpful. Screenshot below:

Can anyone offer me some pointers?


texts: Mapbox GL: Possible to move navigation controls?

texts: Would anyone know if it is possible to move the navigation controls from top-right to top-left? By default, it is top-right, but I would like to add a company's logo at that position and I've seen other Mapbox maps with controls on the left.


texts: Request bbox with QgsRasterLayer WMS in QGIS

texts:
Inside a plugin i'm trying to load a wms layer, with python in QGIS.
Everything works but the size doesn't match: is downloading the entire wms layer not the extent of the window i'm using in QGIS (cropped size).
What i'm missing? Here is a sample of the request:
e = self.iface.mapCanvas().extent()
bbox = str(e.xMinimum())+","+str(e.yMinimum())+","+str(e.xMaximum())+","+str(e.yMaximum())
url = "http://myserver.com:8080/geoserver/ows"
serviceType = "wms"
urlWithParams = f"crs=EPSG:{epsg}&format=image/png&layers={layer_name}&styles&url={url}&bbox={bbox}"
QgsRasterLayer(urlWithParams, temp_layer_name , 'wms')


texts: GRASS dependency issues in Windows QGIS (i.vi for instance)

texts:
I have QGIS (2.14.3) and GRASS (7.0.4) desktop installed on Windows 7 and have successfully installed the "GRASS" plugin in QGIS (using plugin manager). However, some of the GRASS dependencies, such as i.vi or vegetation indexmodules, are not showing up in the QGIS Mapset section.

Next, I opened up the GRASS and noticed that the above dependency is actually available in GRASS under imagery > satellite image products > i.vi menu.

Is this an installation problem? Should I modify the plugin path in Windows?


texts:
Im working with QGIS. When dividing a map with a grid I set in the paramters X=0.1 and Y=0.1. Could someone please explain what does that mean? I mean would the map be divided into 100mx100m squares? How do I know the exact measure in meters of each side of the square? Maybe its a silly question but its not that clear to me. Any help would be appreciated!

Thanks!


texts:
In QGIS 2.6.1, I have merged two polygon shapefiles. After the merge I am facing a problem at boundaries.

So I thought of changing Polygon -> lines, manually deleting one of the lines at the boundary, and finally change line -> Polygon.

The problem is I am not able select a line between two nodes to delete; the entire polygon is getting selected when I try to select. Is there any better approach to do this?


texts: Seeking code behind Fix Geometries tool in QGIS?

texts:
I used the very nice "Fix Geometries" tool in the Processing Toolbox of QGIS 3.6 to repair a vector layer that was giving me trouble. The tool in QGIS worked just fine, but I was hoping to find out the actual algorithm or Python/C++ code used to implement this feature. I was wondering if anyone knew what code was used to actually implement this tool and fix geometries. Or if there is a different an better way to fix the geometries, any such recommendation would work too.

So let me explain a little bit. I did look up the QGIS documentation for the "Fix Geometries" tool in the Processing Toolbox. Unfortunately, that documentation does not explain how it fixes the geometry--meaning the algorithm. Here is a little excerpt from the documentation:

  Fix geometries

  Attempts to create a valid representation of a given
  invalid geometry without losing any of the input vertices. Already
  valid geometries are returned without further intervention. Always
  outputs multi-geometry layer.

I also looked in the Log window of the Fix Geometries tool, but it did not seem to show the GDAL code or such, as in other tools.

In the past I have resolved this issue by wrapping a narrow buffer around the geometry of interest, to overcome invalid geometry errors. Of course this does create some errors in the resulting layers. However, I saw that the Fix Geometry tool did not seem to add a buffer, and hence seemed to solve the invalid geometry problem without introducing some error. Hence I was wondering if anyone knew the code behind that tool.


texts: Fixing encoding error when exporting GeoPandas dataframe to Esri shapefile

texts:
When I try to export a list of geodataframe to esri shapefile like this:

out_path = os.path.join(out_folder, filename)
df.crs= "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
df.to_file(out_path, driver='ESRI Shapefile')

It used to work fine but now it give this error:

UnicodeDecodeError: 'utf-8' codec can't decode byte 0xe9 in position 73: invalid continuation byte
Exception ignored in: 'fiona._env.log_error'
UnicodeDecodeError: 'utf-8' codec can't decode byte 0xe9 in position 73: invalid continuation byte
Traceback (most recent call last):
  File "D:\script.py", line 72, in <module>
    result.to_file(out_path, driver='ESRI Shapefile')
  File "C:\Users\xxx\anaconda3\lib\site-packages\geopandas\geodataframe.py", line 533, in to_file
    to_file(self, filename, driver, schema, index, **kwargs)
  File "C:\Users\xxx\anaconda3\lib\site-packages\geopandas\io\file.py", line 191, in to_file
    filename, "w", driver=driver, crs_wkt=crs_wkt, schema=schema, **kwargs
  File "C:\Users\xxx\anaconda3\lib\site-packages\fiona\env.py", line 398, in wrapper
    return f(*args, **kwargs)
  File "C:\Users\xxx\anaconda3\lib\site-packages\fiona\__init__.py", line 267, in open
    **kwargs)
  File "C:\Users\xxx\anaconda3\lib\site-packages\fiona\collection.py", line 157, in __init__
    self.session.start(self, **kwargs)
  File "fiona\ogrext.pyx", line 1147, in fiona.ogrext.WritingSession.start
  File "fiona\ogrext.pyx", line 619, in fiona.ogrext.Session.get_schema
UnicodeDecodeError: 'utf-8' codec can't decode byte 0xe9 in position 9: unexpected end of data

I don't have any problem if write: dt.to_file(out_path, driver='geojson')


texts: Is it possible to display the length of each segment with the Openlayers measure control?

texts:
I am trying to make a OpenLayers measure tool that displays the length of each segment on top of each segment and the total length of the measurement as shown below:

I have used this example http://openlayers.org/dev/examples/measure.html.
I can see in the example that the geometry created is one big LineString.

Therefore my general thought was to get the measure control to create a new line segment (LineString) for each of the user click except the first,
then add an attribute with the segment length and display it with the label parameter in the style.

But I don't know if it will conflict with the current design of the measure control or if this is the correct approach to my problem?

All comments or guidance would be appreciated.


texts:
I see that it is possible to access Google Fusion Tables in QGIS using OGRs Virtual Format layer. Shown here and details about the Google Fusion Table driver here.

I was wondering whether something similar was possible with cartodb using the Postgres/PostGIS driver? I had tried a few different combinations with out any luck what so ever, I just keep getting an empty layer in QGIS.

Anyone have an idea on how to make this work? or if it is at all possible?

Using dev QGIS 1.9.

Thanks

Ando


texts: Using new Australian Datum GDA2020 in SQL Server 2016?

texts:
We have several hundred datasets/spatial layers within SQL Server 2016 which has the geometry datum/projection of GDA1994 Zone 54 (EPSG:28354) and we wish to reproject/shift the spatial locations to the new Australia GDA2020 datum, Zone 54 (EPSG:7854), how is this done?

Or does the spatial data need to be reprojected in say FME and re-loaded into SQL?


texts: How to pan (rather than zoom) when using a 2-finger drag on a Mac

texts:
Using Mac OSX, dragging the trackpad with 2 fingers performs a scroll of the web page, equivalent to using the mouse wheel.

However, if the cursor is located within an OpenLayers map, the 2-finger drag performs a zoom-in on the map.

Is it possible to change this action to perform a pan instead of a zoom? Eg, see the ArcGIS Server JS API's behaviour on a Mac with a trackpad.

Thanks


texts: Is this error "invalid context...fatal error" due to GRASS GIS bug or the system?

texts:
Copy of error message that appears in the terminal (Macbook Pro, OS X 10.9.5):

  MacBook-Pro.local Python[2653] :CGContextRestoreGState: invalid context 0x0. This is a serious error. This application, or a library it uses, is using an invalid context  and is thereby contributing to an overall degradation of system stability and reliability. This notice is a courtesy: please fix this problem. It will become a fatal error in an upcoming update.

What could cause this? Can I ignore it? I'm importing Lidar data into GRASS, and the application completes the process without an error in the application's command window.


texts: Result Window in ArcGIS 10.1 is not Opening

texts:
I booted up ArcMap 10.1 today and found out the Results Window was not opening up.

Background info:

ArcMap is running through a virtualization program through at my company.
Similar issue has came up in the past but with the search window not showing up.
IT couldn't figure out what caused the issue with the search windows not opening up so they uninstalled the virulization program cleared the cache associated with ArcGIS.
I noticed whenever I click on Geoprocessing->Results (I have also added in the results button through customization), the right edge in the following image blinks.
I noticed a windows update went through last night so that may be the culprit?  Then again, how did the search windows get borked?

I've scoured  stack-exchange/internet and could not find any info relating to this issue. Tools successfully run and I can view the history log as mentioned here but i really need to view the results window so I can rerun various tools again without filling out the parameters again.

Does anybody know of any registry keys/services that I can look into that that may need to be tweaked?

I don't  want to resort to IT since I know for a fact they will uninstall everything, in turn loosing all my defaults e.g. shortcut keys, toolbars, and document customization.


texts: Best SRID for worldwide data in Microsoft Sql Server geometry field?

texts:
I'm using Microsoft Sql Server 2012 to store geographical data related to hydrometeorological stations which are scattered around the globe (point layers and multipolygon layers).
My data is in WGS84 latitude-longitude coordinates, and I initially stored it in geography fields of SqlServer.

Later I had to move to geometry fields, because we started using GeoServer, and the plugin for GeoServer can't read geography fields, only geometry ones.

As a start I saved my WGS84 data in the geometry fields (keeping the SRID 4326), and more or less it works, but it doesn't seem to me the best choice.
WGS84 coordinates are ellipsoidal ones, but in the geometry field they are treated like planar. Moreover, they are in degrees, and functions like STArea() return results in square degrees (?) that are meaningless to me.
I'd need a reference system whose measurement unit is in meters.

Given that the objects I have to manage are scattered around the globe, I can't use a local projection system, like UTM33N for example.

So which could be a better SRID to convert my data to fit in geometry fields?


texts: Automating drawing of polygons (square) around number of points in QGIS?

texts:
In QGIS, how can I automate drawing a 100km x 100km square around a number of points, with the point being dead in the center?

I'd need to loop through all the points and create a polygon, but not sure what would be the automated way of doing so.

I was able to create buffers (circular) around the points, but now looking for a solution to convert those circles to squares -- something like the Envelope approach described here (but need it in QGIS): Creating square buffer around point feature using ArcGIS for Desktop?


texts: Match by location and attribute in QGIS?

texts:
I have 2 different layers.  One layer is the field layer which contains field agronomic data, name, variety, etc.  The other layer is a layer with about 2000 polygons. Each polygon is tied to a variety.  Some contain the same variety and are in different locations.  Also there are some areas that have polygons stacked on top of each other that contain different varieties.

I want to determine if there are any fields inside a polygon that contain the variety inside the polygon.  Another way of saying, I want to know if there are any fields inside the different polygons that have that variety.

I have tried join by location.  (1:Many)Instead of telling me which polygons contain the field, I only get back how many polygons as a count.  since I am not getting the information from the polygon layer, I cannot query if the varieties match.

I have also tried joining by attribute.  Doing this will link up the varieties, but it is not spatially accurate.

I am using QGIS.


texts:
I'm trying to add a vector layer with a WFS layer and it's not working. Can anyone point a solution?

OpenLayers code:

var wfs = new OpenLayers.Layer.Vector(
                "Stavros Features",
                {
                    strategies: [new OpenLayers.Strategy.BBOX()]
                    ,projection: new OpenLayers.Projection("EPSG:4326")
                    ,protocol: new OpenLayers.Protocol.WFS({
                        version: "1.1.0"
                        ,srsName: "EPSG:4326"
                        ,url: "http://localhost:8080/geoserver/wfs"
                        ,featureNS :  "http://postgis.org"
                        ,featureType: "polygon"
                        ,geometryName: "the_geom"
                    })
                });

In the GeoServer configuration, I have a workspace 'myWorkspace',
a store 'myStore' and a layer 'myLayer' and WFS is enabled. Also in the 'publishing' tab of my layer in the 'WMS Settings
Default Style' I have selected 'polygon' as myLayer is multipolygon.

What am I doing wrong?

Am I supposed to specify the workspace name in the 'featureNS' property?

Is there a good tutorial for OpenLayers?
I have read the docs which are very poor and the examples doesn't help me always. I have also read the user manual of GeoServer and I found it very good.


texts: Which raster resolution (DEM) can I get out of lidar data with variable grid width?

texts:
I have ASCII files in *.xyz format and want to generate a DEM out of it. According to the information about the data the *.xyz files are containing variable grid widths - from 1 meter to 4 meter. Following a small section of the data (converted to points):

Is it possible to generate (with interpolation) a raster with a resolution I want (1,2,3 or 4 meter)? If I generate a 1-meter-raster out of it, how is the quality of the resulting raster in your opinion?

I use ArcGIS in general.


texts: Automating Pop up configuration of web app in ArcGIS for Portal?

texts: I have been tasked to customize several hundred layers (pop up windows) in Portal for ArcGIS (web app). I was wondering if there is a faster way to accomplish this, instead of using the on client pop up config (custom attribute display) for each layer over and over which is repetitive, doing it with code instead? I heard there is a config.json file that has the attributes from the layer schema I have and I can populate these pop ups with code instead?


texts: Fill null pixels resulting from Landsat 8 SR image cloud masking (Earth Engine)

texts: I want to collect cloud-free surface reflectance images using a cloud masking function and then compute NDVI. However, the output image has some transparent areas that can't be used in the next step. How can I solve this problem?


texts: Raster Files Sizes from Point Data Extremely Different

texts:
Using ArcGIS 10.1.

I created rasters from several groups of point rasters using the Point to Raster tool.  All of the points groups are similar size and extents.  The first group of points had a file size of 250MB.  The other groups of points had a give size of about 900GB each. All of them are GRID float point files.


texts: Creating TIN raster from point layer from Google Earth in QGIS?

texts:
I created height layer from google earth from https://www.youtube.com/watch?v=bLbY3iMBW-A. The attribute table of the layer is:

The "ALT" field is the height from. My problem is when i try to create a TIN raster there is no option to choose the ALT field. The only options to choose is the xy coordinates. Why is it?

The QGIS 2.14 can recognized the ALT field, because when test the layer and create categorized style- the QGIS create it:


texts:
I am developing a QGIS Plugin in PyCharm. I use the QGIS Plugin "Plugin reloader" to reload the plugin when I want to test the plugin. But changes in the code are not applied on running the plugin most of the time.

I have to reload the plugin a second time or oftener to update the plugin (I am waiting till the reload-info-dialog disappears, before testing the plugin). The likelihood to reload the code, in the plugin-base-path is higher than in a containing directory.

After restarting qgis, the plugin works as expected. But I cannot restart QGIS for every test!

Edit:
I use QGIS 3.8 on Windows 10 (Plugin reloader 0.7.8)

The "plugin-base-path" is the plugin-directory where my_plugin.py is. Here i have subdirectories to handle the plugin-files.


texts: polygonize raster so that every cell becomes a polygon?

texts:
Is there a proper way of polygonizing a continuous raster so that every cell becomes a polygon? What I do now is make the cells differ enough so that the algorithm makes different polygons of them. for example:

raster1 ^ 8 -> raster2
polygonize raster2 -> vector1

This is however not a guaranteed-to-work approach. Anyone has any ideas on a better approach?


texts: Why do I get an OGR error when attempting to save an edited Shapefile?

texts:
I often get an ogr error when I try to save a shapefile after working on it in Qgis ( adding features, removing, modifying attributes.. )

The error says :
"OGR ERROR : Attempt to delete shape with feature id (-28 ) which does not exist"

Does anyone have an idea about the reason why I'm getting this error and/or how to avoid it ?


texts: Calculating common length of two lines using ArcGIS Desktop

texts:
I have a shapefile that consists of around 800 different lines.
Now I have another shapefile which also consists line data.
Is there a method which can help me calculate the shared length in the second shapefile with the first one?
The lines of both feature classes are not necessarily overlapping.
For example both the lines are crossing through the same road but on different corners. As seen in the picture I want to calculate the length, separately of the black line where it passes near the red line and also where it does not.
For one or two lines I used the split tool and calculated the length manually but I want to do it for all the lines.


texts: Using file geodatabase with HTML in attribute column?

texts:
I am working with a file geodatabase (.gdb) where one of the attribute columns contains html code:

My issue is that the attributes that I need are embedded in the html code. The screenshot below shows a "Prj_type" field with a value of "Sediment Diversion" for the highlighted feature. I need a "Prj_Type" column that would have a value of "Sediment Diversion" for this feature:

(The attribute table viewer in ArcGIS only shows the first line, so I exported the table to Excel to see the full html code.)

Is there a way to parse the html in this attribute column to create new columns based of the fields and values in the html?


texts: Unable edit anything in geodatabase using ArcMap?

texts:
I'm using ArcMap 10.5.1. A client sent me a geodatabase (Project.gdb) and asked me to update/add shapefiles and information in the attribute tables. I can start Editing but shapefiles in the geodatabase don't come up in the Create Features box and I can't edit any of the attribute tables.

Is there something I have to do to "unlock" a geodatabase?

I am not accessing the geodatabase through any other programs - just ArcMap. When I look at the geodatabase files through Windows Explorer, my computer is the only one with a lock on it.


texts: Revert QGIS2.6 installation to QGIS 2.4 through OSGeo4W Advanced installer

texts:
It is possible to return to QGIS 2.4 OSGeo4W installation when I have QGIS 2.6?

I want to do this with advanced installer.


texts:
I am using QGIS 3.10.1, and I am trying to fix invalid geometry in a number of polygons.

This polygon (orange) surrounds a road, but also has a sliver going straight across. As a fairly new QGIS user, I don't know how to get rid of it - I have messed around with the Editing tools, but can't find anything that works.


texts: Extracting coordinates of random points using QGIS?

texts:
I have successfully created a layer that contains a variable number of random points per polygon by using Vector > Research Tools > Random points inside polygons (variable). As I needed 500m distance between the points in each polygon, I used a UTM coordinate reference system for the original layer, which was of course continued on to the random points layer.

Now I would like to extract the GPS coordinates of the random points.

I thought I could do this via:

1) Right click on 'Random points' (the name of the points layer)

2) Click 'duplicate'

3) Right click on the duplicated layer and click 'Save As'

4) Choose the CRS to be EPSG 4269 and save.

However, the 'duplicate' option is greyed out, even when I enable editing. I tried this solution and this one.

How can I get the coordinate information from the random points, preferably extracted to a csv or similar?

EDIT: Better way to phrase the question: how do I convert from EPSG:32614, WGS 84 / UTM zone 14N to EPSG: 4326; WGS 84?


texts: Zoomify and Vector Layer Alignment in OpenLayers 3?

texts:
Need help figuring out how to make a vector overlay for a zoomified image using OL3.

Am updating a project from OL2, where what I'm attempting worked fine. Project isn't true map data, the zoomify layer is a scanned document, converted to zoom tiles.

Polygons represent the coordinates of terms/phrases in text found in the OCR text.

This fiddle shows a page with highlighting for 3 hits on the term "texas". Note the position of the terms is correct on one axis, but simply too high on the screen.

I'm sure this is easy, but I'm at a loss as to the problem (if it has to do with extents or projections, or how I'm centering the layer, etc.)


texts: How can I add a featureCollection to another one to export into one sigle .csv file?

texts:
Sorry if this is an silly question, as I am starting to work with GEE.

I am iterating (mapping) over a Sentinel-2 image collection, using different coordinates and dates in each step.

My goal is to download the mean band properties of the available images in each step, but I want to "join" the FeatureCollections obtained in each step to avoid downloading multiple files.

How should I proceed?

Here is some code example:

var S2 = ee.ImageCollection("COPERNICUS/S2");

// List with LAT, LON, and dates

var list = [
[42.454,3.543, '2017-05-23','2017-06-24'],
[42.7,3.8, '2018-04-19','2018-05-20'],
];

var loop = function(point)
{
  var x = point[1];
  var y = point[0];
  var date0 = point[2];
  var date1 = point[3];
  var point = ee.Geometry.Point([x,y]);

  //Filtering
  var S2TemporalFiltered = S2.filterDate(date0,date1);
  var S2SpatialFiltered = S2TemporalFiltered.filterBounds(point);

  //Box coordinates to reduce collection
  var x1 = x - 0.002;
  var x2 = x + 0.002;
  var y1 = y - 0.001;
  var y2 = y + 0.001;
  var box = ee.Geometry.Rectangle(x1,y1,x2,y2);

  //Cutting collection to box
  var S2Clipped = S2SpatialFiltered.map(
    function(img) {return img.clip(box)}
    );

  //Mean values
  var eBox = ee.Feature(box);
  var func = function(image){
    var meanList1 = image.reduceRegion({
      reducer: ee.Reducer.mean(),
      geometry: eBox.geometry(),
  });
    return ee.Feature(null, meanList1);
  };
  var total = S2Clipped.map(func);

//Now I need to add "total" to a "matrix" where
// band properties will be added in each step of the loop.
//Each line relates to an image

//This next part should be out of the loop, dowloading all images properties
//at once (instead of having two files, I would have one file with two rows)

Export.table.toDrive({
collection: total,
description: 'Total',
fileFormat: 'CSV',
});

};

var a = list.map(loop)


texts: Given Lat and Lon identify if point over land or ocean using R

texts:
In R, how do I identify if a point is over land or ocean?
e.g. for:
set.seed(0)
lat <- runif(10, -80, 80)
lon <- runif(10, -180, 180)
points <- expand.grid(lat, lon)

I would like to create a column called "Land" of type "logical".


texts: Does ArcPy work for Engine or Server the same as Desktop?

texts: I am contemplating a project where there is not going to be access to ArcGIS Desktop and I amwondering if ArcGIS Engine and/or ArcGIS Server install ArcPy and Python?


texts: Producing CSV file from .tif file with elevation data

texts:
I have a .tif file that contains all the elevation of a country.

How can I convert that to a CSV file that contains x,y,z values where x and y are the latitude longitude respectively and the z column is the elevation using QGIS?


texts: Are there more events in OpenLayers 3, besides "moveend"?

texts:
I am building an API on top of OpenLayers 3. I can't display here the whole code for professional reasons. This is a bit how my classes look like.

// external file xgis.js

xgis = {};

// xgis.map
xgis.map = function(options) {
    // map definitions ...
};
// Inherits from ol.Map
goog.inherits(xgis.map, ol.Map);

The goal is to build an SDK, which I named here "xgis". We want to implement some callbacks. For example, when I perform a zoomIn, I have a callback named onZoomIn that executes some code. I know that it is done with the event "moveend", like in this fiddle.

Now, I want to implement some callback functions:

onZoomIn --> fires after a zoom in
onZoomOut --> fires after a zoom out
onMapReady --> fires when the map is ready (has loaded all the layers)
onPanEnd --> fires when the drag/pan has ended

Are there other events in OpenLayers 3, like "zoomend" or something?

Thanks


texts: Building Final Web App from QGIS Web App Builder?

texts:
I have used Web App Builder plugin for QGIS from Boundless to build a web application, but the application only shows in a browser when I click Preview.

The application created using Create App does not open up in browser. It says Build final app using Boundless webSDK.

I have read some documentation about the SDK, but I've got no idea how to proceed.


texts: Identify coordinate system used near Kuala Lampur

texts:
I have some coordinates that were given to me, and I do not know much about how they were taken. I thought they were UTM coordinates, but when I try to convert them using some online tools I always obtain a value that is far away from the supposed point. That makes me think that they are not UTM coordinates, so I need some help to identify what system is used.

Here is an example coordinates data straight from the file I was provided, copied exactly as it was written. The positions are located in Kuala Lumpur, Malaysia, and were surveyed around 1997; I used 47N as UTM zone.

            Provided  Converted
easting   106374.634   0.874278 latitude
northing   96819.608  95.464439 longitude

(Converted is even out of the 47N zone!)

   Expected (approx)  converted in UTM
latitude    3.145241  easting  806771.55
longitude 101.759941  northing 348052.77

Any idea on what kind of coordinates are being used here?


texts: Failing to Extract by Mask with Small Groundwater Basins?

texts:
I am using ArcGIS for Desktop.

I am working with 4km resolution precipitation rasters from the PRISM website because I am trying to calculate the total annual precipitation for each groundwater basin I am interested in.

I have been using the "Extract by Mask" tool to try and accomplish this, but the tool has not been able to clip the raster data to the very small groundwater basins in my dataset.

Is there another tool I can use to accomplish this without altering the resolution of the rasters?


texts:
I know that mapserver drop the support for ecw due to licensing issues. But if i understand well, as mapserver rely on GDAL, if GDAL support ECW, mapserver can serve ECW data.

So I've compiled libecwj 3.3 , then compile GDAL. Now GDAL can read ECW.

But in mapserver, ecw data don't work. The error is :

  MsDrawMap(): Image handling error. Failed to draw layer named
  'ortho'.; MsDrawRasterLayerLow(): Unable to access file. Corrupt,
  empty or missing file 'my_file.ecw' for layer 'ortho'.

I tried to recompile mapserver. But with no effect. Does someone know if it's possible and where is the mistake ?

My config is :
- debian 7
- libecwj 3.3
- gdal 1.10.1
- mapserver 6.4.0

Thanks for your help !


texts: Installing SAGA in QGIS 3.12 on Ubuntu 18.04.04

texts:
I have installed QGIS 3.12  through the Ubuntu Store and installed SAGA 7.3.0
through:
sudo apt-get install saga

I tried to install QGIS through apt-get, as instructed in the documentation of QGIS but this resulted in many errors on dependencies which were impossible to install. I did not run into these issues when installing QGIS from the Ubuntu store.
However, I cannot find the SAGA tools in QGIS. I ticked the activate box in the Processing settings but still, no sign of any SAGA tools in the processing toolbox.
Printscreen of Processing Options


texts: How do I interface with bundles using ESRI's Javascript API?

texts:
I've found this example... http://resources.esri.com/help/9.3/arcgisserver/apis/javascript/arcgis/help/jssamples_start.htm#jssamples/layers_custom_tiled.html

Which uses an exploded cache. It notes at the beginning of the article...

"You would most likely use ArcGISTileMapServiceLayer when accessing a service from 9.3 and beyond"

Digging into the ArcGISTiledMapServiceLayer object it seems like it only connects to services via the rest url. I need to accomplish what is achieved in the Portland example (the link above) but with cache bundles instead of the exploded cache.

Potentially Useful Background: I need to make an entirely self contained web map app with the javascript api that would use bundles that are stored with the app (as opposed to out on a server with arcgis installed) to render the map.


texts: Windows script/application to convert GPX to KML?

texts:
OSRM seems to have everything I need to draw bicycle routes, but it only exports routes as GPX while the smartphone application I use only supports KML.

Does someone know if there's a Windows script or application that I could use to easily convert GPX to KML?

I know the Gpx2Kml site, but I'd like a one-click solution.


texts: Can't change data range on r.scale in qgis

texts:
I'm trying to rescale a raster to values of 0-1 in order to do a multicriteria analysis in QGIS. I have a raster with integer values from 0-18. However, when I come to try and change the data input range, it wont let me move it past 1.

Any ideas on why this might be happening?


texts: How to add a raster layer in openlayer 6 and its wmts in Mapproxy

texts:
I want to add a raster layer to my map. I added a TMS layer map using ol6, now I want to add a raster layer to it. How could I implement it. My TMS layer is generated using Mapproxy & mapnik.
My javascript is :
proj4.defs("EPSG:32643","+proj=utm +zone=43 +datum=wgs84 +units=m +no_defs");
register(proj4);

var extent = [291627,904686,958569,1571628];
var projection = new Projection({
  code: 'EPSG:32643',
  extent: extent,
});
var map = new Map({
  target: 'map',
  layers: [
    new TileLayer({
        title: 'District',
        type: 'base',
        source: new XYZ({
           url: "http://127.0.0.1:8080/tms/1.0.0/district/distgrid" + "/{z}/{x}/{-y}.png",
           projection:projection,
           tileGrid: createXYZ({
                  extent: extent,
                  maxResolution: 2605.2421875,
                  tileSize: [256, 256],
                  maxZoom: 13,
                }),
        }),
    }) ],

  view: new View({
    projection: projection,
    units:"metric",
    extent: extent,
    zoom:0,
    maxZoom:13,
    minZoom:0,
    maxResolution:2605.2421875,
    center:[654496.136597752,1155181.26900064],
    numZoomLevels:13,
    constrainOnlyCenter: true,
    constrainResolution: true,
  })
});

My mapnik.yaml is:
services:
  tms:
    use_grid_names: true
    origin: 'sw'
  wmts:
    kvp: false
    restful: true
layers:
  - name: district
    title: District Layer
    sources: [dist_cache]

caches:
  dist_cache:
    grids: [distgrid]
    sources: [dist_tms]
    meta_size: [4, 4]

sources:
  dist_tms:
    type: mapnik
    mapfile: /var/www/Mapfiles/district.xml
    transparent: true

grids:
    distgrid:
        srs:  'EPSG:32643'
        bbox: [291627,904686,958569,1426831]
        bbox_srs: 'EPSG:32643'
        origin: 'sw'

globals:

Can I combine TMS & WMTS in ol6? I want to add as in this example https://openlayers.org/en/latest/examples/shaded-relief.html. How to add wmts layer in mapproxy?


texts: CQL filter causing slowness in Java \ Geoserver app

texts:
I have a Java application that includes an OpenLayers map with WMS tiles pulled from Geoserver. Geoserver and the Java app are looking at the same MySQL database. In the Java application I have a complex filtering system to select data and I would like only the selected data to be shown on the map.

I have tried to do this by passing a CQL filter of the form "ID in (123,3000,6001,5464, etc...)" unfortunately Geoserver 2.3 takes around 12 seconds to render each tile as the list of ID's is so long and I can expect this to get worse as more data is added.

I cannot pass through the SQL query that my Java app is using as it is a complex join that pulls in data from tables not available to geoserver. Also there are security issues as I do not want to expose anything more than the geometry and the feature ID to the public internet via geoserver. Is there a better way of passing the ID's across that is faster?


texts:
I have three feature classes [CircuitSource(A),CircuitBreaker(B),Sub-Station(C)]. In the Geodatabase A and B has a 1:1 relationship based on a primary field namely 'FEEDERID'. C and B has a 1:many relationship based on a field namely 'STATIONID'. I have created a new feature in A(pFeatA) and related this feature to an existing feature(pFeatB) in B feature class so that the value of 'FEEDERID' in pFeatA gets populated in pFeatB. I used the following code.

IRelationshipClass relClass = null;

relClass.CreateRelationship(pFeatA as IObject, pFeatB as IObject);

But the value of 'FEEDERID' in pFeatA is getting populated in all features in B feature class which has 'FEEDERID' = null and the similar value of 'STATIONID' as that of pFeatB.


texts: Problem with calculating the total cell values

texts:
aa= arcpy.RasterToNumPyArray("water_depth")
      sum = numpy.sum(aa)
      print sum
      -2933737

Why I am getting negative values ? In my raster there is no any negative depth.


texts: How to convert geojson to polygon in Leaflet?

texts:
I have this json data:

var states = [{
    "type": "Feature",
    "properties": {"party": "Republican"},
    "geometry": {
        "type": "Polygon",
        "coordinates": [[
            [-104.05, 48.99],
            [-97.22,  48.98],
            [-96.58,  45.94],
            [-104.03, 45.94],
            [-104.05, 48.99]
        ]]
    }
}
}];

It'is easy to visualize this on leaflet map by doing this:

L.geoJSON(states, {
    style: function(feature) {
        switch (feature.properties.party) {
            case 'Republican': return {color: "#ff0000"};
        }
    }
}).addTo(map);

But i need the "type" of states is "polygon" not "geojson", so the result will be like this:

var states_polygon = L.polygon([[
            [48.99,-104.05],
            [48.98,-97.22],
            [45.94,-96.58],
            [45.94,-104.03],
            [48.99,-104.05]
        ]]).addTo(map);

I knew that both variable will be give the same visual on leaflet map but all i need is the object type of that variable as polygon.. How can i do this?


texts: Converting tiff file (DEM) into XYZ file with FME?

texts:
I want to convert a tiff file (DEM) into an XYZ file with FME.  I have tried using the sa - sample tool in ArcGIS and it took 19 hours to accomplish 9% of the task.  I saw a post by Mark Ireland that mentioned using FME to accomplish this task much faster.  Here is an image of the post:

This looks relatively simple, however I am relatively inexperienced with FME and was wondering anyone can help decipher how to accomplish this task by looking at the FME diagram he posted.


texts:
I originated some points with incorrect geometries (ex: 1432.33 , 900.2 for EPSG4326) because lat/lon in data source where also incorrect.

How can I detect in the_geom these incorrect ones?

I think ST_IsValid is not useful, as the ST_Geometry is well formed, even if its data inside is not.


texts: Converting GCS NAD 1983 to CA state plane

texts:
I have a polygon shapefile in GCS_North_American_1983. I'm trying to project it into
NAD 1983 StatePlaneCalifornia III Fips 0403 (US Feet).

This is my original file.

When I run the project tool, this is what happens

What's going wrong? Why isnt it working correctly? What is happening?


texts: Add Rasters to a new GDB/Mosaic Dataset in ArcGIS ModelBuilder

texts:
I'm attempting to add rasters to a new mosaic dataset in a new file geodatabase, all in the same model.  However, the Add Rasters To Mosaic Dataset tool does not recognize the NewMosaicDataset or the NewGDB.gdb because they haven't been created yet.  The first time the model is run, the geodatabase and mosaic dataset are created but no rasters are added, but the Add Rasters tool then recognizes an existing geodatabase and mosaic dataset.  The second time running the model is unsuccessful because the geodatabase already exists and it will not overwrite it.

This is the error produced when I open the Add Rasters To Mosaic Dataset tool.
(Ignore the digit 1 after NewGDB in the error message, I simply entered a new output name in order to take a screenshot of the error message I receive)

I have attempted to set the output mosaic dataset as a precondition as well as exporting the model to a Python script and adjusting variable names, but both attempts were unsuccessful for the same reasons described above.

Would like to know how to do all of this in one model.


texts: Understanding ArcGIS network service areas

texts:
I want help understanding the output from network Service Area analysis in ArcMap 10.8.1.
The example below uses the data here https://sites.tufts.edu/gis/files/2013/11/Powerplant_Proximity_Exercise.zip used for the Tufts Data Lab's Proximity Tutorial, (https://sites.tufts.edu/gis/files/2019/10/Proximity-Tutorial_Nuclear-Plant-Risk-Assessment_10.7.1.pdf.
This image shows Service Areas at 12 and 50 mile breaks and a 3 mile tolerance, as generated with the ArcGIS Network Analyst tools (using the settings specified in the Proximity Tutorial):

What I want to understand is where the 0 to 12 mile break polygons come from.  They aren't enclosed by any of the streets in the network, so how does Network Analyst determine where to locate the boundaries?


texts: What is the coordinate system of this maritime natural gas field map?

texts:
The Israeli Ministry of Natural Infrastructure, Energy and Water Resources published a maritime map of the borders of a new natural gas field.

The map is a part of a full license (PDF), a legal document written in Hebrew.

This map has enormous importance to several environmental and social NGOs. Some of their activists have tried to georeference it, but could not figure out the coordinate system and asked for my assistance.

The coordinates are indeed odd - the first number has 5 digits, while the second has 6 (e.g. 65395/776428). It resembles UTM, but I could not figure how to find the exact zone or why is there a difference in the number of digits.

Could you please help me identify the coordinate system in this map? The location is in the Mediterranean sea, a few dozens kilometers off the Isreali coast.


texts: TrainImagesClassifier (SVM) fails when executed from within QGIS modeler

texts:
I made a small model in QGIS that runs Second Order Image Statistics, the SVM Trainer, and Image Classification, all tools from the Orfeo provider.

The problem is it fails everytime. When I run the tools outside the modeler, using the processing toolbox (not the command line), they run perfectly, but when I put them in the modeler, they fail. More specifically, the trainer does not export the model nor the confusion matrix, which causes the next step, Image Classification to fail.

Is this a problem that has been encountered before? Does anybody have a suggestion on how to solve this?

Update: Below is the log output when the model atempts to execute the Trainer.

2016 May 19 13:47:03 : Application.logger (DEBUG) Entering DoExecute
2016 May 19 13:47:03 : Application.logger (INFO) Elevation management: setting default height above ellipsoid to 0 meters
2016 May 19 13:47:05 : Application.logger (INFO) Number of training samples: 1125
2016 May 19 13:47:05 : Application.logger (INFO) Size of training set: 1125
2016 May 19 13:47:05 : Application.logger (INFO) Size of validation set: 1058
2016 May 19 13:47:05 : Application.logger (INFO) Size of labeled training set: 1125
2016 May 19 13:47:05 : Application.logger (INFO) Size of labeled validation set: 1058
CvSVM::RBF = 2
CvSVM::C_SVC = 100
OpenCV Error: Bad argument (There is only a single class) in cvPreprocessCategoricalResponses, file C:\OSGeo4W64\usr\src\osgeo4w\opencv\opencv-2.4.11\modules\ml\src\inner_functions.cpp, line 729
2016 May 19 13:47:05 : Application.logger (FATAL) The following error occurred during application execution : C:\OSGeo4W64\usr\src\osgeo4w\opencv\opencv-2.4.11\modules\ml\src\inner_functions.cpp:729: error: (-5) There is only a single class in function cvPreprocessCategoricalResponses


texts: Reclassifying based on neighbour polygons using QGIS?

texts:
I have a polygon shapefile and each one has one class attached. I want to reclassify one of them based in the neighbours.

For example, I have a polygon with class street and I want to reclassify to a house if one of the polygons around it (or one that contains it) belong to class urbanization in QGIS.


texts: Exporting ArcMap Feature to KML while retaining Attached Photo?

texts:
I have a set of point features in an ArcGIS file geodatabase. I have attached a photo to one of the points, using the built-in tools 'enable attachments' and 'add attachment'. The photo is stored on my local drive and should be copied to the GDB when 'add attachment' is run.

If I click on an HTML pop-up in ArcMap the image is correctly associated with the point. However, when I export the layer to a KMZ ('layer to KMZ'), the image does not appear in the resulting KMZ file. Oddly, the photo file name is correctly added to the KMZ pop-up.

Note: the KMZ export tool otherwise correctly works. The photo is less than 1 MB so size should not be an issue. I have found tutorials that suggest a similar workflow, but do not encounter this problem. I also manually inserted the full photo filepath in Google Earth by right-clicking on the point, editing its properties, and saving the KMZ again. Unzipping the KMZ shows that the image is not saved with the KMZ though.

Has anyone else run into an error with this ArcGIS tool or Google Earth and found a work around?


texts: How to automatically add identical new polygons to a layer?

texts: I am seeking a way, using QGIS, to automatically add polygon features identical to others on the layer. In this case, they are water valves on a water network. Anyone know of a user-friendly way to add new valves that visually line up with the pipes without needing to manually draw the polygons through editor mode?


texts:
I am trying to split up the world into 2 categories based on elevation (below 2000m and above 2000m) while still including a country identifier.  So in the end, I would like to have 2 types of units (low areas and high areas) for each country in the world, then I will use Spatial Analyst's Zonal Statistics to calculate the mean values of other rasters of interest in my defined areas.  Does anyone have a suggestion for how to do this?

My data:
- 16 DEM tiles for the world
- 1 country borders shapefile

I am struggling with:

1)...working with all the DEM data at once on a global scale because I have been reading the mosaic-king across UTM zones distorts the data.

2)...whether my final units should be in raster form/ file geodatabase grid.  I have heard that I should work in a file geodatabase and forget trying to have my final units be in shapefile form.

I know how to reclassify and have been successful at obtaining a raster with only two categories of elevation for each DEM tiles separately (above and below 2000m).

Any advice would be greatly appreciated,
Liz


texts: Extracting RGB values from symbology of .lyr or polygon file using ArcGIS Desktop?

texts:
I have a polygon shapefile of geological units. Someone has already created a symbology and assigned standard colours to the units, available as a .lyr file.

What I need to do is create a list of the the 40+ units with their assigned RGB values using an ArcMap tool or maybe ArcPy.

I have been able to use polygon to raster, but it doesn't preserve the correct RGB values. That could be a feasible option and allow me to sample the colour if anyone has a suggestion on that path.

I have ArcGIS Desktop 10.5 Advanced License.


texts: Null geometries when importing CityGML data into ArcMap and QGis

texts:
I am looking to convert a CityGML file into a shapefile to perform some spatial analysis on it. I am completely new to the CityGML format and keep running into issues when trying to import them in either ArcMap or QGIS.

In QGIS (Version 2.18.2): I have tried to import the .gml via 'add vector layer'. The validity check states 'unknown geometry type 0'.

In ArcMap (Version 10.3.1): I tried both data interoperability tools > quick import and interoperability connections > add interoperability connection. I still end up with No Geometry.

Is there something wrong with the structure of the CityGML?


texts: Using ArcGIS Online basemap from PyQGIS says connection is blocked (Port:80)?

texts:
I'm trying out some code from the book "QGIS Python Programming Cookbook Second Edition" and would like to add a tiled imagery map to my QGIS Application.
I minimally adapted the code from the book and it now looks like this:

arc_map_server = "http://server.arcgisonline.com/"
arc_map_server += "ArcGIS/rest/services/World_Imagery/MapServer?f=json"
basemap = QgsRasterLayer(arc_map_server, "World Imagery")
QgsMapLayerRegistry.instance().addMapLayer(basemap)

When I use the same URL in QGIS normally (without the ?f=json) and add the layer, it works without problems. Adding shapefiles as vector layers and other raster layers inside the standalone QGIS application works as well without problems.

But when try to run the above code inside the application, I am blocked on Port: 80. Is that on purpose, or am I doing something wrong?

If it is on purpose, is there an easy way to add a tiled imagery map to my QGIS standalone application?

OpenLayers I used when building plugins, but so far I haven't been able to find a way to reliably access it with the application.


texts: How to convert Rainfall Isohyet map into Raster to determine precipitation (mm) at any point?

texts:
Specifically, I have a rainfall map (Isohyets) measured in mm (precipitation)... I wanted to know a particular value (precipitation) at any point between these lines...

How can I do this? Do I have to convert/interpolate these polylines so I could have a Rainfall raster instead? So that I would know the precipitation in any point based on cell values of the Rainfall raster? If so, how do I make a raster map out of my rainfall contours/isohyets? By the way, I have a DEM of 20-m resolution.

Will this be of any help in converting/interpolating my rainfall map?


texts: Convert GDA eastings/northings to WGS84 lat/long in QGIS

texts:
I'm using QGIS to map weed locations and have a mixture of data from GPS, shapefiles of contours, OSM all in a project with CRS WGS84 (lat, long). I am trying to add in a shapefile I have been given which is in GDA eastings and northings but can't get it to mix with the other layers.

Lat for this area is ~-35, long about 149. Can't find an existing CRS definition that works - can get the right lat but the wrong long (~119). What do I do?

Hope this question makes sense as I'm pretty new to this GIS stuff.


texts: Fiona: won't install in Debian7.x-64bits/Python2.7

texts:
I need to install Fiona in my Debian7.x-64bits/Python2.7 machine. I already installed GDAL/OGR, pip, virtualenv, and six. I am using Python 2.7 so argparse is also installed. The installation protocol provided at the Fiona website is not completely clear to me:

$ mkdir fiona_env
$ virtualenv fiona_env
$ source fiona_env/bin/activate
(fiona_env)$ pip install Fiona

(Source: https://pypi.python.org/pypi/Fiona )

I followed the steps from both a non-root- and a root-terminal. In the first case some problems related with permissions arose. In the second the errors were related with an inability to find OGR directories.

I am running the installing protocol from my home directory, e.g. me@mymachine:~($ or #). I am not sure if that is the problem. I already checked around for documentation in the Internet, but it seems that the installation protocol from Toblerity (referenced above) is the most explicit.

I would appreciate the community help with the installation of Fiona in my rather standard Debian machine.


texts: Controlling layer order between cartoDB and Google maps Fusion layer

texts:
I am trying to control the layer order between a cartoDB layer and a Google FusionTables layer, on a google base map. The CartoDB layer always displays underneath, whereas i want it on top. I have found methods to order layers within Google, or within CartoDB, but not how to order them with respect to each other.

I have mocked up an example here:
http://codepen.io/rich_3po/pen/Bjoadr

(the yellow points are from CartoDB, and should display on top of the routes)

How do i go about achieving this?

My question is similar to Ordering Cartodb torque visualization over Google Maps heatmap layer, but the answer wasn't clear and didn't seem to apply to my case


texts: What's the absolute minimum information needed to have a valid coordinate?

texts:
Don't assume anything. You read a specific coordinate value.

What's the minimum information you need for unambiguous interpretation of such coordinate?

Does a coordinate has sense without specifying its Geographic Coordinate System? Should I know its notation? Anything else?


texts:
We are developing a QGIS Plugin.
Considering we can only use Python, is there any way to protect our intellectual property?


texts: Model Builder ArcGIS - Export topology errors (Iteration Tool)

texts:
I was wondering if there’s any option to iterate export topology errors in ArcGIS model builder. My GDB Topology errors has all topologies within the same Topology dataset but model builder seems not have the option to go throughout the same topology dataset and export the errors (line, point and polygon) for each topology file. Here the structure of the topology errors GDB:

My idea is iterate each topology file (e.g. T_RailwaysAlignment3DCentreline_Ln) and export them with same name (e.g. T_RailwaysAlignment3DCentreline_Ln_Line).

The following screenshot shows the test model prepared for this purpose:

Iteration Tool – Iterate Dataset

Export Topology Lines

Any ideas how export the topology errors (within the same topo dataset) using Model Builder? (Python ideas more than welcomed too)


texts: Calculating Focal Statistics for Special Neighborhood?

texts:
I'm looking to calculate focal statistics for each cell of a raster, within a neighborhood of a specified criteria.

Background - I have three binary rasters, each representing a single vegetation type of interest. I'd like to calculate the percent coverage of each vegetation type within (e.g.) 20 km^2 of any cell in my study area (sum/total cells in neighborhood). The problem is that I can't use a simple circle or square neighborhood around each cell because, if I did, the search area used to calculate the sum would incorporate areas outside my study area. This exception is important because the statistics will be used as inputs for a habitat model, and the areas outside of my study area cannot be considered possible habitat - they're urbanized. Including them would give me erroneous statistics. So, what I'm looking to do is write some code in python that will choose a neighborhood representing the n nearest cells (n determined by number of cells required to cover an area equal to my desired neighborhood size) that meet my criteria. The criteria being that they do not fall within an urbanized area. I'm thinking that some form of cellular automata should be used. I've never worked with CA though.

I guess what I'd like is something like starter code, or a point in the right direction.

REPLY TO COMMENT BELOW:

Let's say I'm calculating this statistic for a cell on the boundary of my study site. If I assign all areas outside of my study area to zero (or ignore NoData), then I will get a statistic that represents roughly half of the areal coverage I'm interested in. So, percent coverage in a ~10 km^2 area, instead of 20 km^2 area. Since I'm studying home range sizes this is important. The neighborhood has to change shape, since that is how the animal views/uses the landscape. If they need 20 km^2, they'll change the shape or their home territory. If I do not check ignore NoData, cell output will be NoData - and NoData is no help.

"PROGRESS" AS OF 10/24/2014

Here is the code I've come up with so far using Shapely and Fiona:

import numpy as np
import pprint
import shapely
from shapely.geometry import*
import fiona
from fiona import collection
import math

traps = fiona.open('C:/Users/Curtis/Documents/ArcGIS/GIS_Data/occurrence/ss_occ.shp', 'r')

study_area = fiona.open('C:/Users/Curtis/Documents/ArcGIS/GIS_Data/Study_Area.shp', 'r')
for i in study_area: #for every record in 'study_area'
        sa = shape(i['geometry']) #make a variable called 'sa' that is a polygon

grassland = fiona.open('C:/Users/Curtis/Documents/ArcGIS/GIS_Data/land_cover/polys_for_aa/class3_aa.shp', 'r')
pol = grassland.next()
gl = MultiPolygon([shape(pol['geometry']) for pol in grassland])

areaKM2 = 20
with traps as input:
    r = (math.sqrt(areaKM2/math.pi))*1000
    for point in input:
        pt = shape(point['geometry'])
        pt_buff = pt.buffer(r)
        avail_area = pt_buff.intersection(sa).area
        # works to here
        while avail_area < areaKM2:
            r += 10
            pt_buff = pt.buffer(r)
            avail_area = pt_buff.intersection(sa).area

        perc_cov = pt_buff.intersection(gl).area//areaKM2
        print perc_cov

Unfortunately, it's INCREDIBLY slow.


texts:
Is there any OS library to read/write/process RINEX data?

I've been wanting to start a OS library to read and write RINEX, which is something I need for a project, but not actually post-process them.

Any tips?


texts: Google Earth Engine filterdate ignored by plot?

texts:
I'm starting with Earth Engine (EE) coding. Following instructions from https://developers.google.com/earth-engine/tutorial_api_07 I was able to put together some code and get a plot. But why the range of date on the plot goes from 2016 to 2018 while filterDate('2017-01-01', '2017-12-31')?

var image = ee.Image(
  s2.filterBounds(point)
    .filterDate('2017-01-01', '2017-12-31')
    .sort('CLOUD_COVER')
    .first()
);

var addNDVI = function(image) {
  var ndvi = image.normalizedDifference(['B5', 'B4']).rename('NDVI');
  return image.addBands(ndvi);
};

var ndvi = addNDVI(image).select('NDVI');

var withNDVI = s2.map(addNDVI);

var chart = ui.Chart.image.series({
  imageCollection: withNDVI.select('NDVI'),
  region: point,
  reducer: ee.Reducer.first(),
  scale: 30
}).setOptions({title: 'NDVI over time'});

print(chart);


texts: Use of US Census data to plot house value range

texts:
The house value data from the US Census Bureau is provided as a table differentiating areas into census tracts. My understanding is that they categorize each tract into quantities of homes fitting within a certain value range and provide an estimate of the total number of homes within each tract. For example, tract 123 has an estimate of 200 total houses and 5 houses ranging from 10k to 15k, 6 houses ranging from 15k to 20k, and so on.

I am wanting to derive a map illustrating the tracts fitting a value range. For example, a chloropleth map displaying the tract differences among houses valued within the range of 10k to 15k. In order to do this, wouldn't I first need to take the number of houses within a particular range and divide it by the total number of houses estimated for the same tract? It seems simply plotting the number of homes within a certain value range for a particular tract would mean I was plotting unnormalized values and thus an inaccurate portrayal.

I wanted to make sure I was understanding what the raw census data represented and how to properly use this data.


texts: Integrating R in QGIS with Mac OS El Captain: missing dependencies

texts:
I have a Mac OS (El Captain) and I'm trying to run R scripts within QGIS 2.14. When I open a script (like one from the online list) I receive the error message:

  This algorithm cannot be run :-(
  This algorithm requires R to be run. Unfortunately, it seems that R is not installed in your system, or it is not correctly configured to be used from QGIS
  Click here to know more about how to install and configure R to be used with QGIS

Please note that:

R is correctly installed
I followed instructions from http://osgeo-org.1560.x6.nabble.com/R-in-QGIS-on-OSX-using-Kyngchaos-build-td5254567.html
It is not possible to set the R folder, as I did in Windows (see the figure here below), but only the R script folder.


texts: Shifting a CAD file from ATS77 MTM 5 to NAD83 CSRS UTM 20

texts: Looking to make the coordinate system shift in ArcMAP (for Nova Scotia) between ATS77 MTM 5 and NAD83 CSRS...NTv2 file?


texts: How to remove outliers from spatial polygon

texts:
A source shapefile contains clear outlying points/area. How can they be isolated and removed from the spatial polygon in R?

library(rgeos) #required by maptools
library(maptools)

incheon <- readShapeSpatial("KOR_adm1.shp") #upload Korean provinces
incheon <- incheon[incheon@data$NAME_1=="Incheon",]
plot(incheon) #plot reveals clear outlier in the southeast corner of the map


texts:
I want to install gsutil onto my QGIS's version of python.

Most tutorials say I should be able use the pip/pip3 command.

When I run OSGeo4WShell, it does not have the pip command available. How can I download pip for my QGIS's python?


texts: ASCII Grid to Raster with DX/DY Cell Values

texts:
I have an ASCII grid file that has rectangular cells. The ESRI documentation does not specify how to convert from .asc to raster with a DX/DY value such that:

NCOLS 1900
NROWS 1100
XLLCORNER 115.834717
YLLCORNER 31.638744
XCELLSIZE ...
YCELLSIZE ...
NODATA_VALUE 0

I do not want to read into a NumPy array. So how would one do this, be it manually, in ArcGIS?


texts:
I've noticed that some SVGs are distorted when outputting to a vector PDF from QGIS.  A good example is NorthArrow_07.svg.

The image on the left is how it looks in the composer and the image on the right is the PDF output.

I have tried increasing the composer DPI but it has no effect.  Is there a setting that causes this distortion?  I need vector output - not raster.


texts: Extract 3 colors from JPEG and change it by 3 other

texts:
I would like to import a JPEG image, find and extract 3 colors and change by 3 others to use as ID.

I can find the 3 colors I want by using Raster and I save the 3 images associated to these 3 colors. Then, I open these 3 JPEG images to combine in 1 using EBImage package. The image I got is not nice, the colors are indefinable.

I do not find the way to change these 3 colors in 3 others (by example I would like finally just red, green and blue to identify 3 items)

##################################################
#Plug it into a stack object using package raster

img <- list.files(folder)
img.raster<-stack(img)
names(img.raster) <- c('r','g','b')

## Extract blue
img.blue<-img.raster
img.blue$blue<-0
img.blue$blue[img.blue$r<60 & img.blue$g<80 & img.blue$b<60] <- 1

## Extract grey
img.grey<-img.raster
img.grey$grey<-0
img.grey$grey[img.grey$r>85 & img.grey$g>82 & img.grey$b>77] <- 1

## Extract green
img.green<-img.raster
img.green$green<-0
img.green$green[img.green$r>65 & img.green$r<120
                  & img.green$g>90 & img.green$g<160
                  & img.green$b>45 & img.green$b<70] <- 1

#Export images using JPEG package
setwd("C:/Users/*******/export")
writeJPEG(as.array(img.grey), "grey.jpeg")
writeJPEG(as.array(img.blue), "blue.jpeg")
writeJPEG(as.array(img.green), "green.jpeg")

#Open images using Using EBImage package
x <- readImage(file.choose())
y <- readImage(file.choose())
z <- readImage(file.choose())

#put 3 images in 1
total <- x + y + z

#Export images
writeJPEG(as.array(total), "total.jpeg")

Is it the good way or can it be done easily?


texts: Create a buffer with a pre-defined area?

texts:
I have a road (line) system and need to create a polygon from it to show the right-of-way area. The right-of-way area is calculated in the field from taking measurements from the road centreline. From these measurements a total right-of-way area is calculated (m2.)

I am wondering if it is possible to  create a buffer with a pre-defined area based on the measurements taken in the field. Possibly through a Python script?

Without a tool to do this, it takes some time as I have to guess the average right-of-way with and create a buffer using that measurement. So it takes a lot of trial and error time.

I am using ArcMap 10.2.


texts: Weighting a data frame within another frame for map creation

texts:
I've currently got an inset data frame sitting in a corner of a larger data frame. Some labels are being half covered by the smaller data frame. I was wondering if there are any tips on how to weight it (either in Maplex or Standard) so that the labels that are half covered are either removed, or moved so that they aren't behind the smaller data frame?

Thanks


texts: How do you determine how straight a LineString is?

texts:
I have a bunch of LINESTRINGs and need to classify them based on how strait they are. My idea is to compare the length of the LINESTRING to the distance between the first and last POINTs.

I'm thinking of something along these lines (assume the srs units are in meters):

ST_Length_Spheroid(geom) - ST_Distance_Spheroid(ST_PointN(geom, 1), ST_PointN(geom, St_NumPoints(geom)))

Is there a better way to do this?


texts:
Does anybody know where i can find implementation of feature called "Zoom to Feature".

To use this function from GUI you have to right click on layer then choose "Open Attribute Table" then right click on feature and there should be option called "Zoom to Feature" in context menu.

It looks like this:

[IMPORTANT]

Feature in Qgis is a set of attributes from table.


texts: Finding the nearest distance between points in a grid cell

texts:
I'm trying to find the distance between the nearest points in a grid-cell, such that there is no equivalence relation between these points. Meaning if there are three points {a, b, c}, then a = b; not b = a; though b = c; and c = a.
This of course only works if there are an odd number of points in a cell, otherwise, if there's an even number, or a single point in a cell, then pick any point in the nearest cell. Furthermore, these points cannot overlap between regions, where the dataset has points on a map with boundaries.
The output should return a data from likeso:
Points  cells   region  distance
 pt1     1      canada    2000
 pt2     1      canada    1800
 pt3     2      canada    2100
  .      .        .         .

Code that has points (habitat) on a map in North America:
uk.map <- getData('GADM', country='USA', level=1)

file1 <- "modis_mcd12q1_umd_2010.tif"

dir.create("data", showWarnings = FALSE)
if (!file.exists(file.path("data", file1))) {
  download.file(paste0("https://raw.githubusercontent.com/lime-n/test/main/raster/", file1),
                file.path("data", file1), mode = "wb")
}

r <- raster("data/modis_mcd12q1_umd_2010.tif")

us.map <- st_as_sf(us.map)

rr <- us.map %>%
  st_transform(crs = projection(r)) %>%
  rasterize(r) %>%
  # remove any empty cells at edges
  trim()


texts: how to implement pegman (street view icon) in leaflet map?

texts:
i would like to get street view icon (pegman of google map) in leaflet map. A pegman is a small icon in google map at the right corner. I would like to drag that icon on the map to get the street view..

var button = '
  Street View';

var fenway = new google.maps.LatLng(lat, lng);
    var panoramaOptions = {
        position: fenway,
        pov: {
            heading: 34,
            pitch: 10
        }
    };

$('#myModal').on('shown.bs.modal', function () {
        var panorama = new google.maps.StreetViewPanorama(document.getElementById('streetView'), panoramaOptions);
    });

this works i.e.street view icon is seen on click of button but i want on drag of pegman


texts: Why does a single feature's label not show?

texts:
I am working with a GeoPackage layer, and cannot figure out why one of my features' label refuses to be shown. I am using only data-defined labeling for this layer, not auxiliary storage. What things can I troubleshoot to find the cause?

Here's what I've already ruled out:

Rule-based labeling: labels are displayed by a wordwrap() function; the necessary attributes "name" and "wrap" are correctly populated.
"Show all labels for this layer" is selected in Label properties.
"Only draw labels which fit completely" and "Discourage labels from covering" are unselected.
"Show label" and "Always show": tried setting both to TRUE, without success.
Scale-dependent and pixel size-based visibility are unchecked.
The geometry is valid according to both Check Validity and Topology Checker.
"Simplify geometry" is un-checked in Rendering properties. This is normally the culprit when a single feature is missing its label. If I create a new project containing only this layer, the label WILL be displayed, only if I un-check "Simplify geometry". However, in the project I am using, this does not solve the problem, even after deleting all other layers from the project.

If I split the feature into two new features, the label will show only on the northern of the two. If I split it into two parts of one feature, the label does not show on either part, even if "Label every part…" is selected.

What other settings or conditions will cause a label not to be shown?


texts: GetFeatureInfo template of MapServer not replacing values

texts:
I am using MapServer 7.4.2 and I am trying to setup a WMS GetFeatureInfo using  OUTPUTFORMAT instead of direct TEMPLATE but variables are not replaced

On MapServer file layer definition using:

TEMPLATE 'html.tmpl'

and with a html.tmpl:

:
Value: [value_0]

:

I have the proper value reply on QGIS 2.18

Value:125

When I try to implement  it by replacing the template location with

TEMPLATE 'getfeatureinfo_html'

and creating an OUTPUTFORMAT option as follows

OUTPUTFORMAT
 NAME "getfeatureinfo_html"
 DRIVER "TEMPLATE"
 FORMATOPTION "FILE=html.tmpl"
 MIMETYPE "text/html"
END

The server replies with the correct template BUT the variable is not replaced, resulting in an output like:

Value: [value_0]

Therefore the variable is not being replace and on logs (level 5) I have no indication or an error.

It looks like MapServer is not aware the layer from where to pull the value. Looked at documentation and examples but is not clear how to achieve it (missing layer information).

Later I would like to have a GetFeatureInfo for WMS with multiple mimetype options

Any tips or solutions?


texts: PyQGIS (Spatial Join): Sum features of attribute from point layer and join with vector layer

texts:
I am new to Python and QGIS. I have developed a reasonable idea of where I want to go so and where I am stuck.
I am using QGIS 2.0.1 and Python 2.6.1
What I need: I have a layer consisting of polygons which I want to spatially join with a layer of points. However, I am particularly interested in the sum of features of a particular attribute from the points layer which I want to join with the vector layer to conduct my analysis. Hope that was sufficiently clear.
I am still in the process of understanding what is happening. Thus, would be happy about any advice on how to proceed. What I have understood so far:
#_Set the scene_

qgis.utils.iface
import qgis.core

#_Import stuff_

LayerVector=QgsVectorLayer("x.shp", "x_name", "org")
LayerPoint=QgsVectorLayer("y.shp", "y_name", "org")

#_Define stuff_

VectorI=LayerVector.fieldNameIndex('bearing')
PointI=LayerPoint.fieldNameIndex('bearing')
VectorF=LayerVector.getFeatures()
PointF=LayerPoint.getFeatures()

#_edit_

VectorF.startEditing()

#_loop stuff_

for i in VectorF:
    for j in PointF:
       if j().intersects(i.geometry()) == True:

Some function creating an attribute and accumulating integers from features fulfilling the above statement.
#_edit_

VectorF.commitChanges()


texts: Differentiate between line ending within polygon and line passing all the way through polygon

texts:
I have a line file and a polygon file in QGIS.
Primarily I am looking for any polygon that has a line intersecting it - easy, sorted.
From these intersects, I then want to pull out just the polygons that have a line going all the way through them and ignore those where the line ends within the polygon. See the below diagram for an explanation: I only want polygons 2 and 3.

Does that make sense? I haven't managed to find any predefined tools thus far.


texts: Pause and resume execution of FME workbench

texts:
I started to run a workspace in FME Workbench (desktop, version 2015). The execution takes a long time (already 1 week, lasts 3 weeks).

Is it possible to do this (with my workbench I have already started)?

pause the execution
reboot the pc
resume the execution (without let FME begin at the beginning)

The workspace is quite simple: Reader (ECW) > Reprojector > Writer (ECW).

Besides input and output files there are some temporary files fme has generated:

I know that I am able to pause and resume a workspace with this tool. But as far as I know it does not work if I close FME:

I searched  in FME help and web but I did not found a solution.
I want to reboot the pc because our it department want to install updates to windos and/or to other software.


texts:
I've created a lovely route map in CartoDB, but one of the lines crosses the dateline (from Colorado to Guam, FWIW) and finds its way there by hitting the dateline north of Alaska, bouncing back around past Greenland to the dateline, bouncing again ... you get the idea.

I tried to manually draw a line that crosses the dateline and found that it gets cut off right at the dateline, no matter what I do.

See for yourself:
http://cdb.io/1k1vod6

Is there a way to draw a poly line that crosses the dateline in the Pacific?


texts:
Just trying to find a definition.

I was asked this in an interview and I assumed it was a polygon created by incorrectly modifying a another polygon in a coverage.

According to my best googling resulting in this code: http://www.falconview.org/svn/FalconView/trunk/public/fvw_core/Utility/GeometryUtils.cpp and my limited understanding of it.  A shaved polygon occurs when you eliminate a polygon and your surface no longer forms a convex hull.  I would assume this is important in a situation where keeping area constant (such as city limits).


texts: MapWinGIS Create Buffer Using Data That Has WGS84 Coordinate System

texts:
I am using MapWinGIS with c#. I am trying to create a buffer for a polygon that has WGS84 decimal degree reference system, but it never get displayed.
Here is my code:

// Form load
var gp = new GeoProjection();
gp.SetWgs84Projection(tkWgs84Projection.Wgs84_UTM_zone_40N);
axMap1.GeoProjection = gp; //set map projection to projected mercator
...
if (sfAOI.Open(fileName, null) ) //open shapefile wgs84 lat long RS
{
    axMap1.RemoveLayer(aoiHandler);
    aoiHandler = axMap1.AddLayer(sfAOI, true); //Lat Long

    //till here, all data loaded are in lat long coordinate system

    //create aoi buffer
    sfAOIBuffer = new Shapefile();
    sfAOIBuffer = sfAOI.BufferByDistance(500, 20, true, false);
    aoiBufferHandler = axMap1.AddLayer(sfAOIBuffer, true);
    axMap1.MoveLayerTop(aoiBufferHandler);  //does not appear.
}

Which suppose to create a buffer after successfully loading area of interest.
Having all data are in WGS84, I need to be able to create buffer on the fly in meters regardless to data reference system.

I am not sure if projection will really solve the problem because that deals only with data display.
What am I missing?


texts: How can I get the length of segments from the output of pgr_dijkstra?

texts:
I created a topology from OSM data with osm2pgrouting.

I use a very simple Query

SELECT * FROM pgr_dijkstra('
    SELECT gid AS id,
         source,
         target,
         cost_s AS cost,
         reverse_cost_s AS reverse_cost
        FROM ways',
    13224, 9224, directed := true);

Which works fine. However, I am also interested in the length of the segments.

The ways table has a column length_m so I though I add this columns to the nested SELECT

SELECT * FROM pgr_dijkstra('
    SELECT gid AS id,
         source,
         target,
         length_m,
         cost_s AS cost,
         reverse_cost_s AS reverse_cost
        FROM ways',
    13224, 9224, directed := true);

However, the resulting table does not contain a column with the length of the segment.

Assuming I imported OSM data with osm2pgrouting. Would a join with the ways table gives me the desired output?

SELECT * FROM pgr_dijkstra('
    SELECT gid AS id,
         source,
         target,
         length_m,
         cost_s AS cost,
         reverse_cost_s AS reverse_cost
        FROM ways',
    13224, 9224, directed := true) AS route
INNER JOIN ways w
ON route.node = w.gid;

Two question: Is the join the right way to get the desired information (length of the routing segments).

And: I wonder why I cannot select the stuff within the pgr_dijkstra()?

*edit:
The proper way to join the tables should be

SELECT * FROM pgr_dijkstra('
    SELECT gid AS id,
         source,
         target,
         length_m,
         cost_s AS cost,
         reverse_cost_s AS reverse_cost
        FROM ways',
    13224, 9224, directed := true) AS route
LEFT JOIN ways w
ON route.edge = w.gid;


texts: Python Script create spreadsheet with exception

texts:
In the following Python Script, I would like to create a spreadsheet (Excel would be fine) with the Exception results instead of just printing "There was an error--->..." to be able to keep track of all of the errors that occurred (if any).

Any suggestions?

#import system modules
import arcpy, os

#Set environment options
arcpy.env.workspace = r'Z:\project\....gdb'

# Calculate field with arcpy.da.UpdateCursor
for fc in arcpy.ListFeatureClasses():
    try:
        with arcpy.da.UpdateCursor(fc, ["Depth", "OD"]) as cursor:
            for row in cursor:
                row[0] = row[0]*(-1)
                row[1] = row[1]*(-1)
                cursor.updateRow(row)
                print "Working on {0}".format(fc)

    except Exception as ex:
        print "There was an error --->    {0}".format(ex)

print ("Add Field Finished")


texts: Rearranging columns in table when using arcpy.CopyRows_management in ArcMap

texts:
I am creating a table from a geodatabase using a fieldinfo object. I can create the table, but I need the columns to be in a specified order.  Is there a simple way to rearrange the column order before saving it as a .dbf file? Currently the columns come out as Miles, FType, FName
I need them to be FType, FName, Miles
    #Get the fields from the input

    fields2= arcpy.ListFields(InTable2)

    #Create a fieldinfo object
    fieldinfo2 = arcpy.FieldInfo()

    #Iterate through the fields and set them to fieldinfo2
    for field in fields2:

      #from near table calculation
        if field.name == "IN_FID":
            fieldinfo2.addField(field.name, field.name, "HIDDEN", "")
        elif field.name == "NEAR_FID":
            fieldinfo2.addField(field.name, field.name, "HIDDEN", "")
        elif field.name == "NEAR_DIST":
            fieldinfo2.addField(field.name, field.name, "HIDDEN", "")
        elif field.name == "NEAR_RANK":
            fieldinfo2.addField(field.name, field.name, "HIDDEN", "")
        elif field.name == "NEAR_ANGLE":
            fieldinfo2.addField(field.name, field.name, "HIDDEN", "")
              #From feature class
        elif field.name == "FType":
            fieldinfo2.addField(field.name, field.name, "VISIBLE", "")
        elif field.name == "FName":
            fieldinfo2.addField(field.name, field.name, "VISIBLE", "")
        elif field.name == "Miles":
            fieldinfo2.addField(field.name, field.name, "VISIBLE", "")
        elif field.name == "Direction":
            fieldinfo2.addField(field.name, field.name, "HIDDEN", "")

    #Create the table with the fields listed above

    arcpy.MakeTableView_management(InTable2, "FOC", "", "", fieldinfo2)

    #Append the data from all the other output tables

    arcpy.management.Append(["output1","output2","output3"], "FOC", "NO_TEST")

    #Save table to DBF file

    arcpy.CopyRows_management("FOC", mpars.replace("-","") + "Features_For_GWTool.dbf")


texts: QGIS hides few graduated points after reopening project

texts:
I have a strange problem in QGIS 2.2 (problem is also present in QGIS 2.1).
When I add a point (shapefile) layer, apply graduated style and change point size, everything looks fine. Then I save project, reopen it, and some of the points are missing.

When I open attribute table all points are there, but some of them are not shown on the map. Usually those point are in the middle of the map. I noticed that number of mising points depends on selected "on the fly" CRS (before save). When I remove graduated style then all points show up on the map... Zooming, refreshing view and changing projection doesn't help.


texts: Calculating proportion of area covered by another layer - not just the intersection

texts:
This is a follow up question from $area gives NULL for some polygons in QGIS's Field Calculator and $area is NULL for over-wrapping geometry.
I have a map of Cambodia at district level (dist) and economic land concession (elc). I am trying to calculate the proportion of district that is occupied by the economic land concession. I have tried the following:

created area for elc layer
for dist layer,
100*intersecting_geom_sum('elc','area_elc')/$area

Problem is, there are some districts that have over 100% coverage because the elc layers in these districts have outrageously large area calculated.

I have used the following data, projected them into EPSG:3148, set the output field length as 21 for area calculation as advised. I also ran fix geometry for elc layer afterwards.
How can I correctly calculate the area for these districts.
*** ADDED
I have found that for a small district, area of one big chunk of elc that cross the district and other district was used- instead of part of elc that is within district. Attached as picture below:

When I take the intersection, it works to cut them correctly. How could I possibly correct this problem?
Data using:
GEO2_KH (Cambodia shapefile) https://international.ipums.org/international/gis_harmonized_2nd.shtml
ELCs: All economic land concessions shapefile https://data.opendevelopmentcambodia.net/dataset/economiclandconcessions


texts: Using Extract Data tool without compressing?

texts: I am looking for a way to either use the Extract Data tool within ArcMap without it compressing the files into a .zip at the end, or a way for the Clip tool to have a way to select multiple feature classes (see screenshot). The goal is to put this into ModelBuilder so that I can automate a process. However my knowledge concerning Python and inserting/editing scripts is limited to none and I know of no other way to achieve the results without scripting/adjusting.


texts: Downloading map as shapefile from ArcMap service?

texts:
I am trying to export a polygon file (link is below) available on the ArcMap service as a shapefile.

Which tool I can use to do that?

https://services.arcgis.com/jIL9msH9OI208GCb/arcgis/rest/services/World_Lakes/FeatureServer


texts: Shorter Version of this SpatiaLite Query

texts:
I've been doing some spatial queries using spatialite, and they're large enough that a spatial index becomes a must. An example of using an index is given in the tutorial

spatialite> SELECT NewRegions.Name, COUNT(*) FROM NewTowns, NewRegions
        WHERE NewRegions.Name IN (
        'VALLE D''AOSTA', 'PIEMONTE', 'UMBRIA', 'LOMBARDIA',
        'CALABRIA', 'MOLISE', 'MARCHE', 'BASILICATA') AND
        Within(NewTowns.geom, NewRegions.Geometry) AND
        NewTowns.ROWID IN
        (SELECT pkid FROM idx_NewTowns_geom WHERE
        xmin > MbrMinX(NewRegions.Geometry) AND
        xmax < MbrMaxX(NewRegions.Geometry) AND
        ymin > MbrMinY(NewRegions.Geometry) AND
        ymax < MbrMaxY(NewRegions.Geometry))
        GROUP BY NewRegions.Name;

This is a bit of a mouthful, and for Within queries it's entirely redundant, it will always be xmin > MbrMinX(NewRegions.Geometry)... Is there any way to do this query less verbosely?


texts: Query/find features in ArcPad, hide the rest in map?

texts: Trying to help out another department that does quarterly assessments of field data. Their symbology is all set by the status of the assessment (the only way I could figure out how to hide something was to hide symbology) so that doesn't help much. Their checkouts contain all data for the year but they want to be able to filter/limit what is visible in the map based off the quarter (and potentially other data columns). Been trudging through documentation but can't really find anything. Is this possible?


texts: Ogr2ogr trying to import and convert shapefile to EPSG:3857

texts:
I'm attempting to import a shapefile (along with .prj file) to MySQL using Ogr2ogr and using the t_srs option set to EPSG:3857 (the default for Leaflet and Google Maps) but I'm getting odd results, eg:

POLYGON((-11368056.30134075 7311108.362629404,-11373244.15325952 7311276.412255973,-11372957.24394724 7313966.399114295,-11367933.92414896 7312806.137997068,-11368056.30134075 7311108.362629404)),1

the command I'm using is:

ogr2ogr -f MySQL MySQL:db,host=host,user=user,password=pass sourceDir -nln destTable -t_srs EPSG:3857 -update -overwrite -lco engine=MYISAM

the .prj in the sourceDir is

PROJCS["NAD_1983_UTM_Zone_13N",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",500000.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-105.0],PARAMETER["Scale_Factor",0.9996],PARAMETER["Latitude_Of_Origin",0.0],UNIT["Meter",1.0]]

Any suggestions as to the problem?


texts: Is there a standard for the specification of .prj files?

texts:
We use shape files for data exchange. I've looked at the .prj files: different programs produce different prj files.

Is there a standard how the parameters must be specified in a prj file?

AFAIK the name of the projection (e.g. PROJCS["ETRS89 / UTM zone 32N (zE-N)",…]) can be freely assigned. But I can not believe that the name of the datum (e.g. DATUM["D_ETRS_1989",…]) can be assigned freely too.

I have found no information about this.
Even in ESRI Shapefile Technical Description there is no information about .prj files.
The ESRI KB Article 14056 HowTo: Create projection metadata (.prj) files for shapefiles says: “You have to use the predefined names for map projection and parameter objects”. But the article does not specify where to find these predefine names. Does each GIS-company specify  it own names? Or does EPSG specify these names?

Here are a few differences in a section of .prj files (EPSG 4647) between QGIS, AutoCAD Map and spatialreference.org


texts: Accessing information inside an .sd file using python

texts:
I want to know if there is a way to list all folders/files inside a .sd file using Python.  Below is what I am currently working with.

I am writing a script that will access our ArcGIS Online organization, list all service definitions we have for all our users and download them to my PC.  So now I have a folder with 100+ .sd files and now I want to list out all the mxd names inside the sd file.  Using 7-zip I can open these files and see what the name of the MXD that created the service is.

The only thing I have came close to finding useful is the module libarchive but for whatever reason pip will not install this module.

I currently am using Python 3.6.5 for this particular script but I can switch over to Python 2.7.8 if needed.

EDIT
Here is a screenshot for @enolan

On the left is the .sd file and on the right is the file opened up in 7-zip.


texts: Remove parts of line that fall within a polygon

texts:
I have a polygon and two lines. How can I remove the portion of each line that falls within the polygon?

# polygon
coords <- matrix(c(-1.798123, -1.793072, -1.805767, -1.804129, -1.798123, 55.68066, 55.67369, 55.67508, 55.68139, 55.68066), ncol=2)
myPolygon <- Polygons(list(Polygon(coords)),  "myPolygon")
myPolygon <- SpatialPolygons(list(myPolygon))
proj4string(myPolygon) <- CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0")
plot(myPolygon)

# lines
line1 <- matrix(c(-1.79880, -1.79517, 55.67737, 55.67920), ncol=2)
line2 <- matrix(c(-1.80231, -1.80679, 55.67764, 55.68004), ncol=2)
line1L <- Line(line1)
line2L <- Line(line2)
my.lines <- Lines(list(line1L, line2L), ID="my.lines")
myLines <- SpatialLines(list(my.lines))
proj4string(myLines) <- CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0")
plot(myLines, add=T)