texts: Adding East Sussex ArcGIS Online Services to QGIS?

texts:
Looking for info on how to add web mapping services from the East Sussex  of ARGIS Tithe Map  Resources website to QGIS. I've tried nearly a five different things online and still can't get results. I am asking here but decided to attempt this route for fear that it was impossible (permissions issues, etc.). Link to data I want to add is posted below:

https://escc.maps.arcgis.com/apps/webappviewer/index.html?id=94d9d72603034a85bb68bce793fb5a59

or here :


texts: OD cost matrix not calculating for all origins in ArcGIS Network Analyst?

texts:
I am using ArcGIS Desktop 10.2.

I have created the origins points from a fishnet laid over my area of interest.
Out of 1200, there are some 65 points for which the OD cost matrix cannot be calculated.

I have isolated these 65 points and tried to run the OD cost matrix only for them (just in case they were not computed because the matrix was too big), but to no avail. The warning is still displayed for all 65 locations. Eg: for location 1 "No destinations found for "location 1" in "origins"".

Any ideas on how to solve this problem?

Legend of Picture:
yellow dots are those O that were successfully calculated, the orange are those that were not successfully calculated and the violet lines are each OD cost link.


texts: Showing script tool dialog from button in ArcMap?

texts:
Using Setting Script Tool Parameters Using Python and ArcGIS 9.3.1? I was able to properly set up the parameters of a script tool.

I would like to get the same dialog box from the script tool to appear when used from a button in ArcMap 9.3.1.

I know how to create a button and getting working, but only if it does not contain parameters. I'm assuming that this customization needs to be carried out using VB, or whatever other language.

So, is there some generic code that can be used to create the same dialog box as the one that appears in the script tool?

I'm not a programmer (although I can hack may way through stuff in Python) so I don't know where to start.


texts: Using MMQGIS for geocoding from CSV file?

texts:
I try geocoding locations from a CSV file in QGIS using mmqgis with no success.

is there a log file where I can see the error massages?
let's use one of the records as an example

When I use the GeoCode plugin and geocode Abu Gosh, Israel - I get the geocode but when I use mmqgis I have no luck. I tried putting country in the state field and city in the address field (I have no address).


texts:
I'm trying to get familiar with creating my custom vector tiles. For this purpose I'm creating a simple world map with country boundaries based on natural earth shapefiles. I'm creating the vector tiles using tippecanoe which I'm later on loading into a browser based mapbox map. I'm able render the country boundaries just fine but I'm struggling with the country labels which should be placed at the center of a feature. After some research I found out that it is not possible to simply add a label to a feature based on some property value since a feature might be cut in half by a tile boundary and therefore it is not possible to find the center of the feature. It seems to be necessary to create a separate layer for labels with only point data instead of lines.
Since I was assuming that this is quite a common problem, my question is whether such datasets are available somewhere or whether they have to be created by hand by calculating country centroids or similar. I was hoping to find something like this based on open streetmap data but I wasn't able to find anything.
Any hint is very welcome!


texts: Highest level curve per (Voronoi) polygon?

texts:
I am facing a difficulty in QGIS (2.8.2):
I have three layers:

layer 1: layer with points
layer 2: layer of level curves (contour lines)
layer 3: layer of Voronoi polygons obtained from the points layer (layer 1).

What I need to get is, for each point (layer 1), the highest elevation contained in its respective Voronoi polygon.

How could it be done?


texts: geoserver heatmap tiles doesn't display well (with leaflet)

texts:
I am trying to display a heatmap generated with geoserver on my map using the L.tileLayer.WMS class from Leaflet, but the tiles do not fit well together.

Do you have any idea how to fix it?


texts: Using a different Python IDE/debugger for ArcGIS

texts:
Possible Duplicate:
  How can I effectively debug arcpy Python scripts?

I prefer to use PyScripter over IDLE as my IDE/debugger for writing GP scripts in Python.  In ArcGIS, however, when I right-click on a geoprocessing "Tool" that references a Python script, it always opens in IDLE.

What is the best way to change the default tool that ArcGIS is calling?

Edit: I am currently running 9.3, but am transitioning to 10 soon, so suggestions for both versions is appreciated if they are different.


texts: How to calculate Multispectral Vegetation Indices using Python

texts:
I am trying to calculate Multispectral Vegetation Indices from multi spectral photos, using Python.
So far, I've managed to calculate the NDVI, which is easy to implement as there are many examples online.
#Read the image file
image_file = "my_tif_image.tif"

with rasterio.open(image_file) as src:
    band_red = src.read(3)

with rasterio.open(image_file) as src:
    band_nir = src.read(4)

# Allow division by zero
numpy.seterr(divide='ignore', invalid='ignore')

band_nir = numpy.array(band_nir)
band_red = numpy.array(band_red)

NDVI = (band_nir.astype(float) - band_red.astype(float)) / (band_nir + band_red)

fig1, ax = plt.subplots(figsize=(12, 12))
fig1 = plt.gcf()
ep.plot_bands(NDVI,
              cmap='PiYG',
              scale=False,
              cbar = False,  #Change this to TRUE to display the color bar
              ax=ax,
              vmin=-1, vmax=1,
              title="NDVI")
plt.show()
fig1.savefig('NDVI.png',bbox_inches='tight')

However, apart from the above index, I couldn't manage to find anywhere examples on how to calculate other indexes (at least the ones I am interested in).
For example, instead of the NDVI index, I tried the formula of Burned Area Index:
BIR = 1.0 / (((0.1 + band_red) ** 2) + ((0.06 + band_nir) ** 2))

but I get an empty plot and the calculated array numbers are very small (in the numbers of 10^-9).
My question is, why is it so hard to find examples that calculate and plot indices other than the NDVI? And lastly, if I can't do it with python, is there any other programmatic way to do it?
I've stumbled across the Arcpy, but I wasn't able to import that module. Do I have to have the ArcGIS software in order to use that library?


texts: Run Python 3 scripts on ArcGIS Server 10.5

texts:
Is it possible to run Python 3 scripts on our ArcGIS Server 10.5 running Windows?
I have a script that I want to publish as a geoprocessing tool that exports an .aprx (ArcGIS Pro) file to a PDF using arcpy.mp but this is using Python 3.
I see at Run stand-alone scripts that Python 3 is installed with Pro on the desktop but can't see advice on how to install on your server.
I see The Python 3 runtime for ArcGIS Server on Linux but the Windows option is greyed out.


texts: Modification of Line Density using ArcGIS Desktop?

texts:
I am using ArcGIS Desktop 10.1.

I need to calculate Line Density for polyline vector but the "Search radius" option must be replaced by something like the Neighborhood option in Focal Statistics.

I need to use a rectangle neighborhood instead of a circle neighborhood.


texts:
I need to make image tiles out of a DWG. I need them so I can use them as overlay on google maps. hopefully I can get them in this folder configuration: zoom_level/x_or_column/y_or_row.png or similar cut by 250px by 250px. I tried to follow an example by FME team themselfs and they do not go in details. It was not clear at all.
Here is the dwg file.
here is the FME file.


texts: Best way to visualize Google Maps and Mapinfo (or SHP or SQLite) at the same time?

texts:
I am having a bunch of Mapinfo data of a large area. This data has multiple layers and is also available as SHP, and of course it could be easily converted to SQLite ("SpatiaLite"). Although probably only one layer needs to be visualized at this point.

The problem is: The customer needs to select a rectangle from the data so it can be processed further.

I already though of writing a Google Maps Plugin / Mashup that is fed by some AJAX/JSON. Or writing a Plugin for uDig. (OpenJUMP won't be appropriate because it has problems handling large amounts of data.) Another idea was writing a Java app with Google Maps embedded and maybe some (OSS) overlay that visualizes the SHP/Mapinfo data.

Each solution has its good sides but also its bad sides. So, what do you think is the best way to realize it?

Thanks, Philip

BTW: The solution should not consist of expensive Closed-Source Tools.


texts: Getting wrong scale on georeferenced TIFF in QGIS?

texts:
I just georeferenced a map of Iceland using the ISN93 coordinate system. It looks like it was successful since when I transform some random coordinates to WGS84 and compare them to Google Maps it does match.

For the raster image I used a map 1:2'000'000 that I downloaded from the cadastral office, so it should be accurate.

But when I try to view it on QGIS I need a scale of at least 1:10'000'000'000 to see it.

Why is that?


texts:
In my country, we have a special epsg-code for transformation (#4649). In the definition (http://www.epsg-registry.org/report.htm?type=selection&entity=urn:ogc:def:coordinateOperation:EPSG::4649&reportDetail=short&style=urn:uuid:report-style:default-with-code&style_name=OGP%20Default%20With%20Code&title=), a shTransCom.dll will be listed as a component. I have this.

Now my problem - there are several shapefiles I want to transform by using epsg:4649. One idea is to use QGIS- but how to integrate the dll-file? Normal epsg-code in my knowledge are ony a long list of parameters.

On the other side, my idea is to extract the coordinate - transform in another step - and reinject the new coordinates. One possible packet for perl could be http://search.cpan.org/~jasonk/Geo-ShapeFile-2.51/ShapeFile.pm.

Does anyone have another or a better idea - or could describe how to integrate epsg:4649 into QGIS?


texts:
I am using Python Console and loading a Vector Layer using python commands. When the layer is added, QGIS first asks for the CRS to use, but I don't have all the CRS listed out to select the Reference system from and the "Error Log" shows the error as:

> "Could not open CRS database C:/Program Files/QGISDufour/./resources/srs.db.
Error(14): unable to open database file"

Now what shall I do get this cleared. I cannot load my vector layer because of CRS problem.


texts: How to determine aspect areas for different elevation bands in a basin in QGIS?

texts:
I am working in QGIS and have the dem rater, the vector basin boundary, the hypsometric curve and the aspect raster.

I need to determine the north, south and east-west facing areas in each of the 15 elevations bands already determined.

How can I do this in QGIS?


texts: How to avoir IOError: [Errno 13] Permission denied: when running r.reclass from QGIS not GRASS

texts:
I run r.reclass from QGIS process toolbox (ie it works fine if i run from GRASS) and get this error :

2017-01-26T12:25:14 2   Uncaught error while executing algorithm
            Traceback (most recent call last):
              File "C:/PROGRA~1/QGIS2~1.18/apps/qgis/./python/plugins\processing\core\GeoAlgorithm.py", line 203, in execute
                self.processAlgorithm(progress)
              File "C:/PROGRA~1/QGIS2~1.18/apps/qgis/./python/plugins\processing\algs\grass7\Grass7Algorithm.py", line 274, in processAlgorithm
                func(self)
              File "C:/PROGRA~1/QGIS2~1.18/apps/qgis/./python/plugins\processing\algs\grass7\ext\r_reclass.py", line 45, in processCommand
                with open(tempRulesName, "w") as tempRules:
            IOError: [Errno 13] Permission denied: u'tmp1485429914724

The same process works fine in  GRASS.

I have no admin issue / access in windows and my user is administrator


texts: Is it possible to get the coordinates of the bounding box of a SRTM3 .hgt file using the terminal? Maybe with GDAL?


texts: Tips and advice for Filtering on a large dataset

texts:
I am attempting to build an application using the ArcGIS JavaScript API that will allow my users to research competitor pipelines across the country.  The main layer in this application is comprised of a large feature class containing over 673K features.  One major aspect that I want to give my users is the ability to "filter" features.  For example, show only pipelines of companies X, Y, & Z that are transmission, and greater than 24" in diameter.

Another major aspect is that it has to be simple to use.  I have been experimenting with using ComboBoxes so that the possible options are already populated in the list.  By using more than one of these for each field I can create a simple way to create a SQL query of sorts.  The problem I have with the ComboBox is that only one option can be chosen at a time.  But with nearly 2000 possible options I'm not really sure what else to use.

Any one have any ideas?  Is there something that I am overlooking? Has anyone spent much time using the EnhancedGrid?

I have added the image below that is using an EnhancedGrid.  Could it be possible to use the "Check Boxes" as a selection method.

Thanks for you input.

-Mike


texts:
I have a table with a set of original XY values and a set of new XY values. I want to create new points of the new xy's using the insert cursor. My code updates the ID and Sp row but not the XY row.

for fc in arcpy.ListFeatureClasses('*Clip_1'):

  X = 'X1'
  Y = 'Y1'
  ID ='Sp'
  xy = []

  with arcpy.da.SearchCursor(fc, [X, Y, ID]) as Cursor:
    for row in Cursor:
        xy = [(row[2],(row[0],row[1]))]
        cursor = arcpy.da.InsertCursor(fc,['Sp','SHAPE@XY'])
        for new in xy:
            cursor.insertRow(new)
            print str(new)
  del row, xy

my output table looks like:

ID   Sp    X    Y    X1   Y1
1    21    25   25   21   19
2    22    12   29   12   18
3    23    14   18   22   10
4    21
5    22
6    23

What I want is

ID   Sp    X    Y    X1   Y1
1    21    25   25   21   19
2    22    12   29   12   18
3    23    14   18   22   10
4    21    21   19
5    22    12   18
6    23    22   10

Is there a better way to do this?


texts: What is the meaning of -32768 in gdaldem color config?


texts: Polygon overlap into new layer on CartoDB

texts: New cartodb & PostgreSQL user here. I have two rows in cartodb table and each row has its own polygon. The two polygons have some overlapping area. How can I create a new polygon of just the overlapping area, preferably via cartodb SQL commands?


texts: Categorized symbology inside categorized symbology

texts:
I want to make a complicated symbology for my layer.
Now it's categorized on column use and I want to add a term
"If use=street, different symbology to lamp_type column (H=(x),L=(y)
and others=(z)"

(x) (y) and (z) are symbols

I mean by that, that categorized symbology has another categorized symbology in it

I have two relevant columns one is "use" which contains [road, street, etc.] and the other is "lamp_type" which contains words like [HPS,LED,MH] etc. I need to make special symbology for every item in "use" and if the "use" is "street" I need to make special symbology for every item in "lamp_type"

How can I set up symbology this way?


texts:
Pretty basic question; In this Polyline to Raster tool, there is a field called "Value Field"

I am not sure what to put for this field. My polylines are ZM, so there is Z information for the line, I figured this tool would use that information to build the raster. How can I get the value field to you this z information??


texts: Client Zoom in OpenLayers with Google layer

texts:
I'm displaying some features on a map, but the features are somewhat on the small side and can't be easily viewed at default zoom levels. Using the Client Zoom example, I was able to enable acceptable zoom levels on an OSM layer. However, I'd like the features to be displayed on a satellite/hybrid layer.

Google is my first choice, but I have been unable to get client zoom working properly on their hybrid layer. Attempting to zoom in past the server's maximum zoom level results in the feature itself being zoomed, but not the map tile around it.

Has anyone gotten this to work? Other satellite layers would be acceptable as well.


texts: Plotting bearing and distance from known point in ArcMap

texts:
I have bearing and distance estimates of animals from a survey post. I have the coordinates of the survey post and I already have that in GIS.
How do I plot these bearing and distance estimates as points from the survey post in ArcMap?
Is there a tool for this?


texts: Converting Tab to Shp Using OGR2OGR with Python on linux

texts:
First I downloaded  http://svn.osgeo.org/gdal/trunk/gdal/swig/python/samples/ogr2ogr.py and I am importing this. How can I get ogr2ogr to run in my script?

Nothing seems to be working and I'm not sure how to continue problem solving.

#!/usr/bin/env python
import ogr2ogr
import sys

if __name__=="__main__":
ogr2ogr.main('-f', "Mapinfo File", '/home/firstuser/birmingham/Wood.tab', "/home/firstuser/birmingham/Wood.shp")

With this I get
    main() takes at most 3 arguments (4 given)

and when I do:

if __name__=="__main__":
ogr2ogr.main('-f', "Mapinfo File", '/home/firstuser/birmingham/Wood.tab')

I get:

Usage: ogr2ogr [--help-general] [-skipfailures] [-append] [-update] [-gt n]
           [-select field_list] [-where restricted_where]
           [-progress] [-sql <sql statement>]
           [-spat xmin ymin xmax ymax] [-preserve_fid] [-fid FID]
           [-a_srs srs_def] [-t_srs srs_def] [-s_srs srs_def]
           [-f format_name] [-overwrite] [[-dsco NAME=VALUE] ...]
           [-simplify tolerance]
           [-fieldTypeToString All|(type1[,type2]*)] [-explodecollections]
           dst_datasource_name src_datasource_name
           [-lco NAME=VALUE] [-nln name] [-nlt type] [-dim 2|3] [layer [layer ...]]

 -f format_name: output file format name, possible values are:
 -f "ESRI Shapefile"
 -f "MapInfo File"
 -f "TIGER"
 -f "S57"
 -f "DGN"
 -f "Memory"
 -f "BNA"
 -f "CSV"
 -f "GML"
 -f "GPX"
 -f "KML"
 -f "GeoJSON"
 -f "Interlis 1"
 -f "Interlis 2"
 -f "GMT"
 -f "SQLite"
 -f "ODBC"
 -f "MSSQLSpatial"
 -f "PostgreSQL"
 -f "MySQL"
 -f "PCIDSK"
 -f "DXF"
 -f "Geoconcept"
 -f "GeoRSS"
 -f "GPSTrackMaker"
 -f "PGDump"
 -f "GPSBabel"
 -append: Append to existing layer instead of creating new if it exists
  -overwrite: delete the output layer and recreate it empty
  -update: Open existing output datasource in update mode
  -progress: Display progress on terminal. Only works if input layers have the "fast      feature count" capability
  -select field_list: Comma-delimited list of fields from input layer to
                 copy to the new layer (defaults to all)
  -where restricted_where: Attribute query (like SQL WHERE)
  -sql statement: Execute given SQL statement and save result.
  -skipfailures: skip features or layers that fail to convert
  -gt n: group n features per transaction (default 200)
  -spat xmin ymin xmax ymax: spatial query extents
  -simplify tolerance: distance tolerance for simplification.
  -dsco NAME=VALUE: Dataset creation option (format specific)
  -lco  NAME=VALUE: Layer creation option (format specific)
  -nln name: Assign an alternate name to the new layer
  -nlt type: Force a geometry type for new layer.  One of NONE, GEOMETRY,
  POINT, LINESTRING, POLYGON, GEOMETRYCOLLECTION, MULTIPOINT,
  MULTIPOLYGON, or MULTILINESTRING.  Add "25D" for 3D layers.
  Default is type of source layer.
  -dim dimension: Force the coordinate dimension to the specified value.
  -fieldTypeToString type1,...: Converts fields of specified types to
  fields of type string in the new layer. Valid types are :
  Integer, Real, String, Date, Time, DateTime, Binary, IntegerList, RealList,
  StringList. Special value All can be used to convert all fields to strings.
  -a_srs srs_def: Assign an output SRS
  -t_srs srs_def: Reproject/transform to this SRS on output
  -s_srs srs_def: Override source SRS

  Srs_def can be a full WKT definition (hard to escape properly),
  or a well known definition (ie. EPSG:4326) or a file with a WKT
  definition.


texts: Install of QGIS3.0 on Mac OSX using KyngChaos prebuilt binaries

texts: The installer doesn't work because it requires Python 3.6 to be installed as a framework in /Library/Frameworks (or so I believe). It is insufficient to have Python 3.6 installed by Macports (in /opt/local/Library/Frameworks). I don't know about Homebrew. I am not sanguine about trying to install Python into /Library/Frameworks. It is not clear where a binary from python.org does install, and whether it messes us links elsewhere in the directory structure. Does anyone have information about this? Why is the build from KyngChaos so inflexible?


texts: Enforcing unique values in field when editing in ArcMap?

texts:
I am working on a geodatabase. In some of my features, I would like to add a field, in which I will type the values on my own. On this field no duplicates must be allowed. In other words I want it to be a unique value field. I know ObjectID exists to certain uniqueness of the object, but I am working on a relational db and I want this field to be used as primary key.

I have tried adding a field and through access editing it as a "no duplicates allowed" field, but after that, I cannot edit in arcmap and the message "Warning: Unique index on non-OID field" appears.


texts:
I am trying to build a a library of scripts to call inside the QGis Python console/Python editior. I do not understand why they work independently, but not when being called from another script.
Example:

file c:\scripts\A.py

def create_copy(input,output):
    layer = QgsVectorLayer(input, "input", "ogr")
    error = QgsVectorFileWriter.writeAsVectorFormat(layer, output, "utf-8", None, "ESRI Shapefile")

I can call that function, from within the same script:

create_copy(r"z:\original.shp", r"z:\copy.shp")

works as intended.

Now if I call this from another script, like so:

file c:\scripts\B.py

import sys
sys.path.insert(0, r'c:\scripts')
import A

A.create_copy(r"z:\original.shp", r"z:\copy.shp")

I get

  NameError: global name 'QgsVectorLayer' is not defined.

Even though the call is made from within the QGis Python Editor. What am I missing?


texts:
I am new to GIS and I want to fit a DEM given as an xyz file into an existing map.

My problem is that when I open the xyz in qgis the tile is not at the correct position.

Changing the CRS to the one given in the documentation of the tile results in a white raster image.

Here are the xyz tiles (they are quite large):

https://www.opengeodata.nrw.de/produkte/geobasis/hm/dgm1_xyz/dgm1_05314000_Bonn_EPSG4647_XYZ.zip

For the tiles EPSG:25832 is given in the documentation, here is a link:

https://www.opengeodata.nrw.de/produkte/geobasis/hm/dgm1_xyz/dgm1_meta.zip

Here is the gdalinfo for one of the tiles:

Driver: XYZ/ASCII Gridded XYZ
Files: dgm1_32370_5614_2_nw.xyz
Size is 2000, 2000
Coordinate System is `'
Origin = (32369999.500000000000000,5613999.500000000000000)
Pixel Size = (1.000000000000000,1.000000000000000)
Corner Coordinates:
Upper Left  (32369999.500, 5613999.500)
Lower Left  (32369999.500, 5615999.500)
Upper Right (32371999.500, 5613999.500)
Lower Right (32371999.500, 5615999.500)
Center      (32370999.500, 5614999.500)
Band 1 Block=2000x1 Type=Float32, ColorInterp=Undefined
  Min=46.540 Max=169.160

Ideally I would like to save the tiles all as one georeferenced raster with the correct CRS. The solution does not have to involve QGIS, if there is another alternative.


texts:
I am trying to use the GDAL Lib to read a ShapeFile (.shp). I followed the instructions of http://www.gdal.org/ogr_apitut.html - But it doesn't work!

#include <gdal_priv.h>
#include <ogrsf_frmts.h>

int main()
{
GDALAllRegister();

GDALDataset       *poDS;
poDS = (GDALDataset*) GDALOpenEx( "point.shp", GDAL_OF_VECTOR, NULL, NULL, NULL );
}

The error message is that "GDALOpenEx()" is a undeclared identifier. Also "GDAL_OF_VECTOR" is not declared.

How do I access the GDALOpenEx() method?
Do I have to load some drivers? (Which I thought GDALAllRegister() would do)

I am using xCode 6 on Mac OS 10.10. GDAL 1.11 from  http://www.kyngchaos.com/software:frameworks


texts: Relation between spatial resolution and object size in satellite imaging

texts:
Is there a general rule of thumb that relates the ability of satellite remote sensors to resolve an object and the spatial resolution of the sensors?
For example, a whale that's 15m in length and 7-8 meters in width can not possibly be identified on sentinel-2 imagery which offers 20m spatial resolution but can be identified easily on Worldview 3 images with 10cm spatial resolution.
In general what is the relative spatial resolution required to resolve an object of of a given size? Is it 1:10? 1:5?


texts: How to shade or add grey overlay to non-selected features in ArcGIS API for JavaScript?

texts:
I have an application showing cadastral parcels. When a parcel is clicked, the applications makes zoom to the extent of this feature and prepares a floating pane where it is possible to load some layers specific to that parcel.

The problem is that I think that when the studied feature is clicked, there are lots of lines, points and polygons around that feature, so it is a bit difficult to concentrate in the feature that user is studying.

So I would like to know how can I add a shade or a grey overlay to everything that is NOT the selected feature, so the user has highlighted only what has selected. Then, clicking on a close button or similar, shade has to go out. Yes, it is a pure aesthetic question, but I would like to give to my application a serious and professional look. :-)

Please, consider that at this point, there could be up to 5 or 6 different layers loaded in the map.

Does anyone know how to do this process? I am not sure if it is possible to hide everything except one feature of one of the loaded layers.

PD - Just to provide an example. Check the first item of the "What I would like to do" section in this question.


texts: Centering the map in print composer around a feature (QGIS)

texts:
I have to frequently generate the same 15 maps for work, same data each time (mostly) but at different locations for each project. Is there way that the print composer could actually auto-center itself around a polygon that I created?
I know I can just click "set map extent to match main canvas extent" it does what I want, but I'm hoping to use a plugin that auto-print all my 15 different maps in one go.
Currently by default the composer uses the last extents used, and I have to individually open each layout to re-set the extents to the current extent. Is there any way I can set it so that it loads centered around the data I'm using?


texts: Does Gohlke GDAL break Command-Line GDAL/OGR in Windoiws?

texts:
I recently gave up trying to get the gisinternals versions of GDAL's Python bindings to work, and sought instead to use the 'gohlke' approach - getting the appropriate wheel from the UCI repository.

For me, the 'appropriate' whl was GDAL-2.0.2-cp27-none-win32.whl, as I am running 32-bit Python 2.7.11.

When I was using the gisinternals version, I had used gdal-201-1500-core.msi for the core, and GDAL-2.1.0.win32-py2.7.msi for the Python bindings. My Python is compiled against MCSV1500, as is evident below.

  Python 2.7.11 (v2.7.11:6d1b6a68f775, Dec  5 2015, 20:32:19) [MSC
  v.1500 32 bit (Intel)] on win32

Having installed the gisinternals stuff I did tout ce qu'il faut: i.e.,

added GDAL_DATA to the environment with value C:\Program Files
(x86)\GDAL\gdal-data (having verified that was the directory)
put C:\Program Files (x86)\GDAL at the start of the PATH

GDAL and OGR ran perfectly from the command line, ogrinfo --version gave the right version, and I thought I was off to the races.

Then I tried to run a Python script that tried to import ogr - fail (same for from osgeo import ogr.

I checked double-checked, triple-checked and so forth, and eventually gave up and ditched the gisinternals stuff.

I have set the gisinternals Python bindings up before in such a way that they worked against Python 2.7.3 MCSV1500 and GDAL 1.1, and I was looking var by var at the two machines and there was literally no difference in the PATH or GDAL_DATA variables (since Python and GDAL paths on the two machines were the same: C:\Python27 and C:\Program Files (x86)\GDAL).

OK - long story even longer, I downloaded the whl referred to above, and pip-ed it into place without drama.

Then I fired up Python, typed import ogr, and everything worked.

Enthused almost beyond the human capacity for joy, I opened a ConEMU window and typed ogrlinfo --version, and got the dreaded "OG-what? Never heard of it" message.

I figured that the 'gohlke' package must have put GDAL somewhere other than C:\Program Files (x86)\GDAL - so off I went to find it and amend that PATH and GDAL_DATA variables as appropriate.

And here's the thing: I don't have any GDAL folder anywhere, it seems.

So - finally - the questionS:

Are those who install the Gohlke GDAL, effectively giving up any hope of being able to use GDAL/OGR at the command line? That seems unlike something that's intended, but I've never seen any reference to command-line use of GDAL installed by the 'gohlke' route.
Second question - in the hope that I've got screen-hypnosis and looked straight past teh Gohlke GDAL folder... does Gohlke's whl install GDAL to a non-standard folder? [UPDATE: Yes. Yes, it does. See self-answer below]


texts: How can I use pyqgis to automate a series of field calculations?

texts:
I've got a shapefile with about 30 fields that I'll need to update from time to time by importing a CSV with new data, joining it to the shapefile, and using the field calculator to carry the new values into the old fields one at a time:

However, this is a tedious process, and it seems like something ripe for automation with pyQGIS. The target field names and corresponding join field names will not change over time.

I have a sense of what this would look like in SQL, but the python version escapes me, and this documentation is thin on field calculation.

Is it even possible to automate a batch of field calculations like this?


texts: Intersect points and polygons in shapely

texts:
I have some polygons and points that I know overlay which look like this:

One thing to note is that in the polygon shapefile there can be multiple which overlap and I want to intersect through all of them, thus a single point could potentially return numerous rows in the output.  I have been trying this code:

from shapely.geometry import *

import geopandas as gpd
pts = gpd.read_file("/path.shp")
poly = gpd.read_file("/path.shp")

combined = pd.DataFrame()
for key , grp in pts.groupby('id'):
    inter = poly.intersects(grp['geometry']).to_frame()
    combined = combined.append(inter)

but the output of combined returns False for every single instance.  Is there anyway to intersect points and polygons in geopandas or shapely? id is just the row indicator for each point.


texts: Merging LiDAR tiles without losing data in QGIS?

texts:
I'm trying to merge some LiDAR tiles together but when I do, the clarity of the image drops, so I think that means I am losing some data.

How do I merge the tiles and keep all of the data and the clarity of the LiDAR?

Before merge:

After merge:


texts: I took some ASCII files and created points with attributes. I am trying to join the line features to these but a simple join function will only take the one that is closest. Is there a way to perhaps duplicate a line feature if there are two points that correspond to it? I am trying to accomplish this on ArcGIS.


texts: Compute closest distance between points to polygons with PyQGIS

texts:
I have the following situation:

Some points are outside the polygons and some are in. I want a new field in my point layer with the closest distance to its closest polygon. The points inside a polygon should have the value zero and outside the polygon it should be the distance to the next segment of the polygon. Not the center or the node.

I have tried the "distance to nearest hub" tool, but it gives me the distance to the center of the neighbor polygon.

I have tried the GRASS command "v.distance" but this only makes a join between the two Layers. "v.distance" seems to be the GRASS tool for "QGIS JoinAttributeByNearest". There is no "distance" field in the Attribute table of the new layer.

Same with "JoinAttribteByNearest"... it doesn't give me the distance either.

The Pluging "NNJOIN" works perfectly well for what I need but I want to do this as part of a standalone PyQGIS script and to include "NNJOIN" in it is kind of complicated, thus I'm hoping to find a different solution.

Perhabs PostGIS can help me out here, but do I need to install PostGIS first to use its functionality? And if yes, do the People I'm writing this script for also need to install PostGIS to use my script?

Any other ideas what I could do to get the closest distance from the points to the polygons?

I'm using QGIS3.4


texts: ArcMap export training data for deep learning using RCNN_Masks returns a generic error

texts:
As the title says, ArcMap returns a generic error when I try to export training data for deep learning. The error just says failed to crate raster dateset. I follow the settings from this page  but still, the program fails to run. I can export data using the KITTI Rectangles model format for example but not RCNN Masks. Any suggestions?
EDIT:I was able to run RCNN Masks by trimming down the extent, originally it was probably 1 sqmi so it wasn't too large to begin with.


texts:
I am fighting with the following issue:
I have two labels quite close to each other. When I zoom in everything is alright, but when I zoom out one of them is gone.
See the situation below:

Is there any way to keep these 2 texts even if they could overlap each other?


texts:
I'm trying to put a label over a polygons. My polygons cames from a GeoJSON serialized by GeoDjango. I'm be able to categorize this polygons using the attribute "type" but I wanto to see the height as label. The height is inside the attribute "height", then I've write this:

var strokeEdificiCasalnuovo = new ol.style.Stroke({
    color: 'rgba(0,0,0,1.0)',
    lineDash: null,
    lineCap: 'butt',
    lineJoin: 'miter',
    width: 1,
});

var labelEdifici = function(feature) {
  return new ol.style.Text({
    textAlign: 'center',
    textBaseline: 'middle',
    font: '12px Verdana',
    text: function(feature) {
      var text = feature.get('height');
      return text;
      },
    fill: new ol.style.Fill({color: 'black'}),
    stroke: new ol.style.Stroke({color: 'white', width: 0.5})
  });
};

var baracca = new ol.style.Style({
  stroke: strokeEdificiCasalnuovo,
  fill: new ol.style.Fill({
    color: 'rgba(121,121,125,1.0)',
  }),
  text: labelEdifici,
});
.
.
.
.
var sourceEdificiCasalnuovo = new ol.source.Vector({
    url: '{% url 'edificicasalnuovo_geojson' %}',
    format: new ol.format.GeoJSON(),
  });

var edifici_casalnuovo = new ol.layer.Vector({
  title: '
Edificato
',
  source: sourceEdificiCasalnuovo,
  style: function(feature, resolution) {
            if ( feature.get('type') == 'Baracca') {
              return [baracca];
            } else if (feature.get('type') == 'Cabina di trasformazione') {
              return [cabinaTrasformazione];
            } else if (feature.get('type') == 'Edificio del terziario industriale') {
              return [edificioTerzIndustriale];
            } else if (feature.get('type') == 'Edificio di culto, cappella, campanile') {
              return [edificioCulto];
            } else if (feature.get('type') == 'Edificio diroccato o rudere') {
              return [edificioDiroccato];
            } else if (feature.get('type') == 'Edificio generico') {
              return [edificioGenerico];
            } else if (feature.get('type') == 'Edificio in costruzione') {
              return [edificioCostruzione];
            } else if (feature.get('type') == 'Edificio industriale') {
              return [edificioIndustriale];
            } else if (feature.get('type') == 'Serra') {
              return [serra];
            } else if (feature.get('type') == 'Tettoia o pensilina') {
              return [tettoia];
            }
      }
});

With this in console I see this error:

Uncaught TypeError: t.getFill is not a function
    at e.setTextStyle (VM4026 ol.js:7)
    at MultiPolygon (VM4026 ol.js:7)
    at VM4026 ol.js:7
    at Ua (VM4026 ol.js:7)
    at e.renderFeature (VM4026 ol.js:7)
    at e.<anonymous> (VM4026 ol.js:7)
    at e.prepareFrame (VM4026 ol.js:7)
    at e.renderFrame (VM4026 ol.js:7)
    at e.renderFrame_ (VM4026 ol.js:7)
    at e.<anonymous> (VM4026 ol.js:7)

What is wrong in my code? I'm not an expert of JavaScript...


texts: Howto list EPSG codes as proj.4 strings from the new Proj 6 sqlite3 proj.db

texts: I want to list all EPSG codes from the new PROJ.6 tool in form of proj.4 strings. Is there a tool / SQL script to replicate the former /usr/share/proj/epsg file from the SQLite3 proj.db?


texts: geoserver xml header for cors xmlhttprequest

texts:
I'm running geoserver-2.5.2. From an external server i would like to call it as a wfs source with ol3. For this I would like to allow this external server to access the data. It seams I need a header for cross-origin resource sharing (CORS). Like Zia has tried here, the link wasn't helpful to me as well, as geoserver didn't work anymore after editing the web.xml file the way as described in the link.

I believe I do need the access for calling the wfs data.


texts: Extracting desired LiDAR points with Quick Terrain Modeler

texts:
There may not be a 'better' or 'easier' way to do what I'm asking, but I figured I'd give it a shot and ask around.

I'm using Quick Terrain Modeler v8.0.3.4 (x64) and have used the ZEB1 Handheld LiDAR wand to create a cloud point dataset. As can be seen from the images, the scan is from a storage unit that holds salt. My end-state is to use the "Volume Calculation Analysis" tool in QTM in order to determine how much salt is in this (and any other) particular unit. I have successfully done this, my question however is:

"Is there an easier way to remove LiDAR points that are not pertinent to the volume calculation?"

For example, the way I did it was by using the 'Polygon (Screen)' tool, selecting points that were not salt, and then using the 'Cut' tool, repeat, etc. This is a very tedious process and certainly not an extremely effective way to go about obtaining or achieving the same results every iteration.


texts: How to fix "Couldn't load plugin ..." and "ImportError: No module named ..." in QGIS plugins?

texts:
After installation of Quantum GIS (QGIS) I tried to use Table Manager Module, but I can't do it.

Here is the message:

Nem tudom betölteni a 'tablemanager' modult innen: ['C:/PROGRA~1/Quantum GIS Wroclaw/apps/qgis/./python', 'C:/Users/TamÃ¡s/.qgis/python', 'C:/Users/TamÃ¡s/.qgis/python/plugins', 'C:/PROGRA~1/Quantum GIS Wroclaw/apps/qgis/./python/plugins', 'C:\PROGRA~1\Quantum GIS Wroclaw\bin\python27.zip', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\DLLs', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib\plat-win', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib\lib-tk', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\qgis\bin', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib\site-packages', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib\site-packages\PIL', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib\site-packages\win32', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib\site-packages\win32\lib', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib\site-packages\Pythonwin', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\Python27\lib\site-packages\wx-2.8-msw-unicode', 'C:\PROGRA~1\Quantum GIS Wroclaw\apps\qgis\python\plugins\fTools\tools']

Traceback (most recent call last):
  File "C:/PROGRA~1/Quantum GIS Wroclaw/apps/qgis/./python\qgis\utils.py", line 117, in loadPlugin
    __import__(packageName)
  File "C:/PROGRA~1/Quantum GIS Wroclaw/apps/qgis/./python\qgis\utils.py", line 283, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
ImportError: No module named tablemanager

Python verzió:
2.7.2 (default, Jun 12 2011, 15:08:59) [MSC v.1500 32 bit (Intel)]

QGIS verzió
1.7.4-Wroclaw Wroclaw, 411aff6

Python elérési út: ['C:/PROGRA~1/Quantum GIS Wroclaw/apps/qgis/./python', 'C:/Users/Tam\xc3\xa1s/.qgis/python', 'C:/Users/Tam\xc3\xa1s/.qgis/python/plugins', 'C:/PROGRA~1/Quantum GIS Wroclaw/apps/qgis/./python/plugins', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\bin\\python27.zip', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\DLLs', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib\\plat-win', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib\\lib-tk', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\qgis\\bin', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib\\site-packages', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib\\site-packages\\PIL', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib\\site-packages\\win32', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib\\site-packages\\win32\\lib', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib\\site-packages\\Pythonwin', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\Python27\\lib\\site-packages\\wx-2.8-msw-unicode', 'C:\\PROGRA~1\\Quantum GIS Wroclaw\\apps\\qgis\\python\\plugins\\fTools\\tools']


texts: Installing QGis such that it will work with the Postgres.app implementation of PostgreSQL on Mac (Lion) (open to alternative implementations)

texts:
I just got a new mac and am getting everything installed on it (postgres, postgis, qgis, grass, gdal, brew, ruby, etc) and am about to install QGis.

QGis is a great tool for me, except for that with my other mac I was never able to get it to connect to PostgreSQL.

PGAdmin always worked when connecting, but connecting QGis to the Postgres.app implementation of Postgresql never was able to connect. The error messages were depressingly vague (just recreated): Connection failed - Check settings and try again.

So my question is this. What are the ideal settings that allow qgis to connect easily to postgresql? I've found Postgres.app an easy to use implementation, but should I abandon it and use the normal brew-installed postgresql instead to get qgis to connect?

Or even better, is there a qgis/postgresql co-installation somewhere that would come close to guaranteeing inter-functionality between these two systems?

UPDATE:
(after temporarily removing the Postgres.app PostgreSQL implementation)
When I run this command:

brew install qgis --with-grass --with-postgis

brew starts off by installing postgresql. I'm going with this for now and am hoping the "pipes" between QGis & PostgreSQL will be setup to communicate with one another...


texts: QGIS Print Composer's "Render as HTML" not outputing what is input

texts:
Am I not using the correct version of HTML letter spacing?

In Print Composer using a Label with "Render as HTML" checked, I input

for what looks good in Print Composer only to have the output PDF have the text extending off the page.

I have tried many variations in the settings and have isolated it to  causing the issue.

The full code of the title block here:

FIGURE 4

National Wetland Inventory Map

General Aviation Avenue Tract

Austin, Texas

All of it works fine except

Have also tried
 and see the same over shooting the page on PDF output issue.

Screen shots:


texts:
If we want to get information about a vector layer (shapefile for example) we could use ogrinfo:

ogrinfo -so myshapefile.shp mylayer

This returns a lot of pertinent information about the shapefile.  Is it possible to return whether or not a spatial index exists?  If so, can we access info about the spatial index?

Note: I realize it is trivial to open the folder containing the shapefile and looking for a .sbn or .sbx and so the question is specific to ogrinfo.


texts: Can I Start and End an OSM2PO route in the middle of a road segment?

texts: I need to be able to start and end a route in the middle of a road segment as opposed to a vertex/intersection. I know that you can disable snapping when stating your start/end but the route still shows the start and end snapping to a vertex.


texts: Is there a maximum number of features that QGIS can save into a shapefile?

texts: I have a very large .gdb that I want to reduce to my study area. Unfortunately, the study area is still quite big. I could select features that intersected my study area but when I tried Export > Save Selected Feature As...it crashed and only saved some of the features (4GB). Ideally, I would create another .gdb but it seems that QGIS does not support .gdb creation (read-only) or delete the features I don't need directly in the .gdb but as far as I know it is not possible to modify a layer in a .gdb.


texts: Importing annotation feature class into file geodatabase?

texts:
While using ArcGIS 10.1, I had an issue importing annotation feature classes into a file geodatabase. Using Import Feature Class (Single) or Import Feature Class (Multiple) on an annotation feature class resulted in null values for the font size, angle, text etc. in the annotation feature class, but worked fine when used on any other feature class types.

The Conversion tool Feature Class to Feature Class when used alone or as part of a python script also resulted in null values, but when used with the Iterate Feature Classes iterator in ModelBuilder, the annotation feature classes retained their values for size, angle, text etc., but the annotation was changed into a polyline feature class and appeared as a line on the map where the text should have been.

I have tried to export/import  annotations from file geodatabases, personal geodatabases and coverages and always got the same result – null values, or if used with ModelBuilder, a polyline instead of an annotation.

The only solution that seemed to work is right clicking on the annotation, copying and then pasting it into the other geodatabase.

Last week the computer was upgraded to ArcGIS Desktop 10.2 and the Import Feature Class (Single) / Feature Class to Feature Class tool works perfectly when used as a standalone tool, but when used with the Iterate Feature Classes iterator null values are returned.

Is there a solution other than using the copy/paste method?


texts: Connecting to Folder/File to open mxd document using ArcObjects?

texts:
I have a dialog box written in ArcObjects.

I am programming the Map button on the dialog box.

The user picks a type of station and county.

The program goes to folder containing the mxds and opens the one with the same county name.

I need to know how to code for going to folder and then opening the mxd.

My program is listed below :

Dim pWorkspaceName As WorkspaceName
Set pWorkspaceName = New WorkspaceName
Dim pWorkspaceName2 As WorkspaceName2
Set pWorkspaceName2 = New WorkspaceName2
pWorkspaceName.PathName = "K:\TASS\4_MAPPING_DATA_Maps\2012"
pWorkspaceName.PathName2 = "K:\TASS\4_MAPPING_DATA_SUPPORT\Traffic_Mapping\Urban_Maps"

Dim cboStations As String
If cboStations = "Annual" Then pWorkspaceName.PathName
Else
pWorkspaceName.PathName2

Dim districts As String
Dim cboDistrict As String

If cboDistrict = "Abilene" Then Application.OpenDocument ("K:\TASS\4_MAPPING\Abilene_Base_Map.mxd")
Elif
cboDistrict = "Amarillo" then getmxd "K:\TASS\4_ps\2012\Amarillo"


texts:
I am working with ENC files from noaa, specifically S57 format or ".000" to be as specific as possible.

I am trying to render a layer or even the entire thing in a window but have had no luck. Perhaps you could tell me where I am going wrong.

I am pretty sure the error is happening here:

  layer = QgsVectorLayer(filename, "DEPARE", "ogr")

I get the error:

  ERROR 4: S57 Driver doesn't support update.

Code:

from qgis.core import *
from qgis.gui import *
from PyQt4.QtCore import *
from PyQt4 import QtGui

app = QtGui.QApplication([])

QgsApplication.setPrefixPath("C:/OsGeo4w/apps/qgis", True)
QgsApplication.initQgis()

window = QtGui.QMainWindow()
window_frame = QtGui.QFrame(window)
window.setCentralWidget(window_frame)
frame_layout = QtGui.QGridLayout(window_frame)

canvas = QgsMapCanvas()
frame_layout.addWidget(canvas)

filename = "C:/Users/jacot/Desktop/ENC_ROOT/US1AK90M/US1AK90M.000"
layer = QgsVectorLayer(filename, "DEPARE", "ogr")

if not layer.isValid():
  print "Layer failed to load!"

QgsMapLayerRegistry.instance().addMapLayer(layer)
canvas_layer = QgsMapCanvasLayer(layer)
canvas.setLayerSet([canvas_layer])
canvas.zoomToFullExtent()

window.show()
app.exec_()

Ultimately I want to save this in a raster format like .tif or .png. Maybe it would be easier to do that first? Could you assist me or point me towards documentation for this as well?


texts:
I am having some trouble getting gdalwarp to work. I want to clip my input raster in.img to the boundary of a polygon stored in a shapefile polygon.shp and save it out as out.img. The data type of in.img is 32-bit float. The raster cells external to the polygon delineation have a value of 0.0. I figure gdalwarp ought to be able to clip in.img to my polygon boundary, effectively setting the external 0.0 values to NoData (that is what I want to do, anyway). I am using the following command at the Windows command prompt:

gdalwarp -of HFA -cutline polygon.shp in.img out.img

The problem is, this doesn't seem to achieve anything. out.img is essentially a copy of in.img in that the cells external to the polygon boundary are still 0.0—my understanding was that this would only happen if I used the -crop_to_cutline flag.

What am I doing wrong?


texts: QGIS does not remember paths to previous projects after restarting computer

texts:
I am running QGIS 2.18 on Ubuntu 16.02. After turning on my computer and opening QGIS I see that QGIS will not remember the paths to previously opened projects.  I can see a thumbnail of these projects and the paths to these projects.  To open them I need to track down the project through File - Open rather than select it from the list of recently opened projects.

If I have a QGIS project open, then close QGIS but do not turn off the computer, I can reopen QGIS select a project from the recent projects list.

Why does QGIS not remember the paths to my recent projects and how can I fix this annoying behavior?   I have not notice this behavior with any other software installed.


texts:
I have a segment defined by 2 LatLon points(A and B)  and a third LatLon point (C).
I need an alghorithm, formulae, or some code to find the projection point of C on the segment AB (call this point D) so that CD segment is perpendicular to AB...

i few word I've to find the minimum distance between C and AB segment or in other word the height of the triangle ABC.

I'found this question:
Find the nearest point on a line and this reply: https://gis.stackexchange.com/a/39176/16839

I've translated the algorithm in Java
(GeoPoint is a simple implementation LatLon class and VectorN is a class that mantain an array of doubles and implement basic vector operation like dot product, sum etc...)
but it don't seems to work well...

My source for Geopoints coordinates is googleMap... so they are in WGS84 system... and this may be the source of problems... but in my debug I found that results are very inconsistent... a point that stay at 30meters is given as staying at 0.5meters and a point at 18meters is given as 3m away...

Can you help me?

public static GeoPoint GetOrthographicProjection(GeoPoint start, GeoPoint end, GeoPoint offsetpoint) {
        double EarthRadius=6371;

    // Start node conversion
    double sLat = degToRad(start.getLat());
    double sLon = degToRad(start.getLng());

    double sX = EarthRadius * Math.cos(sLat) * Math.cos(sLon);
    double sY = EarthRadius * Math.cos(sLat) * Math.sin(sLon);
    double sZ = EarthRadius * Math.sin(sLat);

    //End node conversion
    double eLat = degToRad(end.getLat());
    double eLon = degToRad(end.getLng());

    double eX = EarthRadius * Math.cos(eLat) * Math.cos(eLon);
    double eY = EarthRadius * Math.cos(eLat) * Math.sin(eLon);
    double eZ = EarthRadius * Math.sin(eLat);

    //Offsetpoint conversion
    double  oLat = degToRad(offsetpoint.getLat());
    double  oLon = degToRad(offsetpoint.getLng());

    double  oX = EarthRadius * Math.cos(oLat) * Math.cos(oLon);
    double  oY = EarthRadius * Math.cos(oLat) * Math.sin(oLon);
    double  oZ = EarthRadius * Math.sin(oLat);

    // Create vectors
    VectorN p1 = new VectorN(new double[] {sX,sY,sZ});

    VectorN p2 = new VectorN(new double[] {eX,eY,eZ});

    VectorN o = new VectorN(new double[] {oX,oY,oZ});

    // Calculate
    /*Vector<Double> u = p2 - p1;
    Vector<Double> po = o - p1;
    Vector<Double> w2 = po - (u * Vector3.DotProduct(po, u) / Math.Pow(u.Magnitude, 2));
    Vector<Double> point = o - w2;*/

    VectorN u = p2.minus(p1);
    VectorN po = o.minus(p1);
    VectorN w2 = po.minus( (u.times(po.dot(u)).times(1 / Math.pow(u.magnitude(), 2)))  );
    VectorN point = o.minus(w2);

    // Convert to latlon
    double rlat = radToDegrees(Math.asin(point.cartesian(2) / EarthRadius));
    double rlon = radToDegrees(Math.atan2(point.cartesian(1), point.cartesian(0)));

    return new GeoPoint(rlat, rlon);

  }

Example of result:
gp1= 40.8880224061686 , 17.391275949776173
gp2= 40.888020378444516  , 17.39116732031107

output point (that should stay on the segment gp1-gp2)
40.88801433354826 , 17.390843482462408

On gMaps you can see the point (gp1 and gp2 are blue, the output is the red marker... it should be between blue markers... instead is far away...)
see this pic:
http://i.stack.imgur.com/bBMd6.png


texts: Cannot use "CalculateGeometryAttributes_management()" tool to calculate the centroids of features with Arcpy

texts:
I've been trying to run some code in Python with ArcPy library.
Everything is working fine, but in the end of the code there is a part that I want to calculate the centroids, x and y, with CalculateGeometryAttributes_management(), like when we calculate field through ArcMap directly from attribute table, but it returns some error saying that the "module" object has no attribute "CalculateGeometryAttributes_management".
I've searched around and found that maybe the only version of ArcGIS that supports this tool on ArcPy is ArcGIS 10.6, and which I have here is ArcGIS 10.5.
Is there any other way to calculate the centroids?
Here is the part of my code containing the calling to the tool (Everything else is working - already tested):
arcpy.AddField_management(net_clip, 'LAT', 'DOUBLE')
arcpy.AddField_management(net_clip, 'LON', 'DOUBLE')

arcpy.CalculateGeometryAttributes_management(net_clip,
                                             [["LAT", "CENTROID_Y"],
                                              ["LON", "CENTROID_X"]])

And here is the error message I get:
Traceback (most recent call last):
  File "D:/ArcPy_Classes/Scripts/sampling_sig.py", line 90, in <module>
    arcpy.CalculateGeometryAttributes_management(net_clip,
AttributeError: 'module' object has no attribute 'CalculateGeometryAttributes_management'


texts: Surf3D, example from ggRandomForests: Random Forests for Regression

texts:
I have a model and I want to use the surf3D function in R, and produce a plot similar to the following (the image is the example from "ggRandomForests: Random Forests for Regression").
Example data available at
https://www.dropbox.com/s/hho5sgwjhlk4185/data.csv?dl=0
rfsrc_data <- rfsrc(rtp~., data=data)

From this point, I understand that I need to create a partial surface and then use the surf3D similarly to the example from ggRandomForests: Random Forests for Regression, however, I am lost with the partial surfaces.
surf3D(x=srf$x, y=srf$y, z=srf$z, col=topo.colors(10),
 colkey=FALSE, border = "black", bty="b2",
 shade = 0.5, expand = 0.5,
 lighting = TRUE, lphi = -50,
 xlab="Lower Status", ylab="Average Rooms", zlab="Median Value"


texts:
I'm not sure whether to post this here or on Stackoverflow, as it's sort of a crossover. Let me give it a try.

In the proposed HTML5 standard, there's the option of storing web application data in a local cache using a cache manifest. I am looking to use this technique for an offline tiled web map application and did a quick prototype here based on a standard install of Geoserver.

I just ripped the OpenLayers client page that is generated from one of the default layers and changed the <html> tag to

<html xmlns="http://www.w3.org/1999/xhtml" manifest="cache.manifest">

And created the file cache.manifest which simply reads

CACHE MANIFEST
http://lima.schaaltreinen.nl:8080/

Which should, in my limited understanding of the matter, locally cache everything that loads from that URL: all tiles and the OpenLayers JS. In fact, upon first call in Firefox, it does ask me to allow local storage, but when I hit refresh - either online or offline - the page croaks stating that the OpenLayers object is not defined. Looks like it hasn't loaded the OpenLayers JS.

Anyone successfully did a functioning offline OpenLayers-based tiled web map app? Or hints on how to investigate this further?


texts: Open software tool or QGIS plug-in to make land cover classification

texts:
Is there any tool or plugin for QGIS/GRASS to make an automatic land cover classification from satellite imagery --- based on CLC nomenclature for example? Of course, it will not be 100% sure but it would be a first step of the classification process.

I have seen a similar tool expressed into a pdf file on the QGIS blog, but it seems the plugin has not been implemented (at least, not listed into the plugin repository).

I have also the existence of the Molusce plugin. I haven't tested it yet but the approach seems to be based on land use changes.


texts:
BASEMAP

activities_and_attractions = pd.read_csv("data/activities_and_attractions.csv", encoding="ISO-8859-1")

mp = Basemap(resolution='l', llcrnrlon=-150, llcrnrlat=45, urcrnrlon=-110, urcrnrlat=60)

mp.scatter(activities_and_attractions.X, activities_and_attractions.Y, marker='o', color='b', label='Activities and '
                                                                                'Attractions')
mp.readshapefile("data/road_features_inventory", "Road Features Inventory")

activities_and_attractions = pd.read_csv("data/activities_and_attractions.csv", encoding="ISO-8859-1")

plt.legend(loc='lower right')
plt.show()

gives error

/anaconda/bin/python /Users/Chu/Documents/dssg/dssg.py
/anaconda/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3222: MatplotlibDeprecationWarning: The ishold function was deprecated in version 2.0.
  b = ax.ishold()
/anaconda/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3231: MatplotlibDeprecationWarning: axes.hold is deprecated.
    See the API Changes document (http://matplotlib.org/api/api_changes.html)
    for more details.
  ax.hold(b)
Traceback (most recent call last):
  File "/Users/Chu/Documents/dssg/dssg.py", line 19, in <module>
    mp.readshapefile("data/road_features_inventory", "Road Features Inventory")
  File "/anaconda/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py", line 2146, in readshapefile
    for shprec in shf.shapeRecords():
  File "/anaconda/lib/python3.5/site-packages/mpl_toolkits/basemap/shapefile.py", line 543, in shapeRecords
    for rec in zip(self.shapes(), self.records())]
  File "/anaconda/lib/python3.5/site-packages/mpl_toolkits/basemap/shapefile.py", line 515, in records
    r = self.__record()
  File "/anaconda/lib/python3.5/site-packages/mpl_toolkits/basemap/shapefile.py", line 491, in __record
    value = u(value)
  File "/anaconda/lib/python3.5/site-packages/mpl_toolkits/basemap/shapefile.py", line 58, in u
    return v.decode('utf-8')
UnicodeDecodeError: 'utf-8' codec can't decode byte 0x92 in position 6: invalid start byte

and GDAL

os.chdir('data/')
driver = ogr.GetDriverByName('ESRI Shapefile')

data_source = driver.Open('roda_features_inventory.shp', 0)
layer = data_source.GetLayer()
n = layer.GetFeatureCount()
print(n)

gives error

/anaconda/bin/python /Users/Chu/Documents/dssg/dssg.py
Traceback (most recent call last):
  File "/Users/Chu/Documents/dssg/dssg.py", line 32, in <module>
    layer = data_source.GetLayer()
AttributeError: 'NoneType' object has no attribute 'GetLayer'

Process finished with exit code 1

Can somebody help?


texts: Select Random Point Based on Proximity of Other Points

texts: I want to randomly select 50 points from a point layer of 3,000, but no two points can be within 100 miles of each other. Having difficulty moving beyond a theoretical framework into something practical. I have access to ArcGIS and QGIS, if one or the other's better suited for this. I also know some python but haven't figured out what to do.


texts:
I have used phylin package to perform IDW on my data using different approaches since gstat package seems to perform only one type of idw.
rawdata = read_excel("E:\\PhD\\10 years\\R code\\Mapping\\July 2019 - masked.xlsx")
wholearea = read_excel("E:\\PhD\\10 years\\R code\\Mapping\\WA - Jul19.xlsx")
land = read_excel("E:\\PhD\\10 years\\R code\\Mapping\\Land.xlsx")
rw = SpatialPointsDataFrame(rawdata[,c("Longitude","Latitude")], rawdata[,3])
WA = SpatialPointsDataFrame(wholearea[,c("Longitude","Latitude")], wholearea[,3])
land.mask = SpatialPointsDataFrame(land[,c("Longitude","Latitude")], land[,3])
grd = as.data.frame(spsample(WA, "regular", n= 20000))
coord = as.data.frame(rw)
coord1 = coord[,c(2,3)]
coord2 = cbind(coord1, coord[,1])
rw.idw = phylin::idw(coord2, coords = grd, grid = coord1,
                     method = "Shepard", distFUN = geo.dist)

I encounter this error after a while into the run of my code.
Error in diag(values) : 'list' object cannot be coerced to type 'double'
I have used summary function to find out what's wrong with my data and here is the results:
summary(coord2)
   Longitude         Latitude       coord[, 1]
 Min.   :-93.72   Min.   :28.84   Min.   :941.4
 1st Qu.:-92.77   1st Qu.:29.01   1st Qu.:955.9
 Median :-91.04   Median :29.24   Median :957.0
 Mean   :-90.99   Mean   :29.29   Mean   :957.1
 3rd Qu.:-88.92   3rd Qu.:29.49   3rd Qu.:959.6
 Max.   :-88.41   Max.   :30.30   Max.   :978.7

summary(coord1)
   Longitude         Latitude
 Min.   :-93.72   Min.   :28.84
 1st Qu.:-92.77   1st Qu.:29.01
 Median :-91.04   Median :29.24
 Mean   :-90.99   Mean   :29.29
 3rd Qu.:-88.92   3rd Qu.:29.49
 Max.   :-88.41   Max.   :30.30

summary(grd)
       x1               x2
 Min.   :-93.70   Min.   :28.85
 1st Qu.:-92.38   1st Qu.:29.27
 Median :-91.06   Median :29.70
 Mean   :-91.06   Mean   :29.70
 3rd Qu.:-89.74   3rd Qu.:30.12
 Max.   :-88.42   Max.   :30.54

How can I address this matter?


texts: Choosing layers for processing tools in QGIS?

texts: Am I right that if I want to choose layers for QGIS processing tools I always have to choose the shapefiles via the path directory, and that I can't drag and drop them from the layer window on the left like in ArcGIS Desktop?


texts: Sources for population density map of Europe

texts:
I want to make population density map for EU where I can show populations of people aged about 65yrs, where can I get a EU shapefile and population Excel file?

Is there any EU shape file available with country boundaries or will I have to create one myself by collecting each country of EU? I have downloaded NUTS 2013 but its attributes hard to understand.


texts: trying to use mbtiles as a OpenLayers.Layer.TMS

texts:
I am modifying the opengeo MapBoxSource plugin. The plugin is able to render MapBox Mbtiles into a openlayers TMS layer. I am modifying the plugin so that it will render my own tilesets served from tilestream. However I am not able to figure out what tile origin setting to use.

The original working code:

    createStore: function() {

    var options = {
        projection: "EPSG:900913",
        numZoomLevels: 9,
        serverResolutions: [
            156543.03390625, 78271.516953125, 39135.7584765625,
            19567.87923828125, 9783.939619140625, 4891.9698095703125,
            2445.9849047851562, 1222.9924523925781, 611.4962261962891,
            305.74811309814453, 152.87405654907226, 76.43702827453613,
            38.218514137268066, 19.109257068634033, 9.554628534317017,
            4.777314267158508, 2.388657133579254, 1.194328566789627,
            0.5971642833948135
        ],
        buffer: 1
    };

    var configs = [
        {name: "blue-marble-topo-bathy-jul", numZoomLevels: 9}
    ];

    var len = configs.length;
    var layers = new Array(len);
    var config;
    for (var i=0; i<len; ++i) {
        config = configs[i];
        layers[i] = new OpenLayers.Layer.TMS(
            this[OpenLayers.String.camelize(config.name) + "Title"],
            [
                "http://b.tiles.mapbox.com/mapbox/",
            ],
            OpenLayers.Util.applyDefaults({
                attribution: "MapBox | Terms of Service",
                type: "png",
                tileOrigin: new OpenLayers.LonLat(-128 * 156543.03390625, -128 * 156543.03390625),
                layername: config.name,
                serviceVersion: '1.0.0',
                "abstract": '

',
                numZoomLevels: config.numZoomLevels
            }, options)
        );
    }

    this.store = ...

I am then replacing the var configs with my own layer info, I am changing mapbox.com with my tilestream server address, and serviceVersion with 'v2'.
It nearly works, but the tiles are not lining up. I have tried a number of different tileOrigin options without success. It is not clear to me why they are using -128 * 156543.03390625, -128 * 156543.03390625 instead of 20037504?


texts: Merge many objects with equal position while sum up their attribute values in QGIS

texts:
I´m looking for a smart way to merge multiple line objects with the same geometry and same position in one layer while also sum up their attribute values.

The result should be only one geometry (with the sum of all different object attributes)
To merge individual activated objects via expanded digitalizing tools won´t work, because the number of different overlaying objects is much too big to do this step by step. All other merging instruments combine either only two layers or don´t sum up the attribute values (i.e. mmqgis). Also v.dissolve seems not to give usable results. Does anyone have an idea?


texts: Problem thinning out polygon layers to find polygon centerlines

texts:
I have a polygon layer representing road surfacing layers. My task was to "thin out" the polygons in order to come up with a network of 'road centerlines' (dark green).

I generated this new polyline layer as follows in ArcMap:

Rasterized the road surfacing layer to a generously high resolution
Binarized the road surfacing layer where values of 1 indicate presence of road and values of 0 indicate no presence of road.
Thinned out the raster layer
Vectorized the thinned raster layer using "Raster to Polyline" tool

The result is conceptually OK since it does somehow find the centerlines of the polygons. However I am ending up with some really ugly artifacts as follows which originate from the thinned raster:

Zooming in a bit closer:

Now I am doubting whether this approach is appropriate, since it generates so many defects that would take an unreasonable amount of time to fix manually.

How can I either :

automate the fixing of these defects
prevent that these defects occur in the first place? I.e using another method for example...


texts: Use Image in calculation if a certain property is given

texts:
The goal is to select certain images from an image collection of Landsat 5 raw scenes, which possess a property called RADIANCE_MAXIMUM_BAND_6, since not all images contain this property.

if (property === null) {
 var index = 0}
 else {
  index = 1}

does not work, because apparently null is not the right condition, even if print(property) yields null

What's the right statement?

Rodrigo is right, this is neither specific enough nor an API friendly way. I think it would be easier with a filter

var filtered = L5.filterMetadata('RADIANCE_MAXIMUM_BAND_6',"not_contains",null)

and

var filtered = L5.filter(ee.Filter.notNull(ee.List('RADIANCE_MAXIMUM_BAND_6')))

doesnt work.

Which metadata is meant i.e. the properties of the Collection or those of the images?

to be specific: what I want to filter a Collection depending on if an Image has the property 'RADIANCE_MAXIMUM_BAND_6'. Because some images just do not have this property.


texts: Utilize standard visualuzation parameters for ImageCollection in Earth Engine

texts:
I've noticed that for many image collections in EE there's a set of properties describing default vizulatization parameters. For example, for Landsat LANDSAT/LO08/C01/T1_RT it is:

...
visualization_0_bands: B4,B3,B2
visualization_0_max: 30000.0
visualization_0_min: 0.0
visualization_0_name: True Color (432)
visualization_1_bands: B5,B4,B3
visualization_1_max: 30000.0
visualization_1_min: 0.0
visualization_1_name: Near Infrared (543)
visualization_2_bands: B7,B5,B3
visualization_2_max: 30000.0
visualization_2_min: 0.0
visualization_2_name: Shortwave Infrared (753)
...

Is there a quick way of applying those? (other then explicitly getting each parameter at a time with get())


texts:
I am fairly new to Arc Programming and am currently using .net to build my tool.

There is not much back ground information as I am starting from a Blank Map.  The goal of my tool is just to create a short cut (Hotkey) feature in order to place perpendicular lines.

I am attempting to use the function, IConstructPoint2.ConstructPerpendicular() but it requires that I use ISegment as an input.

Prior to this, the only EVENT I am using to capture the points that I place on the map are IEditEvents2.OnVertexAdded(IPoint).  How can I capture the line segment that is constructed immediately after putting down a second point on the map?


texts: I am not able to show attributes table in QGIS,after I converted the file from DXF to shapefile

texts:
I am not able to show attributes table in QGIS, after I converted the file from DXF to shapefile.
the source of file was DGN so i saved it by bentley as DXF file .
to use on QGIS ,is it right way but the source of DGN it has already id for each line so why i do not have at the end

I used another DXF importer or converter to import my map but I found I have the map how it looks like but I am not able to see each line or select each line on the map.

by the way i have FME software if someone can provide me another method to convert from DGN to DWG OR DXF to read on QGIS
Could someone guide me?


texts: Field Calculator: return two or more field value

texts:
I trying to return three field values in same If Statment and conecate them, but i Get NULL value insted of the conecated fields

 getClass (!InDiameter! , !OutDiameter! , !TheridDiameter!, !FREQUENCY! )
     if FREQUENCY== 1:
        return In
     if FREQUENCY== 2:
        return(In &"-"& Out)

     if FREQUENCY == 3:
        return(In &"-"& Out &"-"&Therid)


texts: Handle Curvilinear/rotated grid NetCDF file in R

texts:
I have a NetCDF file with a curvilinear (rotated) grid that contains meteorological data.
I am trying to crop the NetCDF based on a mask shapefile in R, but the rotated grid does project the data in the correct coordinates.
I have tried the following in R:
library(ncdf4) # package for netcdf manipulation
library(raster) # package for raster manipulation
nc_file <- '2000010112.nc'
nc_ras <- brick(nc_file)
plot(nc_ras[[1]])

First, when reading the NetCDF file as brick, I encountered the following error
Warning messages:
1: In .getCRSfromGridMap4(atts) : cannot process these parts of the CRS:
long_name=coordinates of the rotated North Pole
grid_north_pole_latitude=31.7583160400391
grid_north_pole_longitude=87.5970306396484
north_pole_grid_longitude=0
2: In .getCRSfromGridMap4(atts) : cannot create a valid CRS

Second, the resulting figure shows that the data were projected to the wrong location

This is not the correct location of the data as the extents should be lon : -127.0499 to -87.86075 degrees_east and lat : 42.03678 to 62.11168 degrees_north.
R uses the rlon and rlat variables (instead of lon and lat) in the NetCDF file as the coordinates. These extent information were obtained from CDO using the sinfon command and the following was the result:
Grid coordinates :
 1 : curvilinear              : points=42147 (223x189)
                          lon : -127.0499 to -87.86075 degrees_east
                          lat : 42.03678 to 62.11168 degrees_north
                      mapping : rotated_latitude_longitude
                         rlon : 342.1528 to 362.1328 by 0.09000005 degrees
                         rlat : -13.05 to 3.87 by 0.09 degrees

Apparently, I cannot use any shapefile as a mask to crop the data from the file because the shapefile will not be within the same extents as the NetCDF file.
I tried to find the crs(nc_ras) of the NetCDF file but it returned the following: CRS arguments: NA
I have also tried the stars package that handles the curvilinear grid
library(stars)
prec = read_ncdf(nc_file, curvilinear = c("lon", "lat"))

but I encountered the following error
Error in UseMethod("GPFN") :
 no applicable method for 'GPFN' applied to an object of class "rotated_latitude_longitude"
Error in .check_curvilinear(c_v, var, meta$variable, curvilinear) :
 Specified curvilinear coordinates are not 2-dimensional.

I have tried several solutions using command-line tools (CDO, NCO, NCL, which are reported on here, here, and other places) to convert the NetCDF file from the rotated grid to regular lonlat grid, but with no luck. I need to preserve the same values on the new grid (i.e., no interpolation or remapping)
So, my question, is there a way to handle this in R? In other words, handle the rotated curvilinear grid and crop data out of it based on another shapefile.


texts: How to write a CASE WHEN expression for label offset?

texts: I would like to label to feature with id 1 and 2 with offset x of -10 mm and the offset y with 5 mm. But on the GUI there is only one data defind overried. Is it possible to change that for the x and y-value??


texts: ArcGIS Field Calculator String to New String Does Not Remove Colon Formatting

texts:
I imported a CSV from NWS to a file geodatabase and would like to calculate a new field (Time1).
I used the following code. It succeeds with no errors. However, I want to remove the colons. I am not sure why they carry over to the new field? Both time and Time1 are string fields.
Left( [time_],2  ) + Mid( [time_], 3,4 ) + " CST"


texts: Choosing extent and scale while exporting a map

texts:
I'm exporting a map to a PNG file in QGIS 3.10.4 via Project > Import/Export > Export Map to Image. In the dialog that shows up, I can select the extent (exported area), scale, resolution (DPI), and image size (in pixels). Obviously those are related, and changing them completely independently is not possible. QGIS automatically propagates changes made for one quantity to others so that the values stay consistent:

If I change the extent, the image size is adjusted
If I change the scale, the extent is adjusted
If I change the resolution, the image size is adjusted
If I change the image size, the extent is adjusted

My goal is to create an export with a given extent, scale, and resolution (and the image size necessary to achieve that). Given the automatism described above, that's rather cumbersome to achieve: I'd like to take the extent of a given layer (the dialog includes an option for that), but this automatically sets the scale to that of the current view (I think?). When I then set the scale to my desired value, this automatically changes the extent.
My current workaround is to first set the extent to that of the desired layer, writing down the values, setting the scale, and then re-entering the desired extent. Surely there must be an easier way?


texts: Is it possible to create add-ins with VS Express?

texts:
Last year, when version 10 was still beta, I contacted ESRI to recommend a migration path for VBA conversion. I was reassured that the free VS Express will be sufficient.

Now in version 10 and with all VBA customization at "document level" since UI Controls cannot be added to the "Normal.mxt" template created a significant problem with supporting / upgrading VBA customizations.

Seems that I finally need to migrate the tools so i set out to complete the Walkthrough example only to find that I don't seem to be able to complete it.

The Walktrhough is obviously customized for those with full version of VS. I have a very limited knowledge of .Net and so i thought the Add-ins should be simple enough migration path.

Is it even possible to create an Add-in with VS Express?  Does anyone know of a good example using VS Express rather than the full version.


texts: Where can I get a TopoJSON with 3 digit zip codes (zip3) shapes?

texts:
Where can I get a TopoJSON file with 3 digit zip codes (zip3) shapes

I have searched for days and I just can find 5 digit zip codes shapes.


texts: How to delete an image from axMapControl component?

texts:
I have a problem in arcgis c#. I am using Arcgis axMapControl component and adding image to it but I can not delete image when wanted to delete from axMapControl component... I have get a exception such as " The process cannot access the file C:\Users.. because it is being used by another process "

Thanks for your help.

public void Add_map()

  IRasterLayer my_raster = new RasterLayer();
  my_raster.FromFile(path);
  AxMapControl1.AddLayer(my_raster);
--------------------------------------------------
public void Delete_map()

AxMapControl1.DeleteLayer(0);
System.IO.File.delete(path);    //There is a Exception.... This problem...


texts: Field Calculator GUID remove curly brackets from GUID string

texts: I'm attempting to Field Calculate the GlobalID field into a text-equivalent version of it, named GlobalID_T. This new text field will have the exact values for records calculated into it, except for the curly brackets at the beginning and end of the number string. All values the same for all records except the curly brackets. How can I accomplish this?


texts: Getting area covered by intersection or junction?

texts:
I need a way to get following information for each intersection or junction in Berlin:

latitude and longitude
area covered

I could get all Junctions from openstreetmap as nodes, which would give me the lat/lon, but not the area covered.

Any solutions?


texts: Installed plugins disappear after reinstalling 2.14

texts:
I am using Kyng Chaos package QGIS 2.14 on Mac OS 10.11.6, and uninstalled/reinstalled the software recently and each time I re-open a project I lose all my user-specified preferences that I initiated the last time it was open.

My installed plugins need to be re-checked for them appear in the drop-down bars, and the text size/font/zoom ratio changes etc.. Any suggestions for where I should set up a preferences cache? Or where I can look to access my old settings?


texts: Identify adjacent and no-adjacent polygons between multiple shapefiles using R?

texts:
My spatial data are recorded multiple years, providing one polygon for a year. (Polygon is derived from raster data, so I can use  rasters instead if it will be easier to run my calculations).

I start with year t1. Next year, I add my t2 shapefile. I want to know, which of my t2 polygons are adjacent to t1 year and which are not? for t3, which polygon are adjacent to t1 or t2, and which are not?

I think it can be done by some type of proximity analysis, but I can't really figure out which one. My expected results is a table containing the information if the certain polygon is adjacent to any of the polygons from the previous year (TRUE, FALSE) and what is the area of the polygon per year.

here are some dummy data: https://uloz.to/x8o66ow8/t1-poly-zip

I prefer to use R, but I can also use ArcGIS. I'm not python experienced user.

I think, that one of the solutions could be:

calculate area of all year (t1, t2, t3)
merge and dissolve t1 + t2
calculate area of polygons which are not merged with t1 and distract
this area from the total area of t2
merge and dissolve t1 + t2 + t3
calculate area of polygons which are not merged with t1 + t2 and
distract this area from the total area of t3
etc..

but I hope to automate the process..


texts: Does Snapping work for lines as well as polygons?

texts:
I can get snapping working with polygons to form a set with full coverage but I can't get it to work with lines.   What I want to do is have the end of one line segment snap to the end of another.  I.e the two ends have the same geometry.

In my research nearly all the references to snapping involve polygons not lines.

For the record the big picture is that I have a set of line segments representing a network of tracks.  The same line segments also form the boundaries of polygons.  Individual tracks often have more than one segment and I need to have end points snapped so I can generate layers with polygons and tracks that are a single feature automatically.


texts: How to access adjacent rows with cursor?

texts:
In the attached screenshot, the attributes contain two field of interest "a" and "b".  I want to write a script to access the adjacent rows in order to make some calculations.  To access a single row, I would use the following UpdateCursor:

fc = r'C:\path\to\fc'

with arcpy.da.UpdateCursor(fc, ["a", "b"]) as cursor:
     for row in cursor:
          # Do something

For example, with OBJECTID 4, I am interested in calculating the sum of the row values in field "a" adjacent to OBJECTID 4 row (i.e. 1 + 3) and adding that value to the OBJECTID 4 row in the "b" field.  How can I access adjacent rows with the cursor to make these sort of calculations?


texts: Saving band math output as TIFF using ArcPy

texts:
I'm trying to do some band math using ArcPy for ArcGIS Pro. I want the output to be a TIFF file. However, my following code is outputting in GRID format even though the input is TIFF.
band1 = "C:\Geoprocessing\_band_1_rast1"

Raster_Calculator = band1

raster1 = arcpy.Raster(r"C:\image.tif")

band1 = raster1 * 100000

band1.save(Raster_Calculator)

How can I add code to make the output TIFF?


texts: Removing and deleting duplicates in FME?

texts:
I'm analysing a large folder containing huge amounts of data, using FME.

I'm filtering it down as much as possible, and I'm a little stuck at this part.

Right now it is just about creating lists. My current spreadsheet has a list filtered by many things, one of which I set for "filetype" to multiple spatial and non-spatial types. One of the filters was "shp" and "dbf".

Now, I know .dbf can come under one of the files making a full shapefile, but it can also come on it's own as a file itself.

What I'd like to do is to take my list and remove all of the dbf's relating to the shapefile.

I have set a workspace that does the following:

Input Spreadsheet ----> Tester for .dbf and .shp ----> Matcher (on filename) ---->Output from "Matched" into Tester for .shp -----> Output merged with "NotMatched" from Matcher.

This creates a list with dbf's still involved that are not part of a shapefile, and also gives me all of the shapefiles.

I think this works. My issue is now merging this back in with the original spreadsheet, or on the other hand deleting all of the dbf's from the original.

This could be done manually of course but I need it to run in FME.

I need the final output not to include DBF's that relate to the SHP's.


texts:
i have a set of tile which was tiled profile of Mercator and proper to tms scheme mechanism.
unfortunately i don't have any tif or georeference table or anything for this. all was tiled with Mercator profile.

do you know any way for transforming them to geodetic profile?
in gdal re_tile tool can do this but it needs coordinated tif and projection.

what do you advice me for this huge tiled jpeg file?

thanks for your advice...


texts: Adding custom error in published geoprocessing service using ArcPy?

texts:
Is it possible to replace the default arcpy error message with a custom message using arcpy.AddError(msg) or else for a GP Service published into ArcGIS Server?

if I use:

arcpy.AddError("This is my custom error message")

My custom message only appears printed in the details. like this: (sorry its in french):

{
   error : {
       code: 400,
       message: "Unable to complete operation.",
       details: [
                  "Erreur d’exécution de l’outil. AppLogin Job ID:
                      jfe0aa728897644eb9344fa9fc7a48864 : This is my custom error message
                      Échec de l’exécution de (AppLogin).
                      Échec de l’exécution de (AppLogin)."
        ]
}

What I want instead its just:

{
    error : {
        code : 400,
        message : "This is my custom error message"
    }
}


texts: How to measure an irregular road line in QGIS?

texts:
I'm trying to measure the distance between some coordinate points along a road line from an OpenStreetMap source map, I am using QGIS 2.4.0.

My problem is to find a way to measure these distances between points on a irregular line like a road with some precision.

I got some points from a GPS took along some of the roads to verify if the precision is acceptable and help with the measurements. Also I have a table with another points in kilometers without coordinates system reference and I have to create a relation between all these points to relate the distance.

I was able to find some tools that give me a measurement in a straight line and It bring me a lot of imprecisions with kilometer points.


texts: Calculating Reflectance from Parrot Sequoia DN

texts:
I am dealing with the images of a Parrot Sequoia. I could not pre-process the images in Pix4d because I took single shots in different parts (so there are not overlapping images).. for me it is important to have reflectance values

For me still not clear WHAT do I need to read and from WHERE/HOW to calculate reflectance.. I georeferenced them in python and I was thinking that I could also calculate reflectance there.. I saw this approach https://gis.stackexchange.com/a/292426/122737  but I would like to know if somebody knows where/how to find such information


texts: In QGIS, is there a way to have the same layers in multiple groups

texts:
I have a large project with at least 80 layers of various types and growing. These data in these layers changes occasionally and I need to repeat outputs of numerous combinations of layers.

Currently, I have to record each layer I use for each output and turn the layers on as required and sometimes I need to change the rendering order.

I know about the standard groups in QGIS, but for my needs I would have to duplicate dozens of layers and put these into individual groups in the required order to generate my numerous outputs. Another way would be to create different projects for each, but that's just as bad and difficult to ensure the layers keep in sync. There must be a better way of doing this, but I've not found one yet. Anyone any ideas?

I'm using QGIS 3


texts:
NetCDF files downloaded from

could not be converted to geotiff using gdaltranslate command below:

 gdal_translate -ot float32 -unscale -CO COMPRESS=deflate NETCDF:"A2006001170000.L2_LAC_OC.x.nc":chlor_a fulldisk.tif

The error message is the following:

Warning 1: NetCDF driver detected file type=5, but libnetcdf detected type=3
Warning 1: NETCDF:A2006001170000.L2_LAC_OC.x.nc:chlor_a is a netCDF file, but chlor_a is not a variable.
ERROR 4: `NETCDF:A2006001170000.L2_LAC_OC.x.nc:chlor_a' does not exist in the file system, and is not recognized as a supported dataset name.

I got the NetCDF file from NASA ocean color website.

Is there any way of fixing the issue?


texts: Looking for free / open source Lidar web viewer?

texts:
I have lidar data which I want deploy my users using web viewer.

I will probably hold the data in PostGIS (not sure yet).

Is there good open source/free lidar viewer that I can integrate with OpenLayers site?


texts:
Using a geoserver WMS getcapabilities request

http://localhost:8080/geoserver/ows?service=wms&version=1.1.1&request=GetCapabilities&format=text/xml ,

the xml is displayed in the browser, but if we take &format=text/xml the request gives a popup to open/download the xml file.

I am trying to display a geoserver WMTS getcapabilities request,

http://localhost:8080/geoserver/gwc/service/wmts?REQUEST=GetCapabilities ,

in the browser adding the &format=text/xml, but even with the output format ,

http://localhost:8080/geoserver/gwc/service/wmts?REQUEST=GetCapabilities&format=text/xml,

the request still shows a popup to open/download the xml file.

How do I do a wmts getcapabilities showing the xml file in the browser instead of opening a popup?


texts: Using CAST function in ArcMap with ArcSDE?

texts:
I would like to use the CAST function in the definition query property of a layer and/or in the label class sql query window.

This was the most promising link that I've found, including the two links found there.

I haven't tried whether the syntax showed works against a file gdb but it doesn't work on SDE.

Anyone have luck using sql in ArcMap 10 against ArcSDE 10?  Thanks.

This worked on sde when the field is a Double.

CAST(MP AS CHARACTER(12)) LIKE '7%'

It does work but '%5' doesn't work.  Having the percent after works and before doesn't.


texts:
I have a field, "Access" that I want to populate with 4 variables based on the attribute of another field, "Des_Tp".For example, if a feature has "Des_Tp" = 'NP', then I would want to code it as open access, 'OA' in the field "Access".

I've put the lists in a dictionary, but the code runs without any change in the data:

import arcpy
from arcpy import env

arcpy.env.overwriteOutput = True

inFeatures = r"C:\Users\A\AccessTest.gdb\Access_test"
destp_dict ={1:("NP","NM","NCA","NF","NG","PUB","NT","WSR","NRA","NLS","IRA","ACEC","REC","RMA","ACC","OCS","SP","SREC","SLHCA","SRMA","LP","LREC","LHCA","LRMA"),
2:("NWR","WA","WSA","NSBV","RNA","REA","HCA","SDA","SW","SCA","LCA"),
3:("MPA","PROC","FOTH","ND","TRIBL","SO","LO","UNK"),
4:("WPA","MIT","MIL","PCON","PREC","PHCA","PAGR","PRAN","PFOR","PO","PO","CONE","RECE","HCAE","HCAE","AGRE","RANE","FORE","OTHE","UNKE")}

fields = ['Des_Tp', 'Access']

arcpy.MakeFeatureLayer(inFeatures, "TEMP_LYR")

with arcpy.da.UpdateCursor ("TEMP_LYR", fields) as cursor:
    for row in cursor:
        if row[0] == destp_dict[1]:
            row[1] = "OA"
        if row[0] == destp_dict[2]:
            row[1] = "RA"
        if row[0] == destp_dict[3]:
            row[1] = "UK"
        if row[0] == destp_dict[4]:
            row[1] = "XA"
        cursor.updateRow(row)


texts: Creating Dynamic 3D Buffers in ArcGIS Pro

texts: I am currently using ArcGIS Pro to model WiFi router signal distribution and I would like to use buffers to display them in a 3D format. I have used the Buffers 3D tool to create spherical buffers, however in the real world there are a lot of physical external factors that affect the range of the signal (i.e. walls, floors, ceilings, usage, construction materials, etc.). The buffers I have created sit within a constructed 3D building that has walls, ceilings, and floors. I would like to create a buffer that could recognize these physical structures, such as floors, walls, and ceilings, and reduce the buffer by a certain coefficient whenever it passes through one. I realize that this is quite a complex process, but I am looking for some simple solutions that could help make the model more realistic.


texts: QGIS 2.14.1 Processing Plug-In does not work anymore (Python error)

texts:
I did not change anything and started QGIS 2.14.1 today and was greeted with the following error message. So the Processing tab is missing and I cannot activate it.

The error:

2016-04-26T17:58:43 1   Traceback (most recent call last):
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python\qgis\utils.py", line 315, in startPlugin plugins[packageName].initGui()
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\ProcessingPlugin.py", line 57, in initGui
                Processing.initialize()
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\core\Processing.py", line 148, in initialize
                Processing.loadFromProviders()
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\core\Processing.py", line 165, in loadFromProviders
                Processing.loadAlgorithms()
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\core\Processing.py", line 200, in loadAlgorithms
                Processing.updateProviders()
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\core\Processing.py", line 173, in updateProviders
                provider.loadAlgorithms()
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\core\AlgorithmProvider.py", line 54, in loadAlgorithms
                self._loadAlgorithms()
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\algs\otb\OTBAlgorithmProvider.py", line 61, in _loadAlgorithms
                version = OTBUtils.getInstalledVersion(True)
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\algs\otb\OTBUtils.py", line 143, in getInstalledVersion
                out = executeOtb(commands, progress, False)
File "C:/PROGRA~1/QGISES~1/apps/qgis/./python/plugins\processing\algs\otb\OTBUtils.py", line 169, in executeOtb
                os.putenv('ITK_AUTOLOAD_PATH', otbLibPath())
TypeError: must be string, not None

What is wrong suddenly?


texts: Google Earth Engine - remove layer from the layer manager

texts:
Is there a way in Google Earth Engine to remove layers which are in the layer manager?

Example:

//Image selection
var img1 = ee.Image('COPERNICUS/S2/20161231T035142_20161231T040436_T47PNP')

//Location
var center = ee.Geometry.Point(99.256, 12.1682)
Map.addLayer(center,'','Center')
Map.centerObject(center,10);

//color palette
var rgb_vis = {min:0, max:3000, bands:['B4','B3','B2']};

//add images to map
Map.addLayer(img1,rgb_vis,'Image 1',true);

//NDVI Button
var button = ui.Button({
  label: 'Calculate NDVI',
  onClick: function() {
    var NDVIimg1 = img1.normalizedDifference(['B8', 'B4'])
    Map.addLayer(NDVIimg1,'','NDVI',true)
  }
});
print(button);

I would like to get rid of the layer "Image 1" in the layer manager, but the centerpoint should stay. Is there a better way than put a Map.clear() in the button code and re-add the centerpoint again?

//NDVI Button
var button = ui.Button({
  label: 'Calculate NDVI',
  onClick: function() {
    Map.clear()
    var NDVIimg1 = img1.normalizedDifference(['B8', 'B4'])
    Map.addLayer(center,'','Center')
    Map.addLayer(NDVIimg1,'','NDVI',true)
  }
});
print(button);


texts:
My work is finally catching up to the present and buying my group iPads.  My group does extensive field work and currently uses giant laptops to run our GIS software, Google Earth being one of the primary pieces of software.

How can I load a KMZ that contains 3 layers of data that is 600kb in size onto my iPad and be able to view it through the GE iPad app?

Does it have to be hosted or can I load it locally on the iPad?


texts: ImageOverlay or VideoOverlay in ArcGIS JS API?

texts: I would like to lay an image over my map using the ArcGIS JS API.  Whether it be a png, jpg, or gif, I'd like to be able to use a locally hosted image, give the 4 coordinates for the corners, and have that laid over my map.  I know this is super simple to do in leaflet with an L.ImageOverlay or L.VideoOverlay, but I am struggling to see any documentation or examples that can do this with the ArcGIS JS API.  I thought perhaps it could be done with a GraphicsLayer, but I don't see any examples.  Frankly I'm not even sure where to get started with this.  Is anyone familiar with how to do this?


texts: Set NoData values to -999 in ArcGIS 10.3

texts:
I know how to set a value to NoData using SetNull, but how about the opposite?

If I had cloud pixels in a Landsat scene that are currently set to NoData and I want to set them to -999 (an extreme value), how do I do that? All the articles I see discuss how to do the opposite.

Any Arcpy or ArcToolbox solutions?  I'm using ArcGIS 10.3.1 and Python 2.7.8


texts: Create sequential numbers based on field value using Python in the Field Calculator

texts:
I'm trying to modify the autoIncrement() python script to auto number based on another field. For example, I need the autoincrement to restart when the field "STD_ID" changes value. This is my attempt but I'm not getting the results I need.

rec=0
stand = "STD_ID"
stand2 = "STD_ID"
def autoIncrement():
    global rec
    global stand
    global stand2
    pStart = 1 #adjust start value, if req'd
    pInterval = 1 #adjust interval value, if req'd
    stand2 = "STD_ID"
    if (rec == 0):
        rec = pStart
    elif (stand2 == stand):
        rec = rec + pInterval
    else:
        rec = 1
        stand = "STD_ID"
    return rec


texts: How to compare centrality of points in different shaped and sized polygons

texts:
Suppose I have several polygons of varying shape and size, each containing a single point. I want to compare the centrality of each point, and determine which point is more central to its polygon. I'm working in QGIS and R, though for this question I'll provide a generic example.

So far I have found the centroid of the polygon, and measured the distance between the point and the centroid. Initially I thought I could divide the distance by the area of the polygon in order to create a figure for comparison, but this is flawed and gives incorrect results. I then experimentally tried dividing the distance by the square root of the area instead - and this yielded what look like correct rankings. See the below illustration:

I'm having trouble wrapping my head around the mathematics here. Is this a valid way to compare centrality in non-uniform polygon shapes and sizes? Eventually I will need to apply it to more complex shapes like the following: