texts: Getting test data for Lat/Lng -> MGRS/USNG and UTM, and Geotrans rounding errors

texts:
I'm making an MGRS overlay plugin for leaflet. To help with this, I've been using some of the leaflet mosue coordinate plugins. But they have occasional math errors. So I want to write some unit tests for their conversion code.
The first step is having test data.

I ran Geotrans 3.7 on a set of test Lat/Lng data to get UTM, MGRS, and USNG. Following the basic file processing section of the user guide: running the GUI, file open, open the Gold Data map_proj/WGS84/Lat_Lon.csv, and selecting the output coordinate system, and saving.

I ran a diff on the MGRS and USNG data, and it showed a difference.

$ diff mgrs.dat usng.dat
1c1
< COORDINATES: Military Grid Reference System (MGRS)
---
> COORDINATES: United States National Grid (USNG)
100c100
< 40REU0000018785
---
> 40RDU9999918785
208c208
< 31XDP6500529005
---
> ZAA0000033272
258c258
< 41XMP6500529005
---
> ZHD7740366636

In Geotrans, I can convert 40RDU9999918785 USNG -> Lng/Lat(57E, 30N) -> USNG, and it results in 40REU0000018785, which is 10 meters away. (It's also the coordinate given for MGRS)

ZAA0000033272 USNG -> Lng/Lat(0E 84N) -> USNG results in the same thing, but 31XDP6500529005 MGRS -> Lng/Lat is  84 0 0.0N 0 0 0.1W

So GeoTrans seems to be using different algorithms for MGRS and USNG, which is odd, as with WGS 84, they're identical.

So where should I go for example data I can trust, so I can evaluate/fix existing converters, or worst case write my own?

github.com/proj4js/mgrs(Can someone else edit this to be a link, "You need at least 10 reputation to post more than 2 links") looks reasonable, but it doesn't support the UPS part.

Update:

So using GeographicLibs GeoConvert (version 1.45, which was in the debian repos)

$ tail -n +13 Lat_Lon.csv | GeoConvert -m > mgrs.dat
$ tail -n +13 Lat_Lon.csv | GeoConvert -u > utm.dat

Differences between Geotrans and GeographicLib MGRS are as follows:

$ diff mgrs.dat mgrsGEOTRANS.dat
203c203
< ZAA0000033272
---
> 31XDP6500529005
253c253
< ZHD7740366636
---
> 41XMP6500529005

So for the first point, which is 84,0, converting back(with GeoConvert) gives

GeoConvert: 84.00000 0.00004
GeoTrans: 84.00000 0.00001

Second point is 84,60;

GeoConvert: 84.00000 60.00004
GeoTrans: 84.00000 60.00001

GeoConvert will do UTM/UPS at once, while Geotrans just gives Error: Output Universal Transverse Mercator (UTM): Latitude out of range, so I didn't compare those.

So my question becomes, can I just trust the GeoConvert values as perfect and correct, or should I use the two values from GeoTrans?


texts: Calculate 16-bit elevation values from a 32-bit RGBA image?

texts:
Say a single-band 16-bit (UInt) raster with values representing elevation (meters) must be converted to a 32-bit RGBA PNG image (8-bits/band) for use in a web application. The image dimensions must be identical so there is no interpolation etc.

1. Can the original 16bit elevation value be back-calculated (algebraically) from the 32bit RGBA values if the image is rendered using color ramps?

2. How could the original image be converted into a 32-bit RGBA image (using R and G bands in QGIS for example), such that the original 16bit precision persists? How could the original elevation be back calculated algebraically from the utilized bands?

For context - certain web applications built on platforms such as OpenLayers3 have limitations on image format support, as do some browsers. For this reason, though not ideal, it can be necessary to back-calculate values from an alternative image format so that meaningful information can be conveyed to the user.

See this SE question for an example of the OpenLayers3 pixel query methodology.


texts: I was wondering if I can make the polygon with the red arrow bigger than the one with the green arrow based on their difference in size in the attribute table?


texts: How to easily remove loose/isolated edges from pgRouting network?

texts:
I wondered if there's an easy way to 'clean up' a pgRouting network by removing those edges that do not connect with the global network. Otherwise you risk to end up with cost values of -1, meaning there exists no path between the given source and target (pgr_kdijkstra).

I handled the osm.pbf file with osm2po (v4.8.8) and imported the network into a spatial database with the generated SQL script.

Getting back to my question, I found it lacked a bit of context :-) So...

I had generated a set of points as the centroids of a polygon layer. Eventually I want to use these centroids as start/end nodes to calculate shortest paths. As these points do not belong to a network, I cannot use them directly for routing. Therefore I searched for the nodes that are closest to the set of centroids by means of Vector > Analysis Tools > Distance matrix in QGIS. The problem, however, is that it is possible these 'corresponding' nodes are unreachable from the global network. When calculating shortest paths, they are indicated with a cost of -1.

As @Carsten noted, commenting out the following line in the config file of osm2po

postp.1.class = de.cm.osm2po.plugins.PgVertexWriter

will detect all dead ends or 'dangles'. As such you are able to filter them out from the whole set of nodes. But then, the problem indeed is that if an isolated mini-network consists of more than 1 segment, the vertices in the middle are not filtered out as these vertices are no dead ends. The figure from above represents the 'ideal' situation of an isolated mini-network of 1 segment and thus 2 dead ends. Eliminating those nodes results in eliminating the whole mini-network, which is nice.

So in fact, my question should have been 'How to eliminate (or detect) nodes and segments that are not connected to the global network'?
That requires a little more (post)processing, as it is not incorporated as a preprocessing step which @Carsten also mentions in a related question.


texts: How do you retrieve a single attribute from a feature layer by clicking on it using ArcGIS Javascript API v3.2?

texts:
I am building a map that you can hover over and click on. I need to access a feature's project number (ProjectNum) when I click on it, but my click event on the layer doesn't seem to be firing. When I click the map, I get all of the feature layer's project numbers, not the one that I want. How do I retrieve the specific project number of the polygon that I click on?

function clickMap(){

    layer.on("update-end", function(evt){
        console.log(evt);
        map.on("click", function(e){
        console.log(e);
        console.log("project number: " + layer.graphics[0].attributes.ProjectNum);
        });
    })

}


texts: Joining custom shapefile to csv data in QGIS

texts:
I have a custom shapefile which I want to join to geo-coded (Lat and Long) protest data from a .csv file in QGIS. I wanted to know how this was possible without having unique IDs.

In other words, I have a shapefile that doesn't conform to standard boundaries such as administrative levels, districts or municipalities therefore I cannot easily join each dataset by a common attribute.

Is there a way of getting around this, e.g. visualising both the shapefile and point data and then working backwards to spatially join the two?


texts: Extracting point attributes to polygon grid in QGIS

texts:
I have point data of endagered flowers found in a number of locations. Each location may have up to hundreds of flower
registrations.

I want to end up with a 1x1km grid colour graduated from the number of registration (not locations) in each grid cell.

I tried making a vector grid (that worked), but joining attriubutes
by spatial overlap does not work. It runs, but table is blank.

Any idea of a different workflow or ideas on why i wont work?


texts: Adding days in date using QGIS Field Calculator

texts:
I'd like to create a new field adding 10 days in a date (string format) in QGIS; using the Field Calculator.
For example, in the field_1 there is 22/02/2021, in field_new it should become 04/03/2021.
Maybe, I should convert date in number, as you can see in Excel.

Any suggestion?


texts: Using QGIS to create a GeoTIFF raster of Google Earth import

texts: I've managed to link QGIS3 to Google Earth via XYZ tiles.  I would like to convert the Google Earth map shown to a GeoTIFF raster (to save on memory?).  I try and export the Google Earth layer to a GeoTIFF raster and notice a the bottom of the screen I have to choose extent - assume I click mapview extent, because the original Google Earth files must be large? And then I need to select the resolution?  This is the bit where I am struggling, because on the layout composer canvas I want to produce a map at A0 size with a scale of 1:5000, and if this works, at 1:1250.  A previous post I've found noted that the resolution units are metres per pixel - so how do I choose an appropriate resolution, noting the printing specifications of the map to be produced?  Also a previous post queried whether this method of exporting Google Earth data is actually allowed?


texts: I've been searching the documentation for the GUID for a Network Analysis Layer object.  For example, "E156D7E5-22AF-11D3-9F99-00C04F6BC78E" is the Interface identifier for an IGeoFeatureLayer.


texts: How manage user privileges in Portal for ArcGIS (operational dashboard)

texts: I am working on Operations Dashboard using Portal for ArcGIS. Different users have different spatial extent (different district). User name and District Data in PostgreSQL table. Now my issue is how to customize their dashboard based on their required district?


texts: I would like to use height above sea level in the gravity equations. What is the source for getting it? Can i use altitude measurements from Altimeter or GPS 'z' measurement? For your kind information, I am going to use this for aircraft velocity estimation.


texts: Is there an Open Data Cube (ODC) example/documentation on indexing data residing on Google Cloud

texts:
Are there any Open Data Cube (ODC) examples/documentation that cover indexing data residing on Google Cloud.

I know there is documentation on indexing data on AWS S3 at https://datacube-core.readthedocs.io/en/latest/ops/indexing.html#indexing-data-on-amazon-aws-s3

I am mainly exploring Google Cloud for Open Data Cube (ODC), not AWS.


texts: I wanted to find the elevation difference between 2 DEMs. When I got the data, I had to georef one image and had to aggregated it too because their cell values do not match. I also tried to resample but I read that aggregating it is better. After much tweaking, I then subtracted them, but the resulting DEM looks like a hillshade. What does that imply?


texts: QGIS open action slightly changes my url when it tries to open it

texts:
I have an issue with creating active hyperlinks in QGIS. My url (after I deleted the site name) originally looks like this:

?type=canvas_rendering#%7B%22frame%22:%7B%22id%22:481350%7D%7D

When I create an action in QGIS the normal way (Type-open, action-[% "url" %]), the links do not work as expected. It seems that QGIS is confused by the % sign inside the url and changes the exact name of the url by adding the number 25 after each % sign:

?type=canvas_rendering#%257B%2522frame%2522:%257B%2522id%2522:481350%257D%257D%0A

My url comes from a json object encoded to be compatible with JavaScript’s encodeURI() function to obtain a valid  URL. I first assumed that it re-encodes the url but I tried with the original JSON object and didn't work.

I would like to ask if someone knows how to tell QGIS to open the exact url as it is in my attribute without being confused by the % sign


texts:
I am user QGIS and I use composer and atlas to create my maps exports print.
in the atlas settings I use for scale 20% margin but now I need to print in specific scale numbers.for example if my atlas feature scale where i use scale 20% margin is the scale is between 100 and 200 then I need to take scale 200.
if the scale is between 1500-2000 then I need to take 2000 scale.

any idea how to do this ?with some condition ?

summary I want to avoid to print features with scales like this 10006 or 12221 I need to print in scales with specific numbers like 10000 or 15000 or 100 or 200.


texts: How to find objects which are styled twice (in rule based styling) in QGIS 2.18.10

texts:
I have the problem of a complex rule based styling (with 10 different styles and multi-line rules in each style) with several thousands objects. My problem is that the sum of all styled objects in one layer is higher than the actual number of objects in that layer, which means that some objects are included in two different styles and thus are styled twice.

Is there a way to find these objects faster than going through the whole map or through the whole styling rule code?


texts: Creating shapefile containing terrain, buildi ng and road layers from OS UK maps using QGIS?

texts:
I'm new to GIS.

I would like to create a shapefile map using QGIS containing contours, buildings and road information layers which can be exported to other software to build a 3D model.

I have downloaded the OS Terrain 50 map (all of UK) and OS VectorMap District (National Grid Ref Square: SO) from the Ordinance Survey UK site. I am only interested in a specific area of 1 square km within the Grid Ref 50 Map.

Is it possible to extract the Contour information from the Terrain 50 Map and the road and building information from the District Map and overlay them in QGIS and then export as shapefile?


texts:
I'm sure similar answers have been given elsewhere but no matter what I have tried I am unable to optimise this query.

Details: TableA has 3.5mill records TableB Had 57000 records but I have dissolved that down to 170 based off category of data.

I need to return what features in tableB intersect each individual record in tableA in a single string.

Update a
   Set a.[Zone_Codes] = (SELECT SUBSTRING((SELECT '; '+ b.[ZONE_CODE] FROM  [TableB] AS b
                        WHERE a.geometry.STIntersects(b.geometry.MakeValid()) = 1
                        ORDER BY b.[ZONE_CODE] FOR XML PATH ('')), 2, 1000))
FROM [TableA] AS a

Any thoughts on how to optimize this? All tables have spatial indexes. SQL Management Studio 13

--- removed amended code, didnt optimise


texts: Why gdal does not find gcps coordinates on SENTINEL images?

texts:
I use gdalwarp with python to overlay polygons on ENVISAT radar images. Everything goes well for this type of images but when I try with SENTINEL images, gdal does not find the georeferencing. It seems that gdal takes the global coordinates of the pixel at the upper left of the image. When I read the image information with gdalinfo, I realize that ENVISAT images have no gcp while SENTINEL images have. I have this error message when I use gdalwarp with SENTINEL images:

  The transformation is already "north up" or a transformation between pixel/line and georeferenced coordinates cannot be computed for tmp_spill.tif. There is no affine transformation and no GCPs. Specify transformation option SRC_METHOD=NO_GEOTRANSFORM to bypassthis check.

I do not understand why it does not find gcps when they display well with gdalinfo. I use this line for overlay polygon on ENVISAT images : cmd='gdalwarp --config GDALWARP_IGNORE_BAD_CUTLINE YES -q -dstnodata None -cutline tmp.shp %s %s'%('tmp_seep.tif',classout_seep).

Do you know why gdal does not find gcps coordinates ?


texts: How to open none .qgs files in QGIS Wien?

texts:
I have tried tried to open my older qgis files made on a mac on my new PC.
Wien(2.8) only seems to allow .qgs files so all my folders appears empty when I access them from QGIS on the PC.


texts: How to fix broken QGIS project when folders structure has changed

texts: The folder structures for our network was changed the middle of the night without consulting me.  Usually if a shapefile moves location, QGIS will bring up a broken layers dialogue box for me to fix.  Instead, for some projects, I can see the layers in the layer panel, but the status bar pulsates green with nothing happening.  Other projects take a long time to open but eventually they do.  Can someone help me understand what's happening and how to fix?


texts:
I'm trying to sample a line in a bunch of points, but when I use vector/geometry tool/extract summit I get too many points.

I'm working on Switzerland Highway, as the points are generated on every line summit I got 628 000 points for all the lines. I want these points to be separated from 100m.

I would like to set a minimum inner distance between 2 points. What are my options?

I tried to but the result is really random sometimes I get my 100m sample but in the case of very long direct line I got from 1 to 20km without any points even by reducing the simplification

I'm using QGIS 3.8


texts:
I have some raster files of an agricultural region. There are more than 50 small fields in the agricultural region. For each field, a shapefile has been provided. As shown in the image below, each polygon represents a field. In the attribute table of the each shapefile there is the field number along with other information (shown in image).

I also have some satellite raster data and none of them cover the whole agricultural region. Only a portion of the agricultural region is covered by one raster (an example is shown below).

I need to know the shapefiles' (polygon actually) name that falls inside the raster in this image. I can turn the labels on of all shapefiles and find the shapefiles' name one by one. As there are more than 20 rasters, doing that manually will take a lot of time.

Is there any way to know the shapefiles' name that falls inside the shown raster?

I am using ArcMap 10.5.


texts: I used this youtube video: https://www.youtube.com/watch?v=FSblORbDyhc&list=PLWJ2X2wvq78QRGzl8Un5ADF0_4m5t99iD&index=11 to learn how to run flow accumulation and was successful. I realize it calculates where the watershed will accumulate flow (streams/rivers etc.) I just don't understand what changing the symbology to classified (2 values) in order to show it represents. Once flow accumulation is calculated, what can it show you, and how is that done?


texts: Determine correct tolerance/epsilon for line simplification in different zoom-levels when generating vector tiles

texts:
I am writing my own vector-tile generator for a map renderer that I have written, but I am unsure of how to pick the correct epsilon value to use in the PostGIS method ST_Simplify depending on the zoom levels.
I'm using the web Mercator projection, and my tiles widths/heights halve on each zoom level.

Is there a suggested way to pick this epsilon value dependent on the zoom level? Because, ideally, at the highest zoom level, the epsilon value should be 0 such that there's little to no simplification of the lines. I did have some idea relating to getting the horizontal distance covered by a given tile at a particular zoom level, given by the formula:

distance = equatorialCircumference * cos(latitude)/2^zoomLevel

However I feel this is slightly unfounded, any suggestions?


texts: Build your own Street View with Google Street View API

texts: Is there is a way to upload my own panoramic images of our hours/building using Google Street View API, then use Google Street View Engine to navigate through them?


texts:
I have layer of point data of aid projects in Tanzania between 2000 and 2012, and a layer of point data for 42 villages where a survey was conducted in 2003, 2006,2009, and 2013. Within the survey point data my attribute table contains the percentage responses to a number of binary questions for respondents in each village. (There were 30 respondents from each village, so for example each village attribute table contains the % of people in x village who think corruption is a serious problem)

I want to see if the responses (dependent variable) are influenced by the proximity of villages to aid projects/the distance from aid projects (making proximity/distance my explanatory variable).

Should I perform OLS and then GWS? If so, how can I include distance as my explanatory variable? My data is currently geolocated in WGS_1984 datum, so should I change it to a projection?


texts: Error during conversion of shapefile to MySQL using ogr2ogr

texts:
I have installed MS4W in Windows 7. I get an error when I try to convert a shapefile to MySQL database using GDAL/OGR ogr2ogr.

This is the command I'm trying:

C:\ms4w\tools\gdal-ogr>ogr2ogr -f MySQL MySQL:banyuwangigis,user=root,password=
-lco GEOMETRY_NAME=layerkecamatan c:\ms4w\apps\banyuwangigis\petadasar\shp\admin
polygon.shp

This is the resulting error:

MySQL driver does not support data source creation.

What have I missed?


texts: Using Arcpy UpdateCursor to make calculations based on multiple fields?

texts:
I'm new to Python.

I would like to calculate field values  for a specific 'IDF' based on values from another column.  I have something like that:

with arcpy.da.UpdateCursor(fc, ['IDF', 'og70_79_sales_sum','SUM_og70_79_103_sales']) as cursor:
for row in cursor:
    if row[0] == 103:
        row[1] = row[2]
        cursor.updateRow(row)
        print('{} - IDF, {} - value, update og70_79_sales_sum: {}'.format(row[0], row[2], row[1]))

but it does not work.


texts:
I am using the QAD plugin.
The QAD-Dimensioning is working.
But when I want to erase the measurement it should be that way that I have not to select layer text first and erase the text then line-layer erasing line...

Is there a way to erase elements from different layers at the same time?

My layer-tree looks like this

layer-Measurement
    -layer-dim_line
    -layer-dim_symbol
    -layer-dim_text


texts:
I've been having several issues making GRASS algorithms work from the Processing toolbox.

When you fresh install QGIS 3.4 (MacOS High Sierra), first it doesn't find python -- which I solved by creating a custom PATH variable. Then it complains about inconsistent use of tabs and spaces, which I also solved by editing the grass74 file inside the QGIS app. Now I'm getting this error whatever GRASS algorithm I use:

WARNING: Default locale settings are missing. GRASS running with C locale.
Default locale not found, using UTF-8
Traceback (most recent call last):
File "/Applications/QGIS3.app/Contents/MacOS/grass/grass74", line 2029, in <module>
main()
File "/Applications/QGIS3.app/Contents/MacOS/grass/grass74", line 1880, in main
set_paths(grass_config_dir=grass_config_dir)
File "/Applications/QGIS3.app/Contents/MacOS/grass/grass74", line 622, in set_paths
os.environ['MANPATH'] = sys_man_path
File "/Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6/os.py", line 674, in __setitem__
value = self.encodevalue(value)
File "/Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6/os.py", line 744, in encode
raise TypeError("str expected, not %s" % type(value).__name__)
TypeError: str expected, not bytes

Right now I have Python 3.6.6 installed and QGIS 3.4.2 (the most recent 24Nov2018 release), but I've also tried with Python 3.6.5, Python 3.6.7, QGIS 3.0, 3.2, 3.3, 3.4.1... No good, same error.

I found people around here having similar problems (related to encoding) for Windows QGIS, which was solved but setting the Windows region to USA. Although my OS is in Portuguese, I'm using the QGIS in English and tried setting the OS to America region and English language, but the error persists.

Does anyone know anything about this?


texts: Some highways segments ignored by osm2po?

texts:
Where i can read why some roads/segments are disappear from routing? We are trying to comment all "wtr.deny" rules in config file, but nothing changed.
Examples added as links bellow.
Screenshot1
Screenshot2


texts: Extracting X and Y coordinates from .gpx file using QGIS?

texts:
I'm using QGIS 2.18.17 on Linux Mint 19. I have several .gpx files of waypoints recorded using a handheld GPS device. I can load them on to QGIS or Google Earth with no problems. However, I want to be able to access the coordinate data and the values for each GPS point that I have recorded.

When I look at the attribute table of the .gpx files in QGIS, the X and Y coordinates are not displayed, although the elevation, time, name, and comment data are.
When I open the .gpx files in a text editor, the coordinate information is there and I can copy and paste the values for each waypoint.

I assume that there is an automated way to do this, like converting the .gpx to a .csv file.

How do I extract all the X and Y coordinates from a .gpx file using QGIS?


texts: moving the objects in a feature class by a field within the table of that feature class

texts:
I have a feature class and I want to move its objects based on two fields in it I used the following code but get error:

with arcpy.da.UpdateCursor(fc2,["SHAPE@XY","centerpo_4","centerpo_5"]) as cursor:
for row in cursor:
    XM=row[1]
    YM=row[2]

    cursor.updateRow([row[0][0]+ XM,row[0][1]+YM])

The columns "centerpo_4","centerpo_5" are two fields within the feature class.
When I run the code I encounter with following error:

Traceback (most recent call last):

File "C:\Users\Rezars\Desktop\test2.py", line 34, in <module>
cursor.updateRow([row[0][0]+ XM,row[0][1]+YM])
TypeError: sequence size must match size of the row

The number of objects in the table are 2893 (FID starts from 0)which all three columns have the same row numbers(2893).


texts: Does lasfilterdecimate work when homogenizing by pulse?

texts:
In R lidR package, when I run lasfilterdecimate  with the homogenize option excluding use_pulse
 = TRUE it works, but when I include use_pulse = TRUE, it does not work.

For example, the code below does not work:

currenttilethinned <- lasfilterdecimate(currenttilepulse,homogenize(density,res=20,**use_pulse = TRUE**))

The code below works, but I am not sure if it is filtering by pulse (laspulse has already been run on the input LAS file, currenttilepulse):

currenttilethinned <- lasfilterdecimate(currenttilepulse,homogenize(4.25,res=20))

Does homogenize automatically use the pulses when pulseID is present?


texts: QGIS and Raster Mapping in 'Seamless' layer

texts:
I would like to use Ordnance Survey base mapping in QGIS. I would like to set up something like a seamless table in MapInfo. I think that the equivalent in QGIS is a Virtual Raster Table. The files are in tif format as they were when I download them.

I have built a VRT but it seems slow to load when I open the layer. I have also tried building pyramids for the data and this doesn't seem to speed it up that much.

Is there a way to speed this up? Am I missing something? Are there any good resources about steps to take to prepare raster data for use in QGIS?

Would I be best to start to load this information into a PostGIS database with PostGIS Raster installed?


texts: How to smooth NDVI raster (Sentinel 2-A image) in ArcGis?

texts:
My question is regarding the final  display of a NDVI raster.

I currently work with ArcGIS Desktop 10.2.2.

When I work with Landsat images (the image is loaded as composition through the MTL file, after i calculate apparent reflectance function in ArcGis,  and finally NDVI through "image Analisys" module). Using "streched simbology " the image is very visual and useful in a smoothly way, for example:

However, with the Sentinel 2-A images downloaded from ESA recently, I've still not managed the way to do so visual. In this case I downloaded through the Semi-Automatic Classification Plugin (fromgistors) plugin using "DOS 1 atmospheric correction." In this case, I could not (do not know if there form) to load the image as a composition (accompanies an XML archive I guess it's for that, but I have not found a way), so charge the bands individually and calculate the NDVI (even checking the box scientific output in "Image analisys module"). I also use "Streched simbology" to represent the image, but the quality is very low, it appears pixelated. Example:

I don't know the cause why this happens, and that is my question: do you know any way to smooth or blur the final raster (GRID)  in ArcMap?


texts: Difference between python IDLE and eclipse on accessing GRASS?

texts:
I have seen these posts:GRASS and Python, Using Eclipse to develop GRASS Python programs.

I used to use python IDLE to control GRASS, but I can't make python and GRASS connected, as I described in this post:Problem with python script to control GRASS GIS from outside - How to import grass.script under Win 8.1?

So now I am thinking of using eclipse to control GRASS, trying to avoid problems I faced as I used python IDLE. But I don't know whether I will get the same problem after I use eclipse, and I can't tell the difference between using python IDLE alone and using eclipse.

I am working with GRASS 6.4.3 and python 2.7.6 under WIN 8.1.

Can you help me to make decision and furthermore, to solve the connecting problem?


texts: How to add custom CSS class to leaflet tilelayer

texts:
I am working with mapbox.js and leaflet to create a map with a number of overlapping tilelayers. I want to apply some custom CSS to only selected tilelayers, and thus assign these layers to an additional custom CSS class like this:

var layer = L.mapbox.tileLayer(layername).addTo(map);
layer.getContainer().className += ' custom';

Unfortunately, the above did not do the trick, though. Any hints on how to assign a custom CSS class to a given tilelayer in leaflet/mapbox would be very welcome...


texts: How to get rid of artefact rectangles in QGIS Print Composer?

texts:
when inserting a map on the canvas of print composer and moving that map around,
quite often, there is a rectangle that cannot be deleted.

This rectangle is positioned at the initial place, where I draw my map first.

I thought, it is just a "on screen" problem, but in print, the rectangle is there as well.

Then I try to put it into the background, but this doesn't work either.


texts: exporting table to drive from Google Earth Engine returns some NAs where valid data

texts:
My question is similar to exporting table in to a drive from Google earth engine returns blank rows, but instead of getting entirely no data, only a few records come back NA where I have confirmed there is valid data (~20 NAs out of >2500 records).

Briefly, I am calculating the median of moving window and subtracting that value from each Landsat pixel. The PV is photosynthetic vegetation from a spectral unmixing.

Here are 2 screen grabs of the CSV and the console. Note the matching system:index IDs (I have triple checked the dates are the same, etc.)

 //

Here's a snippet of the code used to calculate the normalization and extract the data:

// Adjustment kernel function
var kernel_radius5 = 150
var PV_kernel5 = function(image,kernel_radius5) {
  return (image.select("PV").subtract(image.select("PV").focal_median(kernel_radius5,"circle","meters")));
};
var PV_normalized5 = PVdecimal.map(function(image){return PV_kernel5(image,kernel_radius5)});

// extract values from tree locations
var PVcompute_band5 = function(bandname) {
  var mapfunc = function(feat) {
    var geom = feat.geometry();
    var addProp = function(IMG, f) {
      var newf = ee.Feature(f);
      var index = ee.String(IMG.get('system:index'));
      var date = index.slice(0,20);
      var value = IMG.reduceRegion(ee.Reducer.mean(), geom, 30).get(bandname)
      return ee.Feature(ee.Algorithms.If(value,
                                         newf.set(date, ee.String(value)),
                                         newf.set(date, ee.String('NA'))));
    }
    var newfeat = ee.Feature(PV_normalized5.iterate(addProp, feat));
    return newfeat;
  }
  return treePoints.map(mapfunc);
}

// export to drive
var PVbands = "PV"
var PVtable5 = PVcompute_band5(PVbands)
Export.table.toDrive(PVtable5,
  fileName+"_L8_median5", //job name
  "Ch3_GEE_exports", //folder name
  fileName+"_L8_median5"); //file name

Not sure if the problem is with the approach I've used to extract the data or in the spatial normalization function.

Has anyone else seen this behavior?

I've made quick link with a couple of the specific points exporting as NA
https://code.earthengine.google.com/a5ae71cc6dd0816561b78c6a7acff614


texts: Calculating Flooded area / Inundation in ArcGIS for Desktop?

texts:
I was tasked to create a layer with the flooded areas if the sea / river level was to go up by X feet. There is no requirement to take Flow and Pour Points into account, just a simple raise in water level.
I will be using ArcGIS and have access to all extensions. As an input I have a DTM (tiff) and coastline (shp) layers.
I have found this very helpful video titled How to Model Coastal Flooding in ArcGIS which outlines the following process:

Calculate
Reclassify to have only a single value
Convert to vector
Remove polygons that are not touching the coastline (select by location)

Before creating a model from it, I thought there might be a dedicated tool for that, either in the Hydrology toolbox in the Spatial Analyst extension, in the ArcHydro toolboxes or in the Government toolbox.

Is it indeed the case? Which tool would be the best for that?


texts: Question about the bound box I'm getting from ogr2ogr

texts:
I'm using ogr2ogr to create geojson from a shapefile using the -lco WRITE_BBOX=YES argument. Everything is working great except the bound box (bbox) doesn't seem to be the traditional "left, bottom, right, top" EPSG coordinates. Here is an example:

This should be a district in Alabama -
[-9847815.3009172, 3522127.8720985, 0,-9650071.0310152, 3741251.3025806]

I'm thinking that it's possibly WKT coords.(?)

I've searched all over for more information about this and was unable to find any information on it. How would I go about converting this to EPSG coordinates? Are there any libraries out there that can do this? Preferably in JavaScript or PHP.

Thanks


texts: I installed Qgis 2.16.1 and QuickOSM. Then I activated the plugin and restarted Qgis. But the plugin does not appear under the Web tab.


texts:
I want to query which ward contains the most cycle stations (a cycle station is a 3D geometry i.e. polyhedral) - the wards are 2D polygons - below is what I have worked out so far - The error I get when I run the query is also shown:

 SELECT COUNT(ST_Contains(a.location, b.location))
   FROM wards a, cycle_stations b;

 ERROR:  Unknown geometry type: 13 - PolyhedralSurface

However I am aware that I am trying to query both 3d and 2d coordinates - which is where I am lost:

I guess what im wondering is - whether it is possible to do a within/contains query for 2d/3d geometries - if so, what would it be.


texts: Get a spatial reference object by name from a saved string

texts:
I am creating an ArcMap 10.2 extension is VS 2012 and C#.  I want to give the user the option to select a spatial reference to be applied to a layer - planning to use the Dialog for that.

But then, I want to be able to save the name of the SR the user selected to an custom XML settings file so that when they come back to the same section it will be set by default (and they can use the Dialog to change it if needed). This is where I am stuck.

My plan was to use that string name and recreate the SR object by see no option for it in something like ISpatialReferenceFactory3.

As an example, if I run this code on an mxd using a New Jersey projected coordinate system:

ISpatialReference SR = new UnknownCoordinateSystemClass();
IMxDocument Mxd =  (IMxDocument)WestonGISExtension.ExtensionApplicationHook.Document;
ESRI.ArcGIS.Carto.IMap Map = Mxd.FocusMap;
if (Map.SpatialReference != null)
    SR = Map.SpatialReference;

then, the SR.Name would produce, say, "NAD_1983_StatePlane_New_Jersey_FIPS_2900_Feet".  So, how do I later come back with that string and reproduce the SpatialReferenceClass?

Thanks guys.  Ya, I think you are right in going with the WKID.  Here is what I settled on:

ISpatialReferenceFactory3 SpatialReferenceFactory = new SpatialReferenceEnvironmentClass();
ISpatialReference SpatialReference = SpatialReferenceFactory.CreateSpatialReference(SpatialReferenceCode);

Ernie


texts:
I have a column in my MySQL database table that contains location information (GeoJSON).

I want to import this data and see the points on my QGIS project map but it doesn't seem to work.
What I was trying to do:

Layer -> Add Layer -> Add Vector Layer -> Database -> insert the server configuration -> Add -> import the relevant table.

EDIT:

I missed the query. I used the query builder and it seems to find around 200 relevant rows, but it shows only 1 point on the map. In addition, the point is not in the right spot...


texts: QGIS 3.2 - Forcing column type when importing csv

texts:
I've got this csv file that I'm trying to import to QGIS 3.2. It's just a data table with no geometry attached.

When I just drag it from the Browser on the left-hand side, all the columns get incorrectly read in as string. However, when I use the "Open Data Source Manager" and import the file through there, the column types are all correct (first one is integer and the others are double). But I really didn't specify that the columns should be double instead of string.

So, finally, my question is: how/where can I try to force certain column type at import for each of my columns?

Also, extra brownie points to whoever knows what's going on in the background that makes dragging a CSV from the Browser different from using the "Open Data Source Manager".


texts: Efficiently check for equal adjacent values in raster calculator

texts:
I have a raster consisting of integer values in each cell. I want to check every adjacent cell in a 3x3 grid to see if any one of the adjacent cells has a value EQUAL to the center cell, and if so, give the output raster a value of 1 (else give it a value of 0). I have already accomplished this like so:

Con("raster" == FocalStatistics("raster", NbrIrregular("C:\filepath\kernel_file.txt"), "MAXIMUM"), 1, 0)

Where the kernel file looks like this:

3 3
1 0 0
0 0 0
0 0 0

Then again for another possible combination:

3 3
0 1 0
0 0 0
0 0 0

And so forth:

3 3
0 0 1
0 0 0
0 0 0

Then I have to combine all of the rasters at the end to get one final raster.

Is there a more efficient way to accomplish this? It shouldn't be necessary to loop through the raster this many times to get what I need. Raster calculator solutions will work, but ArcObjects solutions are welcome as well as the end result will likely end up being implemented in ArcObjects.


texts: I'm working with a locally stored SQLite database and have an area field that auto populates via a virtual field using $area. Works fine but there seems to be no way to set the precision in a decimal number when using SQLite. Therefore results are not rounded off to the two decimals places I'd really like. I'm wondering if there is something I could add to the expression to do this? Or is it a constraint in an SQLite dbase?


texts: QGis Print Composer: layer symbology not displayed

texts:
Exactly the same as this problem from way back, only I've had this from 1.8.0 to present (2.1.0 master).

Yay, symbols:

Boo: no symbols, cue MS Paint:

Anyone know why this might be or have the same problem?

Happens on Win7 & Ubuntu (& Xubuntu) 13.04.

This issue raised on hug.qis.org (and can be subscribed to) here.


texts: Ordnance Survey TOIDS: converting from strings to integers

texts:
I'm planning to load a number of OS datasets into Postgres/PostGIS.

In order to save space and speed up queries, I'd like to strip the 'osgb' off the front of the TOID strings and then store the remaining number as a bigint.

Are there risks (or outright problems) with this?

I'm guessing that the numeric part of the TOID string is itself unique, but I know that some TOIDs have 13 numeric digits and some 16, and I'd like to be sure that any leading zeroes that follow 'osgb' aren't material to the number-part. This would be an issue if, e.g., there were two TOIDs in the same dataset as follows, where stripping off the 'osgb' would result in the same number:

osgb0000123412341234
osgb0123412341234

Thank you.


texts: Managing large amounts of geospatial data?

texts:
How do you manage your geospatial data? I have terabytes of data spread out over hundreds of datasets, and have an ad-hoc solution using symbolic links within projects which link back to an domain-name based archive directory for each dataset. This works mostly, but has its own issues.

I'm also keen to hear if anyone manages their geospatial data in a revision control system; I currently use one for my code and small datasets, but not for full datasets.


texts:
I am trying to programmatically access the NoData value in a raster dataset as shown in the screenshot.

Two tools that come to mind to describe raster properties do not seem to be able to describe the NoData value: Describe and Get Raster Properties (Data Management).

How can I programmatically access (i.e. print) the NoData value of a raster?


texts: Using overpass turbo to download all windfarms on admin level 6 in germany

texts:
I managed to write a query which extracts all windfarms:

[out:json][timeout:150];
(
{{geocodeArea:"Lüneburg"}};
{{geocodeArea:"Nienburg/Weser"}};
{{geocodeArea:"Uslar"}};
{{geocodeArea:"Dassel"}};
{{geocodeArea:"Einbeck"}};
{{geocodeArea:"Northeim"}};
{{geocodeArea:"Moringen"}};
{{geocodeArea:"Hardegsen"}};
{{geocodeArea:"Bad Gandersheim"}};
{{geocodeArea:"Oldenburg"}};
{{geocodeArea:"Dötlingen"}};
{{geocodeArea:"Ganderkesee"}};
{{geocodeArea:"Großenkneten"}};
{{geocodeArea:"Hatten"}};
{{geocodeArea:"Hude"}};
{{geocodeArea:"Wardenburg"}};
{{geocodeArea:"Wildeshausen"}};
{{geocodeArea:"Osnabrück"}};
{{geocodeArea:"Osterholz"}};
{{geocodeArea:"Peine"}};
{{geocodeArea:"Rotenburg (Wümme)"}};
{{geocodeArea:"Salzgitter"}};
{{geocodeArea:"Schaumburg"}};
{{geocodeArea:"Stade"}};
{{geocodeArea:"Uelzen"}};
{{geocodeArea:"Vechta"}};
{{geocodeArea:"Verden"}};
{{geocodeArea:"Brake"}};
{{geocodeArea:"Butjadingen"}};
{{geocodeArea:"Wilhelmshaven"}};
{{geocodeArea:"Wittmund"}};
{{geocodeArea:"Wolfenbüttel"}};
{{geocodeArea:"Wolfsburg"}};
)->.myArea;

node["generator:source"="wind"](area.myArea);
out meta;

The Problem is, that some of the administrative regions have the sam e Name as a City in this Region. Unfortunately overpass takes the name of the city instead of the region, so that only the wind farms in the City are selected.

Is there a way I can make clear, that I am searching for the Region and not for a single City?


texts: How can I update position of points with coordinates?

texts:
I have the following use case: I have around 1000 points and their positions in the attribute table. Sometimes I get an email with an updated position (coordinates) of a certain point. The problem is, when I update the position in the attribute table, the point in the map doesn't change. I know, that is just the way how GIS works. I also know I could go into editing mode, use the advanced editing toolbar and copy & paste the new coordinates for the point. But this is a bit too much work because it is a repetitive task.
What would be an easier way to update the position of a point with using coordinates?
From a user perspective it would be great to just copy & paste the new coordinates in the attribute table and the position of the point would update to the new position. I know that doesn't work, but maybe I could realize it with an expression in the attribute formular.
Do you have an idea which expression I could use?


texts: How to resolve this error in google earth engine?

texts:
Error

Too many pixels in the region. Found 30698755252, but maxPixels allows
only 10000000000. Ensure that you are not aggregating at a higher
resolution than you intended; that is a frequent cause of this error.
If not, then you may set the 'maxPixels' argument to a limit suitable
for your computation; set 'bestEffort' to true to aggregate at
whatever scale results in 'maxPixels' total pixels; or both.

The main code I received the error:
 features = features.reproject({crs: aspect.projection().crs(), scale: 30}).float()
 var minMax = features.reduceRegion({
  reducer: ee.Reducer.minMax(),
  scale: 30,
  maxPixels: 10e9})
print(minMax)


texts: How do I assign values of a field, to another empty field when a conditional statement pertaining to a third field is considered?

texts:
When I run the code below I don't get the output I'm seeking. The "max_proportion" field is assigned the text LEVELMEAN instead of the values for each row within the "LEVELMEAN" field when the conditional statement is met, "Highest" = 1. My else statement doesn't seem to be working either. Thanks.

fields ["Highest", "LEVELMEAN", "max_proportion"]
with arcpy.da.Update.Cursor(fc,fields) as cursor:
  for row in cursor:
    if (row[0] == 1):
        row[2] = [1]
    else:
       row [2] = "<NULL>"
    cursor.updateRow(row)


texts: QGIS map tile not printing some given extent in one time why?

texts:
I made a plugin in QGIS using pyqgis and qt designer to plot OSM tile in which I have one button which connects me to function in which I wrote a program to plot raster layer using param and I am setting an extent but when I call function when button triggers once it gives me simple ODM world map and on clicking that button again it plots the given extent.
How to handle this issue?

Here is code of function:

def maptileOSM(self):
        canvas=self.iface.mapCanvas()
        urlWithParams ='type=xyz&url=https://a.tile.openstreetmap.org/{z}/{x}/{y}.png&zmax=19&zmin=0&crs=EPSG3857'
        rlayer = QgsRasterLayer(urlWithParams, 'OSM', 'wms')
        x = float(8667633)
        y = float(3296664)
        scale=1000000
        rect = QgsRectangle(x-scale,y-scale,x+scale,y+scale)'''
        canvas.setExtent(self.rect)
        QgsProject.instance().addMapLayer(rlayer)

self.dlg.pushButton.clicked.connect(self.maptileOSM)


texts: I am using ArcGIS 10.2.2, I am trying to run Kernel Density Estimation on roads to Identify Car Crash hotspots, with ArcGIS 9.3 there is an option "Analysis mask" that one can use so that the analysis focuses strictly on input. I can't find such an option in ArcGIS 10.2.2, that's why the results include even the part were there are no roads, I want the focus to be only on the roads.


texts: QGIS Error while executing Python function in openProject() Macro: "QgsMapLayerRegistry is not defined"

texts:
I have the following script that sets the a feature id automatically after a new feature is digitized.

I want the following function to be executed automatically when project is opened (and layers are loaded) so that users don't need to do anything.

import os

# Set layer by Name

layer=None
for lyr in QgsMapLayerRegistry.instance().mapLayers().values():
    if lyr.name() == "Pozos":
        layer = lyr
        break

# Define function to select added feature and add attribute to field "id"
def update(featureAdded):
    idx = layer.fieldNameIndex('id')
    max = layer.maximumValue(idx) +1
    layer.changeAttributeValue(featureAdded, idx, max)

# Connect "featureAdded" event to "select" function
layer.featureAdded.connect(update)

I tried adding the function inside the project's Macro section inside openProject() function:

But when the project is opened I get:

It looks like QgsMapLayerRegistry is not loaded when the function is triggered?


texts: geoserver manager tutorial add new WFS layer

texts:
GeoServer Manager is a library to interact with GeoServer 2.x. I want to use it. but I didn't find any appropriate documentation for it. Can anyone help me by some tutorial or example?

here is a few examples : geo-solutions examples for geoserver manager
GeoServer-Manager wiki on GitHub: here
For example I want to add a new WFS layer into a certain workspace of geoserver. suppose I have a shafefile(.shp) and I want to add it into geoserver by geoserver-manager. How can I do that? here is the initial code:

String RESTURL  = "http://localhost:8080/geoserver";
String RESTUSER = "admin";
String RESTPW   = "geoserver";

GeoServerRESTReader reader = new GeoServerRESTReader(RESTURL, RESTUSER, RESTPW);
GeoServerRESTPublisher publisher = new GeoServerRESTPublisher(RESTURL, RESTUSER, RESTPW);


texts: Adding granules to GeoServer store doesn't update time dimension of WMS?

texts:
I'm pretty new to using GeoServer, so I'm sorry if this question has been asked already and I just was unable to use the correct verbage to find it. I've been following this tutorial to add a raster time-series to GeoServer using the ImageMosaic plug-in. I've gotten it to work with the data they provided.
In this tutorial I set up a store using a directory named timeseries in the coverages directory of the GeoServer. For my needs, I will need to constantly update this directory by adding a new daily raster.
For example, the tutorial had me add three files, timeseries_20030517.tif, timeseries_20130201.tif, and timeseries_20170803.tif.
To mimic my need of this feature I copied the first file and renamed it timeseries_20170804.tif. However the changes are not reflected when I view the WMS in QGIS. Only the data for 20170803 appears.
I've even tried removing the sample data from the tutorial entirely and inputting my own data with the same timeseries_yyyyMMdd format, but adding additional dates. The actual data clearly changes. However, when I look at the getcapabilities.xml of the WMS the time extent still shows:
<Extent name="time" default="2003-05-17T00:00:00Z">2003-05-17T00:00:00.000Z,2013-02-01T00:00:00.000Z,2017-08-03T00:00:00.000Z</Extent>

Is this a restriction of the ImageMosaic capabilities or is there something happening that I'm not doing correctly?


texts: Error installing GDAL 1.9.0 python package

texts:
I'm trying to install the GDAL python package in a python 2.7.1 environment on RHEL 5.  I am using the command 'pip install GDAL'.  I have no problem installing other packages like numpy, etc. in this manner.

The installation fails with this error:

extensions/gdal_wrap.cpp: In function 'PyObject* _wrap_VSIFTruncateL(PyObject*, PyObject*)':

extensions/gdal_wrap.cpp:6845: error: 'VSIFTruncateL' was not declared in this scope

error: command 'gcc' failed with exit status 1

I also tried compiling GDAL 1.9.0 from source and running 'python setup.py build' (after setting the path to gdal-config in setup.cfg) and 'python setup.py install'.  The GDAL python bindings seem to install okay, but I get the same 'VSIFTruncateL' error when trying to use GDAL...

>>> from osgeo import gdal
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "osgeo/__init__.py", line 21, in <module>
    _gdal = swig_import_helper()
  File "osgeo/__init__.py", line 17, in swig_import_helper
    _mod = imp.load_module('_gdal', fp, pathname, description)
ImportError: osgeo/_gdal.so: undefined symbol: VSIFTruncateL


texts: Using relation tags to nodes in OpenStreetMap?

texts:
I filter an osm planet extract using osmfilter. I am only interested in nodes belonging to route=subway relations. After performing this operation, I use osm2psql to load all the data into my postgres db.

The thing is that I don't know to which relation name belongs each node. I mean, I know all nodes belong to route=subway relations, but I don't know the name of the relations.

So, if I have a relation with these tags: route=subway, name='Line 23', I would like that all the nodes belonging to this relation had also the relation tags. For example: name:'5th street station', relation_route=subway, relation_name='Line 23'.

How can I do this?


texts: Copy-pasting records from one layer to another?

texts:
I am copying some records in a attribute table from layer1 and pasting those as new objects in layer2.

The first time I copied some records/objects from the attribute table of layer1, I pasted them as a new vector layer (it is layer2)
When I copy/paste a new record/object in layer2, it is visible as a new record in the attribute table of layer2, but it is not visible on the map/canvas (although the first copied objects are visible).

See the following screenshots (attribute table with 6 features/polygons and canvas with only 5 polygons)


texts: Combine hover info with onclick info in openlayers?

texts:
I am trying to do some complex thing. I wonna have:

popup with a name of a feature on mouse hover,
popup with more complex info ( name, foto, etc.) after mouseclick on the same feature - the same layer.

It is possible to have those two things at the same time??


texts: Definition Query does not work for flex application against ArcGIS Server 10.5

texts:
I have a flex 3.5 app that has been working for a great number of years. I wrote it from the ground up not basing it on the sample flex viewer.

I had some server problems yesterday and took the downtime as an opportunity to upgrade everything from 10.4.1 to 10.5. For some reason my definition queries on map services on the 10.5 server no longer work. I do not know if this is by design or a configuration that I can change.

All the services draw unless I use a definition query on them. The query I am using works in the service layer query URL (10.5) so I am really at a loss for why it will not work in the flex application.

I pointed the Flex application towards a development 10.4.1 server and all the definitions work as they have in the past. There is not any code to troubleshoot, the only difference between a functional application and an non-functional one is the version of ArcGIS Server.

Can anyone help?


texts:
I have a polyline shapefile with time series data. I want to create an animation map using Carto and show the line width based on the value. I used to create an animation map for point data and it worked. I updated the code to accommodate polyline setting. However, I tried torque this time and it doesn't seem to work. Any thoughts?

Map {
  -torque-frame-count: 64;
  -torque-animation-duration: 30;
  -torque-time-attribute: "datetime";
  -torque-aggregation-function: "CDB_Math_Mode(torque_category)";
  -torque-resolution: 4;
  -torque-data-aggregation: linear;
}

#layer {
  line-color: #4CC8A3;
  line-width: 1.5;
  line-opacity: 1;
  line-comp-op: multiply;
}

#layer[frame-offset=1] {
  line-width: 8;
  line-opacity: 0.45;
}

#layer[frame-offset=2] {
  line-width: 10;
  line-opacity: 0.225;
}

#layer {
  comp-op: source-over;
  line-opacity: 0.9;
  line-color: #FFF;
  line-width: 1.5;
}

#layer[value=1] {
  line-width: 1.5;
  line-color: #ffeda0;
}

#layer[value=2] {
  line-width: 3.0;
  line-color: #f03b20;
}


texts: Sharing data in small remote organization

texts:
My company's GIS "department" currently consists of myself and another analyst, located in offices several states apart. We have been looking for ways to share GIS data between ourselves as seamlessly as possible, instead of constantly compressing shapefiles and emailing them back and forth, as well as having a central GIS data hub so that new hires can have access to everything we do.
My sense from the limited research I've done is that GIS-focused cloud storage seem to be geared towards larger organizations, or at least the costs are. I'm leaning towards trying dropbox, but wanted to see if anyone else has insight on good programs to use. We both use ArcGIS, and I don't think we'd be looking to change platforms.


texts: What is the relationship between Map Service Instances and SDE connections?

texts:
We have a GPS project based on ArcGIS Mobile 10.0. There are two web map services and each were kept as Nonpooled with Minimum instances 15 and maximum instances 60 by original developer. The database is ORACLE and connection is made through SDE. Recently whole enterprise moved to 10.1 except GPS project which has some crashing issues with 10.1. But that led hosting of this GPS project on an old Test server.

Two incidents happened in last two weeks when the connections for GPS projects reached alarming high in few minutes (122 and 58). First time it affected the SDE server and caused total chaos- we had to kill the processes and bounce the SDE server back. Second time, we restarted the services before things out of control. I investigated the first time - only 20 GPS devices were connected that time. I am new to ArcGIS mobile and don't know SDE in detail(It is still under DBA control and I can only investigate the appserver). How will the number of SDE connections being used one time transforms to map services? I want to go to management with the reasons and recommendations.


texts:
The Q&A How to cleanly upgrade QGIS on Windows? outlines steps to upgrade QGIS. An important step recommended is "Make copy of plugins directory located here C:\Users\yourUserName\.qgis2\python\plugins and then delete the .qgis2 directory."

Are there not other elements besides the plugins folder to save for the new installation, that would otherwise be removed when uninstalling the version to be replaced? Is there a recommended process for storing such elements in protected locations outside of QGIS install folders?


texts: Click link inside Leaflet Popup and do javascript

texts:
I have a leaflet map up and running.  It overlays a series of polygons (via GeoJSON) on the map and attaches popups to each polygon. Each of the popups display information about that polygon.

I'd like to have inside the popup a link that, when clicked, runs a javascript function that pulls further smaller polygons via AJAX and shows them.

I can't get the script to catch a click on the link via the normal jQuery/Javascript click events.  Here's what I mean by normal (the following doesn't work):

$('a .smallPolygonLink').click(function(e){
  console.log("One of the many Small Polygon Links were Clicked");
});

The bindPopup part is as follows.  It runs on each polygon when made and it pops up correctly on clicking on a polygon.  It does show the link, just won't run the above code on click.

var popupContent = "Basic Information..." + 'Click here to see the smaller polygons';
layer.bindPopup(popupContent);

Any help is appreciated.


texts:
I want to create an SQL View in my local Geoserver (2.4.4), acccording to the handbook.
http://docs.geoserver.org/stable/en/user/data/database/sqlview.html

They say:

  In order to create a SQL View the administrator invokes the Create new layer page. When a database store is selected, the usual list of tables and views available for publication appears, A link Configure new SQL view... also appears:

My Questions:

Where can I find the "create new layer page"?
Where do I select a database store then.

I simply cannot find any of these.


texts: Exporting CartoDB dataset for create report?

texts:
How could I export the Carto's dataset to then generate a report of that data?

With the idea of giving that report to some architects.

Please excuse my limited English.


texts:
I tried using rgdal::writeOGR in R to overwrite a SHP file that was open in QGIS.

Error in writeOGR(x, "C:/myFiles", "testFile", :
   Creation of output file failed

Not only did it fail to overwrite the file, but it also deleted all the following auxiliary files:

testFile.cpg
testFile.prj
testFile.shx
testFile.sbx

Basically everything except the dbf and shp files. This seems to be standard behavior, and happens every time that I have a file open in QGIS. This also happens whether or not the file is being edited in QGIS.

I'm aware that it's possible to regenerate some of these to save the file, but I'd prefer to avoid this problem altogether by catching a QGIS file lock.

Is there any way (in Python, R or any other) to detect that QGIS has locked a file BEFORE trying to overwrite it? My operating system is Windows 10.


texts: Why is my interpolated map seemingly in the wrong position relative to the points?

texts: I am trying to interpolate a surface based on these points but the output map looks like it's offset. I can't figure out why. They both have the same CRS WGS84 4326.


texts: ArcSDE 10: How to read/write to Oracle with ST_GEOMETRY using Java

texts:
I found the following in the deprecation plan of ArcGIS 10.2:

  ArcGIS 10.2 will be the last major release to support the ArcSDE SDK
  with the ArcSDE C and Java APIs. Today many other options are
  available for developers, including SQL, which is available as a
  result of the widespread adoption of spatial types...

In our recent Java/Oracle project, we want to write/read Features to an Oracle DB, which uses ST_GEOMETRY instead of Oracle's SDO_GEOMETRY.

So my questions is:

Is it a bad idea to use the Java SDE API when starting from scratch?
If SQL is the alternative - how can I do this? Shouldn't there be JDBC Support for ST_GEOMETRY in Oracle? I could not find anything like that.
If there is no JDBC Support - isn't the SDE API the only path to go currently? I don't see "many other options".


texts: Seeking website to do geoprocessing online?

texts:
Is there a website in which I can do geoprocessing online?

I know there are a lot of sites that allow the visualisation of your data online or some basic map-making, but I need something more...
Kind of like a free online GIS.

For example, go on a website, draw a polygon and extract the land cover inside.


texts: Is there any autocomplete filter widget for QGIS connecting to PostGIS?

texts:
In QGIS, there is Value Relation widget for fields, that contain references to data in external tables or views. While it works fine for small number of choices, it can quickly become overwhelming, if there are many of them (e.g. when trying to relate multiple businesses with one map item. It seems, that this is a case for autocomplete widget, similar to http://ivaynberg.github.io/select2/ but I am not aware of such a widget and Google didn't help me.

Do you know about such plugin, or should I turn to custom programming?


texts:
How can I copy "system:time_start" to other formulas from collection?
Like EVI:
var evi2020 = med2020.expression( '2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))', { 'NIR': med2020.select('B8'), 'RED': med2020.select('B3'), 'BLUE': med2020.select('B2') }).rename('EVI2020');

or to other formula:
var formula = image.expression( '(ndvi2019 - ndvi2018) / mean', { 'ndvi2019': image.select('ndvi2019'), 'ndvi2018': image.select('ndvi2018'), 'mean': image.select('mean') }).rename('formula');


texts: How to find coordinates of pixels of a GeoTIFF image with Python

texts:
I have a GeoTIFF image that looks like the following
f = rxr.open_rasterio('myFile.tif')
f5 = f.sel(band=5)
f5.plot.imshow()

I would like to create a Pandas dataframe with the values of each pixel and the coordinates of the centroid of the pixels. How can I find the values of the centroid of each pixel and and two columns lat and lon to the following dataframe?
a = np.ravel(f5)
df = pd.DataFrame({'LSTDay':a})
df.head()

    LSTDay
0   15631.0
1   15647.0
2   15624.0
3   15624.0
4   15590.0

Here the link to the GeoTIFF file


texts: Looking for interactive map visualization software that allows user input

texts:
I am looking for a software tool that allows to create maps which can be customized by the user later on. To be more concrete: the user should be able to choose different indicators, create his own aggregate scores, attach weight to the indicators and so on. One example is the OECD better life index, it contains these features (you can choose which indicators you want to include, you can weight them, etc). Do you know such a visualization program? The well known ones like Geocommons, ArcGIS, Spatialkey, etc. do not have these features.

Thanks a lot for your help,

Stephanie


texts: Applying selected features from output group layer to same existing group layer?

texts:
I am doing a geometric network project, and am running 3 trace process within a model to find connected features to a set location on the network. Each trace function wants to output its own group layer of all participating feature classes in the geometric network, with each group layer containing the selected features of the 5 feature classes from that operation (eg, 4 valves to turn, 2 impacted hydrants, 7 selected orphan junctions, ect.). The first operation's selected features are copied into the second, but the third has a unique section set based on a downstream trace, so I really only need the second and third group layer outputs.

Is there a way to apply these individual selected features from the two group layers to the original group layer or combined them into a single group layer?

I believe I would need them in one group layer because I need to be able to publish the mxd and model to a server, have the model run, and be able to display valves that need to be turned to isolate the selected point, and the mains, hydrants, and buildings that are impacted.

Here is an example of the different selections from the group layers, top being the second trace, middle is from the third, and bottom is the original group layer used for data viewing, editing, ect.

Is there a way to have these apply to the original group layer?

Also, the features from the layers are all referencing the same data source.

Also, i know a bit of arcpy programming, and can find the layers that contain selected features, but am stuck on where to go next. Is there a simpler way then just doing a bunch of select by locations to apply the selections to the original feature classes?


texts:
I want to open a "Raster Dataset" in QGIS. In ArcMAP I can just add it like shown in the Screenshot (only the raster shows up).

If I am trying to do this in QGIS. However there I don't see the raster files itself. I find three files (.ovr, .aux, and .ovr.aux.xml) and a folder (with .adf files inside it). How can I open the raster?

A similar question was already asked here and here, but never really answered.


texts:
I'm developing a GIS solution which uses only geometry data. I designed my spatial DB with PostGIS and now I'm performing some trials connecting the DB with GeoServer map server (using the administration page).

In the layer definition of the administration tool there is the mandatory field SRS(Spatial Reference System) to set with the appropriate SRID

I'd like to know which is the default SRID value to set in order to take into account the geometry type. I read that should be 0, but it doesn't work...


texts:
I am trying to use the distance matrix to measure the distance between a crime (input ID) to the closest subway entrance (target layer). The problem is that sometimes the matrix comes out with results, sometimes it doesn't and instead, the algorithm fails.

I get this error in the distance matrix log

  Traceback (most recent call last): File
  "C:/PROGRA~1/QGIS3~1.0/apps/qgis/./python/plugins\processing\algs\qgis\PointDistance.py", line 123, in processAlgorithm matType, nPoints, feedback) File
  "C:/PROGRA~1/QGIS3~1.0/apps/qgis/./python/plugins\processing\algs\qgis\PointDistance.py", line 155, in linearMatrix fields, out_wkb, source.sourceCrs())
  Exception: unknown

  Execution failed after 1.39 seconds

  Loading resulting layers Algorithm 'Distance matrix' finished


texts:
I am trying to classify a shapefile based on rule set, following this tutorial.

This is the result:

What I want to ask is:
1. How to insert the classification result (water, land etc) into the attribute table? I want to make each feature assign to a certain class in the attribute table.

How to join the vector contain same class? I tried vector-geoprocessing-dissolve, it is not working. I tried to select the polygons using same rule, then using "merge selected feature". However, the attributes value only following the first feature, meanwhile, I need the mean of each feature. Is there any way to merge the same category, and calculate the mean attributes value?


texts:
Does anyone know how to offset a stroke in Geoserver? For example I will color the left side of the street (line feature) red and the other side (same line feature) blue. I believe I will be using 2 LineSymbolizer for this but I don't know how to offset them. Otherwise, they will just be rendered on top of each other.

Any thoughts?


texts: Does cropping and adding an alpha layer to a raster map yield better georeference results?

texts:
I have a collection of nine high-resolution map scans illustrating former borders between two nations. The scale is about 1:13,265 for each of the nine images. I want to georeference these maps using QGIS so that I can trace the historical borders, but I want to know whether processing the images - in GIMP, for example - will influence the accuracy of the results.

Right now, each map is on a rectangular sheet with white background. If I crop everything outside the frontier zone and replace the cropped out space with an alpha layer, would this image produce more accurate, less distorted results than the rectangular originals? Using the same ground control points for any given transformation type (e.g., linear, polynomial-1, thin plate spline), would it yield larger or smaller residuals?


texts: Can't export geometry generator styling to DXF - alternatives?

texts:
I've been using QGIS' wonderful geometry generator to dynamically generate root protection zones and crown areas for some tree inventories I've been working on. It's been a massive time saver.

Unfortunately our clients usually insist on CAD drawings for tree plans, and the robust DXF export in QGIS does not support this level of symbology. Neither does OGR2OGR from what I can gather?

I'm still new to GIS, and our workflow is such that I'm not yet confident migrating our data to PostGIS and using ST_MakeLine etc. to dynamically generate geometries.

Are there any other quick alternatives to dynamically generating geometry, which can be successfully exported to DXF?
(or is there a way to export this sort of styling to DXF?)

Just as a reference, here's what I've been using to generate the geometries from the tree points. RPZ and N/E/S/W are fields in the data, in meters.

RPZ (Red):

Radius line:

make_line(
 $geometry,
 make_point($x + "RPZ"*cos(radians(45)), $y - ("RPZ"*sin(radians(45))))
)

(2019 update: The above doesn't work for multipoint geometries; you need to convert them to single parts by getting a centroid of the buffer like so:)

make_line(
centroid(buffer($geometry,0.1)),
 make_point($x + "rpz"*cos(radians(45)), $y - ("rpz"*sin(radians(45)))))

Circle:

buffer($geometry,"RPZ")

Crown area (blue):

make_polygon(make_line(
     translate($geometry,0,max("N",0.2)),
     translate($geometry,0.6*max("E",0.2),0.6*max("N",0.2)),
     translate($geometry,max("E",0.2),0),
     translate($geometry,0.6*max("E",0.2),-0.6*max("S",0.2)),
     translate($geometry,0,-max("S",0.2)),
     translate($geometry,-0.6*max("W",0.2),-0.6*max("S",0.2)),
     translate($geometry,-max("W",0.2),0),
     translate($geometry,-0.6*max("W",0.2),0.6*max("N",0.2)),
     translate($geometry,0,max("N",0.2))))

And the current output looks like this (pretty fugly crown area, I know :-) - are splines even an option?)


texts: Python script to convert Lat Long to ITM (Irish Transverse Mercator)

texts:
I am currently working on a python script that converts Lat Long (WGS84) to ITM. The formulas used are based on a paper published by OSGB (Appendix B):

http://www.leica-geosystems.co.uk/downloads123/gb/gps/gps1200/other/OS%20Transformations%20and%20OSGM02%20user%20guide_en.pdf.

I get good results in the Eastings but not in the Northings (from few mm to approx. 1.5m) and I can't understand why. Does anyone know the reasons?
Thank for your help.

Andrea

# this script convert Latitude and Longitude to ITM

# import modules
import math

# GRS80 Ellipsoid constants
a = 6378137.0
b = 6356752.3141
e2 = (a**2-b**2)/a**2

# ITM projection constants
F = 0.99982
Lat0 = 0.933751149817
Long0 = -0.13962634016
E0 = 600000.0
N0 = 750000.0

# input Lat degrees, minutes and seconds for latitude
dLat = float(raw_input('Enter latitude degrees: '))
mLat = float(raw_input('Enter latitude minutes: '))
sLat = float(raw_input('Enter latitude seconds: '))

# convert Lat to decimal degrees
ddLat = dLat + (mLat/60) + (sLat/3600)
print 'Latitude dd:', ddLat

# convert Lat to radians
rLat = (ddLat*math.pi)/180.0
Lat = rLat
print 'Latitude r:', Lat
print

# input Long degrees, minutes and seconds for longitude
dLong = float(raw_input('Enter longitude degrees: '))  #The degrees have to be entered with a posite sign
mLong = float(raw_input('Enter longitude minutes: '))
sLong = float(raw_input('Enter longitude seconds: '))

# convert Long to decimal degrees
ddLong = dLong + (mLong/60) + (sLong/3600)
print 'Longitude dd:', ddLong

# convert Long to radians
rLong = -(ddLong*math.pi)/180.0 # the negative sign account for west direction
Long = rLong
print 'Longitude r:', Long
print

# calculate constants for transformation
n = (a - b)/(a + b)
v = a*F*(1-e2*(math.sin(Lat)**2))**-0.5
p = a*F*(1-e2)*(1-e2*(math.sin(Lat)**2))**-1.5
n2 = v/p-1
print 'n:', n
print 'v:', v
print 'p:', p
print 'n2:', n2
print

M1 = (1+n+5.0/4.0*n**2+5.0/4.0*n**3)*(Lat-Lat0)
M2 = (3*n+3*n**2+21.0/8.0*n**3)*math.sin(Lat-Lat0)*math.cos(Lat+Lat0)
M3 = (15.0/8.0*n**2+15.0/8.0*n**3)*math.sin(2*(Lat-Lat0))*math.cos(2*(Lat+Lat0))
M4 = 35.0/24.0*n**3*math.sin(3*(Lat-Lat0))*math.cos(3*(Lat+Lat0))
M = b*F*(M1-M2+M3-M4)
print 'M1:', M1
print 'M2:', M2
print 'M3:', M3
print 'M4:', M4
print 'M:', M
print

I = M+N0
print 'I:', I
II = v/2*math.sin(Lat)*math.cos(Lat)
print 'II:', II
III = v/24*math.sin(Lat)*(math.cos(Lat))**3*(5-(math.tan(Lat)**2)+9*n2)
print 'III:', III
IIIA = v/720*math.sin(Lat)*(math.cos(Lat)**5)*(61-58*(math.tan(Lat)**2)+(math.tan(Lat)**4))
print 'IIIA:', IIIA
IV = v*math.cos(Lat)
print 'IV:', IV
V = v/6*(math.cos(Lat)**3)*(v/p-(math.tan(Lat)**2))
print 'V:', V
VI = v/120*(math.cos(Lat)**5)*(5-18*(math.tan(Lat)**2)+(math.tan(Lat)**4)+14*n2-58*(math.tan(Lat)**2)*n2)
print 'VI:', VI
print

# calculate Eastings and Northings
N = I+II*(Long-Long0)**2+III*(Long-Long0)**4+IIIA*(Long-Long0)**6
E = E0+IV*(Long-Long0)+V*(Long-Long0)**3+VI*(Long-Long0)**5
print 'Easting: ', E
print 'Northing:', N


texts:
I would like to add to the attribute table somehow quickly a column containing the entry and exit kilometre data.
I have the available point layer with kilometre, road axis and polygon as shown:

Usually I did it by changing the training ground to a line and then creating points at intersections. These points were later connected with the help of NNJoin plugin or similar tool with kilometre. The last step was to manually glue the values and organize the tables.

I'd like to automate (or speed up ) this process but unfortunately I don't have enough experience with using Python to write appropriate script. (I'm just learning it).

So does anyone know how to do it in another faster and a better way?


texts: Why are the parameters in a ModelBuilder generated tool disabled?

texts:
I created an interactive tool for generating Euclidean Distance Rasters and Reclassifying them afterwards. For the main input features I checked the:

"A list o values radio button" inside the
"This value contains:" fieldset, from the Input Layers' properties.

When running the model from ArcToolBox I am able to insert the layers normally but only to change the parameters for the first layer. What is causing this problem?

The next image shows how the parameter fields are grayed out:


texts: Finding all polygons inside of a polygon

texts:
How to find all the polygons inside of a polygon? I tried ST_Intersects and ST_Contains but they are booleans variables. I have a list of polygons in WKT format, so I am trying to pull all the polygons given a region in their WKT format.

Does anyone have an idea how to approach this?

select geom from polygon where
ST_Intersects(ST_GeographyFromText
('POLYGON((210000 2400000, 300000 2300000, 330000 3708400, 210000 2400000))'),
polygon.geom);


texts:
I wish to display land seismic lines on a map. Land seismic lines are often not straight so often there is survey data for each shotpoint or every 5 shotpoints. Typically a shotpoint (SP) could be 25m or 50m apart. There are often many hundreds of SPs. The TDL loads as discrete points. I need to be able to identify the seismic line (that is in the TDL file) and also say every 100 SP.

I cannot find a way to label SPs every 100 points, or to just label the line at the start and end.

A typical file looks like:

Rad PROFIL_ID   SP  SWRF99TM_N  SWRF99TM_E
1   MC-81-313E  4   6414971 735192  6414387 1687468
2   MC-81-313E  5   6414988 735192  6414404 1687468
3   MC-81-313E  6   6415004 735192  6414420 1687468
4   MC-81-313E  7   6415020 735192  6414436 1687468

PROFIL_ID = Line number, SP is Shotpoint and the other 2 columns are XY coordinates.


texts: How to add a lat/lon grid to a map using azimuthal projection?

texts:
I have a north azimuthal projection, equal areas.
I want to create a map from it (in QGIS) and display a grid on it for better orientation.

Thats how it looks like:

Why doe my degrees have so many zeros?
Why isnt it changing from South to North when it obviously is in the north?
And why does my grid look like this?
Any ideas on that?
Can QGis not display a good grid on this kind of projection?

My project is in the projection and all my shapefiles have this projection...


texts: Renaming layer in multiple maps in ArcGIS Pro with ArcPy?

texts:
Running ArcGIS Pro 2.2, I have a layer in ~20 maps that needs it's name changed. I pulled this code together from a bunch of snippets, and have done my best to make sure it's updated appropriately, but I've hit a point where it runs without an error, but there are no changes afterwards.

What am I missing here?

import arcpy
aprx = arcpy.mp.ArcGISProject(r"C:\GIS_Local\Working.aprx")
for m in aprx.listMaps():
  layers = m.listLayers()
  for lyr in layers:
    if lyr.name == "old name":
      lyr.name = "new name"

There's a similar Q&A from a couple years ago - Use ArcPy to rename MXD layers, but I haven't been able to get the code there to update to Python 3.0 correctly, so I pulled this from the samples on the ESRI help and this GeoNet post


texts: shp2kml license expired, looking for free solution

texts: So I work for a small African non-profit. I do all the GIS and data work. Most of the layers I create, I convert to kml layers so my boss can then use them quickly and easily in Google Earth. But now in the last months my shp2kml program says its license is expired and to download the newer version. I've downloaded the most up to date version I can find which isn't that new Feb 2011 and I'm still getting the same error. Anyone know another freeware program that does the same thing or know a way around this or to get a more up to date version of shp2kml. My operating system is windows 7 and i'm using QGIS 1.8.0. The non-profit can't afford Google Earth Pro edition. I need free or practically free solutions. Thanks for your help.


texts: ArcGIS Time Slider that controls two maps/dataframes

texts: I have two maps/dataframes on the same webmap with various data throughout the USA. Both maps use a geocoder to be able to zoom to your area of interest in the USA and compare one place with another. I am using the same arcgis server layer for both maps/dataframes, but I can't figure out how to get the timeslider to control both maps/dataframes...any clues or examples out there?


texts:
I'm running MapServer on Windows (with the msplugin_oracle.dll plugin) to create some WMS layers from SDO geometry stored in Oracle 11 database. I take the binaries from this very nice site. The problem I have is that the client who runs the database is planning to swap to using Oracle 12.

Can this cause any problems?

Or will it run smoothly?

Does anyone have any experience with getting SDO geometry out of Oracle 12 using MapServer?

Couldn't find much info online.


texts: Installing LASTools to QGIS on a Mac (or other LAS processing plugins)

texts: I recently saw this link how to install LAStools toolbox to QGIS. The problem is that the instructions are only for Windows users. Now I wonder if it's possible to install LAStools on a Mac with QGIS or if there are any other LAS plugins that could be used on a mac?


texts: QGIS Calculate Field with Function Editor invalid expression

texts:
I have some data which is built like this:

a,bbb,ccc

where bbb and ccc vary in length.

I want to extract the last segment of the string, after the second ",".
Using the Function Editor:

from qgis.core import *
from qgis.gui import *

@qgsfunction(args='auto', group='Custom')
def sp_split(value1):
    return value1.split(",")[2]

I save it as a function called sp_split, but when I call it in the expression with a field for a value, it doesn't recognize it.

Pressing More Info gives me this explanation:

Eval Error:
sp_split() takes exactly 1 argument (3 given)

What might be the reason for that? I'm at a loss here.
I'm using the function because of the split, if there is another way of using the built-in expression to retain the "ccc" part, I'll gladly use it.


texts: Dark Object Subtraction using Top of Atmosphere Reflectance?

texts:
I'm working on a WorldView 3 image in QGIS and have been told to use Dark Object Subtraction (DOS) to atmospherically correct it to surface reflectance using the raster calculator. I have the Top of Atmosphere (TOA) reflectance and the original Digital Number (DN) image.

I have only been able to find guides and papers (e.g. Song et al 2001 Classification and Change Detection Using Landsat TM Data: When and How to Correct Atmospheric Effects? and Chavez 1988 An improved dark-object subtraction technique for atmospheric scattering correction of multispectral data) using DN values to calculate path radiance. These don't use TOA reflectance.

As I have TOA reflectance if I identify dark objects and subtract the mean TOA reflectance value of these from the image would this be a valid DOS method?


texts: Can not load reprojected shapefile - GDAL Python ERROR 4

texts:
I've created a function based on the code snippet of the gdal cookbook to reproject a shapefile to a new crs.
However if I try to load the reprojected shapefile I'll get ERROR 4: 'SamplepointsEPSG_5678/Reprojected.shp' not recognized as a supported file format.
What I've found so far - if I try to load the reprojected shapefile to QGIS no crs is set for it. When I set it there manually the location fits.
def reprojectShape(inputShape, EPSGout, outPathName):

  driver = ogr.GetDriverByName('ESRI Shapefile')

  # input SpatialReference
  inLayer = inputShape.GetLayer()
  proj = inLayer.GetSpatialRef()
  EPSGin = int(proj.GetAttrValue('AUTHORITY',1))

  inSpatialRef = osr.SpatialReference()
  inSpatialRef.ImportFromEPSG(EPSGin)

  # output SpatialReference
  outSpatialRef = osr.SpatialReference()
  outSpatialRef.ImportFromEPSG(EPSGout)

  # create the CoordinateTransformation
  coordTrans = osr.CoordinateTransformation(inSpatialRef, outSpatialRef)

  # get the input layer
  inLayer = inputShape.GetLayer()

  # create the output layer
  outLayerName = "Reprojected"
  if os.path.exists(outPathName):
    driver.DeleteDataSource(outPathName)
  outDataSet = driver.CreateDataSource(outPathName)
  outLayer = outDataSet.CreateLayer(outLayerName, geom_type=ogr.wkbMultiPolygon)

  # add fields
  inLayerDefn = inLayer.GetLayerDefn()
  for i in range(0, inLayerDefn.GetFieldCount()):
    fieldDefn = inLayerDefn.GetFieldDefn(i)
    outLayer.CreateField(fieldDefn)

  # get the output layer's feature definition
  outLayerDefn = outLayer.GetLayerDefn()

  # loop through the input features
  inFeature = inLayer.GetNextFeature()
  while inFeature:
      # get the input geometry
      geom = inFeature.GetGeometryRef()
      # reproject the geometry
      geom.Transform(coordTrans)
      # create a new feature
      outFeature = ogr.Feature(outLayerDefn)
      # set the geometry and attribute
      outFeature.SetGeometry(geom)
      for i in range(0, outLayerDefn.GetFieldCount()):
          outFeature.SetField(outLayerDefn.GetFieldDefn(i).GetNameRef(), inFeature.GetField(i))
      # add the feature to the shapefile
      outLayer.CreateFeature(outFeature)
      # dereference the features and get the next input feature
      outFeature = None
      inFeature = inLayer.GetNextFeature()

  # Save and close the shapefiles
  inDataSet = None
  outDataSet = None

  # return
  outName = os.path.join(outPathName,".".join([outLayerName,"shp"]))
  return(outName)

inName = os.path.splitext(Psamplepoints)[0]
outName = "".join([inName,"EPSG_",str(EPSGout)])
PreprojectedShape = reprojectShape(inputShape = input_samplepoints, EPSGout = targetEPSG, outPathName = outName)
input_samplepoints = gdal.Open(PreprojectedShape)


texts:
I have install leaflet-sidebar-v2 via npm. However when I try to simply initialize the side bar I get this error.

  L.Control.Sidebar is not a function.

here is my index.html import scripts etc

<script src="https://unpkg.com/leaflet@1.0.1/dist/leaflet.js"></script> <script src="https://unpkg.com/leaflet-sidebar-v2"></script>
here is my component html snippet for the sidebar and map

here is my ngOnInit Method
ngOnInit() {
this.isBusy = true;
this.initializeForm();

this.overlayKeys = getOverlayKeys();

// our markers
this.markers = new Array();

// our geoJson urls sources
if (!this.geoJsonUrls) {
  this.geoJsonUrls = new Array();
} else {
  // keep this in for testing
  console.log(this.geoJsonUrls);
}

// the geoJson from each of the geoJsonUrl resources
this.geoJsonResults = new Array();

// let create our layers
this.layers = new Array(this.geoJsonUrls.length);

// get our geo location and if we have it initialize our map
if (navigator.geolocation) {
  this.latitude = 0;
  this.longitude = 0;
  navigator.geolocation.getCurrentPosition(position => {
    this.latitude = position.coords.latitude;
    this.longitude = position.coords.longitude;
    this.map = new L.Map('map', {
      zoomControl: false,
      center: new L.LatLng(this.latitude, this.longitude),
      zoom: 19,
      fullscreenControl: true
    });

    this.map.addControl(L.control.zoom({position: 'topright'}));

    // carto retina layer
    const retina = '@3x';
    const params = new Array();
    params[0] = retina;

    // get our overlay
    this.overlay = getOverlay(this.overlayKeys[1], params);

    // create the title layer
    const wmsLayer =
      L.tileLayer.wms(this.overlay, null).addTo(this.map);
    wmsLayer.addTo(this.map);

    this.currentLocation = L.marker([this.latitude, this.longitude], {icon: this.customIcon}).addTo(this.map);

    // add popup to marker and open automatically
    this.currentLocation.bindPopup('<b>Your location</b>').openPopup();
    // our map click handler
    this.map.on('click', this.mapClick);

    // our mouse move handler
    this.map.on('mousemove', this.mouseOverMap);

    // wait until our map loads
    this.map.on('load', this.openGeoJsonResources());
  });
} else {
  alert('No Geolocation');
}
this.showButton = document.getElementById('showBtn') as HTMLButtonElement;
this.txtLat = document.getElementById('txtLat') as HTMLInputElement;
this.txtLng = document.getElementById('txtLng') as HTMLInputElement;

// initialize our services
this.jsonpClient = new JsonpClientComponent(this.jsonpService);

// add our local storage event listener appears to fire only when storage is deleted
/*    window.addEventListener('storage', (e) => {
      let a: any;
      a = e;
      if (a.data) {
        this.setGeoJson(a.data);
      }
    });*/

// we need to do a deep copy of the passed in geojson urls into the selected items to select them all
for (let i = 0; i < this.geoJsonUrls.length; i++) {
  this.selectedItems[i] = this.geoJsonUrls[i];
}

// lets setup of our dropdown information
this.dropdownSettings = {
  singleSelection: false,
  text: 'Select GeoSpatial Items to View',
  selectAllText: 'Select All',
  unSelectAllText: 'UnSelect All',
  enableSearchFilter: true
};

this.customMarkerHelper = new CustomMarker();

debugger;
this.initializeSideBar();

this.isBusy = false;
}


texts:
I am having issues bringing an AutoCAD drawing into QGIS. I currently use AutoCAD with Carlson (surface mining and civil modules). My current workflow for maps is this...
Import points using Lat/Long.
Use the "Adjust Coordinates" command to convert from Lat/Long to the state plane of that particular job.
Do my drawings I need to in CAD.
I would then like to visualize my drawing in QGIS. I know this is possible because a surveyor sent me his .dwg the other day, I saved it as a .dxf, and dragged and dropped into QGIS and it went to the correct location.
When I try this with my drawings, I cannot get it to show in the correct location. I can get it close (a couple hundred miles) by messing with different CRS.
I don't know if it's not working because of how I am importing the data into CAD, or exporting it from CAD, or importing to QGIS, or what.
I am using QGIS 3.8.1-Zanzibar and AutoCAD 2020.


texts: Searching by attribute using each row with ArcPy cursor?

texts:
I am trying to write a script to select each individual neighborhood using one FC, then use select by location to see the number of lots from another FC that are in the neighborhood and store that number in the neighborhood FC. I am trying to figure out how to write the script of selecting each neighborhood individually using arcpy.da.UpdateCursor.

This is my first time trying to write a real script outside of school.


texts: Are there any AppStudio for ArcGIS tutorials?

texts:
Any other documentation available regarding a detailed tutorial on how to use the new beta 2 app studio.

I have trolled the Internet and found the minimal information docs on esri regarding the app studio but still need assistance.

Does anyone know where to locate some additional resources, or document the baseline process from a Web application to the final product app on an android?


texts: In which way is height information stored in a DEM?

texts:
I'm struggling with different DEM data, because I don't know how to convert the height information in meters.

When I load ASTER and SRTM tifs in QGIS, it shows me the spectrum of the grayband (0 to max). All the time I thought I could just take this value and see it as the height of each Point in m. Later I compared the graybandvalues with the actual height of the objects and realised that it's terribly off (German Zugspitze 500m, Test Area Swellendam 600m, a mountain in Malaysia more than 2000m). I tried it with different SRTM and ASTER DEMs, but no one was correct.

I would be very happy, if someone had an idea, how to solve this problem.


texts: How do I change the font in the QGIS Script Editor?

texts:
How do I change the font in the processing script editor?  I loaded 3.4 today and it works, but the font in the editor is not set to a monospace font for English, making spaces and positions extremely confusing.  I would like to change it to a standard monospace font.  I only see a size control in the editor, and I don't see anything in the main QGIS options list.

Update: I installed the mac version, and it has the same problem.  This is pretty awful.  It wouldn't be a big deal, but the QGIS team doesn't seem to have planned in support for external IDEs.  Here's a pic.


texts: Format of raw LiDAR data (time-of-flight data)?

texts:
I am doing a project where I need to simulate LIDAR data, and I am looking for specifications on how raw LIDAR data looks.

I want the format of the data before it is translated to x, y and z coordinates (using GPS, IMU etc.), i.e., I imagine each data point having a direction and a distance (or time-of-flight).

I have not been able to find this anywhere, every source reference raw data given in x, y and z coordinates in some form.


texts: Removing (defined) sliver polygons using PostGIS?

texts:
How do I write SQL (for PostGIS in PostgreSQL) to remove sliver polygons (including long narrow slivers and some small irregular polygons that were erroneous 'gaps' in the original dataset) from my dataset by merging with an adjacent polygon with the longest edge?

I have a dataset (tablename: "Vegetation" containing several thousand polygons identified as either a sliver or a 'real' polygon in a description column (column name: "desc") as "sliver" or "main". [edit: Prior to importing my dataset into PostGIS I cleaned it up by importing into GRASS (which is also how I converted all the gaps into small polygons, which I subsequently merged with the main dataset.  I found if I increased the snapping tolerance and minimum area past '1' in the GRASS import, I started loosing desired features.]

My understanding is that I may need the following steps:

Identify the polygon (Vegetation.desc="main") adjacent to each sliver (Vegetation.desc="sliver") sharing the longest edge.
Identify the edge between each sliver and it's adjacent main polygon sharing the longest edge identified in Step 1. (Possibly using ST_Relate)
Remove the edge between each sliver and it's adjacent main polygon sharing the longest edge using the edge identifier returned in Step 2. (Possibly using ST_RemEdgeModFace)

I have only a very basic understanding of SQL and the above process is beyond my ability.


texts:
I want to reclassify a raster (variable) in model builder in two classes. First class ist 0-treshold = NoData and second is treshold-max value = 1. It works fine if i use an existing raster file as input (where i can state and know the max value) but not for variables in ModelBuilder.

Any ideas how to do it in ModelBuilder, or with ArcPy?

EDIT:

Here is my edited model, but it's not running yet.


texts: on Arcgis Android SDK Getting codedvalues from ArcGISTiledMapServiceLayer

texts: i have arcgis mapview, and when user click on Mapview, i get point which is clicked after that i search with IdentifyParameters class and i am getting information about layers which is clicked place. But there is a problem, when i am getting information, alias and value is shown but my value has coded value so user see number but exactly value is string. how can i get codedvalues  ?


texts:
I have a polygon layer and every polygon has a unique gml_id. I want to select multiple at the same time. Normaly I use selection by expression and it works very well. I fount the tool select by value and wonder how it works to select multiple features.
Different inputs that doesn't work

I´ve tried different spellings in the input mask, but it doesn't work.
What is the correct synthax?
This is my "select by expression"

I use the QGIS Version 3.4.8 LTR


texts: Is it possible to create road-type dependent 'shields' in QGIS?

texts:
In QGIS, I need to place road numbers on a map (they are point labels) in solid rectangular boxes to make them stand out. Ideally they need to be the same colour as the road they refer to, so a motorway number should be in a blue box because the motorways are blue etc.

The boxes also need to scale with the label size, as some roads are for example 'A9', whereas others are 'A4567'. The roads and their labels have codes defining their type so this could be used to set the colours.

Is this possible?


texts: How to select the vertices from a polygon?

texts:
I am writing some VBA code to extract vertices from a polygon.

I would like to ask, how to select the vertices from a polygon ?


texts: Min and max in a field attribute with ArcPy?

texts:
I have a shapefile with a field attribute called "Value".

I would like to find the min and max in this field and

if max-min > 0.014 return 1 (or Yes)
if max-min < 0.014 then return 2 (or No)

in the field called "YES_NO".


texts: Both-sided transparency in ArcScene for extruded features

texts:
I'd like to extrude my polylines and make them transparent, but from both sides. As you can see below, the transparency works only on one side. I think it's something with the side-orientation.

I can create two halves of circle with inverse orientation and extrude them separately, but it won't work after a change of view.


texts:
I recently upgraded QGIS by building QGIS 2.11 from source against GDAL 2.0 on my Ubuntu 14.04 system. After a few teething troubles I got it to compile using QT creator. I now find that I can't open any spatialite databases e.g. openstreetmap with the error being returned:

obsolete libspatialite: connecting to this DB requires using v.4.0 (or any subsequent)

When I check it appears that QGIS was built against SpatialLite Version 4.1.1 so not really sure what has happened here. I can connect to all of my old spatialite dbs using the DB Manager in QGIS which allows me to view the attribute tables (but no previews). If I try and 'add to canvas' from DB Manager I get an error message telling me to check the logs with the following output:

dbname='/home/dori.osm.db' table="dori_buildings" (geometry) sql= is an invalid layer - not loaded

Equally I can explore the data with the QSpatiaLite plugin but if I try and load the data from the plugin into QGIS I get the following error:

  Impossible to Load table 'dori_buildings' into QGIS

Anyone have any ideas what's happening here - did I miss out a dependency when building perhaps?


texts: Labeling polygons with the number of related table entries in QGIS

texts:
Using QGIS 3.6.1 I have a layer with areas and a table with entries related to the areas:

A relation is defined like this:

What I would like to achieve is to label the areas with the number of related table entries using relation_aggregate function as shown here:

relation_aggregate ('my_relation', 'count', 'text') || ' Table entries'

As can be seen in the screenshot (Output preview), the label expression delivers a valid result as expected (e.g. "4 Table Entries"), but labels are not shown in the map:

Another thing I tried using the aggregate function is the following:

with_variable('my_id', "id", aggregate('table', 'count', 'id', "flaeche_id" = @my_id)) || ' Table Entries'

Which gives a totally weird result depending on which area is shown in the current map extent, i.e.

1) Both areas labelled with "4 Table entries", but the lower has only one related table entry (see table above):

2) If the map extents only to the southernmost area, the label shows correctly "1 Table entries":

So all in all with my first approach I have no labels and with the second approach I have wrong labels.

I think I'm little stuck with the aggregate... functions in QGIS 3, am I using them the wrong way or are there any other approaches to the labeling I want to achieve?