texts: arcpy iterating field names, values to create new feature classes

texts:
I want to create new feature classes based on selections from a list of fields. For now I'm using a user prompt instead of GetParamatersAsText for the selection criteria. I need the script to look for the first field name from the field list and then look for the first value in the value list, make a selection, and then create a new feature class based on the selection with the output name = current field name&""&current selection value. Then the the script should look for the next selection value and apply it to the current field, and so on. When there are no more selection values the script moves on to the next field name and makes a selection based on the first selection value again, and so on until there are no more selection values or field names.

#designate tower generated service location layer
SL = r"J:\workspace\Crown_Geocode\gisdb\layers.gdb\Rt8_15"

#Select route numbers
rtlist = []
while True:
    route = raw_input("Route Number (press enter when finished: ")
    rtlist.append(route)
    if route = "":
        break

#create new feature classes based on route number and day of service
arcpy.MakeFeatureLayer_management("Rt8_15", "Copy_Rt8_15")
fields = ("MonRoute","TueRoute","WedRoute","ThuRoute","FriRoute")
for field in fields, route in routes:
    arcpy.SelectLayerByAttribute_management("Copy_Rt8_15", "NEW_SELECTION", ' field = route ')
    arcpy.CopyFeatures_management("Copy_Rt8_15", ??????)

I should note that I will eventually want to use GetParametersAsText.


texts: Rasterizing polyline data with QGIS/GDAL: custom line width

texts:
I need to make a raster from a polyline shapefile. However, when I use gdal_rasterize (either directly or through QGIS), all my lines are drawn one pixel wide, while I need them to be drawn with the stroke weight of 3. I've examined the documentation at http://www.gdal.org/gdal_rasterize.html, but couldn't find a parameter to do so.

I have in mind some workarounds, but is there a simpler way to do so? My workarounds are:
- render a raster 9 times smaller and blow it up
- write a script that reads the raster and outputs a 3x3 block into another raster whenever it reads a non-zero pixel


texts: Snapping points and returning features in PostGIS?

texts:
I'm creating a point layer and I already have a line layer

There's a column in my point layer called line_id. I need to return automatically the id from the line to which the point is snapped to

I have written the following code, it works but ONLY when a point is snapped to a vertex in the line feature, not when the point is snapped to the rest of the segment

CREATE TRIGGER line_id_update
BEFORE INSERT OR UPDATE
ON pointlayer
FOR EACH ROW
EXECUTE PROCEDURE line_id_update()

CREATE OR REPLACE FUNCTION line_id_update()
RETURNS trigger AS
$BODY$
BEGIN
new.line_id := (select line_id from line_layer b where st_intersects(b.geom, new.geom) limit 1 ) ;
    RETURN new;
END;
$BODY$
LANGUAGE plpgsql VOLATILE
COST 100;

On QGIS I'm snapping to both vertex and segment, how to make it so that it will always return the line when it's connected to the segment


texts:
I have QGIS 2.18.15 and I'm trying to make a raster layer from polylines or points.
every time I start the process, it looks like it's working but won't really do anything (the cursor turns to round loading cursor, not seen in print screen)
left it all night to raster a single point
it seems like an ordinary problem but yet I couldn't find a solution

the software itself is ok and if I press "cancel" everything ok,
its seems like I'm missing something


texts: Need help translating small code block from VBA to Python

texts:
I need help from someone to translate a small piece of code from VBA to Python, due to the fact that the company I am working for has upgraded from ArcGIS 9.3.1 to 10.0 and VBA is no longer supported. The code block is as follows:

static x0 as double, y0 as double
dim pPoint as IPoint
set pPoint = [Shape]
x = pPoint.X
y = pPoint.Y
d = sqr((x-x0)^2 + (y-y0)^2)
x0 = x
y0= y

This is used with the calculate field tool inside model builder and as far as I can tell retrieves info about two points and calculates the difference in distance between them. I have attempted a few translations on my own with little success using !shape.firstpoint.x!, !shape.firstpoint.y!

Unfortunately my knowledge is lacking and I need some assistance. If you are able to help, that would be greatly appreciated!


texts: Choosing GTFS editor for creating and exporting routes?

texts:
I am a beginner and I need information about GTFS editors. I want to create routes,trips,patterns,etc... and export it to GTFS-format.

I have setup locally two GTFS editors:

TransitDataFeeder
conveyal-gtfs-editor

In TransitDataFeeder we need to create each file seperately through UI. I have created a sample agency and created a trip,route and exported it but its not containing the trips and as per my experience TransitDataFeeder does not meet my requirements for creating and editing GTFS.

I have setup conveyal-GTFS-Editor locally and trying to create new route in which its not possible to create trips using timetable option. I don't have experience with the play framework. I am trying to debug it but not able to solve it.

From above two, conveyal-gtfs-editor appears to meet my requirements but it is also in development phase so having problems as I saw in this blog .

So I am looking for best alternative to conveyal-gtfs-editor or need help to resolve conveyal-gtfs-editor with use time-table functionality.


texts: Query processing time for network table creation using pgRouting 1.05

texts:
I am trying to follow these steps for creation of a network table using pgRouting 1.05.

In my case I have two road geometry tables, one is bangalore_road and the other is all_india_roads. Obviously the bangalore_road table is a subset of big all_india_roads. So when I run the network query (after creating ext_view and node table) on bangalore_road_network it took almost around 3 hours to finish it up. In this case the size of my bangalore_roads table is around 8MB.

Now, when I run the same procedure on my all_india_road table (size ~280 MB) it seems like the query is just going on and on. i.e. now it's almost a week i.e 6 days and this query is still running at-least that's what I can see in the pgAdmin.

So just wonder if anyone have previously gone through the same blog post with huge dataset/tables and if yes then is there any way to speedup this process?

Like, running the query part by part (North_India, West_India, East_India, South_India etc.) on several machines and then combining all network tables together at end.not sure though as start_id and end_id may collide.

My system configs are:

PostgreSQL: 8.4
Processor: Intel Xeon @3.10 GHz
RAM: 4GB
OS: 64 bit win 7.


texts:
I have a main addresses point layer and a areas point layer. Both layers have XYZ coordinates in the attribute table. I need to create a table based from the addresses layer which calculates the 3D distance to every area point. The output I'm after is similar to the table below

ADDRESSES | AREA1 | AREA2 | AREA3 | ...
addressA  | ..... | ..... | ..... | ...
addressB  | ..... | ..... | ..... | ...
addressC  | ..... | ..... | ..... | ...

I've written a custom script which does this but it gets very slow since I have around 15,000 address points and 6,000 areas. I've read here on GSE that PostGIS is good for large datasets so installed it. But I do not know to how to do the query to achieve what I'm after.


texts: Translating and renaming attributes from one layer to another

texts:
I have two shapefiles, one that contains data from our Corporate GIS system and another that I want to import into a software model.

I want to select a subset of the Corporate GIS data (for instance all reservoirs in a given area) and create a new SHP for use in a model. I want to utilise only a subset of the attributes and rename them to match the schema for the model. I am happy to recreate the model SHP every time.

At present, I export it into Excel, and manually manipulate the data. However, I would like to ultimately do this step in a Process Model in QGIS so I can automate it. Is this possible to do? I tried joins but to no success.

For example, I would use an expression to Filter the Corporate SHP to the records I want. Then translate the fields and selection into a new SHP

Corporate SHP (AssetID) -> Skip
Corporate SHP (Name) - > Model SHP (ID)
Corporate SHP (Value) -> Skip
Corporate SHP (Level) -> Model SHP (RL)
Corporate SHP (Easting) -> Model SHP (X)
Coprorate SHP (Northing) -> Model SHP (Y)


texts:
Introduction:

I'm still fairly new to GIS and have only had about 3 weeks experience. I set up an OSM server locally with .osm.pbf data and it worked just fine. When I pulled the map in the browser, which implemented leaflet, all the names of streets, buildings, and other objects were already included in the data.

The Problem:

To increase my knowledge I went on to create a TinyOWS server with WFS-T standardization. To get the data for this server I pulled multiple states in the form of shapefiles from http://download.geofabrik.de/. When I load these shapefiles into QGIS I get all the streets, buildings, etc.., but no labels and names for any of the information. I can't find much helpful information on the internet for this, which means I'm thinking of model all wrong, or I'm visiting the wrong resources.

Question:

How do I get labels for streets, building, etc..., for shapefiles?


texts:
I have an Excel file with coordinates and a measured value belonging to this coordinates. In ArcMap I can display the points and color it concerning the measured values. Since most of the values close to each other have the same value, I want to build polygons out of them.

I thought the easiest way to do so would be to interpolate via kriging. However, now I don't find a way to create a polygon out of the resulting kriging grid. By using Conversion Tools --> from Raster --> Raster to Polygon the Input Field for Raster tells me:

  The resulted Kriging Grid is not in the correct data format.


texts: QGIS: Delete all records/rows with NULL values

texts:
I have joined a CSV file to a shapefile by a common key in QGIS. Not all records in the shapefile and CSV had corresponding keys.

This left some records in the shapefile as NULL.

i.e.

 City     School_ID   ID_no   Total
London    1234       1234     12
Glasgow   4321       4321     10
Bristol   4444       NULL     NULL

This is just a small slice of the data. The actual data has about 12000 records.

Is there any way to delete all records/rows with NULL values?


texts: Generating PDF for each feature in feature class based on attribute with ArcPy?

texts:
I have a feature class of points that represent storm water inlets.  I am trying to automate a process that does the following.

finds all of the inlets that have a field called "Marker" populated with the string "No"
generate a PDF for each of these features that is a centered view at 1:50 scale.

I'm pretty new to Python.  I've been working on this for awhile and keep getting an output of a single PDF (not 800, which is the aprox. # of inlets with "No") and the ArcMap crashes as the script is finishing.  I'm sure my code is pretty whack!


texts: Is every cell in flowdirection draining towards an outlet

texts:
I have created a Flowdirection map in a defined catchment.

I want to find out if every cell in my map is flowing towards the creek/outlet.

Is there an easy way to do this?

I have created the catchment using watershed.


texts: How to query subdistrict with multilinestring

texts:
I use postgis database and I download subdistrict shapefile of thailand from:

https://data.humdata.org/dataset/thailand-administrative-boundaries

and I have my point latitude and longitude (13.670524, 100.634237).

I try to query geom column as below query:

SELECT ST_ASTEXT(geom) FROM subdistrict limit 1

The result of geom column are multilinestring not multipolygon and data in multilinestring are coordinates; not lat, long.

I try to open shapefile in QGIS:

How do I query in subdistrict table with my point (13.670524, 100.634237)?


texts: Updated features from GDB after Check-out Replica created from ArcSDE GDB

texts:
I have created Check-out Replica from the main ArcSDE GDB. I have below questions.

If GDB Admin update few things (Attribute or Spatial) after version was created, will it be reflected in the child versions.


texts: What is the easiest way to find the ground altitude of an image taken by a drone

texts:
I understand that the Exif Altitude metadata that is included in images taken by DJI Drones is based on sea level. It is not the Drone altitude from the actual ground.

What is the easiest way to find the ground altitude of the drone? (means the vertical distance between drone and ground)


texts: Set Vector Layer line-width property to be dependent on attribute with PyQGIS

texts:
I am trying to use the QGIS Python API to get the same effect as manually going through the GUI:

Layer Properties>Symbology>Width>Assistant:

Input:

Source = float attribute
Hit "Fetch" button

Output:

Size from = 0.05
to = 3

Opening Layer Properties>Symbology>Width>Edit shows coalesce(scale_exp("<attribute>", <minValue>, <maxValue>, 0.5, 3, 0.57), 0) so I tried using it to set the data-defined attribute:
project = QgsProject.instance()
layer = QgsVectorLayer("path/to/layer", "VectorLayer", "ogr")
project.addMapLayer(layer)
f_n_index = layer.fields().indexOf('<float attribute>')
max_n_f = network.maximumValue(f_n_index)
min_n_f = network.minimumValue(f_n_index)
n_expression = "coalesce(scale_exp('<float attribute>', " + str(min_n_f) + ", " + str(max_n_f) + ", 0.5, 3, 0.57), 0)"
n_width = QgsProperty.fromExpression(n_expression)
n_prop = network.renderer().symbol().symbolLayers()[0].PropertyWidth
QgsVectorLayer.renderer().symbol().symbolLayers()[0].setDataDefinedProperty(n_prop, n_width))

But it doesn't seem to work. However, if I do the manual setup through the GUI, import another layer through the console and then do
GUI_layer.renderer().symbol().symbolLayers()[0].copyDataDefinedProperty(console_layer)

It works. What should I be doing differently?


texts: problem in read_shp() to generate a networkx.DiGraph from point and line shapefiles

texts:
I tried to use the Python module networkx and specifically the read-shp () to ensure the transformation of a graph to a file sharpe but I always find this error noting that this is the first time I developed Python (level zero).

       >>> import networkx as nx

       >>>G=nx.DiGraph()

       >>> G=nx.read_shp('/home/test.shp')

Traceback (most recent call last): File "", line 1, in AttributeError: 'module' object has no attribute 'read_shp'

thank you in advance for helping me


texts: How to access GDAL utilities from command line using Rgdal installation?

texts:
I have just installed rgdal from within R install.packages("rgdal").  I can now successfully do library(rgdal) with the following message:

Loading required package: sp
rgdal: version: 1.3-9, (SVN revision 794)
Geospatial Data Abstraction Library extensions to R successfully loaded
Loaded GDAL runtime: GDAL 2.2.3, released 2017/11/20
Path to GDAL shared files: C:/Program Files/R/R-3.4.4/library/rgdal/gdal
GDAL binary built with GEOS: TRUE
Loaded PROJ.4 runtime: Rel. 4.9.3, 15 August 2016, [PJ_VERSION: 493]
Path to PROJ.4 shared files: C:/Program Files/R/R-3.4.4/library/rgdal/proj
Linking to sp version: 1.3-1

Based on the startup message, I hoped to find tools like gdalbuildvrt in the directory C:/Program Files/R/R-3.4.4/library/rgdal/gdal but they aren't there, and aren't in the path on the command line.  How do I access these binary tools?

My ultimate goal is to use the system() command from my R script to call out to gdal directory to do some raster processing.  I don't want to read these big rasters into R because of performance issues.


texts:
I'm trying to join some downloaded CANVEC data to associated style file, but I am unable to join on the code field.

The field used to join the stylefile for ArcMap has been removed from the geometry, and updated with a new code, but the stylefile that comes with the data, has not been updated with the same code.

Would there be a simple solution, besides manually inputing the new code into the stylefile?


texts: Samples of 5 cm per pixel satellite imagery?

texts:
Are there samples of 5 cm per pixel or better resolution satellite imagery especially for a city or urban area ( street view level)?

Also are there providers of 5 cm per pixel satellite imagery for cities in Asia and Middle East that I can buy from and more important if they have such imagery available for selling to general public?


texts: ArcObjects custom GP tool: Multiple extensions using IGPFileDomain

texts:
I am using IGPFileDomain to make the user select between different file extensions, for example .type1, .type2, .type3. These are added using domain.AddType("type1") etc.

The dialog box then lists all these filters. In addition, it also adds the option to use all 3 extensions/filters at the same time. Example:

"File (*.type1)"
"File (*.type2)"
"File (*.type3)"
"File (*.type1;*.type2;*.type3;)"

Is it possible to tell/force this dialog to use the last option as default ?

Oyvind


texts:
I am working with c++/QT and c++api gqis 2.8.2

I want to change the geometry of a feature,
I have this code :

QgsVectorLayer * lLayer = mMapLayer.value(TRACKS);
QgsVectorDataProvider * lDataProvider;
lDataProvider = lLayer->dataProvider();
//create point
gsGeometry * MyPoint;
MyPoint = QgsGeometry::fromPoint(QgsPoint(PosX,PosY));
QgsGeometryMap MyGeometryMap;
QgsFeatureId MyFeatureId;

MyGeometryMap.insert(MyFeatureId,MyPoint);

lDataProvider->changeGeometryValues(MyGeometryMap);

I don't understand why i have an error on the insert :

error: no matching function for call to
 ‘QMap::insert(QgsFeatureId&, QgsGeometry*&)’
  MyGeometryMap.insert(MyFeatureId,MyPoint);
                                          ^
src/IHM/mainwindow.cpp:599:42: note: candidate is:
In file included from /usr/local/Trolltech/Qt-4.8.5/include/QtCore/qvariant.h:49:0,
                 from /usr/local/Trolltech/Qt-4.8.5/include/QtCore/QVariant:1,
                 from /usr/local/include/qgis/qgis.h:25,
                 from src/Core/common.h:13,
                 from src/Core/libcarto.h:12,
                 from src/IHM/mainwindow.h:24,
                 from src/IHM/mainwindow.cpp:20:

Any idea ?


texts: Change projection didn't change position in a vector shape in qgis

texts: I got a shape of points and when i save it with another datum and open it again, it stays in the same position. I was hoping for a slight change in position when save from epsg:5361 to epsg:32719.


texts: Print out the properties of the attribute table in qgis?

texts:
I have to merge a lot of shapefiles into one great shapefile.
So I want to compare the properties of the attribute tables to see if there are differences between the tableproperties of each shape (type, length, precision, name, etc). With the intention to standardize all the shapefiles

eg:

Is there a way to print out a list of the table properties as shown in the tab 'fields' or to export the properties into an .xls or .dbf?

Or is there an automated way to do the standardisation? (Perhaps this is a subject for a new thread...)

I am using QGIS 1.8 (QGIS 2.2 is also possible).

Thanks


texts: Writing GeoPackage with FME - no geometry

texts:
I have a large amount (400 shapefiles, 11 gigs total) of shapefiles which I would like to store into one file (or DB). I am also using PostGIS, but just out of pure interest I would like try test out how the files would work in a GeoPackage.

I am trying to write the files with FME into a GeoPackage and the features do get transferred and the file gets created, but when I try to open the resulting GeoPackage file on QGIS, it doesn't recognize any geometry in the file.

I was trying to go through the FME help and tried to force the file to create a spatial column, but it still doesn't seem to work as it should. The files contain only polygons and no Z-values. Any ideas what I'm doing wrong?


texts:
What is ArcGIS Pro .NET SDK ?

Can I code with ArcObjects?


texts:
I need to define an layer in OpenLayers just showing world map. The most simpliest variant is to create OSM layer. When I do like this:

var layer = new OpenLayers.Layer.OSM("OSM");
map.addLayer(layer);

How can I define restrictions/bounds to show the world map?


texts: Regrid a GRIB2 file and convert to netCDF

texts: I have a GRIB2 file on a Gaussian grid. I'd like to convert it to a netCDF file with to an equally-spaced 0.5-degree grid. How can I do this?


texts: What's the correct request to get the geometry from Overpass reverse geocoding with GeoJSON?

texts:
I found a query results a list of all the boundaries contains a point.

[out:json];
is_in(28.2,-16.6);
out;

But it lacks the shape information, preferably in GeoJSON format.

What should be included in this query to request such information?


texts: Update a table from list or np.array using arcpy.da.UpdateCursor?


texts: How to use MapInfo PEN patterns in MapServer?

texts:
I want to use tab files in MapServer WMS. My layer config:

 LAYER
   NAME Background
   TYPE LINE
   STATUS DEFAULT
   CONNECTIONTYPE OGR
   CONNECTION "D:\dev\Mapserver\ms4w\apps\example\tab\road.tab"
   STYLEITEM "AUTO"
   PROCESSING "GETSHAPE_STYLE_ITEMS=all"
   CLASS
     NAME "road"
   END
END # Layer

In MapInfo road.tab looks like:

But my WMS shows just default black lines. After I had added file with symbology to MAP section:

SYMBOLSET      "D:\dev\Mapserver\ms4w\apps\etc\symbols_mapinfo.txt"

I found this file in Mapserver docs , and it's contains mapinfo symbols like this:

# 67 - tent
Symbol
Name 'mapinfo-sym-57'
Type VECTOR
Filled TRUE
Points
0 1
.5 0
1 1
.75 1
.5 .5
.25 1
0 1
END
END

in ogrinfo I find that my layer use this symbology:

PEN(w:1px,c:#000000,id:"mapinfo-pen-73,ogr-pen-0")

But what next? How can I declare this style in symbology file?


texts: Creating shapefiles using textfiles in QGIS

texts:
I have a shapefile in QGIS which has polygons in it . There are several other attributes for each feature .

I want to add another attribute to each feature . The value for this attribute is available to me in a text file . How do i add this column programatically to the existing shapefile's attribute table .

It is not feasible to do it manually as i will have to change the value of that attribute quite often and the number of features are in thousands.

I am open to suggestions to use other GIS software if that fits my requirements.


texts: How to set up a Revision Control for ESRI File Geodatabase?

texts:
I'm interested in doing revision control (Mercurial) on some file geodatabases. Does anyone have an idea of which files are important (need to be saved) and which can be ignored? Note that I'm not interested in ESRI database versioning. These are essentially empty (no data), and I want to make sure they are in sync with code changes.

For example, if I add a new feature to the database, I add some things in code that uses that feature. If I want to go back to an old revision of the code, I want the database that matches. I don't even need to be able to compare versions, just get back to a particular one if needed.

Some of the files include: gdb, timestamp, *.freelist, *.gdbindexes, *.gdbtable,  *.gdbtablx, *.spx, *.atx


texts: Converting geocentric coordinates to geodetic coordinates?

texts:
I would like to start by saying that I am very new to this.
I have 160 different points which all have similar layout as below:

at point | X= 153346.357 | Y=354912.111 | Z=0 |

Doing some research I believe that these are geocentric coordinates. What is the formula which uses the above information (and the Datum: WGS 84) in order to convert these values to geodetic coordinates (i.e. GPS coordinates in the forms of DMS)?


texts: Is there any way to do queries in postgis using multiple cores in a transparent way, other than splitting my queries among several client connections ?


texts: When selecting input don't minimize window - pyqgis

texts:
I am making a plugin and I have a window in which you have to select an input file. There is a pushButton and lineEdit. I am typing a path and when clicking OK button to go back to window it minimalize. How can I avoid it?

Here is my code:

def select_input_file(self):
        filename = QFileDialog.getOpenFileName(self.iface.mainWindow(), "Choose input file ", os.getcwd(), "*.qml")
        self.d.lineEdit.setText(filename)


texts: Geometry column read as Numeric after operation in SpatiaLite (QGIS)

texts:
Pretty weird problem with SpatiaLite used within QGIS (2.13) with the DB Manager.

I have 2 layers (result of the same algorithm) and I use this function to merge them in a new table:

CREATE TABLE new_layer AS
SELECT * FROM lay1
UNION ALL
SELECT * FROM lay2;

The function works fine except for the fact that the new Geometry column is recognized as NUM Type and not as MULTILINESTRING (the 2 starting layers have a correct Geometry column type).

When I recover the geometry column, it is correctly recovered (so the layer is transformed into a geometric object) but still the Geometry column is a NUM type:

SELECT RecoverGeometryColumn('new_layer', 'Geometry',  3003, 'MULTILINESTRING', 'XY')

I think that this misunderstanding of SL can cause future incoherence with spatial operations.

Some idea?

Am I missing something?


texts: Logical judgement of Google Earth Engine

texts:
I am using logical judgement in GEE. Here is my toy code:

    var aa = ee.Number(0);
    print(aa.neq(0) );
    print(aa.eq(0) );
    if (aa.eq(0)){
        print('01 is zeros');
    }
    if (aa.neq(0)){
        print('02 is not zeros');
    }

Rather simple, but I got the following feedback:

Is there any recommended and STABLE logical judgement form in GEE?

By the way, this code doesn't work either; it returns the same feedback as shown in the figure.

var aa = ee.Number(0);

print(aa.neq(0) );
print(aa.eq(0) );

if (aa.eq(ee.Number(0))){
  print('01 is zeros');
}
if (aa.neq(ee.Number(0))){
  print('02 is not zeros');
}


texts: How to Add Basemap from ArcGIS Online to ArcMap 9.3?

texts:
I want to load a basemap into ArcMap 9.3 using "Add Data from ArcGIS Online" tools.

Does anyone know the URL needed to load it?


texts: Calculating mean value with field calculator for selected records using ArcPy

texts:
I want to calculate the mean value in the field z  which has Type Zero and write it in Ground field.
How can I do that with arcpy?

arcpy.MakeFeatureLayer_management(Shapefile, "layer")
arcpy.SelectLayerByAttribute_management("layer", "NEW_SELECTION","Type = 0")
arcpy.CalculateField_management("layer", "Ground", "mean(!{0}!)".format("Z"), "PYTHON_9.3")

I have got error for this code.


texts: Filtering (lasfilter on) huge LiDAR dataset

texts:
My dataset is about 8.3 GB after loaded into R environment by using lidR package. Some points which have return number greater than number of returns need to be removed.

With a small dataset, I could just use:

l <- lasfilter(las, ReturnNumber <= NumberOfReturns)

which is not possible with my PC.

What is the right way to filter (to use 'lasfilter' on) huge LiDAR data?
I am not sure which one to chose from the list of filter expression for argument 'filter' in the 'readLAS' function.  I modified a sample dataset from lidR package so it has points that return number are greater than their number of returns:

library(lidR)
LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
las = readLAS(LASfile)

# make rn > nr
las@data$ReturnNumber[4] <- 3L
las@data$ReturnNumber[8] <- 3L

writeLAS(las, "/TEMP/rnGTnr.laz")

# using "-keep_last", guessing from its source code (and it's wrong)
las2 <- readLAS("/TEMP/rnGTnr.laz", filter = "-keep_last")

# > Warning message:
# > Invalid data: 2 points with a 'return number' greater than the 'number of returns'.

So, the problem persists.

Can it be done using LASCatalog object (by dividing the data into chunks)?  I don't know how to filter using LASCatalog object. There is a filter option for LASCatalog, and again I don't know which one is appropriate for my case.


texts:
I'm trying to figure out (and eventually convert to WGS) projection system of coordinates I get from one web API service. As an example - lat, long = 71804125, 180232666 is a point in Cracow, Poland.
Any ideas what it can be?


texts: How to project WGS84 coordinates to calculate the distance in meters

texts:
I have geo coordinates from Google Maps which are WGS84 in decimal degrees for India. I want to calculate the distance between them in meters. Which projection should I use?

I use ArcMap 10.4.1.


texts: Current status of open source 3D GIS software and libraries

texts:
I'm an open source advocate and been asked to present to a top management in government regarding the availability of open source 3D tools and software that leads to SmartCity/SmartNation initiatives.

It is very difficult to keep track of what's happening around the world in 3D GIS space. Hope you can enlighten. Looking for open source initiatives, as a priority.

I have broken down into categories, if you have anything new other this list, that's what I'm interested.

3D GIS is tightly related to Point Clouds so let's start from there

LiDAR/Point Cloud

LASTools - LiDAR processing
PgPointCloud - PostgreSQL extension for storing Point Clouds
PDAL - Point Cloud Data Abstraction Library like GDAL

CityGML

3D City Database The CityGML 3D City DB

JS Libraries

CesiumJS JavaScript library for world-class 3D globes and maps
Cesium 3D Tiles 3D tiles spec for CesiumJS
ThreeJS 3D javascript library. QGIS has a plugin.

Projects

iTowns JS/WEBGL FRAMEWORK FOR 3D GEOSPATIAL DATA VISUALIZATION
Horao PostGIS 3D

Non-GIS, but of interest

Unity 3D engine The most popular 3D gaming engine
Blender Open source 3D model creation

Obsolete

Glob3 open source 3D GIS multiplatform framework written in java

I know GRASS GIS has 3D capabilities.

Can you let me know what else is out there?


texts: Running iteration at end of another in ModelBuilder?

texts:
I would like to run two iterations. One assigning the respective attributes from the polygon file containing the cols and rows. And after that iteration another iteration that selects all of those attributes with the same name and saves them in a separate SHP file. However I have not been able to nest one iteration at the end of another iteration process?

Supposed parent iteration model:

The supposed sub model that should use the result of the parent model:

I haven't been able yet to find a way to connect them together, or to transfer the name of the main model's resulting file to the sub model's iteration.


texts: Which coordinate system is this data in? Kilometers from the South Pole (90° S) and the 180th meridian

texts:
I'm using some data which, after some hours trying to figure out what its numbers mean, I found that its points are represented by two quantities: kilometers from the South Pole (90° S) and kilometers from the 180th meridian.

For example, these are some points I have with their equivalent lat-long degrees representations (I've rounded them):

16384,20642   41.6° N, 0.0° E
16384,20203   38.8° N, 0.0° E
15759,20203   38.8° N, 6.9° W

I'm not sure that the second quantity represents the kilometers from the 180th meridian. That's because as you go North the arc in the sphere should become smaller, but as you see in the first and second point in the example, although they're both in the same meridian both have the same value.

The software that is using this representation is SigPac (a well-known GIS for agriculture in Spain). They use this coordinate system internally to get parcels data.

Is there any projection that uses that reference system or something similar to what I've described?


texts: Parcel shapefile not lining up with basemap

texts:
I have a parcel shapefile for a municipality I'm working with. When I bring the shapefile into ArcMap, I change the fill color to clear and add a basemap. The parcel lines match with the the imagery from the basemap a little, but not 100%. For example, notice the parcel lines for the cul-de-sac in the upper right portion of the picture below don't line up.
The parcels are in the correct state plane coordinate system, although I just added the parcels into ArcMap, projecting on the fly rather than using any projection tool. I don't know if I need to do any formatting to the basemap.


texts: Plotting Smith and Sandwell global seafloor topography grids in ArcGIS Desktop?

texts:
Does anybody have any experience with plotting Smith & Sandwell 2-minute Mercator-projected grids in ArcGIS?

I've downloaded one in two formats: ers and img. I've found tools for processing .img data called Generic Mapping Tools (GMT).

Within this library there are commands which allow manipulating img grids, which is just what I need.

This example uses the img2grd command which I can't find within the GMT library.

The only tool I found so far was img2mercgrd but it wasn't successful.

How can I process this dataset to GIS format?


texts: Mapping raster values to nearest grid cell on a polyline

texts:
I'm not sure if the following task is possible in ArcGIS, or whether I need to use something like MatLab - I wondered whether the community could possibly advise if it is possible, and if so, what tools would be best to use.

I have a grid set up (see image below), with raster values (orange) to extract to a buffer zone [polygon] (light green) along the coastline [a polyline] (where the coastline acts as source inputs into the ocean). Not all of the grid cells in the buffer are on the coastline however, so I was wondering if there was a way to map their values to the nearest coastline grid cell.

I'm happy for any output of the raster into the grid. All I require at the end is that each grid cell on the coast has a co-ordinate (already defined) and a total population value (the raster) that is on the coastline AND that derived from other cells in the buffer being allocated to their nearest coastal cell.

There are datasets for the coastline, the buffer, and the raster.

Desktop Version 10.2
Co-ordinate system: WGS84

I am happy to provide the datasets if it would be of assistance.


texts: How to determine the appropriate size of a 3 dimesional raster block/chunk

texts:
I am processing large MODIS time series (12 years) dataset with an array dimension of 4800 x 4800 with a data type of Int16-Sixteen bit signed integer. I want to build a 3 dimensional array (row,col,time(10 years)) using np.dstack. This is because I need to access the elements of each index in the 3D array and perform interpolation to those values flag as clouds. I know I cannot build my 3D array with the a dimension of 4800 x 4800 x 644 array my computer simply crashes after processing 10 rasters/array, hence I plan to process it by chunks/blocks. How to set the optimal block size (1200 X 1200 x 644 or 2400 x 2400 x 644) so I can process it faster without my computer crashing?
I am running a 32-bit python with 12GB ram and 64-bit operating system.


texts:
I have the following table :

          Table "osm_france.france_2po_4pgr"
Column     |           Type            | Modifiers
---------------+---------------------------+-----------
id            | integer                   | not null
osm_id        | bigint                    |
osm_name      | character varying         |
osm_meta      | character varying         |
osm_source_id | bigint                    |
osm_target_id | bigint                    |
clazz         | integer                   |
flags         | integer                   |
source        | integer                   |
target        | integer                   |
km            | double precision          |
kmh           | integer                   |
cost          | double precision          |
reverse_cost  | double precision          |
x1            | double precision          |
y1            | double precision          |
x2            | double precision          |
y2            | double precision          |
geom_way      | geometry(LineString,4326) |
Indexes:
   "pkey_france_2po_4pgr" PRIMARY KEY, btree (id)
   "france_2po_4pgr_gix" gist (geom_way)
   "idx_france_2po_4pgr_source" btree (source)
   "idx_france_2po_4pgr_target" btree (target)

As you can see, a spatial index is available for the geom_way column.
I have performed a VACUUM ANALYZE on the table for safety.

Now I run this very simple query :

select * from osm_france.france_2po_4pgr where
ST_DWITHIN(St_TRANSFORM(geom_way,3035),ST_transform(ST_SetSRID(ST_MakePoint(2.3, 48.8),4326),3035),100);

So I'm looking for all the geometries within 100 meters from a given point (I transform the SRIDs from 4326 to 3035 to specify a range in meters).

The explain analyze yields :

Seq Scan on france_2po_4pgr  (cost=0.00..1630779.00 rows=69644 width=293) (actual time=35.044..9407.356 rows=3 loops=1)
Filter: ((st_transform(geom_way, 3035) && '0103000020DB0B0000010000000500000046D8C5565EA84C41A950F47CDFFF454146D8C5565EA84C41A950F47C4300464146D8C556C2A84C41A950F47C4300464146D8C556C2A84C41A950F47CDFFF454146D8C5565EA84C41A950F47CDFFF4541'::geometry) AND ('0101000020DB0B000046D8C55690A84C41A950F47C11004641'::geometry && st_expand(st_transform(geom_way, 3035), 100::double precision)) AND _st_dwithin(st_transform(geom_way, 3035), '0101000020DB0B000046D8C55690A84C41A950F47C11004641'::geometry, 100::double precision))
Rows Removed by Filter: 5225288
Planning time: 0.517 ms
Execution time: 9407.390 ms

Why is this query not using the available spatial index ?


texts: Perpendicular lines from point to another feature

texts:
Scenario: Measurements made along road centerline, location of which shown by dots.

Objective: Compare measurements with building properties which are perpendicular to it (on both sides).

Workflow: Draw Perpendicular lines from points to building line. Spatial Join  these perpendicular lines with buildings and points.

Question: How to draw perpendicular lines from these points till it encounter first building polygon with road centerline acting as base for angle.

Edit: ArcGIS 10.5 Advance with Spatial and Network Analyst extensions.


texts:
I have two equal point layers (76300 features) with different attributes as a result of queries with the "point sampling tool". Now I want to merge the attribute tables using the "join attribute by location" tool. Whatever "geometric predicate" i use (although i think it should be "equals") the result is just one point.

What am I doing wrong?


texts: How to deal with overlapping polygons from GeoJSON in OpenLayers 4

texts:
I have a series of GeoJSON files (generated by QGIS) loaded on to a map in OpenLayers 4 as different layers. Each file contains multiple polygons. When mouse-clicking each polygon a pop-up window appears with info on the polygon.
Layers (each layer corresponding to a single GeoJSON file) are loaded as:
//adding layers
map.addLayer(ocean_map);
map.addLayer(layer_mediterranean_aoi);
map.addLayer(layer_mediterranean_cimma);
map.addLayer(layer_mediterranean_imma);

While the function to trigger the pop-up is written in this way:
     var info = document.getElementById('info');
     var target = document.getElementById('map');
     function displayFeatureInfo(pixel) {
            info.style.left = '30%';
            info.style.top = '20%';
            info.style.height = 300 + 'px';
            var feature = map.forEachFeatureAtPixel(pixel, function(feature, layer) {
        return feature;
      });
      if (feature) {
        var geometry = feature.getGeometry();
        var coord = geometry.getCoordinates();
        var text = '<button class="close-button" id="close-button">X</button>

';
        //if (feature.get('AOI')) {text += '

 ' + feature.get('AOI') + '

';} else {text += '';}
        text += '

 ' + feature.get('Title') + '

';
        text += '
Summary: ' + feature.get('Summary') + '

';
        text += '
Region: ' + feature.get('Region') + '
';
        if (feature.get('AOI')) {text += '';} else {
        text += '
 Criteria: ' + feature.get('Criteria') + '
';
        }
        if (feature.get('AOI')) {text += '';} else {
        text += '
 Species: ' + feature.get('Species') + '
';
        }
        if (feature.get('Status')) {text += '
 Status: ' + feature.get('Status') + '
';} else {text += '';}
        if (feature.get('URL')) {text += '
  MORE INFO
';} else {text += '

';}
        text += '

';
        info.style.display = 'none';
        info.innerHTML = text;
        info.style.display = 'block';
        target.style.cursor = "pointer";

        //handles the closing "x" button to shut down the popup
        var closer = document.getElementById('close-button');
        closer.addEventListener('click', function () {
              info.style.display = 'none';

        });

      } else {
        info.style.display = 'none';
        target.style.cursor = "";
      }
    }

  //handles the opening of the popup
  map.on('click', function(evt) {
      if (evt.dragging) {
        //it doesn't show the info if I just drag the map
        info.style.display = 'none';
        return;
      }
      var pixel = map.getEventPixel(evt.originalEvent);
      displayFeatureInfo(pixel);
  });

Few of these polygons may overlap, and obviously the one that is in the background (if smaller) is not selectable. How can I fix the problem?
At the moment I figure out that one way is to prepare the geojson files with the polygons ordered in a way that the smaller areas appear after the bigger ones where they overlap, so that when the  layer is loaded on the map the smaller polygons are loaded last (therefore on top).
I was wondering if there is a better way to solve the problem on OpenLayers side.


texts: How to project a GeoJSON on the fly in OpenLayers 3?

texts:
I have been trying to load a GeoJSON without success. I have a map with a background layer from OSM. I have some data in 4326 I want to overlay on that map. I do not want to transform coordinates one by one on the server-side and I'd rather have OpenLayers transform them on the fly ... if possible at all.

I have tried using the "projection" attribute under the vector Layer, vector Source, main background layer, ... everywhere I have thought of ... to no avail: the sample point keeps showing near 0,0.

I set up a fiddle with my test code, any idea? Is it simply not possible?

http://jsfiddle.net/zzahmbff/2/

The data here is just a sample. Once I get it to work I'll load the data with AJaX (Asyncronous Javascript XML).


texts: QGIS, how to create a grid that does not exceed the extent, constrained within the extent

texts:
I have to create a grid of 15mx15m over an area much bigger, 25000mx25000m. Of course it means that the grid will never properly fit within the extent, but it will exceed a bit (25000/15=1666.6667). I tried with both the "research tool, vector grid" and the "create grid tool" with MMQGIS. The first create a grid that is bigger than the given extent, the second smaller.

The question has somehow already been asked, but solved poorely with the clip tool... In this case it will take ages by using this tool... Anybody with a better idea?


texts: Select the most distant vertex from polygon's centroid using ArcPy

texts:
I have a polygon and the polygon's centroid. I want to select the most distant vertex from centroid (picture: point A). I used arcpy.FeatureVerticesToPoints_management to create a feature class containing points generated from polygon, but I don't know what to do next. Could you help me?

Also I want to find intersection point (point B) of the outline of polygon and line created by centroid and the most distant point (line OA).

I would be very grateful for any help.


texts: problem creating and editing features at extent edges

texts:
I'm in the process of creating a polyline feature layer of tectonic plates. I'm working with WGS-84 reference system, which means the extent is +-180 in lon, and +-90 in lat.
In order to create the features in a correct way, I imported a CSV file with the 4 "corners" and displayed them as points I then snap my features to when I create them.

The problems arises when I want to snap the line to the edge point. It simply doesn't want to create the line all the way, cutting it shy of the actual point.

Here's what I see:

The coordinate is that of my cursor. The green point is (-180,-90).

The green background is another layer I used (continents) to see if there's some kind of error or a general problem, but as you can see, it sits perfectly in place.

I have tried different editor tools (line, polygon, autocomplete), but the problem persists.

I have also tried to extend that segment the the point, as well as moving a node to the point, but it doesn't edit it, and jumps back to it's previous position.

Does anyone know why this happens?


texts:
How can I import anaconda-installed modules from arcgis's python shell? I have tried the method described in https://sites.google.com/site/advancedgisprog/other-resources/arcgis-anaconda-packages.

Steps I followed:

navigate to the site-packages directory in the ArcGIS python install (C:\Python27\ArcGIS10.3\Lib\site-packages). In that folder, I created a new text file, and named it AnacondaPath.pth. Opened this file in Notepad.
navigate to the Anaconda site-packages folder (C:Users\anwar\Anaconda3\Lib\site-packages). Once there, Paste this path in the Notepad window of the AnacondaPath.pth file, and save the change. Now, open the ArcGIS python interpreter, and try to import a dependency known only to be available to Anaconda. For example, I want to import pandas which is installed into my anaconda. But it showed me the error:

import pandas

Runtime error

Traceback (most recent call last):

   File "<string>", line 1, in <module>

ImportError: No module named pandas

Also I followed the same process I did for importing arcpy into spyder:

Importing arcpy into spyderI copied the arcpy folder from 'C:\Program Flies(x86)\ArcGIS\Desktop10.3' and pasted that into my 'C:\Users\anwar\Anaconda3\Lib\site-packages'. Now I can import arcpy into spyder of Anaconda3.

Importing pandas into arcgis's python I copied the pandas folder from 'C:\Users\anwar\Anaconda3\Lib\site-packages' and pasted that into 'C:\Program Flies(x86)\ArcGIS\Desktop10.3'. Also added this path into Desktop10.3 from 'C:\Python27\ArcGIS10.3\Lib\site-packages'. Still I am unable to import pandas from arcgis's python.


texts:
I'm trying to write a Python script in QGIS 2.18.16 to use zonal statistics of a multi-band raster (GeoTIFF) using a polygon shapefile.

QGIS crashes turning off so I can't get any error message to correct my code.

I'm working in Ubuntu 16.04 64 bit.

Am I doing something wrong or is there a bug in QGIS?

import processing
vectorlayer = qgis.utils.iface.mapCanvas().layer(0)
rasterfile = qgis.utils.iface.mapCanvas().layer(1).source()

processing.runalg("qgis:zonalstatistics", rasterfile, 3, vectorlayer, "Zonal_", True, "/media/giacomo/4A48-5E67/prova/out.shp")


texts: Export Query Layer does not return all values

texts:
When I create a query layer, the attribute table has 1004 records. When I attempt to export this data to a file geodatabase, it only has 523 records. I have tried exporting to a personal geodatabase, a shapefile, and also using Copy Features to export the data, each time only 523 records are exported.

What is going on?

When exporting to a geodatabase or using copy features tool, there is no error message. Only when I export to a shapefile, do I get three separate pop up windows:

Error exporting data
Error reading some rows from the source. These rows have not been converted.
Warning: Some records may not have been exported correctly.
There was an error exporting the layer.

I have managed to export many features successfully, but there appears to be a few tables that struggle with this.

I'm new to this type of work.

As a note, there are no joins or selections on the query layer when I attempt to perform the export. I am using ArcMap 10.4.1.5686, Advanced License.


texts: Incomprehensible speed lowered situation

texts:
I recently started to learn QGIS system and used quite well to develop a project.
I wanted to have faster outputs so I bought a new desktop.

The old one specification is like this:

CPU intel i5-7400 CPU @ 3Ghz
8GB RAM
No extentional VGA

and new one is like this:

CPU AMD Ryzen 9 3900 (12 core, 24 threads)
64GB RAM
VGA GeForce RTX2080 Ti

Yes, I put much money to see higher performance but it is totally slower than before for making even simple joining by location.

Is there any hint of reason?

Other environments are same like
Windows 10 / QGIS 3.8 / sets all the settings and options between old and new one.

I can't understand...


texts: Removing points outside of standard distance circle in ArcGIS for Desktop?

texts:
I have a feature class of companies where each may have multiple locations (for example, company A may have a cluster of points associated with it). For each company, I used the Standard Distance tool to create a feature class of circles, which surrounds the group of points for that company that is within 1 standard distance of the mean point. See the image below:

I highlighted the points that belong to a company and the standard distance circle for those group of points (the big green circle) for that company.

What I need to do is to remove all the points that fall outside a company's standard distance circle. So in the image, I need to remove the lone highlighted point in the top right corner among the cluster of red points. I need to do this for every company in the dataset (around 33,000 of them).

Here is what I have tried:

Use the Select Layer by Attribute tool two times to select company A, for example, in

the point feature class containing the points belonging to company A
the standard distance circle for company A

Use the Clip tool to keep the points for company A that fall inside company A's standard distance circle
Do 1 and 2 for every company in the dataset

I can write a Python script to execute this, but given that there are 33,000 companies, my sense is that this will be slow. Any suggestions for a more direct approach?


texts: Is there any number of rows limitation or data limitation for CartoDB "Mercator Plan"?

texts:
We are using Google MAP API V3 and CartoDB with "Mercator Plan" and having around 200 MB empty space. I am trying to visualize whole UK wide properties on the map with points and I have added about 1.4 million properties geo location into CartoDB table. Its working fine till now but I have added more data into CartoDB table and its going very slow as well as stuck some images of the layer on the map and giving the error "NetworkError: 400 Bad Request" in the console for these images.

Is there any number of rows limitation in table or data limitation for "Mercator Plan"?


texts: Georeferenced Image Overlay and export to new image

texts:
I have two Geotiffs (one satellite image, one clipped image). I want to put the clipped image over the satellite image, because it shows flow velocity from this shown glacier (see image below). I tried qgis, Saga and Matlab with no problems to load and visualize that on screen. But I'm not able to export this overlay to a new geotiff. I've tried Raster Merge in qgis with layer stacking, but I can't visualize like same before. I need geotiff to put more information on it.

Does anyone have an idea for how can I resolve that?


texts: Google Tracks convert kml to shape and preserve attributes

texts:
I've been attempting to convert a Google Tracks kmz file to an ESRI shapefile using the GDAL and libkml(http://www.gdal.org/drv_libkml.html) which I built on my Ubuntu 14.04 (64-bit) system. The kml is a line vector but when converting with GDAL I end up with point data and all attributes concatenated into one field 'Description'. I've read a ton of articles which suggest that this may be the case with kml as they often have non-standard schema, but I wondered if anyone has had any jy converting Google Tracks to shapefile or other GIS format.

Here's how I'm attempting:

Unarchive the kmz to a kml file
Use the command line ogr2ogr to convert kml to shapefile as follows: ogr2ogr -f "ESRI Shapefile" myfile.shp myfile.kml

As I said the shapefile ends up as point vector layer with two points as opposed to a line vector showing the actual path that I walked and recorded using Google Tracks. If I try and load the kml file directly into QGIS I also get point data not a line (note: QGIS is built against an older version of GDAL without libkml support - I built another version of GDAL with libkml support to test this conversion).

Alternatively, if I try and convert the kml to csv the attributes are preserved but no coordinates.

I have also tried converting to a geojson using the link I found here: gpsies.com/convert.do from a previous stackexchange question (How to convert KML with tracks to Shapefiles?) and loading straight in QGIS as a geojson. This preserves the line vector (track) but again loses all attributes.

link to an example kmz file here :https://drive.google.com/file/d/0B-TLTa7JtM_FSk9JY0wwdmkwcHc/view?usp=sharing


texts: Create a combined Landsat and Sentinel mosaic using SEPAL

texts:
I'm new to the SEPAL platform and I cannot find lot of information on the web so far.

I would like to create a Landsat and Sentinel based mosaic on a AOI that I stored as an asset in GEE : users/bornToBeAlive/aoi_sandan
How should I proceed ?


texts:
I want to change the color to a layers after the user has selected the color from the interface
what have I done so far:
 var imobile = new ol.layer.Vector({ source: vectorSource.teren });

    $('#color-picker-imobile').on('change', function(){
                        var imobile_color = $('#color-picker-imobile').val();
                        console.log(imobile_color);

                        $.post( "ser_edit_pdf.php",
                             {
                                  save_setari_imobile: 1,
                                  imobile_color: imobile_color
                             },
                             function( data ) {
                                console.log('save setup');
                                console.log(data);
                                var color = JSON.parse(data);

                                    //imobile.setStyle(new ol.style.Style({fill: color.color_imobile, opacity: '0,9'}));
                                    imobile.getSource().forEachFeature(function(feature){

                                        //console.log(feature.getProperties());

                                        console.log(color.color_imobile);

                                          var style = new ol.style.Style({
                                                fill: new ol.style.Fill({color: color.color_imobile}),
                                                stroke: new ol.style.Stroke({
                                                  color: [0,0,230], width: 1
                                                })
                                            });
                                            feature.setStyle(style);
                                        });

                             });

I want to do this without reloading the map with the new setting
I store the layers as follows:
var vectorSource = {
                            teren: new ol.source.Vector(),
                            constructii: new ol.source.Vector(),
                            uat: new ol.source.Vector()
                        };

I could not find any concrete example for my situation :(


texts: OS VectorMap Local not styling RailCLine and RoadCLine

texts:
I am running QGIS 3.14.1 (64 bit version) on a Windows 10 PC. I have downloaded two 5km tiles of the VectorMap Local (3 July version) in Geopackage format. My Parish Council holds the relevant OS Licence. I have downloaded the relevant styles from Github, in QML format. When I load them, Text, Line and Area seem to work but RailCLine and RoadCLine do not - see this screengrab:

The two A roads should be coloured and also the East Coast main railway line should be shown. When I look at the symbology and the symbol levels the colours are there as the various fills.
Clearly I am missing something, but what? I have looked at the other posts here, but most of them are nearly seven years old and relate to different levels of QGIS and map. There is one fairly recent question but that does not seem to have been resolved.

Thanks for the comments. I've looked at the references but they don't help me (though I may be overlooking something, of course.) There are differences in that I am using QGIS 3.14 (not 2) and the map is VectorMap Local not District. This may or may not be significant. I have looked at the Layer properties, clicked to edit the rule (but not actually changed it) then clicked Test. It's finding the features, just not displaying them properly:

Any thoughts?


texts: Embeding ESRI compatible ArcGIS projection using GDAL?

texts:
Is there any way to embed ESRI compatible ArcGIS projection using GDAL -tools?

Currently I am using gdalwarp to create raster file which translate raster from WGS84 (EPSG:4326) to Spherical Mercator (EPSG:3857). It also create Geotiff with tags and tfw -file.

gdalwarp -q -overwrite -s_srs '+init=epsg:4326' -t_srs '+init=epsg:3857' -r bilinear -of GTiff -co COMPRESS=LZW -co TILED=YES -co TFW=YES...

It still not open correctly when opening in ESRI ArcGIS Desktop 10.4 but complining missing projection. Is it possible to create raster file which open directly and correctly in ArcGIS environment and can be produced by using GDAL tools?


texts:
When I click on Start Editing button following error occurs:

  Could not edit the data in the folder or database you selected. The featureclass is currently in load-only mode.

What do I need to do to start editing?

Currently I am not able to do editing in ArcMap, ArcGIS 9.2 and data is an enterprise geodatabase stored in Oracle 10gR2.

I do data migration on electric dataset while loading data into existing dataset on database server. When I imported data in ArcMap, after clicking on Start Editing, the following error highlighted:

  Could not edit the data in the folder or database you selected. The featureclass is currently in load-only mode.


texts:
With FME Desktop I'd like to write a text file with a dynamic filename, consisting of the date of today and a guid.
I know how to create both parameters (date and guid), but I don't understand how to integrate it into to filename, e.g. 2020-09-02_ABCDEFG123id.txt
The problem about it is that only feature type names can be created dynamically. But the filename has to be chosen when adding the writer, so it cannot be changed depending on the input data.
I tried the TXT-, the XML- and the Generic-writer.


texts: Shapefile PRJ to PostGIS SRID lookup table?

texts:
I was wondering if there is such a thing as an shapefile PRJ to PostGIS SRID lookup table?  Something that can translate the most standard shapefile PRJ definitions into the likely SRID.

When using PostGIS and pgAdminIII, if you use the postgisgui to import your shapefiles, the SRID is left as "-1".  It seems like the tool should be able to parse the Esri PRJ and determine the correct (or at least a couple of options) that are the likely SRID, rather than just leave the default.

Or does the importer have the capability to reproject on the fly if you choose another SRID?

It may seem lazy on my part, but to me it seems curious that this function hasn't already been put in place.  Does anyone know if this concept is in the works, or good reason's why it has been left out?


texts:
I would like to use Geoserver SLD to style points using one attribute for the size in mapunits and the second attribute for the colour.

This works well as separate rules styling either the colour or the size of the symbol.

This is for Geoserver 2.8

I don't get any result by combining the rules like so:

        <ogc:And>
          <ogc:And>
            <ogc:And>
              <ogc:PropertyIsGreaterThanOrEqualTo>
                <ogc:PropertyName>tt_grid</ogc:PropertyName>
                <ogc:Literal>818</ogc:Literal>
              </ogc:PropertyIsGreaterThanOrEqualTo>
              <ogc:PropertyIsLessThanOrEqualTo>
                <ogc:PropertyName>tt_grid</ogc:PropertyName>
                <ogc:Literal>159174</ogc:Literal>
              </ogc:PropertyIsLessThanOrEqualTo>
            </ogc:And>
            <ogc:PropertyIsGreaterThanOrEqualTo>
              <ogc:PropertyName>y15plus_t_2014</ogc:PropertyName>
              <ogc:Literal>2.2</ogc:Literal>
            </ogc:PropertyIsGreaterThanOrEqualTo>
          </ogc:And>
          <ogc:PropertyIsLessThanOrEqualTo>
            <ogc:PropertyName>y15plus_t_2014</ogc:PropertyName>
            <ogc:Literal>7.2238</ogc:Literal>
          </ogc:PropertyIsLessThanOrEqualTo>
        </ogc:And>

Here is a screenshot directly off Ian's solution.


texts:
I am trying to migrate my ArcMAP model over to ArcGIS Pro but experiencing some issues. The model is using an iterator to run the model on multiple files. In ArcMap the model would stop once it has gone through all the files in the specified folder, but ArcGIS Pro just keep looping the model indefinitely.

Please advice


texts: Combine tool: Does it perform a union or intersect operation?

texts: Can I use combine with different shaped rasters? Does it perform a union or intersect algebra operations with the fields?


texts: Problem merging/clipping MODIS rasters in QGIS

texts: Im pretty new to QGIS and am having some trouble making a NDVI map for Mali. To cover the whole country I need 5 MODIS NDVI tiles, however every time I try to merge them I get a "the following layers were not correctly generated" error message. I get the same message when I try to clip the rasters with my underlying polygon (state borders). I have tried to read through past threads but nothing seems to work and I am getting a bit desperate.


texts: Do you use a proven workflow to edit and document ArcGIS ModelBuilder models?

texts:
Do you have any reference to best practice examples , or can you share your expertise on how to keep track of edits to ModelBuilder models , track changes between model/tool versions and maintain  model stability, robustness and function?

I am aware of the ESRI guidance to managing, copying and naming models here.

I would be interested if anyone has established practices, or specific expertise in the use of version tracking, archiving, editing of ModelBuilder tools that they successfully use to track model development and checking and to maintain models in a stable condition.


texts:
I have made a map from QGIS by using RT MapServer Exporter plugin, and I want to browse it to MapServer, just like this:
https://trac.osgeo.org/openlayers/wiki/MappingYourData#no1
and the tutorial on the youtube: https://www.youtube.com/watch?v=WJt0rAyD2xU

However, when I browse it to MapServer, there is always a error:

msLoadFontset(): Unable to access file. Error opening fontset ./fonts/fonts.txt.

Someone said "This means that that the map file is requesting for fonts or symbol definition files but that MapServer cannot find those." But where I can download the font file and how I can put it into MapServer?


texts: Advance text style of feature in openlayers

texts: Is there any way to create table or any html element in text style for features? i tried other ways like putting overlay on feature then in overlay i used table and html elements but had some problems with it also control of overlays to disable and enable on map for each feature was hard.


texts: Adding X and Y text as pre string to Grid Label Tags using Dynamic Text - ArcPro

texts:
Is there any text command to add the "X" and "Y" value when creating Grid Label Tags in ArcPro with Dynamic text?
By using the following I can tag correctly except for the X and Y values:

dyn type="grid" units="meters" showDirections = "False" showPlusSign
= "True" separator="True"/
= "+2 660 000m"

I am aiming for "X +2 660 000m" and Y where Y is relevant


texts: Network Analyst: arcpy's 'Solve' running out of memory

texts:
I've run OD cost matrix analysis using GUI tools n ArcGIS 10.

My goal is relatively simple - find closest 50 neighbors for each origin point.

The dataset I'm working with is relatively large:

network dataset consists of ~ 1.16 million elements
there are ~ 100k origins
and ~ 1.3 million destinations.

Origins are part of the destinations dataset.

The limits I set for the geoprocessing tool are: default_number_destinations_to_find set to 51 (since first 'neighbor' will be the same, overlapping point) and default_cutoff = 50000 (since I know from previous analyses that this value cannot be larger than that).

In order to save memory (I think) I go with output_path_shape set to NO_LINES - I do not need graphical representation of OD matrix - just tabular data (correct me if I'm wrong here).

It takes a while to process this task (~ 24h), but nonetheless I get the results.
So far so good.

In order to automate this task I tried putting all this steps into python script and came up with:

try:

    import arcpy
    from arcpy import env

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

    # set environment settings
    env.workspace = "D:\\GIS\\04-buildings.gdb"
    env.overwriteOutput = True

    # Local variables:
    NW_FIX_ND = "D:\\GIS\\02-network-restrict.gdb\\NW_FIX\\NW_FIX_ND"
    DESTINATIONS = "buildings\\destinations"
    ORIGINS = "buildings\\origins_1"
    OD_COST_MATRIX = "OD Cost Matrix"
    LINES = "OD Cost Matrix\\Lines"
    bn_1 = "results\\bn_1"

    # Make OD Cost Matrix Layer
    arcpy.MakeODCostMatrixLayer_na(NW_FIX_ND, "OD Cost Matrix", "Length", "50000", "51", "Length", "ALLOW_UTURNS", "restriction", "NO_HIERARCHY", "", "NO_LINES")

    # Add Destinations
    arcpy.AddLocations_na(OD_COST_MATRIX, "Destinations", DESTINATIONS, "Name h_gebaeude_id #;CurbApproach # 0", "5000 Meters", "", "TLM_STRASSE_RESTR_FIX SHAPE;NW_FIX_ND_Junctions NONE", "MATCH_TO_CLOSEST", "APPEND", "NO_SNAP", "5 Meters", "EXCLUDE", "TLM_STRASSE_RESTR_FIX #;NW_FIX_ND_Junctions #")

    # Add Origins
    arcpy.AddLocations_na(OD_COST_MATRIX, "Origins", ORIGINS, "Name h_gebaeude_id #;TargetDestinationCount # 51;CurbApproach # 0;Cutoff_Length # #", "5000 Meters", "", "TLM_STRASSE_RESTR_FIX SHAPE;NW_FIX_ND_Junctions NONE", "MATCH_TO_CLOSEST", "CLEAR", "NO_SNAP", "5 Meters", "EXCLUDE", "TLM_STRASSE_RESTR_FIX #;NW_FIX_ND_Junctions #")

    # Solve
    arcpy.Solve_na(OD_COST_MATRIX, "SKIP", "TERMINATE")

    # Select OD Lines
    arcpy.SelectData_management(OD_COST_MATRIX, "Lines")

    # Select
    arcpy.Select_analysis(LINES, bn_1, "")

    print "Completed successfully"

except Exception as e:

    # If an error occurred, print line number and error message
    import traceback, sys
    tb = sys.exc_info()[2]
    print "An error occured on line %i" % tb.tb_lineno
    print str(e)

However, when I attempt to run this script I get the following error:

An error occured on line 39
ERROR 030024: Solve returned a failure.
Out of memory.
Failed to execute (Solve).

Why does arcpy fail when GUI could do the job?

My reasoning was that by running the task in python, and avoiding the graphical overhead, it should even work faster.


texts: How to get QGIS canvas in plugin window?

texts:
I'm working on project where I need insert QGIS canvas as a widget in QGIS plugin. When user click on button it should be loaded into widget. Something like this, this is example where satellite map should be loaded in widget:

How to implement this in python code?


texts: Getting error converting .tif to .map format?

texts:
I am getting this error when I try to convert .tif to .map format:

  ERROR 6: PCRaster driver: Cannot convert cells: Usetype is not version
  2 cell representation, VS_LDD or VS_BOOLEAN  Input file size is 2240,
  2160.

How do I deal with it?


texts: How do I rename rasters in an mxd by using my csv?

texts:
I have the following data in a csv:

    RasterName  First name
    Raster 1    Aaron
    Raster 1    Paul
    Raster 2    Bryan
    Raster 3    Cranston
    Raster 3    Vince

I have added these rasters to my mxd as many times as necessary (twice for Raster 1, once for Raster 2 and twice for Raster 3), so that my TOC looks like the RasterName column.  Now I need to rename the raster in my map document to the appropriate First name.  I read the csv into csvdict using reader and when I use the following:

csvdict = []
for row in reader:
        csvdict.append(row)

for layer in arcpy.mapping.ListLayers(mxd):
    for i,j in csvdict:
        if layer.name == i:
            layer.name = j

I get this:

    Aaron
    Aaron
    Bryan
    Cranston
    Cranston

When I really want:

    Aaron
    Paul
    Bryan
    Cranston
    Vince

I'm sure I'm missing something simple, but I haven't been able to figure it out.


texts: ST_Union Topology exception despite of valid geometries

texts:
I have an isochrones layer and I want to calculate the population lays on the overlapping regions.

My first step is to create an intersection between all the isochrones. After that, I try to create a union with the intersections.

Unfortunately, ST_Union throws GEOSUnaryUnion: TopologyException.

I've checked that all the geometries used for the union are valid.

Union query:

select st_union(the_geom) from isochrones_intersection;

  ERROR:  GEOSUnaryUnion: TopologyException: Input geom 0 is invalid:
  Self-intersection at or near point 2.3739051999999998
  48.544292458332933 at 2.3739051999999998 48.544292458332933

Check geometries:

select count(*) from isochrones_intersection where st_isvalid(the_geom)=false;

count
-------
     0
(1 row)

Here a pg_dump of isochrones_intersectionlayer

I've tried several workarounds without success:

Reproject the layer to 3857.
ST_SnapToGrid to reduce precision.

Something similar appears at this post, but it looks they implemented a workaround which it's not useful for me.

I agree it's not an easy union.

EDIT

Adding versions:

PostgreSQL 10.1 on x86_64-pc-linux-gnu, compiled by gcc (GCC) 4.9.2 20150212 (Red Hat 4.9.2-6), 64-bit

POSTGIS="2.4.4 r16869" PGSQL="100" GEOS="3.5.0-CAPI-1.9.0 r4084" PROJ="Rel. 4.9.1, 04 March 2015" GDAL="GDAL 2.2.2, released 2017/09/15" LIBXML="2.9.4" LIBJSON="0.11" TOPOLOGY RASTER

I've also created a smaller version of the layer with the issue

Any ideas?


texts: Features shift with perspective change in Ol3-cesium Globe view

texts:
http://i.imgur.com/G7YWY1p.gif

In the gif above, the feature I'm talking about is a line string on the river in center of the image. As you can see, as I'm changing the perspective, the line string appears to be shifting from its place.
I'm uploading a geojson for the said feature:

    {
  "type": "FeatureCollection",
  "features": [
    {
      "type": "Feature",
      "geometry": {
        "type": "LineString",
        "coordinates": [
          [
            -87.46780363235473,
            41.64921533718473
          ],
          [
            -87.46604410324096,
            41.650626292892184
          ],
          [
            -87.46428457412719,
            41.65168448939255
          ],
          [
            -87.4596068016052,
            41.65473071559933
          ],
          [
            -87.45874849472045,
            41.65556439450464
          ],
          [
            -87.45282617721557,
            41.65970056487012
          ]
        ]
      }
    }
  ]
}

I'm using a very basic geojson here. Currently, there is no 3rd coordinate (or z-coordinate) for the feature. But the same issue occurs even when I supply a z-coordinate.

I'm not getting why the feature is shifting. Also, the feature seems to move below the surface of map. I'm kind of new to this thing so not quite understanding why this is happening.


texts: ST_Area of a buffer differs by coordinate

texts:
Why are these results different? The only difference is the coordinate. Seem like the area should be the same.

 select
   st_area(st_buffer(st_setsrid(st_point(-70, 40), 4326)::geography, 10)) AS ex1,
   st_area(st_buffer(st_setsrid(st_point(-71, 40), 4326)::geography, 10)) AS ex2;
       ex1        |       ex2
------------------+-----------------
 312.338323043769 | 312.17021877505
(1 row)


texts: Automating a dissolve and merge function using a lookup table with a postcode shapefile

texts:
I'm not any kind of QGIS power user.
I have some CSV data which includes postcodes, local authority names, and a reference label. I have access to the postcode shapefile for the UK.
What I need is a shapefile with one polygon for each of the reference labels, made up of the polygons of each postcode area that it is recorded against in my lookup table. Some postcodes could have two or more reference labels, so might appear multiple times in the table, so I expect to have overlapping polygons in my final output.
As a manual process, I can select all the postcode polygons for a particular reference label, save them as a new shapefile, dissolve them into a single shape, rinse and repeat for every label, and then save all those shapes as a single shapefile layer.
But sadly that is a huge task, way beyond the time or resources I have available.
Any ideas as to how I might achieve this with some clever piece of code?
(Also asking in QGIS subreddit)


texts:
I use bands 654 for vegetation classification and bands 432 to prepare the map layout.

Is there a way to pan-sharpen the landsat 8 image using the bands 65432 ? So I swicth the GRB to classify and to Layout.


texts: Data defined Style gives not same color in toc as on the canvas

texts:
QGIS 2.12.0 Lyon

I have an object with the following attributes:

shape [QString] ...circle, square
size  [integer] ...8, 10
color [QString]

With the widget value map I fill these fields according to csv files.
My color file looks like:

 rgb;description
'215,12,16';red
'51,36,215';blue

I did a rule-based style and on the map canvas it just looks like it should. But the colors in legend (TOC) are the predefined colors from the page "rule properties". What am I doing wrong?


texts: Georeferenced overlay not being displayed properly on Leaflet

texts:
I have an image of Europe and Asia with an unknown projection that I want to display on Leaflet. So I georeferenced it on QGIS with EPSG:3857.

The result looks pretty good. For example I get the following points/coordinates:

Northernmost point of France: 51.07, 2.53
Northernmost point of Hokkaido: 45.46, 141.89
Cape Town: -33.92, 18.46
Southernmost point of India: 8.09, 77.46

The result matches the world shapefile I have downloaded.

So, still from QGIS, I saved the output in a tif file. Which I then converted to PNG with:

gdal_translate -of PNG  -co worldfile=yes file.tif file.png

Here is what gdalinfo tells me about file.tif:

Size is 1269, 637
Coordinate System is:
PROJCS["unnamed",
    GEOGCS["unnamed ellipse",
        DATUM["unknown",
            SPHEROID["unnamed",6378137,0]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    PROJECTION["Mercator_2SP"],
    PARAMETER["standard_parallel_1",0],
    PARAMETER["central_meridian",0],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]]]
Origin = (-27.056867121684334,81.529534762914096)
Pixel Size = (0.217559535863945,-0.217559535863945)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  ( -27.0568671,  81.5295348) (  0d 0' 0.88"W,  0d 0' 2.64"N)
Lower Left  ( -27.0568671, -57.0558896) (  0d 0' 0.88"W,  0d 0' 1.85"S)
Upper Right (     249.026,      81.530) (  0d 0' 8.05"E,  0d 0' 2.64"N)
Lower Right (     249.026,     -57.056) (  0d 0' 8.05"E,  0d 0' 1.85"S)
Center      ( 110.9846584,  12.2368226) (  0d 0' 3.59"E,  0d 0' 0.40"N)
Band 1 Block=1269x2 Type=Byte, ColorInterp=Red
Band 2 Block=1269x2 Type=Byte, ColorInterp=Green
Band 3 Block=1269x2 Type=Byte, ColorInterp=Blue

Everything is looking good. So I added it on Leaflet with:

var bounds= L.latLngBounds({lat: 81.5295348, lng: -27.0568671}, {lat: -57.056, lng: 249.026});
L.imageOverlay('file.png', bounds, {}).addTo(map);

On the map, the corners are all at the right place, but between them, it's not projected as I was expecting and the overlay not match the EPSG:3857 projection of Leaflet.

Which step did I fail?

By the way here is what I have on Leaflet: http://i.imgur.com/2O3k4sh.jpg

And here is what I have on QGIS: https://i.imgur.com/nRXC9ft.png


texts: I need to create a proportional symbol map for a GIS class. I am wondering, how can I manually change the symbology to make the dots grow in size depending on their magnitude represented, and also how can I manually change the color scheme?


texts: Arranging elements horizontally in ArcGIS tool

texts: Is it possible to arrange multiple check boxes horizontally in a tool and changing the parameters dynamically (adding/removing) based on user selection of check boxes?


texts: Projecting as USA_Contiguous_Albers_Equal_Area_Conic_USGS in QGIS?

texts:
Is it possible to project EPSG: 6350 (NAD83 (2011)) as WKID: 102039 (USA_Contiguous_Albers_Equal_Area_Conic_USGS) in QGIS?

There is no EPSG code and I'm unsure how to create a custom CRS.


texts: Unable to export or edit XY data using ArcMap

texts:
I am using ArcMap.
I have a polygon shapefile. I added columns for the x and y of the centroid of each polygon using the Add Geometry Attributes tool. I then used the Make XY Event Layer tool to create a new point file of the centroid x,y coordinates. I also didn't need all of the records in that file so I selected the ones I needed and then used Create Layer from Selected Features.
I basically can't do anything further with this file and I don't know why. I realized there was a point that I could delete from my selected layer but I am unable to delete it. I need to move some of the points to make the map look better but keep getting an error that this layer is not editable. I thought maybe I needed to export it to a new shapefile, but when I try that, I get the following errors:

Error exporting data.  The Fields collection contained multiple
geometry fields.  The table was not found. (this line six times)

My layer has 41 records. There are no null values. There are no weird characters in the column names. I am attempting to export to the same default geodatabase that I export everything to. I get the same result whether I choose to use the same coordinate system as the source layer or the data frame.
How do I export and/or edit this file?


texts: Geocode area/neigborhood into polygon vertices

texts:
I have been successfully using the Google Maps API using the R package ggmap to find location coordinates for e.g. addresses or city names.

However, Google Maps now shows areas as well, indicating their extent with a transparent red polygon, e.g. when searching for "South Extension I, New Delhi", it shows the extent of this area/neigborhood (figure below).

The Google Maps API only seems to return the coordinates of the geometric center, and the viewport (rectangular) boundary coordinates (after zooming to the area extent): https://developers.google.com/maps/documentation/geocoding/intro

How can I obtain the coordinates of the vertices of the polygon which indicates the area/neighborhood?


texts: Variations in postal codes: canonical forms?

texts:
I was wondering, are there many variations that postal codes can take around the world?

For example, in Switzerland, every zip code is a four-digit number, but it can be written as XXXX or CH-XXXX.

I would basically like to have a function that makes a postal code canonical, given the country. For example f(XXXX, Switzerland) == f(CH-XXXX, Switzerland) == CH-XXXX.


texts:
I have to draw all the rows connecting the 2 edges with a line. It's not so hard to it manually, but it takes a lot of time. Is there a way to connect the two edges of the rows automatically?

For example, if i call the edge on the right side X and the left edge AX (X will be the row number) can i give some kind of command to connect them automatically?


texts: How can I pass user input from an HTML form to an arcpy script?

texts: I currently have an arcpy script that selects stores that are open after the user inputs their date and time. I'm looking to embed the results on a website using leaflet. If I have the user input their date and time using an HTML form, how can I pass the results to my arcpy script?


texts: Can I define a transformation between 2D and 3D in Geotools?

texts:
I'm using Geotools, and I need to convert coordinates between 2D and 3D coordinate reference systems.  I recognize that there's no universally-correct way to do this:  2D -> 3D means "making up" a Z coordinate, and 3D -> 2D means losing information.  But, is there some way to teach Geotools that, in my particular context, I want to use a specific transformation.  E.g.:

+-----+                 +-----+
|     | --- Z = 15 ---> |     |
| 2D  |                 |  3D |
|     | <-- drop Z ---- |     |
+-----+                 +-----+

I know I can transform any particular coordinate manually this way, but I'd like to  integrate it so this Just Works.  e.g. I can do CRS.findMathTranform(a, b) and if a is convertible to my 2D CRS, and my 3D CRS is convertible to b, it works.


texts: Tool to convert mapinfo wor to qgis project file with symbology?

texts:
Is there a tool to convert mapinfo wor to qgis project file with symbology?

It should be possible to create one as both wor/qgs are text based.


texts: Filtering out spatial-outliers from cluster of points using PyQGIS?

texts:
I have a cluster of points as seen below from a shapefile layer.

I need an automated procedure (in PyQGIS or any) of filtering out those points that are outside or further-away from the center of the clustered points, by creating a new attribute field for the shapefile layer.

The output would be the same shapefile but with an added field (named "Outlier"), where all the outliers points will be marked with a "1" and the rest with a "0".

That is the points in red circle will have value of "0" while those outside the red mark will have the value of "1".

Code Edit:

from PyQt4.QtCore import QVariant

point_layer = iface.activeLayer()
data_provider = point_layer.dataProvider()

count_ft = data_provider.featureCount()

# Create Outlier attribute field
data_provider.addAttributes([QgsField("Outliers",  QVariant.Int)])
point_layer.updateFields()

# Get points average coordinates
features = point_layer.getFeatures()
for f in features:
    #calculate average xy coordinate
    pass

# Buffer from average xy coordinate to get outliers

# update "outlier" attribute field
## if point in buffer, update with 0 else update 1


texts: How to convert xml into readable gis format


texts:
Is there an ArcGIS function to sample x% where x lies between 0 and 100, of a raster?

I do not want point sampling and ArcGIS functions seem to focus on point sampling.

For ArcGIS, I use version 10.2

--EDIT: Desired output: This output raster will have x% of cells selected randomly


texts: Problem using field calculator to add numbers in new field with preceding zeros

texts:
I wish to add a field in an attribute table that should work as a georeferencing key. It needs to take into consideration values from a series of fields. The last field is a housenumber and it contains from one to five digits (e.g. 2, 20, 20A, 200A, 200AB etc.) What shall I write in Field calculator to ensure that the above housenumbers becomes the following values in my new field:
002
020
020A
200A
200AB

So that the numbers are always represented with 3 digits (containing preceding zeros if necessary), and possible subsequent letters are kept?


texts: I am having a shape file which consists of various fields for example: shape area ,shape length, locality, sub-locality,district etc.In my company they are demanding seprate shapefiles say for example sub locality shape file must have only the attributes of sub locality not any thing else like shape area, shape length etc. this is the actual problem kindly help me in solving this?


texts: failed to load dom raster via raster2pgsql

texts:
I want to load a raster file (dom) to postgis via raster2pgsql. I used C:\Programme\postgresql\11\bin>raster2pgsql -s 25832 -I -C –b -t auto –R D:\gis\dom.tif public.dom_test| psql -h IP -d database -p 5432 -U AdminUser
It takes a bit of time but the raster get loaded and can be seen in the database. However, when I try to load it via QGIS browser (3.12) per drag and drop a get following message:

CRITICAL    Invalid Layer : Raster layer Cannot instantiate the 'postgresraster' data provider

When I use DB Manager I get the message:

No primary key defined for this table!

But in the pgadmin I see the primary key is set for the column "rid".
Any ideas what is the reason?


texts: In an ESRI Shapefile, what is the "Record Number" field of each record header in the .shp file?

texts:
The shapefile specification (http://www.esri.com/library/whitepapers/pdfs/shapefile.pdf) only says "Record numbers begin at 1." (page 5), which doesn't specify what the field actually is. As per page 25, each shape feature must correlate to exactly 1 record in the DBF, and their order must be identical. So the "Record Number" field is not for any purpose related to the DBF. My interpretation is that the "Record Number" field is required to count up from 1 to the number of shape features in the shape file. But this being completely redundant, I am unsure whether I missed something?

Some clarifications after reading your answers and comments:

I am writing a reader/writer for shapefiles, so I need to parse them correctly as well as produce valid output.
Shapefiles as well as DBF files have an "intrinsic" record number (let's call them A and B). The first shape / record after the header in the shapefile / DBF has A=1 / B=1, the next has A=2 / B=2, and so on. The shapefile specification clearly states that the shapes and records have a one-to-one relationship based on these numbers, so the attributes for shape A=1 are contained in record B=1, and so on. Neither A nor B is the "Record Number" field I am asking about.
Each shape in the shapefile also has an (explicit) INT32 field in its shape header, which is called "Record Number" in the shapefile specification. This is where my confusion starts. Let's call this value C. The specification requires that this field starts at 1 (so for shape A=1, it must be C=1). But it doesn't give any other rules (such as "increment C by 1 for each shape"), nor does it give any semantic meaning to this field. So it seems to me that C is indeed meaningless. But if that is true, the shapefile is bloated by 4 unneeded bytes for each shape, which seems weird to me. So I want to be sure I am not missing anything...
I have no idea where things like "FID, OID, OBJID, or OBJECTID" come from. I am not using any ESRI products, but glancing at http://support.esri.com/es/knowledgebase/techarticles/detail/37480 I conclude that those are properties / attributes some ESRI products use and store in the DBF file associated with the shapefile. If that is so, this has nothing to do with my question.


texts: Create kmz file out of Illustrator maps?

texts:
I'm migrating my ERP from an old-style mapping system which relied on jpg files and coordinates on the file to locate my assets, to a vectored system which uses KMZ files.

I'm new to this - I'm here sitting with my illustrator files and I'm not sure where to start.


texts: Using a different projection to calculate areas using the Field Calculator

texts:
I have a dataset that is in a WGS84_UTM projected coordinate system, can I force the Field Calculator to use another coordinate system when I use it to calculate areas using pythons !shape.area@squareyards! ?

My question is similar to checking "Use coordinate system of the data frame" in calculate geometry when the data frames projection is different; I would use this method but it does not have the option to calculate Square Yards. I could reproject the data to a new feature class, then calculate it and perform a join, or calculate the area in square feet using the data frames projection and dividing by 9, but since I will be doing this often the field calculator route would be the best way to go.


texts:
I'm newbie to Android Offline Map Integration, I have a problem in offline map integration, I take a month to analysis the offline map integration.  Here I referred to Osmdroid Osmand and few other programs, I have an idea about map data in vector format, but till now I don't have initial level progress.

Please can you to give some reference code for Map Integration using vector format and POI data Integration.


texts:
I have a shapefile with lines (these are administrative borders in India). I want to split each line into a max length of 10 kilometres.

I am using the GRASS tool, v.split.length, following advice here: How to split a line into a specific number of parts?

However, I get a really large number of segments; whatever change I make to the 'Maximum segment length' option seems to have no effect.

One mistake could be using the wrong projection: I was earlier using the coordinate projection EPSG::4326. Is that it, if I use an equidistant projection for India, will that allow me to create 10 km segments?

Which equidistant projection can I use?


texts: Data Security for Geoserver - Users errantly seeing Workspaces not assigned to them

texts:
I'm having trouble with Data Security in my Geoserver 2.5 application.  Let me preface my situation. I am using Geoexplorer composer which is making GeoServer API calls through ows and wms.  I have multiple clients which should not be allowed to see each other's wms services. The implementation of my data security rules are scatter shot when ows and wms are used through GeoExplorer.  CompanyAUser can see all layers in the CompanyAWS and 1 layer from CompanyB's ws.  CompanyBUser can see 10 out of the 60 layers in CompanyA's ws.  A further twist to this problem is that if I sign in as each company's user in Geoserver's webapp and I navigate to Layer Preview each user can only see the layers which I have specified that they are allowed to see in Data Security.  I have made some edits to the Geoexplorer Composer code but nothing major.  I don't think that this should effect the security in any way since security is implemented by the Geoserver application.  I have included my security parameters below.

Groups:
group, role
CompanyA, COMPANY_A
CompanyB, COMPANY_B
CompanyC, COMPANY_C

Users:
user, group
CompanyAUser, CompanyA
CompanyBUser, CompanyB
CompanyCUser, CompanyC

Data rules:
Roles, rule path
Company_A, COMPANYAWS.*.r
Company_B, COMPANYBWS.*.r
Company_C, COMPANYCWS.*.r

Service rules:
*.*     ADMIN,COMPANY_A,COMPANY_B,COMPANY_C

Roles are independent and do not have parent/child relationships.

Are my roles, groups, data rules and service rules set correctly?

EDIT1
I've noticed that the only layers that are errantly shared to users belonging to different company groups are 'Group Layers'.  This leads me to my next question. What causes a regular layer to be hidden and a group layer to be shown to roles that do not have permission to do so?


texts: Adding two vector shapefiles together (one on top of the other) in QGIS

texts:
I understand how to merge two vector shapefiles together but I don't know how to layer one layer on top of another one and export it to make it permanent?

I would like to add the GreaterMan shapefile on top of the GreatBritain shapefile to highlight the location of GreaterManchester (as shown in the attached picture).

I know this is a quick fix but I can't work out how to do this.

I have merged the two shapefiles together using the union tool so that only the Greater Manchester in the attribute table is displayed and I have run a simple select features using an expression "DISCRIPTIO" LIKE 'Metropolitan District' to highlight all of the Greater Manchester districts.

The purpose is to export the shapefile, so I can import the shapefile in Inkscape.


texts: Dynamically Change Leaflet PopUp Background

texts:
Is there a way to set the background image of a leaflet pop up dynamically?

I have multiple markers with popups, which all need a different background image. The filepaths for the images are being pulled from a database.

I have tried using jQuery to change the background property of the ".custom-popup .leaflet-popup-content-wrapper" classes, but the change does not occur. When I type $(".leaflet-popup-content-wrapper") into the console, it returns an empty array which tells me that the current page cannot even see this class.

How can I dynamically apply a background image to each of the popups?


texts: ArcGIS: Assigning values to missing attribute data for line layers

texts:
I have a large road network layer with an AGE attribute associated per segment. In some occasions, I would have segments A -> B -> C with segment B with no value for the AGE field. In other situations I would have segments D-> E with either segments D or E with no AGE value.

In the first case, the age of B could be assumed as either the average or max of A and C and in the other case assign D or E the age of the other segment. Is this possible in ArcGIS?

Note: Thought must be given to the plausible situation, which I haven't seen in the data yet, of the following case G <-> H <-> ... <-> I <-> K where H, I, and values in between have no AGE values.


texts:
I have over 140 layers that need to appear as overlays on basemaps(Google and OSM). My users can pick and choose which layers they want to see and store them in their preferences. While we have built the solution, the performance is very slow. These layers get served with openlayers making WMS calls to GeoServer and Postgres DB.

I tried using GeoWebCache and Tilecache but then realized that they are not all that helpful when I have dynamic layers. The XML definitions of what layers to cache are static.

What are my options here? Is there a way to cache such combinations, or do I take a different route and investigate on HTTP caching with tools like Varnish?


texts:
I'm wondering how can I get lat and lon out of the function in a variable with global scale in Leaflet. The code to get the coordinates per click in the map works fine but then I can't write my lat and lon in global variables. My code which doesn't work look like this:

var lat;
var lon;

function clickEvent(e) {

        if (typeof(newMarker)==='undefined'){
            newMarker = new L.marker(e.latlng, { draggable: true });
            newMarker.addTo(satMap);}
        else { newMarker.setLatLng(e.latlng);}

        lat = e.latlng.lat;
        lon = e.latlng.lng;
}

satMap.on('click', clickEvent);
document.write(lat);


texts: gdal2tiles format of argument for the '--srcnodata' flag

texts:
When you run gdal2tiles, by default it sets pixels outside your map area to black-transparent.
I would like to set them to white opaque.
Looking at the documentation...
http://www.gdal.org/gdal2tiles.html

it looks like the setting --srcnodata=NODATA might do it,
but it doesnt explain what the allowed values for the 'NODATA' argument are.
I have made various attempts but they had no effect e.g...

gdal2tiles -r near --srcnodata=255 -z12-12 in.tiff outDirectory


texts:
I'm trying to add the attribute $area. I made a intersection and Tot_aera isn't the good one, so I want to add a new.

I made that :

from PyQt4.QtCore import *
try:
    if ind==-1:
        if caps & QgsVectorDataProvider.AddAttributes:
            res=provider.addAttributes( [ QgsField(nameattr,QVariant.Double) ] )
except:
    print "FALSE"

So I have a new attribute but NULL. So my question : How put "$area" in nameattr?