texts: GMT Pen name c not recognized - Trying to fill polygons

texts:
I am trying to fill my polygons like this person tried:
Fill color in polygons using attribute variable in GMT (Generic Mapping Tool) with some hints here:
https://gmt.soest.hawaii.edu/boards/1/topics/7553
They suggest to:

  Try what the man page for psxy -W says:
  +c Controls how pens and fills are affected if a CPT is specified via -C:
  Append l to let pen colors follow the CPT setting.
  Append f to let fill/font colors follow the CPT setting.

I tried but I get the error:

  psxy: Pen name c not recognized!

when I try the solution suggested.

Code drawing the borders of the geologic layers with the palette colors:

gmt psxy "$repo/geology.gmt" $reg $proj -aZ=Code -C$colornames
-Zzval  -O -K

Code supposed to fill in those layers with the palette colors, but which gives an error:

gmt psxy "$repo/geology.gmt" $reg $proj  -W+c -aZ=Code -C$colornames -Zzval  -O -K

What do I do wrong?


texts: QGIS 2.16 Nødebo installation package for OSX

texts:
Im using qgis on my mac, few days ago the notification of the new available 2.16 version appeared in user interface. I went to the official webpage of the qgis but did not succeed. qgis.org is redirecting you to the kyngchaos for download http://www.kyngchaos.com/software/qgis but there is still  only old 2.14 version for download.

Does anyone knows when it will be eventually available? or maybe someone can suggest where can i get fresh dmg.


texts: Cost surface in R: topographic distance betwen cells' centers in place of planar distance

texts:
Aim
Calculating slope-dependent walking distance taking into account topographic distance NOT planar distance between cell centers.
Tool
'gdistance' R package
Issue
I do not really know how to code to consider topographic distance in place of planar distance between cells' centers
Executable code
What follows is after https://www.jstatsoft.org/article/view/v076i13/0
r <- raster(system.file("external/maungawhau.grd", package = "gdistance"))

altDiff <- function(x) x[2] - x[1]

hd <- transition(r, altDiff, 8, symm= FALSE)

slope <- geoCorrection(hd)

adj <- adjacent(r, cells = 1:ncell(r), pairs = TRUE, directions = 8)

speed <- slope

speed[adj] <- 6 * exp(-3.5 * abs(slope[adj] + 0.05))

Conductance <- geoCorrection(speed)

In the previous code, the Tobler's hiking function has been applied and the speed of movement has been calculated when moving between adjacent cells.
The very last line of code takes into account the distance between cell centers because moving with the same speed along a diagonal connection between cells takes longer than moving along a straight connection.
What I am not sure of is about how to proceed to take into account not the "planar" distance between cell centers, but the topographic distance. I do know that the Pythagorean theorem should be used in the calculation, but I really cant wrap my head around how to implement that in the context of the above code.


texts: How to create a circular sector in Mapinfo?

texts: I want to build a tool that creates a cell sector from data like that below. The sector diagram should be like the second picture. The first picture shows a Mapinfo tab file with the required fields Antenna beam, Azimuth, Band, (the radius), Latitude, and Longitude. How could I achieve this? I had a look at a link in mapinfo-L but the solution was not conclusive!


texts: Map Series Extent with Rotation not working properly

texts: In ArcGIS Pro I have enable map series with % extents and the map series work properly. Then when I try and rotate the layers or dataframe and refresh/re-enable map series the map extents no longer work correctly. I do not have this issue when using arcmap. Any suggestions?


texts: ArcMap: Copying features class to Oracle DB fails: Shape integrity error

texts:
I have a feature class in an ESRI file geodatabase where I have isolated one feature. When I copy the feature class to a oracle db (ST_Geometry as defaults; by copy and paste; also fails with gp tools or by python cursor) I get the error "Shape integrity error". Copying the feature or the feature class to other file gdb feature classes works. It only fails on my oracle db.

Check geometry returns no error. What could be wrong? What can I do?

The shape looks like this:

The domain of the feature class is X=(-20037700 to 900699887774,099) Y=(-30241100 to 900689684374,099). The resolution 0,0001 meters and the tolerance 0,001 meters.

EDIT: It also occurs when I copy it to a SQL Server database.


texts: Buffers of public transportation transit times

texts:
I'd like to be able to create buffers/zones/rings of transit times using Denver's public transportation system.  For an example of what I'm talking about, check out triptropnyc.com - I don't need it to look nearly as nice, and I only have one origin point of interest.

My thought so far was to distribute random points around the area of interest and then pass those points to the Google Maps API a bit at a time.  Once I had the travel times, I could then interpolate between the points to get the buffers - going back from more points if necessary to smooth out the estimates.

Would that be a valid approach?  Are there other ways of estimating transit times that I should consider?  Perhaps other sources of data that would be more useful or easier to get?

Let me know what you think.


texts: Is it safe to combine data for zipcodes across different demographic sets?

texts:
My team is working with demographics data across different data sources (some paid sources and some free data sources available online).

Each of these data sources comes with a shapefile and some attributes associated with each demographic area and could be defined across different cuts of time.

However, when we display these attributes to our end users, we wanted to abstract the multiple datasources concept and show zip codes as a single demographic unit.

We were planning to combine the attributes of all the datasources into one single data source and point that to one of the shapefiles (For the time being, we are willing to look past issues related to granularity or precision in the definition of polygon across these data sources.

However, should we be concerned that the ZIP code to actual geographic area might not be consistent across demographics datasources taken across different cuts in time? E.g. Zip code 12345 used to Map to an area in State A till 2010 but points to an area in State B for all datasets after 2010?


texts: Convert from EPSG 3857 to EPSG 4326 QGIS2.18

texts:
I need to convert EPSG 3857 to EPSG 4326 in a plugin.

            #example EPSG 3857
            x = 5062443.00656
            y = -952576.546977
            #convert here to EPSG 4326
            #where x2 =41.3387949319,
            #where y2 = -8.55714071443
            return x2,y2 #return in EPSG 4326

something like the code above. Is there a way to do that?


texts:
I have a legacy program that attempts to read USGS 1/3 arcsec data.  When I compare results to GDAL, I find that the lat/lon to pixel conversion is off by 1 line.  The reason (I think) is because the legacy code calculates the latitude of line 0 as:

lat0 = yllcorner + (nrows - 1) * cellsize

where yllcorner, nrows and cellsize come from the .HDR file.

I see in the GDAL code in frmts/raw/ehdrdataset.cpp that there is no "subtract 1" term:

 dfULYMap = dfYLLCorner + nRows * dfYDim;

I'm assuming that GDAL is right.  Can anyone explain why so that a noob can understand?


texts:
I am stuck with placing Label on the Grid lines. Below is the data grid on the map. I created the grid using data driven pages.
After I turned on the Label, the labelling appeared inside the grid. I want it to be placed on the line itself with a little buffer around it. I have explored the placement properties of the grid layer and it looks as if there is no option there

I wish I could placed the Labeling on the gridLines/borders instead of inside the grid.


texts: How to find wmsURL? And how to use QueryTask method?

texts:
I have couple of questions.

var wmsLayer = new esri.layers.WMSLayer(wmsURL); with this line of code i can display wms layer over any base Layer. Now From where I get wmsURL ? I have published my shapefile on GeoServer, but i dont know how i get the URL (wmsURL) to use as parameter in above line of code ?
Now, Let's suppose I am successful to display wmsLayer. Now i want to query data against this wmsLayer. Like code below:

queryTask = new esri.tasks.QueryTask("wmsURL"); Can I use wmsURL as a parameter in QueryTask method ? From ArcGIS API reference, I come to know that the URL parameter of QueryTask method is the URL to the ArcGIS Server REST resource. But as u know i am using GeoServer, so how it could be fixed. Need some explanation ?


texts:
I'm trying to create a button that sets the map scale to 24000. I'm working in ArcMap 10.1 and my Python programming skill is rudimentary.

I followed the Python Add-In instructions and the creation was successful - I have a new toolbar and a new button - but nothing happens when I click the button. Here's the code:

import arcpy
import pythonaddins

class ButtonClass13(object):
    """Implementation for ScaleBtn_addin.Scale24 (Button)"""
    def __init__(self):
        self.enabled = True
        self.checked = False
    def onClick(self):
        mxd = arcpy.mapping.MapDocument('current')
        df = arcpy.mapping.ListDataFrames(mxd)[0]
        df.scale = 24000

What's missing?


texts:
I want to store a large LIDAR dataset as a VRT so I can quickly query for elevation values in the entire state of Massachussets. However, I am having trouble actually creating the mosaic. The naive way of performing this operation, with gdal_merge.py is impractical since the dataset is 21GB of memory. Because of the size of the dataset, it also spans multiple UTM zones so I need to translate at least some of the geotiffs. I wrote a brief python script to facilitate this:

from glob import glob
from os import path, system
from tqdm import tqdm
from osgeo import gdal, osr
from gdalconst import *
import subprocess

#proj4 = '+proj=utm +zone=18 +ellps=GRS80 +units=m +no_defs'

#tgt_proj4 = osr.SpatialReference()
#tgt_proj4.ImportFromProj4(proj4)
ds = gdal.Open("LIDAR/19_03344719.tif", GA_ReadOnly)
prj = ds.GetProjection()
tgt_srs = osr.SpatialReference(wkt=prj)
files = list(glob("LIDAR/*"))
files.sort(reverse=True)
print(tgt_srs.ExportToWkt())
for filename in tqdm(files):e
    newfilename = path.join("LIDAR_TRANS", path.splitext(path.basename(filename))[0] + ".tif")
    subprocess.call(['gdalwarp', filename, newfilename, '-t_srs', tgt_srs.ExportToProj4(), '-multi', '-q'])

```

I then simply merge all the items in LIDAR_TRANS into a giant VRT. However, I had issues doing this as I keeping large gaps in my mosaic VRT as shown below

As such, I was wondering how my pipeline needs to be changed to avoid these gaps so I can properly use the dataset?


texts:
Is there a way in qgis 2.0 to prevent overlapping labels?

I have a layer with a large geographical extent, which means a hell of a lot of labels. Many of them collide and it looks really bad. Is there a way for me to get qgis to automatically place the labels on an alternative 'location'?

I know it is possible in MapInfo but I'd rather use qgis :)


texts: Selecting features inside a polygon - PyQGIS

texts:
I have just started working with PyQGIS and I would like to ask this question. I have a rather big dataset (cell phone towers) and I want to narrow it in the towers that exist inside a certain  country and I want to do it with PyQGIS. I used that piece of code:

lyrPts = QgsVectorLayer("/path/Cell-tower.shp", "CellTowers", "ogr")

lyrPoly =QgsVectorLayer("/path/GRC_adm0.shp", "GRC", "ogr")

QgsProject.instance().addMapLayers([lyrPts,lyrPoly])

ftsPoly = lyrPoly.getFeatures()

for feat in ftsPoly:
    geomPoly = feat.geometry()
    bbox = geomPoly.boundingBox()
    req = QgsFeatureRequest()
    filterRect = req.setFilterRect(bbox)
    featsPnt = lyrPts.getFeatures(filterRect)
    for featPnt in featsPnt:
        if featPnt.geometry().within(geomPoly):
            lyrPts.select(featPnt.id())

iface.setActiveLayer(lyrPoly)
iface.zoomToActiveLayer()

Although when I run it I don't get any errors I also don't get any results even when I left it running all night.

Can anyone help me and tell me if the code is good?

For reference the cell tower data I got them from mylnikov.org as a .csv file and made it to .shp file and the polygon from diva-gis.org.


texts: Problem With the Edges When Clipping Raster

texts:
The following code clips the TIFF images based on the GeoJson file. However, when clipping the tiles that has half image and half nothing (white area in QGIS), the nothing half changes to black, which is what I want but the half that has the real image turns to white. How can I get the half that has the real image to have the same thing (data) while the half that has nothing to turn to the black color?
import os
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
from shapely.geometry import mapping
import geopandas as gpd
import rasterio as rio
from rasterio.plot import plotting_extent
from rasterio.mask import mask
import earthpy as et
import earthpy.spatial as es
import earthpy.plot as ep

shape_file = '/.../abc.geojson'
crop_extent = gpd.read_file(shape_file)

lidar_chm_path = '/.../reprojected.tif'
output = '/.../Data/'
if not os.path.exists(output):
    os.mkdir(output)
i, j = 0, 0
while i < crop_extent.shape[0]:
    with rio.open(lidar_chm_path) as lidar_chm:
        try:
            lidar_chm_crop, lidar_chm_crop_meta = es.crop_image(lidar_chm,[crop_extent.iloc[i]['geometry']])
            lidar_chm_crop_affine = lidar_chm_crop_meta["transform"]
            j += 1
            lidar_chm_crop_meta.update({'transform': lidar_chm_crop_affine,
                       'height': lidar_chm_crop.shape[1],
                       'width': lidar_chm_crop.shape[2],
                       'nodata': -999.99})
            name = output + str(crop_extent.iloc[1][2]) + '_image_' + str(j) + '.tif'
            with rio.open(name, 'w', **lidar_chm_crop_meta) as ff:
                ff.write(lidar_chm_crop[0], 1)
        except:
            None
    i += 1


texts: How to render PostGIS raster with Mapnik XML?

texts:
I'm trying to manage raster files from PostGIS, which to be rendered using Mapnik. I used raster2pgsql to load those GeoTIFF images to PostGIS, and I think the loading process turned out OK. The problem is to render from Mapnik XML using the PostGIS datasource: the rasters rendered out blank.. The raster files rendered perfectly with Mapnik XML using GDAL datasource, so I suspect my PostGIS raster configuration is at fault here..
Is there any documentation out there that resembles what I am trying to do?

Below is a snippet of my XML file to be used by Mapnik:

<Datasource>
            <Parameter name="type">postgis</Parameter>
            <Parameter name="user">postgres</Parameter>
            <Parameter name="dbname">gis</Parameter>
            <Parameter name="geometry_field">rast</Parameter>
            <Parameter name="estimate_extent">false</Parameter>
            <Parameter name="table">staging.hillshade</Parameter>
            <Parameter name="extent">101,3,103,4</Parameter>
</Datasource>

I appreciate any help, thanks!


texts: Replacing null field in generated ArcMap Report?

texts:
I am creating a report from ArcMap and I would like for the null value fields to be replaced with either a blank or a 0 , but only in the generated report not in the actual attribute field.

Is this possible and, if so, how would you do it?


texts: How to split circles in 12 sections in PostGIS?

texts:
I have a point layer in my PostgreSQL db and I would like to split a large number of circles into 12 sections each. Sample scenario is:

The circle corresponds to a 50 meter buffer of a point layer, and for each point, I need to split the buffer in 12 sections (pie-like wedges). The desired scenario is shown in following figure:

Can anyone suggest how to split circle into 12 sections?


texts:
When I want to upload a .zip-file of an ESRI-shapefile, I always get the encoding error 2006.
When I got this error the first time, I encoded the shapefile to UTF-8 in QGIS. Furthermore I uploaded the encoded shapes but the encoding error 2006 came again.

Does anyone know an answer?

I want to upload the zip file on the online mapping platform carto.com. For uploading georeferenced data, they need to be zipped.


texts:
I am currently working on a project in which spatial data are stored in a Microsoft SQL Server 2008 database. The tables which hold the spatial table contain 3 rows and are:

RoomId - BIGINT,
Description - VARCHAR(40),
Geometry - geometry.

I use QuantimGIS 1.8.0 for Windows and I use the tool for adding an MSSQL Spatial Layer. I succeed in connecting to the Server and I get a list of the tables which is as follows:

Schema: dbo
Table: tbl_lv1
Type: Polygon
Column: Geometry
SRID: 0
Primary key column: (empty)
Select at id: checked
Sql: (empty)

And when I click add, I get the following error:

dbname='gis' host='CAPDEV' user='qgisuser' password='*********' srid=0 type=POLYGON table="dbo"."tbl_lv1"(Geometry) sql= is an invalid layer and cannot be loaded

Why is this happening? Am I missing something?

Thank you

Solved: The solution was to add an auto-increment identifier to use as a primary key, because QuantumGIS needs a primary key attribute in order to import data.


texts: Editing of features with ArcGIS Pro SDK for .NET

texts:
I have programmed several ArcMap AddIns by using C# and ArcObjects the last 10 years. Now there is the idea to migrate all the AddIns to ArcGIS Pro, to support the 64 bit versions and be prepared for the end of life cycle for ArcGIS Desktop. I started with this great Getting Started with the ArcGIS Pro SDK for .NET by Rob Burke and had a quick look at the ArcGIS Pro SDK for .NET.
One of our most used issues (using ArcObjects) is the following. By activating an ESRI.ArcGIS.Desktop.AddIns.Tool the user can draw a rectangle in the map. All features of a defined feature class inside this rectangle are listed in a custom dialog. Double-clicking one feature opens another custom dialog which allows the user to edit all attributes. After pressing the OK button the feature is updated in the geodatabase by using the ESRI.ArcGIS.Geodatabase.ICursor.
Am I right, that we don’t have that deep access to ArcGIS Pro with the Pro SDK like we have had with ArcObjects to ArcMap?
This would be very disappointing.
Or am I overlooking something in the API how to read and write features?


texts: Displaying Wales, England & N.Ireland instead of UK in CartoDB?

texts:
I'm a journalist, and I would like to create a map to explain some things (about the countries that are playing the Euro 2016). I'm using CartoDB to create the map.

But for some reason, I can't find a way to display Wales, England and N.Ireland instead of UK. I have no idea how to fix this thing. Does anyone have an idea?


texts: I have CAD files from AutoCAD and from MicroStation. I need to know if it's 2D or 3D. I have FME at my disposal. How can I know if it's 2D or 3D in FME?


texts: ArcGIS Export Map Operation is returning imageData string, not the expected png Url

texts:
I am executing an ArcGIS server Export Map operation /service. The JSON response returns a property called imageData, and not the expected href property containing the Url to the generated png file. I tried a different ArcServer with the exact request and get back the expected result. The MapServer endpoint lists the following Supported Image Format Types: PNG32,PNG24,PNG.

What needs to be done to the ArcServer to return the expected PNG Url?

Current Response:

{
"imageData": "iVBORw0KGgoAAAANSUhEUgAABRQA...
"contentType": "image/png",
"width": 1300,
"height": 385,

Expected Response:

{
"href": "http://myserver/arcgis/rest/directories/arcgisoutput/Dir/MyMapServer/abcd.png",
"width": 1300,
"height": 385,
....


texts: Increasing marker size based on number of points at same location in QGIS

texts:
I have a dataset that contains many points (representing specimens) with GPS locations taken at the same point. When I map these points as a vector layer, only the uppermost point shows - so, at one given location, for example, there could be 50 data points yet it looks exactly the same as a location that may have just 1 data point.

I am wondering if there is a way to increase the marker size relatively to show this abundance of points at the same location using QGIS? I have tried to use the point displacement tool but I am not looking for the result that this gives with a centroid and a ring of points surrounding. I basically want to create a buffer that has a size in accordance with the number of points at that same lat long to show the relative abundance of specimens at that point.


texts: Writing Python script in Field Calculator in ArcGIS Desktop?

texts:
I'm a programmer, I have programmed in C#, Java, but not with Python. And I'm not a geographic worker.
How can I write a python script into Field Calculator to be able to check the other column value. If that other column value is "a certain string", the value in this column should be 0, otherwise it should be computed with an equation.
The value I want to look is a Text type, and I want to count in a Double type column.


texts: Obtaining Centroid from Polygon in SQL Server 2008 R2

texts:
I have several spatial layers created in QGIS that are held as spatial data in SQL Server 2008 R2. What I want to do is to find the centroid (OS Easting and Northing) of the polygon from the geometry. But I want to do this in a view that I have of the data. Each table I want to do this on, has a Geometry field. I've looked at examples of using ST.Centroid, but none of them seem to use a geometry field from a table.

Help on this will be invaluable.


texts: Extract By Mask Output is Empty Using ArcGIS Pro

texts:
I am trying to extract part of a large raster using ArcGIS Pro 2.2. Each time I try I get a message saying it was successful but the outputs do not show up. The attributes are all copied over correctly and the count values are the same as those in the original raster. I've used a raster and a vector as mask data and made sure both input and output are in the same projection. I've tried exporting to a .gdb as well as a .img and .tif all with the exact same results. The other characteristic of this issue is that the processing takes only a few seconds which doesn't seem correct.

In the screen shot below the Cop2_ss.img is the output from the extract by mask tool.

Does anyone have experience with this happening?

~Natasha


texts: How can I separate consecutive points from overlapping points in ArcGIS?

texts: I have a point shapefile with tracks of a machine. I'm looking to peel apart the overlapping point on a single track. The idea is that I'm looking at machine compaction on the soil and I have data for this at every point, but I'd like to see if with every pass of the machine, I can separate (i.e. peel back like an onion), the different machine passes to show them one at a time. I attached an image of what the points look like. In some rows, the machine only did one pass, but in others it went back and ran over it's previous tracks. Is there a tool that I can use, or a method to cluster all the points from the first pass into one shapefile, and every consecutive pass on top would be separated out? I have date/time information for each point as well.


texts: How to process GRIB2 Data for Visualization?

texts:
My goal is to extract wind information from data provided to me in the GRIB2 format as it has been done in the attached image.

I was curious if anybody had any experience with this format. I'm not sure if I am going in the right direction or what format I should be exporting it to. My end goal is to display this information through openlayers. So far, I have used gdal_translate to extract specific information ex. Temperature from band 25 (see command below).

gdal_translate -of Gtiff -b 25 weather.grib2 band25.tiff

I also have access to FME Desktop professional edition.


texts: Analysing many MODIS AOD tiles on the same map

texts:
I wish to analyse several days' worth of Aerosol Optical Depth retrievals from (MYD04_L2) data over the Atlantic Ocean. I have limited experience of using GDAL, but only analysing a single MODIS tile at a time. I'm still relatively new to using remote sensing data.

Could anybody point me in the direction of how to visualise several MODIS tiles simultaneously on a map?


texts: Unknown error concurring when developing QGIS 2.0.1 plugin in Windows 7

texts:
I am following the tutorial found here for developing a simple QGIS plugin within QGIS 2.0.1 installed on Windows 7.  However, as I have began to try developing a link between the GUI and each the mouse click of the map canvas, the main plugin code does not seem to recognize the fields in the GUI.

To illustrate with the example in the tutorial.  Within in my plugin.py code I have the function:

class vector_selectbypoint:
     def __init__(self, iface):
         # Save reference to the QGIS interface
         self.iface = iface
         # refernce to map canvas
         self.canvas = self.iface.mapCanvas()
         # out click tool will emit a QgsPoint on every click
         self.clickTool = QgsMapToolEmitPoint(self.canvas)
         # create our GUI dialog
         self.dlg = vector_selectbypointDialog()

    def handleMouseDown(self, point, button):
         self.dlg.clearTextBrowser()
         self.dlg.setTextBrowser( str(point.x()) + " , " +str(point.y()) )

that simply writes the lat and long into a TextBrowser with the GUI by calling the clearTextBrowser and setTextBrowser functions in a seperate pluginDialog.py script where the following code is ran:

class vector_selectbypointDialog(QtGui.QDialog, Ui_vector_selectbypoint):
    def __init__(self):
        QtGui.QDialog.__init__(self)
        self.ui = Ui_vector_selectbypoint()
        self.setupUi(self)

    def setTextBrowser(self, output):
        self.ui.txtFeedback.setText(output)

    def clearTextBrowser(self):
        self.ui.txtFeedback.clear()

which interacts with the ui_plugin_ui.py where txtFeedback is the TextBrowser in the GUI.  Below is a snippet of the ui_plugin_ui.py file:

class Ui_vector_selectbypoint(object):
    def setupUi(self, vector_selectbypoint):
        vector_selectbypoint.setObjectName(_fromUtf8("vector_selectbypoint"))
        vector_selectbypoint.resize(400, 300)
        self.buttonBox = QtGui.QDialogButtonBox(vector_selectbypoint)
        self.buttonBox.setGeometry(QtCore.QRect(30, 240, 341, 32))
        self.buttonBox.setOrientation(QtCore.Qt.Horizontal)
        self.buttonBox.setStandardButtons(QtGui.QDialogButtonBox.Cancel|QtGui.QDialogButtonBox.Ok)
        self.buttonBox.setObjectName(_fromUtf8("buttonBox"))
        self.chkActivate = QtGui.QCheckBox(vector_selectbypoint)
        self.chkActivate.setGeometry(QtCore.QRect(20, 250, 81, 20))
        self.chkActivate.setObjectName(_fromUtf8("chkActivate"))
        self.txtFeedback = QtGui.QTextBrowser(vector_selectbypoint)
        self.txtFeedback.setGeometry(QtCore.QRect(80, 20, 221, 151))
        self.txtFeedback.setObjectName(_fromUtf8("txtFeedback"))

The issue is that when I run the plugin in QGIS and click on the map canvas the following error is thrown:

AttributeError: 'ui_plugin_ui' object has no attribute 'txtFeedback'

I have verified that there are no syntax issues but the issue persists. If anyone has any insight I would eternally grateful!


texts: ArcGIS Feature Service - Restrict Editing to Specific Feature Classes

texts: In ArcGIS Server 10.1, is there a way to restrict editing rights on a Feature Service so that only specific feature classes can be edited?


texts: Loading multiple shapefiles in GeoDjango?

texts:
I have 3 three shapefiles and would like to upload all of them and project the result in GeoDjango admin.

a) state.shp (for the shape of a state)
b) boundary.shp (to identify a specific portion in the state.shp)
c) beat.shp (to further sub divide the region obtained from boundary.shp file)

What I want to do is the following:
 Uploading state.shp gives the following result.

Uploading boundary.shp gives the following result:

Please note new operations are performed on results obtained by previous operations. shapefile can be of different types (i.e poly, lines or points).

And similarly any other shapefile file which is uploaded should be projected. After following the official documentation I can project single shapefile in GeoDjango admin.

If I create separate models for each shapefile I will not get combined result.

How can I display result of multiple shapefile in GeoDjango?
EDIT: I want my end result to look something like this (load all shapefiles on openstreetmap in GeoDjango)

But however I can load single shapefile only (one shapefile at a time). Something like this
state.shp

boundary.shp

I would like to combine both of my results in a single result as shown in QGIS output.


texts:
I am looking for a tutorial of how to create a plugin for QGIS.

I want a tutorial that:

Explains about the latest version of QGIS
Working in windows environment.
adequate for absolute beginners and explains step by step.

A search on the web brought up a few tutorials but none of them answers all requirements.

p.s.
I have knowledge in python programming.

Can anyone help with finding such a tutorial?


texts:
I’m hoping to find a solution. I’m running Windows 10 on a Dell XPS, and I have a single install license of ArcMap 10.5. The Windows 10 fall 2017 Creators Update has killed ArcGIS Desktop. It can’t open after installing this update. However, ArcCatalog runs without any trouble.

I’d love any troubleshooting suggestions. I can run Desktop for a while after I roll back the Windows update, but since I can’t disable Windows updates it’s short lived, because ArcGIS Desktop starts crashing as soon as the update downloads.

Update: I checked in my license, uninstalled Arc, rebooted, reinstalled Arc, re-checked out my license, and Desktop still can't open.  ArcCatalog still runs with no trouble.


texts:
I'd like my QGIS map canvas scale to be accurate, e.g. so if the scale is 1:10000, then 100m on the ground should show in 1 cm of screen.
My Windows laptop screen has 157 ppi (dpi) based on its physical dimensions and pixel resolution. That, of course, varies on the laptop and/or external screen.
QGIS however behaves as if the ppi is 96. This is if the overall system scaling (under Windows' Settings / Display / Scale) is 100%. If it is 125% (which improves font legibility), then QGIS uses a ppi of 120. 150% (fonts too large on all Windows applications) yields a QGIS canvas ppi of 144. These figures are confirmed by running iface.mapCanvas().mapSettings().outputDpi() from the Python console.
How can I adjust the canvas dpi? I've tried iface.mapCanvas().mapSettings().setOutputDpi(157.0) but that has no visible effect on the canvas, and .outputDpi() continues to return 120 (or 96).
Either a Python console solution or a startup.py solution would be fine. Or, some magic incantation in QGIS settings. Of course, something that magically identified the actual physical ppi (the 157 in my case) would be even better, but I doubt Windows supports that.
I've also tried fiddling with High-DPI compatibility settings in the Windows program options for the QGIS app. As expected, those affect the scaling of QGIS' UI fonts and icons, but don't change the canvas dpi.
This is not about changing the DPI of a print layout, which is easily customizable; it's about the canvas on-screen.
Editing to add: I've followed some breadcrumbs to https://stackoverflow.com/questions/16561879/what-is-the-difference-between-logicaldpix-and-physicaldpix-in-qt, and note in PyQGIS, iface.mainWindow().logicalDpiX() returns 120 while iface.mainWindow().physicalDpiX() is correctly 158. So I'm trying to find how to force the QGIS map canvas to set its DPI from .physicalDpi, or a custom figure, rather than from .logicalDpi. Ideally, it should keep UI font and icons scaled according to .logicalDpi (as is usual behaviour in all applications), and shouldn't force me to adjust .logicalDpi via operating system settings, which has side effects (on all open applications!) and allows only DPIs of 96,120,144,168 or 192.


texts:
I'm working with a SqlGeometry SqlServer DB and am a novice to GIS.  I'm tasked with getting the centroid of a series of geometries from C#'s SqlGeometry object.  This is fine as it has a GetCentroid method, but it blows up when I try to get the centroid of lines and multilines.

Is averaging the points for lines and multilines a valid method of getting the centroid?
If not, why not?


texts: Trimble GeoXH Geoexplorer 2008 COM1 problem

texts:
I am using TerraSync and getting an error:

  GPS: The COM port is already configured for use by another device.

Does anyone know how to correct this so I can start using my bluetooth connection again for logging information?


texts: How To Turn Off Map Drawing When Using Python in ArcMap?

texts:
I have a Python script that creates layers in ArcMap.  The script works fine, but when a new layer is created a map is always  automatically drawn, which makes the script take longer to process.  I really don't need the layers to be drawn on the map.  Is there any way to set my script to not draw the layers after they have been created?

During a search I found a similar question:  How to disable and enable drawing in ArcMap using Python?

Attempting that solution did not work; "No module named SendKeys"

I'm using ArcGIS 10.2 for Desktop.


texts: Reclasifing raster with raster calculator

texts:
I am using QGIS 3.8 and I have a raster named "Ri" with one band with values ranging from 0.00299858 to 0.198849 and second raster named "IP_Q" with one band with values ranging from 0.00100708 to 14.5811. I am trying to create one raster out of these two rasters using Raster calculator.

The output should have three cathegories based on this table:

cathegory 1: Ri < 0.01
cathegory 2: 0.01 =< Ri < 0.1
cathegory 3: Ri >= 0.1 or IP_Q >= 2

The raster calculator expression I've used is as follows:

( "Ri" < 0.01 )  * 1 +  (  ( 0.01 <= "Ri" )  AND  ( "Ri" < 0.1 )  )  * 2 +  (  ( "Ri" >= 0.1 )  OR  ( "IP_Q" >= 2 )  )  * 3

However, the result is a raster with values ranging from 1 to 7. I am trying to get a raster with values ranging only from 1 to 3 so it corresponds to the table above.

It seems I am missing something.

Any ideas?


texts: Using Python to iterate through folders and join like shapefiles, based on name

texts:
I'm trying to automated a process of searching through folders and merging shapefiles of the same "type" (Roads, Water, Buildings)

On the surface, what I want to do seems pretty simple. I have several data sets in individual folders, each folder containing point, poly and line data for my study area.

The files in each folder are listed by map sheet, then description of the shape, then shape type. For example:

Each file folder contains a similar structure, only changing the map sheet name, but keeping the rest of the description. I want to join the files together based on the last 7 characters that make up their names. For example, each of the folders I am working with contain a building point file. (BL_POINT)

I've set up my script so it searches through all of my folders and lists the names of the files.

import arcpy, sys, os
from arcpy import env

arcpy.env.workspace = sys.path[0]
workspace = env.workspace

feature_class = []

walk = arcpy.da.Walk(workspace, datatype="FeatureClass", type="All")

for dirpath, dirnames, filenames in walk:
    for filename in filenames:
        feature_class += filename + ";"
        print ("Found Feature Class: " + filename)

Are there any suggestions how I can automate the process of picking out files with identical final 7 characters once I add all my files to the feature_class list?

Ideally, I would like to have it so that it is automated, and does not require me to input the final 7 characters for each type of feature (Roads, Buildings, Water, etc)

Or, if I am going about this backwards, any suggestion of a better process would be swell.

Again, all I want to do is:

Search through folders
Merge files that are similar, based on their names.


texts: LAScatalog and shapefile fail to line up in lidR

texts:
I am currently trying to read in a catalog of las files (statewide 10TB) and extract around 500 plots from a shapefile of central points. The shapefile reads in correctly, and I have checked it in ArcGIS to be sure. However, when I attempt to display the catalog of las tiles, they appear disjointed and fail to line up with any of the points from the shapefile (which should fall well within the spatial extent). It appears that the northing may be correct, but the easting could be off. I attempted to change the projection of the las files to match the shapefile, but this did not remedy the issue.

catalog:

class       : LAScatalog
extent      : 70000 , 3195000 , 940000 , 2380000 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
area        : 997518.4 kunits²
points      : 379.02 billion points
density     : 0.4 points/units²
num. files  : 41353

Shapefile:

Object of class SpatialPointsDataFrame

Coordinates:
                min       max
coords.x1 -87.97372 -84.84951
coords.x2  37.80275  41.72992
Is projected: FALSE
proj4string :
[+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0]
Number of points: 530
Data attributes:
      PLOT          LAT             LON
 1005   :  1   Min.   :37.80   Min.   :-87.97
 1024   :  1   1st Qu.:38.25   1st Qu.:-86.71
 1031   :  1   Median :38.70   Median :-86.34
 1035   :  1   Mean   :38.82   Mean   :-86.39
 1042   :  1   3rd Qu.:39.22   3rd Qu.:-86.21
 1046   :  1   Max.   :41.73   Max.   :-84.85
 (Other):524

Here is the code I used.

shp <- readOGR("D:/All LAS/Shapefile/plot4.shp")
ctg <- catalog(folder = "D:/All LAS")
plot(ctg)
EABplots <- lasclip(ctg, shp, radius = 7.3)


texts: Display Legends on esri map using javascript

texts: I have different module like distance calculator,radius search and more. Each pages I need to show Legends div so that user can identify places by images on map(or identify the layers).


texts: How to know if shapefile is georeferenced?

texts:
I have a shapefile which is located in province in indonesia. it contains many sub-province

if i removed sub regions inside a shapefile, and it leaves a small region, and my question is that small region's scale is as same as a previous scale into a map how to reset a new scale that i want?


texts: Color ramp resources other than Color Brewer?

texts:
I am looking for more resources on color theory in cartography. When to de-saturate a color? Or if you start with one color, how do you create a ramp from it? RGB or HSV? I can hack at it and come up with something, but I'm wondering if there's a methodology to it.

I say this because my work requires I constantly come up with new color ramps, and we are limited in certain colors that we cannot use because of said colors' associations (red is only used for hotspots, green is for positive topics only, etc). I constantly riff off of Color Brewer ramps, defaults, and make my own. I'm wondering if there are additional resources available to help me as a cartographer who primarily does thematic population maps, to come up with new ramps. I know there is a lot of information out there for designers, but I am looking for information more specific to cartography.


texts:
Is there a way to define multiple points to snap in the construction mode in advanced digitizing?
Let's suppose I want the value -26- from point 3 to 4 but also want a value -31- from point 1 to 4. How would I do that?


texts: Open a local mbtiles map in Mobile atlas creator "Mobac"

texts:
I'm trying to add local mbtile into mobac.
I have tried to create an XML file in the map sources folder with this information
<?xml version="1.0" encoding="UTF-8"?>
<localTileSQLite>
   <name>Custom SQLite atlas</name>
   <sourceFile>D:\myplace\newmap.mbtiles</sourceFile>
   <atlasType>MBTiles</atlasType>
   <backgroundColor>#000000</backgroundColor>
   <tileImageType>PNG</tileImageType> <!-- optional -->
</localTileSQLite>

but I can't find anything new when I try to add the map inside Mobac.


texts: Road graph tool crashes, but works with sample shapefile

texts:
When I try to calculate the distance between two points using the "Road graph" plugin in QGIS (1.8.0-2, on OSX 10.8.2 ML), QGIS just crashes. It does work with the sample shape file from

http://manual.linfiniti.com/en/vector_analysis/network_analysis.html

I'm sure it's a problem with my shapefile, but I can't see why it won't work. Even if I extract one single line as a new shapefile, it still crashes (program closes itself).

EDIT:

"QGIS quit unexpectedly
Click Reopen to open the application again. Click Report to see more detailed information and send a report to Apple.
Ignore / Report... / Reopen"

EDIT2:

Info from the .prj file:

  PROJCS["ETRS_1989_UTM_Zone_32N",GEOGCS["GCS_ETRS_1989",DATUM["D_ETRS_1989",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",500000.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",9.0],PARAMETER["Scale_Factor",0.9996],PARAMETER["Latitude_Of_Origin",0.0],UNIT["Meter",1.0]]

Filenames:
gsk3c_ruhr_gew_line.sbn
gsk3c_ruhr_gew_line.sbx
gsk3c_ruhr_gew_line.shp
gsk3c_ruhr_gew_line.shx

Plugin settings are the default (do work with sample shape file)


texts: Extracting non-duplicate features between layers in ArcGIS Desktop?

texts:
Just making an example: Layer B has 5 points. There's a layer(A) contains 10 points and 5 of which  duplicate(totoally same and overlap) with layer B.

How to output the other 5 non-duplicate points in layer A as a new feature?

Or just delete the 5 duplicate points from it?

I don't think CLIP tool can solve this.

Any suggestion ?


texts: Creating Python Scripts from within ArcGIS?

texts:
I used to use SAGA GIS and automated simple tasks by dragging and dropping (or clicking something and copy and paste - don't recall sorry) to create a simple script that I could then run.  I've used this for example to mosaik a bunch of files and save as a new file.  I can then just update the list of files I want to mosaik and run the script.

I am working with a team who uses ArcGIS and was hoping that there was a simple way to automate simple tasks.  I know that you can use Python with ArcGIS, and I'm looking for an easy way to help them get started along the Python path.

I don't have access to ArcGIS itself yet, so can't just play and see.  I've checked out quite a few sites on this, but they focus on more advanced tasks.

I am somewhat guessing at the version being used - but it's a new system so I presume the most recent.


texts: Custom Projection and Alignment of Large-Scale RGB Time Series Mosaics

texts:
We are running an experiment to observe time to germination and growth in potato plants in a greenhouse. To that end we have built a gantry system equipped with an RGB camera to capture the plants which are placed on a bench 1.8m wide by 6.5m long (distance from bench to camera approximately 2m). In Pix4D I mosaic the images, but because this experiment is so large-scale, the geolocation tags are discarded and the resulting mosaics have no projection, are therefore not aligned in space, and are not measureable. It is possible to add scale in Pix using "scale constraints", however, this resamples and therefore lowers the quality of the mosaic, and still doesn't solve the projection issue.

In ArcMap, my initial solution was to define a custom projection for the first mosaic, and georeference/resample/snap the rest so the pixels of each aligned as they should (pixel size should all have been the same from the beginning). But, my projection was simply a Transverse Mercator and I did not fiddle with the parameters. Although this aligned the mosaics, they are still not measureable and the metadata does not truly reflect pixel size.

I also tried co-registration in ENVI, but this does not precisely align pixels nor does it allow measurements.

Is there another, more precise solution out there so that I have a time series of perfectly aligned mosaics that I can measure on the mm-cm scale, and that exist in the space that they should?


texts: QGIS 2.6.1: "no enhancement" option not working on DEMs

texts:
When I import DEM raster layers and try to change the contrast to "no enhancement" I get a blank white image in the map window. Changing "black to white" to "white to black", the map remains a blank white image.

I want to import DEMs without stretching the histogram at all as some of them may be 8-bit images.

How can I get an un-stretched image in the map window?


texts: Identify points contained inside a linestring

texts:
This question is, I think, very simple but I can't resolve it by myself.

I have a table 'cable' like this :

CREATE TABLE cable (
  id uuid primary_key,
  the_geom geometry(LineString, 2154)
);

And another table let's call it partpoint

CREATE TABLE partpoint (
  id uuid primary key,
  the_geom geometry(Point, 2154)
);

I want to know if, for a specific partpoint, this partpoint is part of linestring in the cable table. The first request was :

SELECT id
FROM cable
WHERE ST_CONTAINS(the_geom,
                  (SELECT the_geom
                   FROM partpoint
                   WHERE id=<anUuid>))

But this does not work :-( Do I need to loop over all the table ? (with a for loop).


texts: pgrouting calculate length of shortest travel time route

texts:
I currently use the pgr_kdijkstraCost() function to calculate the travel times of a few origins and destinations. But I would also like to return the distance length of these trips. Is this possible? If so, how would you go about it?

If it is not possible for a k matrix, is it possible for a single origin and destination?


texts:
I've got a block-level census dataset that I am attempting to aggregate using a 2-mile grid. Basically, I want each grid cell polygon to have a value equal to the sum of the total population of each census block inside it, plus the fraction of the population of each census block polygon intersecting it that is inside the grid cell (so that if a block is 50% in one cell and 50% in the other, each cell is assigned 50% of its population).

Is this something I should be thinking about in terms of a raster conversion? At a bit of a loss here, this type of geoprocessing is outside of my normal QGIS use.

Thanks!


texts: Export QGIS Atlas frames showing ALL the features on a given date

texts:
Apologies if this has been covered elsewhere, but I'm trying to generate a bunch of frames (to animate later with ffmpeg) using QGIS Atlas + NASA FIRMS data. I'm using the date column from the FIRMS data, in YYYY-MM-DD format. QGIS version is 3.10.1

Using this technique from this other post to show cumulative points.

And here are my Atlas settings before I export them as pngs.

Unfortunately this exports a separate png for each features from of the data, rather than grouping all the points under a given date and exporting a single image for that date.

What could I try to export all the points for each date?


texts: Multiply a stack of floating point rasters( composite raster ) with a constant

texts: I am trying to multiply a floating point stack of rasters (10 layers) with a constant and round up to the next integer. I have the formula RoundUp(10*composite_raster) in raster calculator, but for some reason the output raster has only 1 band while i expect 10 bands? Can you please suggest any way to multiply the composite raster by a constant ?


texts:
I'm trying to get RemapRange to accept a variable as input for the remap statement. The reason I need to do this is the input is variable based on user input.
The reclassify statement in the code block below does not accept this input.

When I check DirRange it outputs [[157,166,1],[166,175,2],[175,184,3],[184,193,2],[193,202,1]]. This assumes Direction = 180

However when I check DirRemap I get [; [; 1; 5; 7; ,; 1; 6; 6; ,; 1; ]; ,; [; 1; 6; 6; ,; 1; 7; 5; ,; 2; ]; ,; [; 1; 7; 5; ,; 1; 8; 4; ,; 3; ]; ,; [; 1; 8; 4; ,; 1; 9; 3; ,; 2; ]; ,; [; 1; 9; 3; ,; 2; 0; 2; ,; 1; ]; ]

How do I get Remap to give the correct value or get Reclassify to accept the input?

import arcpy
from arcpy.sa import *

Direction = arcpy.GetParameterAsText(0)
Dirout = arcpy.GetParameterAsText(1)

DirMin = Direction - 23
DirMax = Direction + 22

if DirMin <= 0:
 DirMin = DirMin + 360
else:
 DirMin = DirMin

if DirMax >= 360:
 DirMax = DirMax - 360
else:
 DirMax = DirMax

Part1 = DirMin + 9
if Part1 >= 360:
 Part1 = Part1 - 360

Part2 = Part1 + 9
if Part2 >= 360:
 Part2 = Part2 - 360

Part3 = Part2 + 9
if Part3 >= 360:
 Part3 = Part3 - 360

Part4 = Part3 + 9
if Part4 >= 360:
 Part4 = Part4 - 360

Cat1 = "[" + str(DirMin) + "," + str(Part1) + ",1],"
Cat2 = "[" + str(Part1) + "," + str(Part2) + ",2],"
Cat3 = "[" + str(Part2) + "," + str(Part3) + ",3],"
Cat4 = "[" + str(Part3) + "," + str(Part4) + ",2],"
Cat5 = "[" + str(Part4) + "," + str(DirMax) + ",1]"
DirRange = "[" + Cat1 + Cat2 + Cat3 + Cat4 + Cat5 + "]"

DirRemap = RemapRange(DirRange)
DRaster = arcpy.sa.Reclassify(Dirout, "Value", DirRemap, "NODATA")


texts:
The problem is that the Chukotka Peninsula is located in the Western Hemisphere. So it is on the left side of Natural Earth Datasets. It needs to be moved to the right. I tried to use following srs(Coordinate system in which the map is rendered) - centered on 11°E :

+proj=mill +lat_0=0 +lon_0=11 +x_0=0 +y_0=0 +R_A +a=6371000 +b=6371000 +units=m +no_defs

And get this result (notice horizontal straight lines):

Obviously I did something wrong.

This is similar unanswered question:
http://lists.berlios.de/pipermail/mapnik-users/2009-January/001497.html

I have tried some actions from there as well.

How can I get the world to wrap around a different meridian?

Do I need to process the data through org2org or something?

More images with +lon_0=+-15 :

Images from linked question from mailing list :

What's the problem with this pictures? How do I make it right?


texts: What is the fastest/proper way to see if a FeatureCollection contains no Features

texts:
Take this code example, where mylabel is a ui.Label:

var mypt = ee.FeatureCollection(ee.Geometry.Point(Map.getCenter().coordinates()));

var states = ee.FeatureCollection('TIGER/2016/States');

var selectedstate = states.filterBounds(mypt);
var selstatename = selectedstate.first().get('NAME');
selstatename.evaluate(function(result) {
    if (result) {
       mylabel.setValue('State:  ' + result.toString());
    }
});

In the above simplified example, if the user has their map centered, say, over the ocean, then selectedstate will be a FeatureCollection with no Features in it.  In that case, I can skip the following selstatename.first().get('NAME'); and selstatename.evaluate

I'm pretty sure there's got to be an easy way to check that, hopefully combined in the .filterBounds step so I don't have to send requests to the server twice, if that makes sense?

THIS IS ANSWERED HERE: How to get several property values at once so there is only one call to .`evaluate`


texts: CartoDB local installation - how to create users other than development

texts:
I have a local install of CartoDB. It all works fine following the instructions at https://github.com/CartoDB/cartodb. I manage to create a user called development, and then log in successfully via my browser at http://development.localhost.lan:3000/login as development/development. I can insert data and plot beautiful visualisations.

Now I'd like other people to access this marvellous thing and have their own accounts and data stored there. So for that, I did this:

sh cartodb/script/create_dev_user newusername newusername newusername newusername@email.com

This went fine, and created a new empty DB I can see with PgAdmin.

Now if I go to http://development.localhost.lan:3000/login and try to log in as newusername/newusername, I'm redirected to http://newusername.localhost.lan:3000/login, which doesn't exist. I'm not sure why this is happening.

Is there some other configuration I'm missing for the server?

I'm confused about the instructions in the website naming development to both the user and the server.


texts:
I have a table with 4 column. There are thousands of data. The column names are

serial || name || X || Y

From 500 I made a mistake and put the X data into Y field and Y data into X field. So can I change the data making a range condition where serial no in greater than 500 or serial no equals to 500 to 1000? I think this is possible. I am new to postgres.

Thanks


texts: Fetching field values from uncommitted features on Qgis

texts:
I'm trying to fetch a feature's field value on Qgis, I use  QgsFeature.attributeMap().value(cont int & key).toDouble()
This works fine for committed features but it always returns " 0 " if the feature is not yet committed to the vector layer !!


texts: QGIS python Console Error - list index out of range

texts:
I have a text file which has the details of many .tif files which have to be converted to a single merged .vrt file in QGIS Python Console.

text = "W:/04 Tools/07 MapMaker/workingfiles/text.txt"
with open(text, 'r') as f:
    list_text = [line.strip() for line in f]
    list_no = len(list_text)

I created a list called list_text which contains all the .tif files.

vrtSublist = []
for sllist in range(0,list_no):
    output_vrt = "W:/04 Tools/07 MapMaker/workingfiles/tif" + "_tile" + str(sllist) + ".vrt"
    if os.path.isfile(list_text[sllist]):
        general.runalg("gdalogr:buildvirtualraster",list_text[sllist],0,False,False,output_vrt, progress=None)
        vrtSublist.append(str(output_vrt))
    else:
        print ("The file " + list_text[sllist] + " is Not Found. Please add this file.")

So, all the .vrt files are now in vrtSublist. Now, I want to merge 7 .vrt files each into a single .vrt file; so that I will group all the .vrt files in another list called as sublists (where sublists[0] = 7, sublists[1] = 7,.... ).

y_count = 7
sublists = [vrtSublist[x:x+y_count] for x in xrange(0, len(vrtSublist), y_count)]
if len(vrtSublist) % y_count == 0:
    sublist_no = len(vrtSublist)/y_count
else:
    sublist_no = (len(vrtSublist)/y_count) + 1
Latlist = []
for slsub in range(0,sublist_no):
    output_vrtsublist = s"W:/04 Tools/07 MapMaker/workingfiles/tif" + "_Lat" + str(slsub) + ".vrt"
    general.runalg("gdalogr:buildvirtualraster",sublists[slsub],0,False,False,output_vrtsublist, progress=None)
    Latlist.append(str(output_vrtsublist))

Now, I have merged each of the 7 .vrt files into single .vrt files; and now, I can see the list of output .vrt files in another list Latlist.

Latlist_no = len(Latlist)
output_vrtfinal = "W:/04 Tools/07 MapMaker/workingfiles/tif" + "_merged.vrt"
general.runalg("gdalogr:buildvirtualraster",Latlist,0,False,False,output_vrtfinal, progress=None)
output_vrtclipped = "W:/04 Tools/07 MapMaker/workingfiles/tif" + ".vrt"
extents_vrtclip = "%f,%f,%f,%f"% (lon, longmax, lat, latmax)
general.runalg("gdalogr:cliprasterbyextent",output_vrtfinal,"",extents_vrtclip,5,4,75,6,1,False,0,False,"",output_vrtclipped, progress=None)

Unfortunately, I am getting an error

Traceback (most recent call last):
IndexError: list index out of range

I am getting this error in the line :

sublists = [vrtSublist[x:x+y_count] for x in xrange(0, len(vrtSublist), y_count)]

How can I solve this issue?


texts: How can I change the name of a coordinate system in SAFE FME?

texts:
Every time I use SAFE FME to write data in the coordinate system EPSG:25832 FME names it NUTM32. I them have problems reading the data in ArcGIS Desktop because it doesn´t recognise the name NUTM32.

Is it possible to change the Name of NUTM32 to "ETRS 1989 UTM Zone 32N" in Safe FME?


texts: Setting DockableWindow Initial Size ArcObjects .NET SDK

texts:
In an ArcMap custom tool I used DockableWindow for user interaction. I want to set a size for the window while initial opening. I changed the height and width of InitialPlacement. After changing it looks like this.
  <DockableWindows>
    <DockableWindow id="First_DockableWindow" class="FirstComponent+AddinImpl" caption="First Component" image="Images\connected-component.png">
      <InitialPlacement height="600" width="500" state="unpinned" position="left" />
    </DockableWindow>
  </DockableWindows>

I experienced that it ignores the size whatever I set. Did I miss anything?


texts: Using list of numbers to Select By Attribute in ArcPy?

texts:
I am trying to use Python to select all records whose attribute (which is a number) in a certain field is within a list.  I cannot get the syntax correct.
my list is:

upstream = [9841978, 9842008, 9841966]

I have tried several expressions without success, including:

expression = '"COMID" in ' + str(tuple(upstream))

I am receiving the following error:

Runtime error  Traceback (most recent call last):   File "<string>", line     2, in <module>
File "c:\program files (x86)\arcgis\desktop10.4\arcpy\arcpy\management.py", line 7182, in SelectLayerByAttribute
raise e ExecuteError: ERROR 000358: Invalid expression Failed to execute (SelectLayerByAttribute).
arcpy.SelectLayerByAttribute_management(catchments, "NEW_SELECTION", expression)


texts: ArcPy Create Point Along Line A Specified Distance From A Specified Point

texts:
In ArcGIS Basic I am trying to take a road polyline(one segment for the entire length of road) and get a point placed on the map a specified distance from an intersection. I have a layer of intersection points and the layer of roads. I need to tell it to place a point on the polyline where the starting intersection point is. Then travel the specified distance down the line and place a point for a sign. Any suggestions?

I am currently trying to write a script to select a road and points and create a temporary layer of each. Then from that temporary road I CreateRoutes_lr. Then take the intersection point and LocateFeaturesAlongRoutes_lr. From locating the point I would like to take the "fmp" from distance as a variable. Is there any way to extract it?


texts: Which time format does ESA use for Sentinel-3 file names?

texts:
regarding the ESA Sentinel-3 file name convention I have got a question; which time zone (GMT, CEST, ...) does ESA use? I could not find anything yet, even it should be a critical point for many users...

There is no information given on the official website of ESA for the Sentinel-3 SLSTR-product as you can see here URL:https://sentinel.esa.int/web/sentinel/user-guides/sentinel-3-slstr/naming-convention


texts: Wide borders look weird for small polygons in Mapnik

texts:
I am experiencing this kind of visual defects for small polygons in Mapnik. These, in particular, are the small islands in the Pacific. As you see, instead of being uniform and nice, the borders have turned into some sort of dried bubble marks.

The .mss file used to generate this style for mapnik with kosmtik is very simple:

Map { background-color: #bbffff; }
#countries { polygon-fill: #ffffba; line-color: #8dd; line-join: round; line-cap: round; line-width: 14; }

I have played with the line-join and line-cap properties, but this is the best result I have got and it is not enough.

Can this problem be solved at all? Am I missing something? Must something be done differently?


texts: Cannot run SAGA commands due to PATH problems

texts:
I have tried to upgrade SAGA on MacOS 10.15.7 to use a version more recent than 7.3 in QGIS 3.16 since it has been unable to correctly export files using the Ordinary Kriging and Semivariogram Surface tools (while I am new to geostatistics, I understand the principles and half of the tool parameters, but I don't see how that would affect the ability to execute outputs). According to the SAGA-GIS wiki and QGIS's training manual for "Running an external algorithm", it's only a matter of adding the SAGA directory to my PATH variables to use it in QGIS, but I can't even do that. I've installed osgeo-saga (Saga 7.7.0) with Homebrew and have folders for Saga For Processing and Saga 7.8.2 in my home folder, but I've only tried to add osgeo-saga to my PATHs since it is the only folder that has the saga_cmd executable necessary to use SAGA (I probably have to compile Saga 7.8.2 to get an executable but my Unix skills end there for now). I've edited .bash_profile and .bashrc (which I had to touch/create) in my home folder by adding /usr/local/Cellar to the former and "export SAGA_MLB=/usr/local/Cellar" to the latter. echo $PATH outputs the following:
/opt/local/bin:/opt/local/sbin:/Library/Frameworks/Python.framework/Versions/3.8/bin:/usr/local/Cellar:/usr/local/bin:/usr/local/bin:/usr/bin:/bin:/usr/sbin:/sbin

Using bash as my default shell, this still doesn't allow me to run any commands for saga or any osgeo brew formulae, while all my other brew programs seem to work. What am I doing wrong? I can add one of the SAGA tool logs if you think that's more important to my original issue.


texts: Combining python (parallel) multiprocessing toolbox with arcpy

texts:
I'm trying to cut down on arcpy processing time by using the python multiprocessing toolbox. Using the example from http://wltrimbl.github.io/2014-06-10-spelman/intermediate/python/04-multiprocessing.html as a guideline.

Basically, I'm trying to select all tracks (point features) with a specific trackID (in field 'track'), give it a buffer and extract all road features within that buffer as preprocessing for another script. But the data set is so large it would take me about a week to run it.

I now have the following code, but I'm having a problem with the program executing normally, but not outputting the files in Tracks.gdb. Does anyone have advice?

import arcpy
from arcpy import env
import multiprocessing
def MyFunction(i):
    TrackFeatureOut = "C:\tryingOut\Tracks.gdb\Track_" + str(i)
    RoadFeatureOut = "C:\tryingOut\Tracks.gdb\Road_" + str(i)
    if not (RoadFeatureOut in featureclasses):
        BufferFeatureOut = "C:\tryingOut\Temporary.gdb\Buffer_" + str(i)
        Query = "track =" + str(i)

        arcpy.Select_analysis("C:\tryingOut\Input.gdb\Tracks", TrackFeatureOut, Query)

        arcpy.Buffer_analysis(TrackFeatureOut, BufferFeatureOut, "100 Meters")

        arcpy.Clip_analysis("C:\tryingOut\Input.gdb\Roads", BufferFeatureOut, RoadFeatureOut)

if __name__ == '__main__':
    global featureclasses
    env.workspace = "C:\tryingOut\Tracks.gdb"
    featureclasses = arcpy.ListFeatureClasses()

    # Getting all unique track numbers
    Track = set([])
    cursor = arcpy.da.SearchCursor("C:\tryingOut\Input.gdb\Tracks",["track"],sql_clause=(None, 'ORDER BY track' ))
    for row in cursor:
        Track.add(row[0])

    # Trying to apply multiprocessing
    pool = multiprocessing.Pool(multiprocessing.cpu_count())
    for i in Track:
        if Counter % 1000 == 0:
            print "Counter = ", Counter, "Current Track = ",i
        Counter += 1
        pool.apply_async(MyFunction, i)


texts:
I've seen other related questions concerning this issue, but I can't figure it out. QGIS doesn't show me Google Layers from the OpenLayers Plugin, and the following message comes out:

This page was unable to display a Google Maps element. Please contact the site administrator. If you are the administrator of this site, please check the JavaScript console or check the following page for troubleshooting: http://g.co/mapsJSApiErrors

My knowledge on QGIS is pretty basic.


texts: ArcGIS Pro .pdf export - Colors are wrong

texts:
I am exporting a .pdf map from ArcGIS Pro (2.5) and the coloring is coming out completely different when I open it in in Adobe Acrobat. (Pictures below)

It seems like this is only a problem on my machine, as others can open the .pdf and the colors render correctly.There is no rasterization, all exported layers are vectors. Any ideas on what could be causing this?


texts: Find points within X distance with a specific variable

texts:
I am trying to use QGIS to find points within 5 meters of each other that all have a value of <50% for a given parameter that has been loaded through a CSV file. I looked in to doing it with a spatial query, but I only have one vector layer.

The CSV file I have loaded in has about 10 parameters such as the percentage value, lat, long, name, etc. Is it possible to search through these points based on distance and the percentage identifier with only using one vector layer (that is, not with spatial queries)?

Edit: Here is an example with a screenshot.

The colours are coordinated with the percentage values I mentioned above. I'm looking to do a search for groups of points that are below 50% and within 5 meters of each other. The black circle wouldn't have any but the cluster of red with the one green to the right should be caught within the search.


texts: How can I create a multi layer exportable file (.xcf .psd) from multiple raster image in QGIS?

texts: I have many raster file of the same zone but in differents spectral bands. I would like to create an exportable file of a small portion that can handle different layers. So in practice I would like to implement a sort of clipper (Raster->extraction->clipper) with two main features: be able to have multiple file input (maybe with an option: "all the raster in the ROI") and create an output like an .xcf file or .psd file with many layers as the number of raster. How can I do that? I know a little bit of Python can I use the code write for the clipper function?


texts: Process finished with exit code 139 (interrupted by signal 11: SIGSEGV) - PyQGIS in PyCharm

texts:
When I try to execute this code:

 from qgis.core import *

 QgsApplication.setPrefixPath("/usr/lib/qgis", True)
 app = QgsApplication([], True)
 app.initQgis()

 layer = QgsVectorLayer('/home/city.shp',"p", "ogr")
 QgsProject.instance().addMapLayers([layer])

 shp_path = '/home/points.shp'
 csv_path = '/home/a.csv'
 tag_save_as_shp = '/home/sep.shp'

 layerInput = QgsVectorLayer(shp_path, "cells", "ogr")
 csvInput = QgsVectorLayer(csv_path, "Table_csv", "ogr")
 QgsProject.instance().addMapLayer(layerInput, False)
 QgsProject.instance().addMapLayer(csvInput, False)
 # joining csv with shp
 csvField = 'Square-ID'
 shpField = 'cellId'
 joinObject = QgsVectorLayerJoinInfo()
 joinObject.setJoinFieldName(csvField)
 joinObject.setTargetFieldName(shpField)
 joinObject.setJoinLayerId(csvInput.id())
 joinObject.setUsingMemoryCache(True)
 joinObject.setJoinLayer(csvInput)
 joinObject.setPrefix('')
 layerInput.addJoin(joinObject)
   QgsVectorFileWriter.writeAsVectorFormat(layerInput,tag_save_as_shp,"CP120",layerInput.crs(),"ESRI Shapefile")

I get an error in PyCharm:

Process finished with exit code 139 (interrupted by signal 11: SIGSEGV)

But in the python console in QGIS everything works fine. Also, this code worked before two weeks from PyCharm also. I read that is about some memory violation, but I don't understand it. :/

UPDATE: If I deleted the QgsVectorFileWriter line, everything is ok. So I have problem when I try to make output..


texts:
I'm trying to update a column to a value in Mapinfo/Mapbasic using

Add Column "MMT_PARCEL_MP" (Limited_StratumTF )From Limited_Stratum Set To "T" Where intersects

and just getting an 'Unable to update table' error message
The table is editable and I can update it manually with no problems.
The intersect statement works with no problems in a query.


texts: How to use the Oracle Spatial functions inside ms4w MAPFILE DATA string?

texts:
How to use Oracle Spatial functions in the DATA queries of MapServer(ms4w 3.0.3)? I am getting the following error while trying to highlight selected features which lie inside a particular boundary :

  msDrawMap(): Image handling error. Failed to draw layer named
  'thematic_layer'. msOracleSpatialLayerWhichShapes(): OracleSpatial
  error. Error parsing OracleSpatial DATA variable.

Following is the query I am using in the DATA string :

DATA "geom from (select S.XFM_ID as XFM_ID,S.LOCATION as LOCATION,S.OGC_GEOMETRY as geom from STRAND S,NODEBOUNDARY NB where SDO_RELATE(S.OGC_GEOMETRY,NB.OGC_GEOMETRY,'mask=inside')='TRUE' and NB.NODENAME='PUGA2003' and S.LOCATION in ('A','2')) as FOO USING UNIQUE XFM_ID SRID 2622 VERSION 10g"


texts: Joining table in ArcGIS Desktop but join field is greyed out?

texts:
I'm trying to join a shapefile with an Excel file I converted into a CSV. I had done this before with the same shapefile but with a different CSV but this time when I try to select the column it's removed from the join list. I have tried exporting the CSV as DBF then turning it back into a CSV as well which did not work. I have included two screenshots below. I want to use the column PCODE3 because it shares the same data as MUN_P_6. Screenshots:

Full list of fields:

Limited list of fields:


texts: Does SDO_WITHIN_DISTANCE also return points inside a polygon?

texts:
Using Oracle Spatial Database, just wondering if the SDO_WITHIN_DISTANCE returns points inside a polygon, not just within a certain number of metres from the boundary.

I've already checked the documentation, it trivially states:

  Uses the spatial index to identify the set of spatial objects that are
  within some specified distance of a given object, such as an area of
  interest or point of interest.

I would have thought that it should as being inside a polygon would count as being a distance of 0m from a polygon...


texts:
I am trying to read a shapefile and plot it using matplotlib. Here is the code:

import matplotlib.pyplot as plt
import shapefile

shpFilePath = "D:\test.shp"
listx=[]
listy=[]
test = shapefile.Reader(shpFilePath)
for sr in test.shapeRecords():
    for xNew,yNew in sr.shape.points:
        listx.append(xNew)
        listy.append(yNew)
plt.plot(listx,listy)
plt.show()

However, i get lines connecting my polygons. How can I draw the polygons such that they are the way in the shapefile. Here are screenshots of the plot and the shapefile when it is opened with ArcGIS.


texts: Carto legend colours not showing in embed?

texts: Carto legends are showing no colours in the bullets after manually adjusting colours in the HTML window -- but only in IE. No issue in other browsers. Has anyone else seen this?


texts:
I am on the ModelBuilder of ArcGis10.3.
I have a While iteration and I want to add a counter which indicates the number of the run at each run. This value can be also useful for in line variable substitution to save/load some data using the number of the run.

What I tried to do is to link the Continue with Calculate Value (2) which is my counter whose value is stored in the Main Value.
Calculate Value (2) contains this:

  Expression:
  x(%Continue%)

  Codeblock: def x(n):
  if n==True:   x=x+1

The code should be syntactically ok, because it allows me to run it, but it doesn't recognize "True" as a boolean value and it think it is a variable.

I was using a For iteration before and the counter was already part of the tool. Now with the While cycle I am not even sure if what I am trying to do is the easiest way to create a counter.


texts: OpenLayers 4: Load baselayer outside projection extent

texts:
I have a map that uses EPSG:3912. Then I have base layer from CartoDB. When I zoom out I see:

My question is, if it's possible to have CartoDB layer loaded outside map extent (so that it fills white border) and with keeping projection as is?


texts:
I have installed r.denoise, put mdenoise into QGIS2.18/bin and installed "pip install pyproj".
However when I try to run r.denoise on STRM DEM, a get this error:

ERROR: Unable to open file <tmp_xyz_5660.xyz> for writing
Removing temporary files...
Traceback (most recent call last):
  File "C:\Users\admin\AppData\Roaming\GRASS7\addons/scripts
/r.denoise.py", line 231, in <module>
    main()
  File "C:\Users\admin\AppData\Roaming\GRASS7\addons/scripts
/r.denoise.py", line 166, in main
    grass.run_command('r.out.xyz', input=in_raster,
output=tmp_xyz, separator='space', overwrite=True)
  File "C:\PROGRA~1\QGIS2~1.18\apps\grass\grass-7.2.0\etc\py
thon\grass\script\core.py", line 410, in run_command
    return handle_errors(returncode, returncode, args,
kwargs)
  File "C:\PROGRA~1\QGIS2~1.18\apps\grass\grass-7.2.0\etc\py
thon\grass\script\core.py", line 329, in handle_errors
    returncode=returncode)
grass.exceptions.CalledModuleError: Module run None
['r.out.xyz', '--o', 'input=DEM_sliv_STRM_30x30@PERMANENT',
'separator=space', 'output=tmp_xyz_8292.xyz'] ended with
error
Process ended with non-zero return code 1. See errors in the
(error) output.

For some time I cant figure out whats the problem and cant find anywhere similar problem.
I am using GRASS GIS 7.2.0 with QGIS 2.18.7 on Windows 10.


texts:
I'm very inexperienced in using R for GIS. I looked at other posts but couldn't find an answer.

I'm trying to create a map that shows a municipality, roads through it and shop locations. I've installed and added packages sf, raster, rgdal and ggplot2.

I have two questions.
First, I've found two ways to load shapefiles to R:

munboundary<-shapefile("Filename.shp") and
munboundary<-readOGR("Filedirectory", "Filename")

Are there differences between both? Both seem to work.
Once I've loaded the three shapefiles I'm trying to plot all three in one map.
I've looked up that I can us add=TRUE to overlay shapefiles but it doesn't work in my case.

plot(munboundary, col="#FFCC99")
plot(roads, add=TRUE)
plot(supermarkets, add=TRUE)
title("Supermarket access in municipality X")

The first one (munboundary) is plotted correctly but the others don't appear. There is no error message, in the console it seems to have worked but the shapefile isn't visible on the map.

I tried to find an alternative online that works and ran into this method:

ggplot() +
  geom_sf(data = munboundary, fill="#CCFFCC") +
  geom_sf(data = roads) +
  geom_sf(data = supermarkets) +
  ggtitle("Supermarket access in municipality X")

but this one doesn't work either. Again looking at the console there seem to be no problems but now nothing shows up in plots.

I've just been stuck for multiple hours and feel like the ultimate noob right now.

Summary information:

summary(munboundary)
Object of class SpatialPolygonsDataFrame
Coordinates: min max x 146.17037 146.51399 y -39.15919 -38.77096
Is projected: FALSE
proj4string : [+proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs]

summary(roads)
Object of class SpatialLinesDataFrame
Coordinates: min max x 440884.8 454233.8 y 5668710.8 5681023.7
Is projected: TRUE
proj4string : [+proj=utm +zone=55 +south +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs]


texts: How to add layerSwitcher to a specific div

texts:
I'm trying to add a control in a specific panel (<div>) which will groups others controls. I tried many things but nothing with succes.
It seems I didn't catch the option target.

Maybe somebody could help me with this ?


texts: Negative values in Landsat 8 surface reflectance higher level product

texts: I have downloaded Landsat 8 higher level surface reflectance product from USGS to genetare chlorophyll map in lake. But the individual bands contain some negative surface reflectance values. So, when I am applying the regression model, the logarithmic scale turns those values to NA. Does someone know what is causing those negative values in the surface reflectance bands? Is it possible to neglect those values ? I am processing my data in R.


texts:
I'm a little stumped with this issue that I'm having.. I have a large zip code map (the whole U.S.) with a Geocoder. When I do a search for a particular zip code, and press enter or select it, the map is programmatically clicked on the longitude/latitude of the resulting zip code area. This works initially. The map flies down to the particular zip code, and highlights (just like it's supposed to do upon click). If the user enters a different zip code that is within the viewport, it will activate again (like it's supposed to).. HOWEVER, if the user selects a zip code OUTSIDE the view port, the view port "flies" over to the new place but does not fire the click event, and so therefore, nothing is selected.

I have a JS Fiddle documenting the behavior.

https://jsfiddle.net/wbeelergofish/t5heL2vz/

EDIT: I'm not sure that the 'fire' event is fired or not, but it does not work like it's supposed when the fly action is finished taking place outside the current view port.


texts:
I have a vector layer which is a grid. In the table attribute I calculated, using Zonal Statistics, the mean value of two raster layers, for each cell of the grid.

Now, I would like to use GroupStats to have some means value of a column compared to the second one.

For example :

R1   R2
10   15
10   25
11   14

I want to obtain :

R1    R2
10    20
11    14

So that is easy :

But as I have to do this a lot of time, I would like to write a Python script for this or to use Graphical Modeler.

But it seems that there is no module for Python, and no possible use in Graphical Modeler as it's a plugin.

Does someone would know how to manage with this ?

I use Qgis 2.18.9


texts: Does a Settlement Hierarchy plugin exists for QGIS?

texts: Is there a Settlement Hierarchy plugin for QGIS ? When I say settlement hierarchy I mean like a plugin that would designate where a village area begin and end ? Same with towns. So you could separate towns from villages from hamlets and isolated dwellings. Here is more information about settlement hierarchy : https://en.wikipedia.org/wiki/Settlement_hierarchy


texts:
How do I reclassify many rasters in ModelBuilder when each raster has its own unique values?

For example, when you create a model in ModelBuilder and you add the Reclassify tool you may want to reclass several rasters with the Iterate Raster tool as an input. However, when you open the Reclass toolbox you choose a value for reclassification and then name it as an output raster. The problem is how can you, when iterating the raster, reclass with several reclassification values for each raster? This is what I am trying to do.

@Michael Miles-Stimson Here is the python I copied and paste here.

# -*- coding: utf-8 -*-
# ---------------------------------------------------------------------------
# Reclassify2006.py
# Created on: 2014-06-23 16:22:40.00000
#   (generated by ArcGIS/ModelBuilder)
# Description:
# ---------------------------------------------------------------------------

# Import arcpy module
import arcpy

# Check out any necessary licenses
arcpy.CheckOutExtension("spatial")

# Load required toolboxes
arcpy.ImportToolbox("Model Functions")

# Local variables:
Final_NDVI = "E:\\FY2013_projects\\PMC_project\\USGS EarthExplorer\\2006\\Final_NDVI"
v_outFGDB___Value__new = "E:\\FY2013_projects\\PMC_project\\USGS EarthExplorer\\2006\\Vectors\\NDVI_reclassify.gdb"
v_Name_ = "E:\\FY2013_projects\\PMC_project\\USGS EarthExplorer\\2006\\reclassify\\%Name%"
v_Value__new = "E:\\FY2013_projects\\PMC_project\\USGS EarthExplorer\\2006\\Vectors\\NDVI_reclassify.gdb\\%Value%.shp"
Name = "ast_l1b_00302032006175500_20101030104009_12125_ndvi"
v_Value__new__2_ = "E:\\FY2013_projects\\PMC_project\\USGS EarthExplorer\\2006\\Vectors\\NDVI_reclassify.gdb\\%Value%.shp"
ast_l1b_00302032006175500_20101030104009_12125_ndvi_img = "E:\\FY2013_projects\\PMC_project\\USGS EarthExplorer\\2006\\Final_NDVI\\ast_l1b_00302032006175500_20101030104009_12125_ndvi.img"

# Process: Iterate Rasters
arcpy.IterateRasters_mb(Final_NDVI, "a*", "", "NOT_RECURSIVE")

# Process: Reclassify
arcpy.gp.Reclassify_sa(ast_l1b_00302032006175500_20101030104009_12125_ndvi_img, "Value", "-0.17647059261798859 0.022727273404598236 1;0.022727273404598236 0.10313901305198669 2;0.10313901305198669 0.1428571492433548 3;0.1428571492433548 0.16883116960525513 4;0.16883116960525513 0.19117647409439087 5;0.19117647409439087 0.21538461744785309 6;0.21538461744785309 0.24324324727058411 7;0.24324324727058411 0.29870128631591797 8;0.29870128631591797 0.56451612710952759 9", v_Name_, "DATA")

# Process: Parse Path
arcpy.ParsePath_mb(Name, "NAME")

# Process: Raster to Polygon
tempEnvironment0 = arcpy.env.parallelProcessingFactor
arcpy.env.parallelProcessingFactor = "ast_l1b_00302032006175500_20101030104009_12125_ndvi"
arcpy.RasterToPolygon_conversion(v_Name_, v_Value__new, "SIMPLIFY", "VALUE")
arcpy.env.parallelProcessingFactor = tempEnvironment0

# Process: Add Field
arcpy.AddField_management(v_Value__new, "Acres", "DOUBLE", "", "", "", "", "NULLABLE", "NON_REQUIRED", "")

# Process: Calculate Field
arcpy.CalculateField_management(v_Value__new__2_, "Acres", "[Shape_Area]", "VB", "")


texts:
I have problem with atlas generation in QGIS 2.2.
I'd like to export pdf map with 2 maps window and use coverage layer function.
Atlas doesn't iterate across the features in coverage layer.

First problem is output preview which doesn't work at all.

Anyway when i export map as pdf i have same map for all features.
It seems that atlas can count my features in coverage layer but it can't iterate across it.

I use shp file with points and lines + shp with polygons for coverage layer.

When I try it with same file in QGIS 2.0 everything works fine.
I have installed both version QGIS 2.0 and QGIS 2.2.

Have you ever seen a similar problem as I have or do you have any idea for a solution?


texts: Why does my query not run since SRID are the same?

texts:
I am running this query:

SELECT source,ST_AsText(ST_ClosestPoint(roads.geom,ST_GeomFromText('POINT(22.256  39.569)')))
FROM roads
ORDER BY ST_Distance(roads.geom,ST_GeomFromText('POINT(22.256 39.569)')'))
LIMIT 1;

I have followed this answer and reprojected my table     geom_ways in my database to 4326 which is the SRID of geometry points but the query will not run because there are not matching geometries?


texts: Fiona and Pyproj incorrect transformation (feet to meter)

texts:
I am trying to reproject a shapefile from EPSG:2229 (NAD83 / California zone 5 (ftUS)) to EPSG:2163 (US National Atlas Equal Area (meters)).

Below is the code that does this transformation for a point in the shapefile:

import fiona
from fiona.crs import from_epsg
from pyproj import Proj, transform

# Open the shapefile we want to reproject and initiate original projection
shape = fiona.open(inputfile)
oriproj = Proj(shape.crs)

# Initiate output projection
outcrs = fiona.crs.from_epsg(2163)
outproj = Proj(outcrs)

# Create a new file where the transformed point shapefile will be written
with fiona.open(outputfile,'w',shape.driver,shape.schema.copy(),crs=outcrs) as output:
    for feat in shape:
        point = feat['geometry']['coordinates']
        x,y = point
        new_x,new_y = transform(oriproj,outproj,x,y)
        feat['geometry']['coordinates'] = (new_x,new_y)
        output.write(feat)

This process works well with lat/long to meter transformations and transformations that involve meters to meters (EPSG:2163 to EPSG:32620) but for some reason it is breaking when trying to transform from feet to meters.

After this transformation, the points maintain their original distribution, but are located off of the coast of Delaware (when they should be in Los Angeles County). I have a feeling this might be an origin issue but I'm not sure.

Has anyone had this issue before? And if so, how did you fix it?


texts: As I understand, DTED0 should be available for free. Yet, I am not able to find any download link. Everything that Google search brings up (eg. http://geoengine.nima.mil) is dead and not working.


texts: Is WFS:native mandatory to get OGC compliance?

texts: I am a newbie to this domain. I am using WFS 1.1.0 and found out that is has a wfs: native operation. I wanted to know whether I can get OGC compliance without implementing WFS: native or is it mandatory to implement it to get OGC compliance.


texts: apply gdal to calculate raster that is already in qgis raster calculator

texts:
I am using the following code in qgis to calculate some raster.It is not giving the correct result as it is in arcpy.I want to apply same in gdal.Can anyone help me how to use the gdal api.Please give me code

from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry
layer= iface.activeLayer()
entries = []
# Define band1
ent1= QgsRasterCalculatorEntry()
ent1.ref = 'boh@1'
ent1.raster = layer
ent1.bandNumber = 1
entries.append( boh1 )

# Process calculation with input extent and resolution
calc = QgsRasterCalculator( 'boh@1 * 0.03705882+ 3.2',
                            'D:\test\testtif',
                            'GTiff',
                            layer.extent(),
                            layer.width(),
                            layer.height(),
                            entries )

calc.processCalculation()


texts: Horizontally slicing point cloud into a cylinder

texts:
I have some point clouds from mobile terrestrial laser scanning of circular sampling plots in woodland environments. I am interested in studying the structural complexity of the woodland understorey. I have used the clip_circle() function to crop the point clouds to the 15 metre radius of the sampling plot.
Now, I would like to horizontally slice the point cloud so that I can obtain the first two metres of points from the woodland floor. Therefore, I would like to extract a 15-metre radius cylinder of two metres depth, positioned on the floor.
I managed to do this in CloudCompare and reimport the cloud back into R like this:

However I would like to do the whole pipeline in R to save time when I have more clouds to process.


texts: Geocoding with county level .Locator Style

texts:
I'm using ArcGIS 10.1 locators to geocode some addresses on a state by state level,
using the dual ranges locator style.

Since I have the administration levels city/county/state on both datasets (reference data and addresses to geocode), I would like to add them when creating the locator so that it geolocates on groups: city/county/state. The issue I'm finding is that on locator creation of Dual Ranges style, I can't find a way to add the county level, I only see the city left and right, and state left and right along with other reference attributes.


texts: Using R, why am I unable to select a column from an sp Spatial Lines Dataframe?

texts:
With utm_RDS being the spatial lines in question, derived from the following shapefile:

library(rgdal)
library(sp)
    RDSYS <- readOGR(dsn = "C:/************")

        utm_RDS <- spTransform(RDSYS, CRS("+proj=utm +zone=14 +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0"))

We can then see

    > class(utm_RDS)
[1] "SpatialLinesDataFrame"
attr(,"package")
[1] "sp"

> head(utm_RDS)
An object of class "SpatialLinesDataFrame"
Slot "data":
   ID   NAME TYPE
0 ******  <NA> <NA>
1 ******  ***   ST
2 ******  <NA> <NA>
3 ******  <NA> <NA>
4 ******  <NA> <NA>

Among some omitted columns.

I am now interested in being able to select lines to plot and analyze given their type, but always seem to be met with error. For example,

utm_RDS@lines[utm_RDS$TYPE == "HWY"]
length(utm_RDS@lines[utm_RDS$TYPE == "HWY"])

[1] 46925

Which gives me the number of distinct lines in the shapefile.

But,

plot((utm_RDS@lines[utm_RDS$TYPE == "HWY"])) will not work, returns:

Error in xy.coords(x, y, xlabel, ylabel, log) :
  'x' is a list, but does not have components 'x' and 'y'

and is of class "list".

Further,

> plot((utm_RDS[utm_RDS$TYPE == "HWY"]))
Error in `[.data.frame`(x@data, i, j, ..., drop = FALSE) :
  undefined columns selected

Does anyone have a workaround for this?


texts: Getting field type information from FieldInfo object in ArcPy?

texts:
I'm looking at the FieldInfo object documentation for ArcGIS 10.2 and am not seeing access for field type (e.g. short, string, double).

Is it possible to get field type information from the FieldInfo object in ArcPy?

I know I can access the information from calling arcpy.ListFields(), but I want it through feature layer FieldInfo.


texts: Getting wrong lat/lang after transforming from "EPSG:900913" projection to "EPSG:4326"

texts:
In OpenLayer3 Map I'm getting lat/lang from nominatim and the map is in "EPSG:900913" projection and then transforming it to "EPSG:4326" projection .

if lat = 12.9791198 and lang = 77.5912997 then i am transforming like this,

var point = new ol.geom.Point([lan,lat]);
point.transform("EPSG:4326","EPSG:900913");

After transformation lat/lang should be [8643375.747091185,1460888.081802196] but i am getting [8637423.97259228,-13864431.524951909].


texts: Google Maps MarkerClusterer - Remove all clusters

texts:
So, I have a PostGreSQL DB with 317 datapoints - all spatially enabled through PostGIS. I have queried this using PHP as GeoJSON and overlaid it on Google Maps (using a converting routine). I have also used the MarkerClusterer extension available in google Maps api to cluster all the markers.

Now, I want to remove all the clusters from the map, to display other data. How do I do this?

My markers are stored in an array called markers = []; and the marker clusterer is called like this markerCluster = new MarkerClusterer(map, markers);. I'm not able to find the right function to do this. MarkerClusterer link - http://google-maps-utility-library-v3.googlecode.com/svn/trunk/markerclusterer/docs/examples.html


texts: How/Where to get reverse geocodes up to 5 decimal points precision for offline use?

texts:
I would like to have a simple text file(s) with following reverse geo-codes information for all the countries (i.e. world):

  Coordinates -- Address

e.g.

  12.96779, 77.65306 -- #7, Jeevan Bima Nagar Main Rd, 10th Main, DOS Colony, New Thippasandra, Bengaluru, Karnataka 560075, India

First I tried to call Google maps API, but it has its limitation of number of API calls per day. In case of OSM, there is no such limitation but at .00001 precision, it will take several years to get just 1 square degree!

We downloaded planet.osm (550 GB after extration) and planet.osm.pbf for offline database/server creation. However, we are unable to find a way out of it to re-create the server locally and make the API calls which will be fast enough.

Is it possible to get/purchase such reverse geo-codes of near Google-Maps like accuracy for our offline usage?

Otherwise, is there a full proof link which teaches, how to create a OSM like local server out of .osm or .osm.pbf? (Have tried many of them so far, but no success)

[Note: Previously I had asked slightly related question]


texts:
I was using lines to polygons tool and I encountered the following error:
"One or more line ignored due to geometry not having a minimum of three vertices."

I'm really confused because it worked perfectly in December. Did anyone encounter that error too and has a solution for this?

@EDIT.
It actually does work in QGIS 2.18. I assume it's because the contours of the smallest "chunks" I created in QGIS 3.4 are actually triangles (see above), whereas when I carried out the same operations in QGIS 2.18 I got squares, so QGIS was able to create polygons from them.


texts: Changing selection color in exporting ArcMap layout with ArcPy?

texts:
I use ArcMap and am looking for a way to change the selection color that would be respected when exporting map layouts programmatically using arcpy.mapping.ExportToPNG(). Currently, it is always the light green color (the default selection color in ArcMap).

Changing a selection color in ArcMap (Selection menu > Selection Options) does change the selection color but it is only respected when exporting the layout manually from the File menu or when exporting  the layout in the Python window (it respects the selection color set in ArcMap since I am inside a live ArcMap session). Exporting a layout using arcpy though always uses the default light green color for selection.

I have found no other settings in the ArcMap menus that would help and AdvancedArcMapSettings utility does not have any either.

Is there any way to make arcpy export layouts using any other selection color than the default light green one?


texts: Finding distance to closest road (polygon to line) in ArcGIS for Desktop?

texts:
For my project I need to find the distance to the nearest road from forest polygons. I'm trying to do a weighted overlay analysis (eventually) for a habitat selection project so I think I need to end up with a raster that shows the distance to the nearest road for each forest cell (among other factors I'm including).

Can anyone give me any advice on how to do this? At the moment my layers are in vector format but I can convert to raster.


texts:
This code adds a QgsPoint to a layer. But how to increase the size of feature ?

layer =  QgsVectorLayer('Point', 'points' , "memory")
pr = layer.dataProvider()
pt = QgsFeature()
point1 = QgsPoint(-4080198.52,1492348.027)
pt.setGeometry(QgsGeometry.fromPoint(point1))
pr.addFeatures([pt])
layer.updateExtents()
QgsMapLayerRegistry.instance().addMapLayers([layer])


texts:
I'm developing a Python Addin for ArcMap.  This will consist of a button that will run another python script:

import arcpy
from arcpy import env
import pythonaddins
import os

class ButtonClass1(object):
    """Implementation for Schaltflaeche_UID_addin.button (Button)"""
    def __init__(self):
        self.enabled = True
        self.checked = False
    def onClick(self):
        os.system(r'C:\Users\admin\Desktop\ToeB 10.2\Addins\Schaltflaeche_UID\Install\Insert_UID.py')

I've tried F_Kelner's answer posted in ArcGIS python add in to call Python script for importing mdb data, but I'm getting the following error, even though I have imported the os module:

      Traceback (most recent call last):
        File "C:\Users\admin\AppData\Local\ESRI\Desktop10.2\AssemblyCache{D9E1A7BE-8FE5-4EE5-9C9E-93C2B867EA3E}\Schaltflaeche_UID_addin.py", line 17, in onClick
          workspace = os.path.dirname(workspace)
      NameError: global name 'os' is not defined

Anyone have any ideas?


texts:
I've published some layers in GeoServer which I want to consume in Mapbox.I've tried several examples like:

Add simple image layers to MapBox / Leaflet?
Mapbox Vector Tiles from Geoserver 2.11 in Mapbox GL JS
Mapbox Docs | Add a vector tile source

But I'm not been able to get any good result, I've got this as the best result till now

 using the code

map.addLayer({
    'id': 'wms-test-layer',
    'type': 'raster',
    'source': {
        'type': 'raster',
        'tiles': [
            'http://localhost:8080/geoserver/gwc/demo/topp:states?gridSet=EPSG:4326&format=image/png'

        ],
        'tileSize': 256
    },
    'paint': {}
}, 'aeroway-line');

What can I do?


texts: Getting an error when I try to import a loaded vector file to the QGIS GRASS Plugin

texts:
I am getting an error code when I try to import a contour line vector file GRASS via QGIS.  I have been trying to create a DEM from the contour lines following this tutorial:
http://linfiniti.com/2010/12/3d-visualisation-and-dem-creation-in-qgis-with-the-grass-plugin/

When I try to import the loaded vector (v.in.ogr.qgis) I get the following error message:

  v.in.ogr "dsn=C:\GIS\Elevation\Florida Counties\Bradford\Bradford.shp" output=BradfordGrass2 snap=-1 min_area=0.0001 -o

  Datum not recognised by GRASS and no parameters found

  Over-riding projection check

  Layer: Bradford

  DBMI-DBF driver error:

  Unable to create table: 'create table BradfordGrass2 (cat integer, FNODE_ double precision, TNODE_ double precision, LPOLY_ double precision, RPOLY_ double precision, LENGTH double precision, TOPO5_ double precision, TOPO5_ID double precision, CONTZ integer, CARRYZ integer, ELEV integer, WMD varchar ( 6 ), ADD integer)'
  Finished with error

I am using five-foot contour lines downloaded from downloaded the data from the FGDL.

http://www.fgdl.org/metadataexplorer/explorer.jsp

Any suggestions as to how to fix this would be greatly appreciated.


texts: Create tiles given images and locations; convert markers to tiles

texts:
What I need

I need to convert a list of markers to a tile layer

I've been googleing without success how to do this. Anyone has faced this problem and has a nice easy solution?

Why I need it

I've developed a HTML application that displays a map. This map has many fixed markers (they do not change in time). Due to that:

I want to speed up the application
I want to hide the marker position from malicious users (so they can't use the information I've gathered for they're own porpose)

I want to convert the markers to a tile layer.

Thanks!


texts:
Is there a way to find out all the octagonal shaped elements within a dataset?

I'm trying to query OSM to give me all the octagonal shaped buildings in a certain area.


texts: Can a Fusion Tables Heatmap be used in Google Earth?

texts:
From the Fusion Tables web interface there's an option to view a heatmap, which draws on top of Google Maps.

From the Fusion Tables map view there's an option to Export To KML, which allows the Fusion Table points to be displayed in Google Earth.

Is there a way to add the heatmap to Google Earth?

Thanks


texts: QGIS - How to update Python Init Function

texts:
In layer properties I can set a "Python Init Function" and it works perfectly.
But I can't find a way to update it when developing (without closing QGIS).


texts:
I want to use PyQGIS function crs.isValid() in normal Python file outside of QGIS like this:
wkt = 'GEOGCS["WGS84", DATUM["WGS84", SPHEROID["WGS84", 6378137.0, 298.257223563]],' \
      'PRIMEM["Greenwich", 0.0], UNIT["degree",0.017453292519943295],' \
      'AXIS["Longitude",EAST], AXIS["Latitude",NORTH]]'
crs = QgsCoordinateReferenceSystem(wkt)
assert crs.isValid()

In the documentation it says that:

If you’re outside the PyQGIS console, the code snippets on this page
need the following imports:

from qgis.core import (QgsCoordinateReferenceSystem,
                       QgsCoordinateTransform,
                       QgsProject,
                       QgsPointXY
                       )

But it doesn't mention how to install qgis.core to my Python. If I do something like pip install qgis I get this error message:

ERROR: Could not find a version that satisfies the requirement qgis
(from versions: none) ERROR: No matching distribution found for qgis

Is there a way to use PyQGIS functions in a standalone Python file?


texts: Create Reduced Resolution Dataset (.rrd) file with modern software

texts:
I need to create, very specifically, a set of reduced resolution dataset (.rrd) pyramid files for a large set of TIFF imagery. I'm using ArcGIS 10.2, and have access to Global Mapper and QGIS. ArcGIS only creates .ovr pyramid files. I'm at a loss on how to do this with software that no longer supports nor creates these files for most image formats.

Is there a way to create these files with my software? Does anyone know an easily accessible software that will do this?