texts: Link QGIS CSV / Excel with raster layers

texts:
I have the following problem. I want to link a value from an Excel file to a raster layer. Uploading the Excel file and linking to other vector layers via properties -> linking is no problem.

However, the ID values from the Excel file should be linked to the same ID values from the grid layer in order to have a region number in the grid, depending on the ID values. Is there a simple tool for this or is it only possible via detours?


texts: Splitting a point layer in equal parts in QGIS

texts:
I need to split a point layer in equal parts, let´s say I have 3252 and want to split in 15 equal groups of points (average is almost 217 points).
I am using QGIS 3.10 in Windows.


texts: How to visualize this country data in Qgis?

texts:
I have a csv table, with data in this format: "name","country". Country code is in the same format as the world map SHP file attribute table (both are in the ISO standard). I learned how to join the tables together, but I want to colorize countries based on how many people are from each country.

Example:

"john","LVA"
"peter","OMN"
"martin","PRI"
"janis","LVA"

In this case Latvia should be darker than the two other countries (because two people are from LVA). I know how to colorize things based on a number that's bigger than other numbers, but how to achieve it without the numbers in the table? I only have raw list of names. Is this possible with Qgis only?


texts: Can GeoNetwork store metadata for a PostGIS table?

texts:
Can GeoNetwork be used to store Metadata for a PostGIS table? Or is it only for WMS/WFS services?

If it can be used to create a store metadata for a PostGIS table, how is it done? I can't seem to find this in the GeoNetwork documentation.

Thanks.

DBlack


texts:
I have spatial data that has an existing related table with photo info.
I have data from another separate table that I want to join and then field calculate to add to my spatial data, then I would remove the joined table.

I checked Esri's site for information and here is advice they offered.

"Relate tableB to layerC, then join tableA to layerC: This scenario does not work. This is different from the previous scenario in that layerC owns the relate. Therefore, when tableA is joined to layerC, the relate is removed.

Relates owned by the join table are not affected by the join. These relates can't be accessed by the target table or layer, however.
 Relationship classes are never removed as a result of a join or unjoin operation. They can be used regardless of whether the data has been involved in joins or relates."

Will my relate be removed? I do not know if my relate is owned by "layerC" or TableB.


texts:
I have a point layer with singular points set with high accuracy to the position of certain sites, each point representing one site. There are 8 different types of site. When representing these points on a smaller scale map, the image becomes confusing, however.

I am looking for a way of clustering points of the same type into bigger points for representation on small scale maps. Therefore the size of the new clusterpoints should correspond to the number of points cluster within them.

There should be a workaround using buffers and spatial joins, but I assume, that there is a more elegant sollution?

I am using QGIS 2.18 on a Windows 7 PC.


texts: Training for GIS professional working in Telecom?

texts:
This question is not about any code or bugs. This is about GIS profession. I am GIS professional working with MapInfo, OpenLayers, SQL Server, Oracle, and System integration.

I learned everything from stackoverflow.com and gis.stackexchange and google.com. Now I want to have some training in this field related to telcom and GIS.

What sort of training should I take that can be applicable to development of telecom and help use vital data to be presented as useful information?

As for any technical jobs the requirement changes and need for update is essential. I want to be expert in the field of GIS related to Telecom with current simple knowledge and without training I don't think it's possible.


texts: Several errors compiling QField App with Qt Creator

texts:
I tried compiling QField on Windows for a week now. But I still have several compiler errors and now I don't know how to fix them. I followed the build manual of the Github repository. So I installed Qt 5.15.2 with Qt Creator. Next I downloaded the latest nightly build qgis-dev with OSGeo4W. In Qt Creator I opened the CMakeLists.txt and added the variables for finding GDAL, QGIS and Proj includes and libs. I use the Desktop Qt 5.15.2 MSVC2019 64 bit kit. With this setting the project configures correctly. When I hit build I get the following errors:

Is this caused by false Include-Files? Or is it maybe a false compiler?
Do you guys have any tips for me?


texts: Converting encoded polylines into GeoJSON and KML files?

texts:
There are plenty of resources available online to help users encode polylines. I am looking for something that does the opposite but my extensive research yields minimal results. I have a massive set of encoded polylines that I'd like to convert into geojson files, and eventually into kml files so I can import into qgis for further analysis. I'm not a coder (at all) so will not be able to complete this tasks using python or javascript.

An example of an encoded polyline:

_{reFxmjVpAeBn@{@w@aAc@m@yAsBJuAC}BIeAO_Ak@sBCQOg@g@aAk@}@qAwAs@s@}I{JgCqCsC_D}DeEyDgEsA}AwB{BcCwCoGyGqBwByEwFgHqHgDqDkIaJyFuGeEsEeEcEcGkGkEaFyDyDgBqBcAgAm@y@w@gAuCcCmFeHk@i@eAgAaDaDcC_CoBaCeBkCwAiCkCuFoPy]eAkCy@eC{@eDq@kDk@cEwAyQiBsVkAsN}Egn@q@sJ_@iGWgFg@mMcAuWQcCg@yFkBkRc@{EiAmKuBiT]sDWiEy@kPg@_Gg@_E{AuOkAgKo@gFe@oDq@oEe@qBmAeEWkAOaAKoAC}ADiBPeBFc@\wA~@oC@sAh@CVgBJoA@c@EWBkAAcBGgB[yC}AuNu@iHG}B@uBJcBJmAjEeZbAuGz@uFh@}DV]H[RkAfAqGJe@l@mBl@eAv@w@n@]TIj@Ml@In@?Df@n@Fh@@|H|AzDj@zDl@l@F^KhAZhCp@pAXnCr@bBh@nBbAfCdB|@z@vCfDPTjAvAhArAjBjBlBdBlEnC\TCzAtEpCtDtBrEdCdExB~DpBnAj@z@XdAPtBLr@CH?pAQZKfBs@t@e@x@s@|@cAt@mAh@gA@mAlAmEFELSpAyD|A_GfAeC|B}HT_At@gDjCaJZq@dB{F~AeFmCuAgCyAv@mCd@}AfA}C~@oBr@iA|BcDn@gA~@iBTo@~@qCh@gCXoB??P{AP_BU_@e@w@d@v@nArBl@v@ClBlBtAdATnAXLCFAXiALc@Ng@X}A^eCvB{T\qDZsBSeATmBr@mE?@ASVkAv@iDlBqHvA_FhAcDv@qBtBaFbIyQpBsDrAoBdAiAp@k@r@g@v@a@lAi@dBk@dG{AfD}@v@YzAq@|BsAfA{@[|A_BvAiBv@oAAiB@A@cAf@_BzBaJrM}j@x@eCz@mBbAgBlAeBhAmArAmAhCkB|MeJnD_CjByAxE}DtG{FlNsL~WaUjMmK|GyFjR}OdI_HhYcV~SoQfCwBnGeF~I}HbHcGnDyCDmCbJyHhA}@~FeF~GwFnNmLzCkCtJkJ|MeNDeDjGkGxGmGtAyApHcHrBiBlAsAdF_GRyZhI{MlJcPBeC~PaYdPaXjOqVbLcRdXmc@dDqFhLgRtGmKj@aAlB}Cz@oAR]tAeC~EqHtWwb@fRc[bDgGzB}Df@{@fIuLxGsKzBwCtCyCbB{AfAw@bC}AlGgDbEsBdIgEl@]fNkH~VcN~f@wWzP_Jr@c@zIeFAo@xEwCbHmEdGyDnIsF~MyIrJmGrUaOlSmMzDgC~k@e@bf@{Z~n@ca@|QqLtRcMdp@ab@naAkn@fBeADaBfDwANBJ?nAWhCm@^Gf@@XHZPp@f@xDfCpF}NpAcDpA_DjF}MtGyPLY]UyDkCoFqDaAC_@z@aApBoAzC]W??\Vt@iBd@gAdAuBNa@r@gBz@h@~B~AdBhAhChBZTLYnCeHtAmDnHgRjHcRtCqH~FcOvBuFdEyKpLoZvL_[v@mBjH}QXaAfAkCx@sBXaAlEwKpA}CfDgIx@mB^y@FE@Ap@uAfCaGrIsSA_CxAeDnGyOlFiMdBeECcGhG_O??AaCnDrC~DvCrNtKz_@~YfUfQtFhE~DD|@p@\Pn@RbATjBBpABzANv@ThAh@|@n@@^ZXf@p@NfA@t@KXSV_@LWA[KU]Ie@?WJk@Le@hA}CBK@g@bAwBlCkGzFoN~KqYpOk@zF_OzAaDfB_DpCmEzJwMfPqTT}XjLoOvJwMh@q@vK}NjEwFbDsErKiOvWi_@d\{d@V_@r@cA~@qAb@u@xB{CrIyLhIoLxKsOnG_JjB}Bp@u@fC_CxCyBnC{A|D}AfFuBrOeGG}Bfh@mShc@}PdIeDhG{BbE{AjD_AzDw@pFq@zJiA|@MxAKbNcBfJaAtOkBbKiA~TwCzFo@nC]bKaAhGm@~Hi@d@AbPaA~AOpCe@pEkAdA]zEkBj@UfB}@|Aw@xBs@|Ac@xJ_Cf@Dr@GrBa@fBYdBQvDQv@MlA]pG{BrAc@xCu@~Cm@zCg@|@WSqAOeA_@oDMqDT}`@JaC^qBZcAZw@xD}HvAsCpBiE~BkGJYg@c@eCaBO^i@tAc@dAo@~A

Is there an online converter somewhere that will help accomplish this task?


texts: Outputting to multiple attribute fields after Spatial join?

texts:
I'm totally newb to GIS.

I'm trying to perform a spatial join of two polygon layers, and want the results to propagate into two separate attribute fields.  Specifically, I have one polygon layer of municipalities in New York state, and another of electric utility service areas.  Some of those municipalities are bisected by two utilities.  I want to create an attribute table (exportable to csv) that lists each municipality, with one attribute field listing one utility, and the next attribute field listing the second one (if there is one).  (I.e. Bedford / NYSGas&Electric / ConEdison).

When I perform a spatial join now (using QGIS on a Mac), using the spatialjoin plug-in and doing an intersect join type, the output only lists one of the utilities when a municipality is bisected (i.e. Bedford / NYStateGas&Electric).

Any idea how I'd get the results I'm hoping for?


texts: Creating vector gradiants showing magnitude and direction by taking input of two rasters in QGIS

texts:
I have two binary rasters and I am focusing on pixels having value 1. Now, I want to know show direction and intensity from which pixels (having value 1) of one raster is shifted other raster. Basically it is land change study and I want to present it with arrows.
I have tried a tool named 'Gradient Vector from Direction and Length', which is not giving desired output because I don't have particular rasters to give input as length and direction.
Attaching raster images here

Here, I want find red pixels which is one raster, shifts to blue pixels which is another raster.


texts: Translate custom MKT projection to proj4 string?

texts:
I need to use a specific custom projection that not exist as SRID into standard spatial_ref_sys of my GIS,

PROJCS["Conica_Equivalente_de_Albers_Brasil",
    GEOGCS["GCS_SIRGAS2000",
        DATUM["D_SIRGAS2000",
            SPHEROID["Geodetic_Reference_System_of_1980",6378137,298.2572221009113]],
        PRIMEM["Greenwich",0],
        UNIT["Degree",0.017453292519943295]],
    PROJECTION["Albers"],
    PARAMETER["standard_parallel_1",-2],
    PARAMETER["standard_parallel_2",-22],
    PARAMETER["latitude_of_origin",-12],
    PARAMETER["central_meridian",-54],
    PARAMETER["false_easting",5000000],
    PARAMETER["false_northing",10000000],
    UNIT["Meter",1]]

so I created new arbitrary SRID for it,

INSERT into spatial_ref_sys (srid, auth_name, auth_srid, proj4text, srtext)
values ( 952019, 'BR:IBGE', 52019, '', 'PROJCS["Conica_Equivalente_de_Albers_Brasil",      GEOGCS["GCS_SIRGAS2000",          DATUM["D_SIRGAS2000",              SPHEROID["Geodetic_Reference_System_of_1980",6378137,298.2572221009113]],          PRIMEM["Greenwich",0],          UNIT["Degree",0.017453292519943295]],      PROJECTION["Albers"],      PARAMETER["standard_parallel_1",-2],      PARAMETER["standard_parallel_2",-22],      PARAMETER["latitude_of_origin",-12],      PARAMETER["central_meridian",-54],      PARAMETER["false_easting",5000000],      PARAMETER["false_northing",10000000],      UNIT["Meter",1]]');

but ST_Transform(ST_SetSRID(geom,4674),952019)  is an error:

 ERROR:  AddToPROJ4SRSCache: could not parse proj4 string ''

So, seems I need to insert also the proj4text column... I need to parse  the MKT definition of srtext, is it?  How to do it?

PS: I am using PostGIS 2 but the question seems generic, about simple parsing procedure.


texts: Set NoData values to cells with a value and vice versa?

texts:
I am working on a ModelBuilder model that calculates natural habitat based on different attributes such as proximity to water, landcovers, elevation (all raster data). I also want to take into account "human disturbance" - i.e. exclude areas that are too close to a developed area.

To do that, I took the landcover raster, used the Con tool to assign all "developed" cells value of 1 and a NoData value to all other cells. Then, I used the Euclidean Distance tool to create a "buffer" around those "developed" cells. So, now I have a raster that has cells with regular values in the "buffered" zone around developed areas and NoData further from developed areas and I need the exact opposite - NoData around the developed areas and regular values further from them.

I am using ArcMap 10.4.1


texts: What is the difference between `Area_2D`, `Area_3D`, and `Volume` in surface volume calculation?

texts:
I want to know the meaning of output text files of surface volume dialog box. For example, the file contains headers and numbers as below:

Dataset, Plane_Height, Reference, Z_Factor, Area_2D, Area_3D, Volume.
raster_to_tin, 50,00, BELOW, 1,000000, 385387,1505728, 414338,45311207, 4544086,3828428

What is the difference between Area_2D, Area_3D, and Volume?


texts: For normalizing las files in a LAScatalog, is opt_filter() only used in the processing?

texts:
I'm using

opt_filter(ctg) <- "-drop_overlap -drop_z_below 0"`
ctg_norm <- lasnormalize(ctg, algorithm = knnidw(k = 10, p = 2), na.rm = FALSE, use_class = c(2L, 9L))

Are overlapping and negative values disregarded only for normalizing purposes, but still written to output las files? (i.e. I am getting negative values in my normalized output, but I am also getting fewer points. Wondering where those missing points went is what lead me to the question in the first place).


texts: USGS Bulk Downloader cannot locate install4j

texts:
Trying to run the the USGS' Bulk Download Application (https://www.usgs.gov/media/images/earthexplorer-bulk-download-application-bda) and keep getting an error message that it could not find the JRE. I choose the options presented to locate it and  navigate to my JAVA folder and the java.exe file. Error persists. Java is installed.
any ideas?
previous to this I set the JAVA_HOME and another related environment variable to the java jre folder.


texts: Problem when using Simplify Tool in QGIS

texts:
I want to simplify district borders. As I got bad results when using polygons I converted those polygons into lines but I still have the same problem that I get new vertexes or corners I don't want.

Neither the v.clean nor the v.generalize tool will work for me to solve the problem.

This is my original border:

The output looks like that :


texts:
What does this error message mean:

  AttributeError: 'QLineEdit' object has no attribute 'toPlainText'

I get this every time I try to run the QPackage Plugin.


texts: Finding problematic rows - spatial join in GeoPandas?

texts:
I have an sjoin function from geopandas that is behaving erratically: it works on some version of the "points" geodataframe but not others.

merged=sjoin(points,polygons, how='left',op='within')

The error I get is always:

  rtree.core.RTreeError: Coordinates must be in the form (minx, miny, maxx, maxy) or (x, y) for 2D indexes

The "polygons" geodataframe never changes. The size of the "points" geodataframe depends on how much data I want to include (in a parameter). Generally the join fails when I include more data (e.g. 100,000 rows), and succeeds on smaller datasets (e.g. 2,000 rows). I assume this is because some rows contain invalid data. However on visual inspection I cannot find anything wrong with any row.

Is there a way to quickly find out which rows are blocking the join, or to automatically ignore them?

I can't easily share the full code and data.


texts: Conversion of files to Mapsforge *.map format

texts:
I am a complete newbie and not a cartographer. I organize off-road events such as 4X4 trail rides, biking events and running events.  Many of these trails/event maps have been created in different file types. To track trail leaders and participants I am currently using APSdroid This software uses the Mapsforge *.map file format.

What are my options and what file formats can be converted to the *.map file format?


texts: Finding the percentage of change of 30 meter resolution raster within buffer area

texts:
How could I find the percentage of change for each buffer area in the image below using ArcGIS Pro 2.5.0?
I want to be able to show the percentage of change for each INDIVIDUAL and "unique" buffer/circle so that I can create a histogram or some type of similar graph. The blue pins are the "Point Locations", the general circles are 250-meter buffers from the point locations, the yellow cells are equal to change, and the dark black cells are equal to no change:


texts: When I add a feature, the feature attribute box comes up but it has nothing in it...completely blank, no fields. It used to come up but stopped. I have re-installed it but the same thing occurs.


texts: Selecting largest feature within each feature using ArcGIS Desktop?

texts:
I need to select the largest feature within another feature using Arcmap 10.3.

I have two sets of data, building footprints and land parcels. Generally there is a many to one relationship where there are many building footprints within one land parcel.

I would like to select the building footprint with the largest floor area in each land parcel. Resulting in a one to one relationship.

As I am working with 7000 building footprints and 5000 land parcels, I would like to find a way to automate this however I do not have a lot of knowledge on scripting. I can use ModelBuilder.

Below is an example of the data. Building footprint pink and land parcels black:


texts: Creating symbol with size and color based on two variables in ArcMap?

texts:
In ArcMap 10.2, I have a layer that has two attributes: depth and magnitude.

Depth should be represented by color: deeper points being lighter and shallower being darker.

Magnitude represented by symbol's size: points with a higher magnitude have a bigger symbol and points with a lower magnitude have a smaller one.

How do i do that?

So, many people has told me to use "multiple attributes" in the symbology dialog. I wish i could upload a screenshot but my firewall seems to have the image service blocked. Let's try with words.

"Multiple attributes" offer these input boxes: "value fields" and "variation by".
in "value fields" i select "depth" and "magnitude"
in "variation by" we have two buttons: "color ramp" and "symbol size"
in "color ramp" i choose "depth"
i open "symbol size". There, "depth" is selected. I change it to be "magnitude".
open "color ramp" again. No longer "depth" is selected. Now "magnitude" is.

Basically, whenever i select something on "symbol size" or "color ramp", the same selection is duplicated on the other button.


texts:
I use next code to get tiles that were rendered by Cesium engine:

    // cesium initialization
    var viewer = new Cesium.Viewer('cesiumContainer', {
            sceneMode : Cesium.SceneMode.COLUMBUS_VIEW,
            imageryProvider: new Cesium.TileCoordinatesImageryProvider({
                tilingScheme : new Cesium.WebMercatorTilingScheme()
            }),
            baseLayerPicker : false,
            terrainProvider: new Cesium.EllipsoidTerrainProvider({
                tilingScheme : new Cesium.WebMercatorTilingScheme()
            })
        });

    // get rendered tiles
    viewer.camera.moveEnd.addEventListener(function() {
        var tilesToRender = viewer.scene.globe._surface._tilesToRender;
    });

the problem is that tiles that I get are not the same as tiles I see on the Cesium globe. Cesium globe picture looks like this:

but tilesToRender variable contains only level 2 tiles and also tiles that are not on the globe (for instance L2X0Y0).

tilesToRender[0]: L2X1Y1
tilesToRender[1]: L2X1Y2
tilesToRender[2]: L2X2Y1
tilesToRender[3]: L2X2Y2
tilesToRender[4]: L2X1Y0
tilesToRender[5]: L2X0Y1
tilesToRender[6]: L2X2Y0
tilesToRender[7]: L2X0Y2
tilesToRender[8]: L2X1Y3
tilesToRender[9]: L2X0Y0
tilesToRender[10]: L2X2Y3
tilesToRender[11]: L2X0Y3
tilesToRender[12]: L2X3Y1
tilesToRender[13]: L2X3Y2
tilesToRender[14]: L2X3Y0
tilesToRender[15]: L2X3Y3

Is it possible to get exactly same tiles that are visible on the Cesium globe?


texts:
Looking at GXP examples here they do not list the dependencies within the JavaScript file, so I was wondering how to make the files work. I would assume just copy the JavaScript file, e.g. viewer into my SDK app.js file and build that. However it just loads a white screen when I load it.

Any hints on how to sue these examples?


texts:
I am trying to clip my Raster file according to a shapefile geometry. I am using the following code

gdalwarp -cutline INPUT.shp INPUT.tif OUTPUT.tif

But it is resulting a black color outside the shapefile geometry extent. I am giving the example here. In the first Image that I want to clip. The second image is the resultant raster but I want the 3rd Raster where the outside of the geometry will be null.

1st Image: INPUT Raster

2nd Image: OUTPUT Raster

3rd Image: I want to get output like this


texts: Difference between GWR in GTWR addin and GWR tool in ArcGIS

texts: I implemented geographically-weighted-regression (GWR) by GWR tool in ArcGIS and bandwidth result is 20307 and next again implement GWR by GTWR addin in ArcGIS and bandwidth result is 0.188 what is difference between them?


texts: Clustering spatial lines in R by direction

texts:
There are tens of ways to cluster spatial points in R but I couldn't find any for spatial lines. What I would like to get is the average trajectory directions that each of these cyclones (spatial lines in my case) has over the territory of France similar to the possibility that trajCluster function from openAir library offers.

library(dplyr)
library(maptools)
library(geosphere)
library(rgeos)

#load the world boundaries polygon
data(wrld_simpl)
contour <- wrld_simpl[wrld_simpl$NAME =='France',] #only the polygon for France

#download the data
x <- read.table("https://forms.naturwissenschaften.ch/imilast/_ERAinterim_1.5_1979_MTEX/ERAinterim_1.5_NH_M02_19790101_20121231_MTEX.txt?_ga=2.18919096.1825595846.1546710263-1112023567.1546710263", sep="", fill = T, nrows = 20000,
                header = F, skip = 2) # download only the first 20000 rows
colnames(x) <- c("Code","CycloneNo","StepNo","DateI10","Year","Month","Day","Time","LongE","LatN","Intensity1","Intensity2","Intensity3")

x <- x %>% filter(Code!=90) ##Just cleaning to remove those 'NA' rows
x <- x[x$LongE > -120 & x$LongE <= 120,] #subset only the traj between -120 and 120 longE
x$DateI10 <-as.POSIXct(paste(x$Year, x$Month, x$Day, x$Time, sep="-"), "%Y-%m-%d-%H", tz = 'GMT')

# create an empty Spatial Line to be filled later
sline <- SpatialLines(LinesList = list(Lines(Line(matrix(0, ncol = 2)), ID = NA)))
proj4string(sline) <- "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
sline$time <- NA
sline$start <-NA
sline$distance <- NA
sline$mean <- NA

#calculate the lifetime for each track; keep only trajectories for France
for (k in unique(x$CycloneNo)){
  message(which(unique(x$CycloneNo) == k), " from ", length(unique(x$CycloneNo)))
  txt2 <- x[x$CycloneNo==k,]
  if (nrow(txt2) > 1){
    d <- SpatialLines(list(Lines(Line(cbind(txt2$LongE, txt2$LatN)), ID= k)))
    proj4string(d) <- "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
    d$time <- difftime(max(txt2$DateI10), min(txt2$DateI10), units = 'hours')
    d$start <-  as.character(min(txt2$DateI10)) # get the formation date
    d$distance <- lengthLine(d)/1000 # calculate distance
    d$mean <- mean(txt2$Intensity1) # average the sea level surface pressure
    if (gIntersects(d, contour, byid=TRUE) == TRUE){
      sline <- rbind(sline, d) #keep only trajectories crossing the France's territory
    }
  }
}

spplot(sline['time'],sp.layout = list(contour, lwd = 3),
scales  = list(draw=T), main = 'Mean lifetime (h)') #plot the tracks

### what I would like to get is something like
library(openair)

traj <- importTraj(site = "paris", year = 2009)

## calculate clusters
traj <- trajCluster(traj, method = "Angle", type = "season",n.clusters = 4)

Is it possible to get something like this using my data? is clustering method appropriate here or I should some else (like PCA)?


texts: Using arcpy.Describe Method on Annotation layers?

texts:
I'm testing layers in the current document using the Describe function. My code is:

for lyr in arcpy.mapping.ListLayers(df):
     if(lyr.isGroupLayer == False and lyr.isRasterLayer == False):
         x = arcpy.Describe(str(lyr))

It's throwing an error on the annotation layer:

  "LSDSectionBoundaryAnno\Default" does not exist

I'm thinking that the Describe function doesn't like annotation layers or the subclasses? Can I use Describe on annotation layers and if not, how would I filter out annotation layers? There doesn't seem to be a 'isAnnotationLayer' property.


texts: Assigning attributes from hydrologic sinks to upstream features?

texts:
I have a large hydrology geometric network with more than 300 sinks and roughly 150 barriers.  I'd like to have the SinkID added to the attributes of the barriers (some individual networks have a single barrier, some have multiple) for later analysis but haven't seen any obvious way of doing so.  I could assign the SinkID to the first stream segment to which it's attached but some barriers may be 10 or more stream segments upstream.

Thoughts?


texts: Clean way to convert from US feet to meters with lidR - LiDAR point cloud

texts:
I would like to convert this point cloud from US feet to meters. Metadata of the point cloud here.
I firstly though to convert it properly by re-projecting to a metric CRS but this caused me a lot of trouble. See the following posts (in their order of publication) :

Non finite transformation detected - lidR spTransform
"non finite transformation detected" using sp R
Reproject using WKT for just the horizontal CRS

But basically, what I want is just metric point clouds! So I came to the idea of simply converting the XYZ coordinates of the point cloud using las_rescale() from lidR R package:
library(lidR)
tmp_las <- readLAS('./USGS_LPC_CA_FEMA_R9_Russian_2017_LAS_2018.laz')
tmp_las_RS <- las_rescale(tmp_las,xscale = 0.3048 , yscale = 0.3048, zscale = 0.3048)
#> X coordinates were moved by 0.0762 on average
#> Y coordinates were moved by 0.0762 on average
#> Z coordinates were moved by 0.0762 on average

But the headers don't take it into account
tmp_las_RS
#> class        : LAS (v1.2 format 3)
#> memory       : 618.7 Mb
# extent       : 6268117, 6270217, 1998802, 2000902 (xmin, xmax, ymin, ymax)
#> coord. ref.  : +proj=lcc +lat_0=37.6666666666667 +lon_0=-122 +lat_1=39.8333333333333 +lat_2=38.3333333333333 +x_0=2000000.0001016 +y_0=500000.0001016 +ellps=GRS80 +units=us-ft +no_defs
#> area         : 4.41 kus-ft²
#> points       : 7.05 million points
#> density      : 1.6 points/us-ft²

Is there a better way to do it through lidR?


texts: How to classify/style rasters correctly?

texts:
I am using a GeoTiff raster data and trying to classify it. I was wondering what are the techniques that I can use in QGIS to do this?

The raster I have looks like

but I tried the reclassify tool in ArcMap and it looks like this

Now I wish do the similar classification using QGIS and I was able to acheive it by adding 9 entries in the Generate new Color map and clicking classify and applying it to the map. The image looks like this,

I again tried going to the Transparent tab and found that the No data value is -3.4e+38 and then clicked on Apply. Nothing changed. I then selected the Band1 in the Transparent band and then applied and the entire map disappears. Is it possible if I can display it in the same manner as the ArcMap's output?


texts: Visualising Multiple Routes from pgrouting in QGIS

texts:
I am able to Visualise the Single Route from Source A to Target B by the following code

in QGIS SQL

"id" in (
SELECT id2 as id
 FROM pgr_dijkstra(
   'SELECT id, source, target, cost FROM ways',
   12746, 12714, false, false
 ) as a left join ways b on a.id2=b.id)

But when I try to load multiple paths (4 in this case)  using K-shortest path, Results are segments are scattered all over in different places and I do not see single proper route.

"id" in (
SELECT id2 as id
 FROM pgr_ksp(
   'SELECT id, source, target, cost FROM ways',
   12746, 12714, 4, false
 ) a LEFT JOIN ways b ON (a.id2 = b.id)

I am very beginner to pgrouting, so help me if i am doing something wrong here?


texts: Specifying Font Size/Style in Label Expression?

texts:
I want to create a label expression with two lines of text, this should be simple but I'm having a lot of trouble. Basically, I need the top line to be 10 pt font, and the next line to be 8pt font that is italicized.

This is what I have:

    def FindLabel ( [LABEL1], [LABEL2] ):
      return "<FNT size = '10'>" + [LABEL1] + "</FNT>" + "\n" + "<FNT size = '8'><ITA>" + [LABEL2] + "</FNT><ITA>"

This is interpreting the formatting as a string, so it's just spitting out exactly what I put in.

The piece of code below will work fine, but as soon as I add the new line everything is messed up.

    "<FNT size = '10'>" + [LABEL1] + "</FNT>"

I would like to do this using Python.


texts: Programmatically block Toolbar choosing in ArcMap 10.1

texts:
With an extension I can block ArcMap's Customize menu using following code line:

ArcMap.Application.LockCustomization("password");

But, it makes no sense since by right-click on Toolbar docking area, user still can add any toolbar s/he wants. Is there a way to block that right click menu too?


texts: Converting shapefile to kml with good projection

texts:
I'm trying to work on to convert from a shapefile to KML however, I could not get projection correctly. I have had success converted it, but I lost the files and cannot reproduce them.
The file i'm trying to convert is the shapefile coming from here: http://www.cartelec.net/?page_id=3609

I tried to do a few times with the CRS with no success, but the only one with something viewable in google earth was WGS84 but it showed my map in France in the Netherlands.

Can someone please help me to how do I get it correct?

EDIT  - Content of PRJ files:

PROJCS[

    "_MI_0",

    GEOGCS[

        "NAD83",
            DATUM[
                "D_NAD83",
                SPHEROID["Geodetic_Reference_System_of_1980",6378137.0,298.2572221009113]
            ],
            PRIMEM["Greenwich",0.0],
            UNIT["Degree",0.0174532925199433]

    ],

    PROJECTION["Lambert_Conformal_Conic"],

    PARAMETER["false_easting",700000.0],
    PARAMETER["false_northing",6600000.0],
    PARAMETER["central_meridian",3.0],

    PARAMETER["standard_parallel_1",44.0],
    PARAMETER["standard_parallel_2",49.00000000001],
    PARAMETER["latitude_of_origin",46.5],

    UNIT["Meter",1.0]

]


texts: Converting colored .Tiff images to Shapefiles using ArcScan?

texts: I have a digital map of water pipelines (water lines in various colors, parcel data in black). I want to convert the .tiff image's water pipelines (not parcels) to a shapefile. I have the .tiff geo referenced and but cannot seem to figure out how to get it to work in Arcscan.


texts: transform longitude and latitude to postgis geometry failed: exceeded limits(-14)

texts:
I'm trying to transform my longitude and latitude column in PostgreSQL,  into PostGIS geometry point. (address in LA area)

This is the sample:(-118.373258481016618, 34.712161771785425) // longitude, latitude

When I run this Command:

 SET geometry =
 ST_Transform(ST_SetSRID(ST_MakePoint(longitude, latitude), 4326),
 2227);

I would have this
Terminal Error :

  Executing (default): UPDATE myTable SET geometry =
  ST_Transform(ST_SetSRID(ST_MakePoint(longitude, latitude), 4326),
  2227); Unhandled rejection SequelizeDatabaseError: transform: couldn't
  project point (-118.001 -1.4873e+008 0): latitude or longitude
  exceeded limits (-14)

I know the command works because, if I hard code in latitude like this:

UPDATE crime_la_counties SET geometry =
ST_Transform(ST_SetSRID(ST_MakePoint(longitude, 34.095360134448367), 4326),
2227);

It's sort of working, and I can see my points in QGIS like this:

(It's a line because I hard coded in latitude).

Anyone knows what went wrong with latitude column?
The column type currently is numeric.
I also tried to trunc + round latitude:

SELECT trunc( CAST("latitude" as numeric), 8) FROM myTable;

SELECT round( CAST("latitude" as numeric), 8) FROM myTable;

It didn't work either...

version:

  QGIS 2.3.3

  PG Admin4

  psql (PostgreSQL) 9.6.5

UPDATE:

After I type in:

select min(latitude), max(latitude) from crime_la_counties

I got this:

so trying to filter the raw data now... in progress


texts:
I have a shapefile which contains landuse types. The Landuse types are written in a field called "Type" and contain also the area in m² and a landuse code. Some polygons contain more than one landuse type separated with a tab

An example or such a field would be

(43001-1200) Brachland 1185 m²
(43002-1300) Laub- und Nadelholz 14 m²

What I would like to do is make a list of the landuse types that interest me, and then search the field for these types. If the word is contained within the field, that polygon should be selected.

I then need to remove the code and the area leaving just the Landuse type. IE int he above example the resulting polygon would simply have Brachland, Laub- und Nadelholz in the field.

I have ArcGIS (Editor License) and FME Desktop. Does anyone have an idea which approach would work best?

Thanks for any replies,

Rob


texts: MapServer Raster WMTS From PostGIS Column

texts:
I'm going to guess this isn't possible, since I haven't seen documentation supporting this in either MapServer, PostGIS, or GDAL.

I have a PostGIS table with two columns -- timestamp and rast -- timestamp is of type timestamp without time zone and rast is of type raster.  I'm trying to set up a WMTS layer in MapServer that uses timestamp as the wms_timeitem and the associated rast as the DATA.

I've tried a variety of approaches -- CONNECTIONTYPE and CONNECTION of course although the MapServer documentation says that isn't supported for raster layers.  I've also used a GDAL PG connection in DATA instead (e.g. PG:host=...) but I see no way to select a specific column let alone tie that to the time request.  I obviously get errors with either of those approaches and no raster is produced.

I have no control over the structure of these tables so if this isn't possible with MapServer I will have to try something else (maybe just serving raw output via GDAL) and I am open to FOSS4G suggestions in that regard.


texts: Unable to add layer to GeoServer - Origin does not correspond to request

texts:
I am unable to add a layer to GeoServer, whenever I click to add layer and select a store from dropdown, I get
HTTP ERROR 400
Problem accessing /geoserver/web/wicket/bookmarkable/org.geoserver.web.data.layer.NewLayerPage.
Reason: Origin does not correspond to request

I have tried setting GEOSERVER_CSRF_WHITELIST to allow-list of proxy with the CSRF filter but no luck, even tried setting GEOSERVER_CSRF_DISABLED property to true but that also didn't work. Recently we shifted out GeoServer from http to https
Can someone tell me how to resolve this or possibly the reason for it to occur


texts: Is there an online tile-based elevation service?

texts:
I am working on an application that requires the elevation at a given lat/long coordinate for lots of closely spaced coordinates along any given irregular path. (The path is a trace of a glider/hang-glider flight.)

I could use SRTM/GTOPO30/GLOBE data but I only need small (but arbitrary) sections of the world at any one time so I'd like to avoid storing a DEM database locally.

There is also Google's elevation API (for use with their maps) but it has a usage limit of 2500 requests a day whereas my data can be tens of thousands of points. Plus querying a server so many times for closely spaced locations seems like a waste a bandwidth.

So, ideally, I'd like a tile server for small DEM tiles. I found this http://www.arcgis.com/home/item.html?id=1b48cd3b6276416784fe90a68c580a89 but the ArcGIS API guide doesn't contain any useful info about it. Has anyone actually used it? Are there any other similar services?


texts: GeoServer WFS 2.0 CSV output multiple features/typeNames

texts:
I'm setting up an export function so that users may download data from a GeoServer WFS 2.0 service in a variety of the supported formats.

The service has multiple FeatureTypes which should be available for export in a single query (along with BBOX and attribute filtering).

For simplicity we are using a get request at this stage (more due to applications requirements).

I have a query like this:

http://www.myurl.com/wfs?
service=WFS&
request=GetFeature&
version=2.0.0&
count=9999999999&
typeName=a_prefix_here:feature_type_one,a_prefix_here:feature_type_two&
outputFormat=csv&

The problem is: the data from feature_type_two is not included in the exported CSV file (only exports the first feature)?

Confusingly though; feature_type_two's data is included when exporting to other formats, e.g. JSON, KML, XML?

Has anybody else had this issue?


texts: Automated matching between the QGIS custom CRS and a custom SRID from Postgis

texts:
According to this answer we know that QGIS does not get the SRS from the postgis spatial_ref_sys table. We have to define the custom SRID under QGIS.

Is there a way that QGIS can make the matching automatically between his custom CRS in QGIS and the Custom SRID in postgis?


texts: Remove/dissolve parts in ESRI line feature class

texts:
I am using ESRI ArcGIS 10.1.

I have a line which has been dissolved, having initially been composed of 4 features. I want to remove these parts, properly dissolving it into one continuous feature. I tried using the Unsplit Lines tool but no luck.

Additionally, once dissolved it places the final vertex at a mid point, not the end of the line:

This proves very difficult to manually correct. Curiously, the position of the end point is the end point of the right most line from the four that this was initially dissolved from.

If anyone knows how to get ridd of these parts, I'd appreciate it.

I am wondering if it is because there is an infinitessimally small distance between parts, but it seems they're snapped together.


texts:
How to set attributes for writers generic objects, for instance, lines?!

Writer generic objects lines are created, but without attribute because I do not know how to set attributes and values.


texts: Importing NetCDF file to ArcGIS Desktop as point feature layer?

texts:
I have a netcdf file of my computational grid. It is a structured curvilinear grid. When I make a feature layer from a it with the Multidimension Tool toolbox, I cannot get the point locations to plot properly. If I don't tell arcgis what netcdf dimensions to use (i.e. row dimensions option empty) all I get is a single point at the 0,0 origin. If I try to use the xi and yi dimensions, which are the indices of of the netcdf variable that I want to plot (i.e. with row dimensions as xi and yi) I get them to plot with the coordinates as the xi and yi indices instead of latitudes and longitudes.

I have the set of lon - lat variables in my netcdf files, but I cannot select them to be used as the dimensions or the dimension values (the tool box will not see them at all, although they are in the netcdf file).

How can I fix this issue?

Using ArcGIs 10.0.


texts: How to Register a Spatial View when ArcSDE is direct connect only?

texts:
I have an ArcSDE 10.0 set up with MS SQL 2008 R2 Database, on direct connect. There is no service running, nor installed.

I have created a spatial view built on a couple of tables. I have followed the steps given in
Creating a spatial view in SQL Server using SQL and registering it

The View can now been seen in MS SQL Server management studio. It is also seen in ArcCatalog as a table.

Now I need to register it with ArcSDE, so that it is seen as a featureClass.
This requires that I run the sdelayer -o register command.
The problem is that since there is no service installed, I cannot run the command.

What is the solution in this case? How do I register this spatial view so that it appears as a layer to ArcGIS?


texts: How to easily edit attribute data using Regular Expressions?

texts:
I am writing a tutorial for MSc students that uses QGIS.
One of the tasks requires the attribute data to be modified
(see pdf draft of tutorial here if interested ).

The question is how to get them to edit the attribute data.
I know this question has been asked before, and the answer used to be simple: just use LibreOffice calc to modify it and then save it, ensuring the row order and column titles are kept the same. I have used this approach before to great effect.

Unfortunately, the solution no longer works: saving the .dbf file, in LibreOffice (5.4.2) or OpenOffice (3.4.1) Calc scrambles the data from the perspective of QGIS (only random integer numbers appear in the attribute table), and creates an additional .dbt file instead.

So to get to the question: how do I tell a class of MSc students to do this next week? I've considered the following options:

Use R instead (great for me, but scary for students)
Edit the file as a .csv, and then join based on timestamp (fine, but seems overly complex, and could introduce more error in the 2 hour practical)
Use an alternative program: I've tried Gnumeric and Excel, neither of which generated good results quickly.

To help answer this question, I've created a replicable example.
Open the file "points.shp" in QGIS, observer the attribute table then close.
Edit a couple of cells in the attribute table "points.dbf" in LibreOffice - save.
Re-open "points.shp" in QGIS and observe the scrambled attribute table. The other shapefiles in the .zip folder have already been corrupted by LibreOffice and OpenOffice to illustrate what the problem looks like for me.

In short, how do I best solve this problem?


texts: How can I change the display name of the "Scale" control?

texts:
Here is the thing. I added the Scale control to my web and is displayed in a div at the bottom. I need to translate the "Scale" text to other language but I can't find how to do it. I searched in the openlayers website for examples or code to do it without results. None of the properties of Scale control tell me how to do this (in openlayers website).

So I ask to you, people, if you know some method to customize this feature.

Here is the code I use in my script:

map.addControl(new OpenLayers.Control.Scale($('scale')));

and this is the div where I displayed:

I tried with writing the label in the div but the scale control overwrite the label. Hope you can help me with this.


texts:
I use ogr2ogr (bundled with QGIS 2.14.3) on Windows to convert an ArcSDE MS SQL Server 2012 database with SQL Geometry to PostGIS (2.1.8). I use the following syntax:

  ogr2ogr --config PGCLIENTCODING ISO-8859-2 -f "PostgreSQL"
  PG:"host=localhost user=postgres dbname=atlas password=secret"
  "MSSQL:server=srv-sql;database=test;trusted_connection=yes"
  -a_srs "EPSG:23700" -lco GEOMETRY_NAME=geom

It seems to do the conversion, all tables are copied, but I get this error:

  ERROR 1: Error initializing the metadata tables :

From other post it seems that it is related to the MS SQL driver or the geometry_columns table, but I can't seem to find the problem as the error message does not continue...

Do I need to worry about this? Or how to get rid of this?


texts: Merging Polygons using ST_Multi ST_Union group by column value

texts:
I'm trying to merge polygons  together to form a larger 1 and I'm getting some odd results
In the  picture the red is giving the correct outline but there are 11 identical polygons I was expecting only 1 (this is related to the number of smaller polygons to be merged) eg the number of polygons with the same gd.gid_4 value
Can anyone say  what I'm doing wrong?  I'm running this query in PostgreSQL and viewing in qgis
Please note using "*" for clarity here
INSERT INTO new_table (*)
SELECT g.*,gd.geometry
FROM gadm36 g
JOIN(
    SELECT gd.gid_4 AS gid_4,
    ST_Multi(ST_Union(gd.geom) ) As geometry
    FROM gadm36 AS gd
WHERE  gd.gid_0 = 'FRA'
    GROUP BY gd.gid_4

) gd ON g.gid_4 = gd.gid_4;


texts: Route hatch labeling as stations (0+00) with vba

texts:
I have a route feature which I have set up to label every 500' with the distance along the line. I'd like to use typical stationing format (000+00). I've been using the following VBA script:

Function FindLabel ( esri__measure )
if esri_measure=0 then
    FindLabel = "0+00"
elseif   esri_measure<1000 then
     varleft = left(esri__measure, 1)
     varright = right(esri__measure, 2)
     FindLabel = varleft & "+" & varright
else
     varleft = left(esri__measure, 2)
     varright = right(esri__measure, 2)
     FindLabel = varleft & "+" & varright
end if
End Function

each "else" works properly on it's own, but when all the script is entered, the first "else" to be satisfied is applied to all labels. For instance, with the script as copied above, every distance label reads exactly "0+00". To see what happens I try changing the script to:

Function FindLabel ( esri__measure )
if   esri_measure<1000 then
     varleft = left(esri__measure, 1)
     varright = right(esri__measure, 2)
     FindLabel = varleft & "+" & varright
else
     varleft = left(esri__measure, 2)
     varright = right(esri__measure, 2)
     FindLabel = varleft & "+" & varright
end if
End Function

With this script every label is formatted "X+00". When the value is greater than or equal to 1000, I need it to read "X0+00". The problem seems to be that the program tests only the first label (distance 0) against the "if" statement, then applies the outcome to every labeled location. How do I get each label to run through the "if" statement independently rather than only testing the first point?


texts: vb.net arcobjects mapping different field values in attribute table

texts:
I need to mapping different field values in an attribute table. The attribute table as follows:

For example, I want to get the "DrainName" value when "BQ no." = "22/22" and "Segment ID" = "22/22-1". The result should return 2 values.

Is there a way vb.net to do this? Thanks.


texts: Is it possible to access fields from a relationship class to use as a kriging value

texts: I have a point dataset of well locations related to a table of calculated production figures by year through a simple relationship class.  I am trying to use the production figures as the z value for a kriging surface.  The kriging tool does not seem to be able to recognize the related table attributes.  The kriging tool acts like the z value must come fom the original input table.  Is it possible to access table fields from a relationship class so I can use the kriging tool?


texts: Creating color table file in PCI Geomatica

texts: What is the simplest way to create a stand-alone color-table file in PCI Geomatica (PCT or STR)? I have some classified imagery with values ranging between 0-255 and I would like a way to automatically have them consistently colored the same way...


texts: Qualitative Color Palette - photocopy safe (grayscale ready)

texts:
I am looking for a 5-class qualitative color palette which can be printed in black and white printers.

Below image shows a map in color (top: qualitative, bottom: diverging) and how each looks in grayscale (right hand side). Both color palettes (Set3 and Spectral) belongs to ColorBrewer.

[UL] Qualitative Set3
[UR] Qualitative Set3  in grayscale
[LL] Diverging Spectral
[LR] Diverging Spectral in grayscale

I especially suffer to discern groups in the top-right image, while the bottom (spectral) is much better, if not perfect (i cannot say which was blue or red). And according to ColorBrewer the Spectral color scheme (the one I used here) is the only photocopy safe for 5-class palette.

Reducing classes (to 3 or 4) may be an option, but I am afraid it was not intended purpose of this map.

This is probably not related to GIS softwares, but I basically use QGIS which offers ColorBrewer, cpt-city, and other palettes and color ramps.

@obrl_soil suggested this is Viridis color set

@underdark suggested this is two types of 4 colors + white approach.

LEFT: Orange-Red + white (sequential), RIGHT: Purple-Orange + white (diverging)


texts:
I am trying to No_Test append a polygon feature class to another polygon feature class. Both fc's are located in the same geodatabase. I am using the Append tool with ArcGIS 10.2.

I'm getting two errors, but to focus the question on one--I keep getting Error 000594: Input Feature 1L falls outside of output geometry domains.

I recall having this issue before and being able to fix it with environmental settings, but the error still occurs with the following set:

  Output Coordinates: Same as input

  Processing Extent: Union of inputs

  M and Z values have  been disabled (although there are no M or Z values for either fc. I've read that if the fc was created in ArcGIS 9.x or earlier, it may have issues with M or Z values.)

  Geodatabase Advanced: Maintain Spacial Index has been checked. Output M and Z are same as input.

Any suggestions on what else I should be troubleshooting?


texts: Extract seasonality and trend values for each time in a MODIS time series

texts:
I'm trying to construct a MODIS temporal series as shown in this thread:

Harmonic MODIS trend in Google Earth Engine?

But I'd like to plot, if possible, the values for trend and seasonality for each time in the time series, and the regression returns just one image as result.

I'd like to plot something like presented in the BFAST-R webpage (http://bfast.r-forge.r-project.org/), which the figure below was extracted.

Can I extract this kind of information from resultant image given in GEE?


texts: What is status of OpenSceneGraph (OSG) plugin for QGIS?

texts:
Does anyone know if the OSG plugin for QGIS is working?
If so does anyone knows what the future plans is for OSG plugin in QGIS?
I am using QGIS version 1.8.0-Lisboa.

I keep getting the following errors when selecting the "launch Globe" menu option, with the WorldWind elevation type. UnsupportedOperationException: GeometryGraph::add(Geometry &): unknown geometry type: N4geos4geom7PolygonE


texts:
I am using an example fracnsites to make a similar map but I don't under stand the ajax code?
The example has js code:

var map,
fields = ["gid", "createdby", "featname", "feattype", "status", "acres"],
autocomplete = [];
function getData(){
$.ajax("php/getData.php", {
    data: {
        table: "fracsandsites",
        fields: fields
    },
    success: function(data){
        mapData(data);
    }
})

};

function mapData(data){
    //remove existing map layers
    map.eachLayer(function(layer){
        //if not the tile layer
        if (typeof layer._url === "undefined"){
            map.removeLayer(layer);
        }
    });

Here is the getdata php code:

<?php
ini_set('display_errors', 1);

//database login info
$host = 'localhost';
$port = '5432';
$dbname = 'sandwikimap2';
$user = 'postgres';
$password = 'postgis';

$conn = pg_connect("host=$host port=$port dbname=$dbname user=$user password=$password");
if (!$conn) {
    echo "Not connected : " . pg_error();
    exit;
}

//get the table and fields data
$table = $_GET['table'];
$fields = $_GET['fields'];

//turn fields array into formatted string
$fieldstr = "";
foreach ($fields as $i => $field){
    $fieldstr = $fieldstr . "l.$field, ";
}

//get the geometry as geojson in WGS84
$fieldstr = $fieldstr . "ST_AsGeoJSON(ST_Transform(l.geom,4326))";

//create basic sql statement
$sql = "SELECT $fieldstr FROM $table l";

//if a query, add those to the sql statement
if (isset($_GET['featname'])){
    $featname = $_GET['featname'];
    $distance = $_GET['distance'] * 1000; //change km to meters

    //join for spatial query - table geom is in EPSG:26916
    $sql = $sql . " LEFT JOIN $table r ON ST_DWithin(l.geom, r.geom, $distance) WHERE r.featname = '$featname';";
}

Here is my code:

    <?php

   //database login info
$host = 'localhost';
$port = '5432';
$dbname = 'thndr';
$user = 'postgres';
$password = 'somepassword';

// connect to database
$conn = pg_connect("host=$host port=$port dbname=$dbname user=$user password=$password");
if (!$conn) {
    echo "Not connected : " . pg_error();
    exit;
}
//get the data
$propBoundaries = $_GET['tm4326'];
$flagPoints = $_GET['tmxy4326'];

//turn propertyBoundaries array into formatted string
$propBoundaries = "";
foreach ($propBoundaries as $i => $propBoundary){
    $propBoundaries = $propBoundaries . "l.$propBoundary, ";
}

//get the geometry as geojson in WGS84
$propBoundaryDatastr = $propBoundarystr . "ST_AsGeoJSON(ST_Transform(l.geom,4326))";

//turn PropertyPoints array into formatted string
$flagPoints = "";
foreach ($flagPoints as $i => $flagPoint){
    $flagPoints = $flagPoints . "l.$flagPoint, ";
}

//get the geometry as geojson in WGS84
$flagPoints = $flagPoints . "ST_AsGeoJSON(ST_Transform(l.geom,4326))";

//create basic sql statement
$sql = "SELECT *
from $propBoundaries.l
inner join $flagPoints.l
on $propBoundaries.gid=$flagPoints.gid";

echo $sql;

//send the query
if (!$response = pg_query($conn, $sql)) {
    echo "A query error occured.\n";
    exit;
}

//echo the data back to the DOM
while ($row = pg_fetch_row($response)) {
    foreach ($row as $i => $attr){
        echo $attr.", ";
    }
    echo ";";
}

?>

I don't understand what the ajax call in the example is actually doing?  Therefore, I can't write a get data php file to connect to my database.  When the ajax in the example says data: {
            table: "fracsandsites",
            fields: fields
What exactly is table:?  is that some preset thing in ajax?  I am trying to convert my map which runs from a geojson file to coming from a postgis database instead.  Here is what I have so far, but honestly I don't understand what I am doing.


texts: Multi-field labeling of shapefile with PyQGIS?

texts:
I wrote the following code in order to draw the line labels with multiple fields. However, the label do not appear automatically in QGIS and I have to switch the "Layer Properties/Labels", and toggle OK.

def labelLine(lyr, champNom1, champNom2):
    label = QgsPalLayerSettings()
    label.readFromLayer(lyr)
    label.enabled = True
    label.fieldName = '''concat('BASE: ' + to_string("''' + champNom1 + '''"),'\nL = ' + to_string("''' + champNom2 + '''") + 'm')'''
    label.placement= QgsPalLayerSettings.Line
    label.bufferDraw= True
    label.bufferSize= 1
    label.setDataDefinedProperty(QgsPalLayerSettings.Size,True,True,'8','')
    label.writeToLayer(lyr)

def saveCoucheFRAC():
    layer=iface.activeLayer()
    nom=nom2+"_FRAC"
    lien=path2+"_FRAC"+".shp"
    QgsVectorFileWriter.writeAsVectorFormat(layer,lien,"utf-8",QgsCoordinateReferenceSystem(22332, QgsCoordinateReferenceSystem.EpsgCrsId),"ESRI Shapefile")
    #Charger cette couche dans QGis
    vlayer=QgsVectorLayer(lien, nom, "ogr")
    if not vlayer.isValid():
        print " Verifier la validite du fichier"
    else:
        # appliquer les couleurs selon le niveau de tension avec la fonction appCouleur ci-dessous
        appCouleur(vlayer)
        # changer l'epaisseur de la couche ligne
        epaiLyr(vlayer, "0.7")
        # ajouter l'etiquette de longueur sur chaque segment
        labelLine(vlayer,"BASE", "Longueur_m")
        # Definir la visibilite de la couche entre ech min 0 et max 100 000
        visLyr(vlayer, 0, 100000)
        #Charger dans le Canvas
        QgsMapLayerRegistry.instance().addMapLayer(vlayer)
saveCoucheFRAC()

I try to add the following to enable the label draw but it does not work:

layer.triggerRepaint()

or

layer.commitChanges()
layer.updateExtents()

A screenshot of the map (After toggling OK, line labels appears in QGIS)


texts:
I am using QGIS 1.8.0 on a windows 7 machine.  Trying to upload a project onto QGIS Cloud.

I have 2 issues. Firstly when i go to publish a map it directs me to upload the data.  As i am using a free account i appear to have the option of 2 databases. fwjdl_kjshzi and eqgmke_bqmbtb.  I don't think i created either database but they are the options.

I go through the process of uploading and when all is uploaded the first problem presents.

I am asked to save the project as a new project. I do so. then it directs me back to uploading again. It is a never ending cycle.  I have the box ticked to replace local data but it seems not to.

Second problem is on some projects which successfully uploaded (some layers did not upload). I would like the layers to be displayed in the groups which i had them in originally. Such as roads in one folder, water data in another, etc.  Is this possible.

I hope this is clear.


texts: Exporting Selected Raster Graphics using ArcGIS Desktop?

texts:
I'm having a problem exporting raster graphics using ArcGIS 10.  It seems simple enough, but for some reason I cannot get ArcMap to work with me. Here are my steps:

Add a raster to ArcMap
Open the raster's Attribute table and select some values
Right click on the raster in the ArcMap table of contents --> Data --> Export Data
I am given 3 options for the data extent (a) Data Frame, (b) Raster Dataset, (c) Selected Graphics.  The last option (which is the one that I want) is grayed out.

I have not been able to activate the "Selected Graphics" option in the Export Raster Data dialog.  The strange thing is that my colleague has the exact same version of ArcMap and when repeating the steps above, the option to export Selected Graphics is not grayed out.   I'm assuming that there is some setting that I have to adjust.  Does anyone know what can be causing this?


texts: How to make use of indexes with inner join in Postgis?

texts:
I got 2 sets of points in 2 separate tables. Table_a got 100k points and table_b got 300k of points. I trying to find nearest points in relation find me any point from table_b that is within 50 meters from tabla_a. After that calculate fall column, group them by table_a a_id column and return highest value.

I wrote a following query that meet this criteira

SELECT DISTINCT ON (a_id) *
FROM (
       SELECT
         table_b.b_id,
         table_b.height - st_3ddistance(table_b.geom, table_a.geom) fall,
         table_b.geom,
         table_a.a_id
       FROM table_a
         INNER JOIN table_b ON _st_3ddwithin(table_a.geom, table_b.geom, 50)) a
WHERE fall >= 0
ORDER BY a_id, fall DESC;

I added 3d geometry indexes:

CREATE INDEX table_a_geom ON table_a USING GIST (geom gist_geometry_ops_nd);
CREATE INDEX table_b_geom ON table_b USING GIST (geom gist_geometry_ops_nd);

However my problem is that i can't make query to use them. Query planer is keep choosing sequence scan that is slow. I run some test changing _st_3ddwithin with st_3ddwithin, <<->> < 50 , creating 50 m buffer and intersect, st_3ddistance < 50 but everytime planner is choosing sequence scan. Is there a way to use indexes with higher performance or changing the query to use indexes?

My query plan:

Unique  (cost=10462593.70..10473018.43 rows=1 width=144)
  ->  Sort  (cost=10462593.70..10467806.06 rows=2084945 width=144)
        Sort Key: table_a.nmbayuid, ((table_b.height - st_3ddistance(table_b.geomgr, table_a.geom))) DESC
        ->  Nested Loop  (cost=0.00..10243762.28 rows=2084945 width=144)
              Join Filter: (_st_dwithin(table_a.geom, table_b.geomgr, '50'::double precision) AND ((table_b.height - st_3ddistance(table_b.geomgr, table_a.geom)) >= '0'::double precision))
              ->  Seq Scan on table_b  (cost=0.00..1459.47 rows=47147 width=96)
              ->  Materialize  (cost=0.00..10.97 rows=398 width=56)
                    ->  Seq Scan on table_a  (cost=0.00..8.98 rows=398 width=56)


texts: False easting/westing with MGI/Austria M28

texts:
After extracting the XY coordinates from geometry data via "save as csv" I get coordinates which are about 150 km more in the west than they should be. North/south seems right. I set the CRS to MGI/Austria M28.

Any suggestions?


texts:
We have the go2streetview plugin installed in QGIS 2.8.6 Wien and when performing a search using this plugin, the following error appeared:

An error has occured while executing Python code:

Traceback (most recent call last):
  File "C:\Program Files (x86)\QGIS Wien\Shared Plugins\go2streetview\go2streetview.py", line 548, in canvasReleaseEvent
    self.openSVDialog()
  File "C:\Program Files (x86)\QGIS Wien\Shared Plugins\go2streetview\go2streetview.py", line 571, in openSVDialog
    print self.bbeUrl
IOError: [Errno 9] Bad file descriptor

Python version:
2.7.4 (default, Apr  6 2013, 19:54:46) [MSC v.1500 32 bit (Intel)]

QGIS version:
2.8.6-Wien Wien, e72e2ab

Python path: ['C:\Program Files (x86)\QGIS Wien\Shared Plugins\processing', 'C:/PROGRA~2/QGISWI~1/apps/qgis-ltr/./python', u'C:/Users/adowning/.qgis2/python', u'C:\Program Files (x86)\QGIS Wien\Shared Plugins', u'C:/Users/adowning/.qgis2/python/plugins', 'C:/PROGRA~2/QGISWI~1/apps/qgis-ltr/./python/plugins', 'C:\PROGRA~2\QGISWI~1\bin\python27.zip', 'C:\PROGRA~2\QGISWI~1\apps\Python27\DLLs', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\plat-win', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\lib-tk', 'C:\PROGRA~2\QGISWI~1\bin', 'C:\PROGRA~2\QGISWI~1\apps\Python27', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\PIL', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\jinja2-2.7.2-py2.7.egg', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\markupsafe-0.23-py2.7-win32.egg', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\python_dateutil-2.2-py2.7.egg', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\pytz-2014.2-py2.7.egg', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\win32', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\win32\lib', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\Pythonwin', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\Shapely-1.2.18-py2.7-win32.egg', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\six-1.6.1-py2.7.egg', 'C:\PROGRA~2\QGISWI~1\apps\Python27\lib\site-packages\wx-2.8-msw-unicode', u'C:/Users/adowning/.qgis2//python', 'C:\PROGRA~2\QGISWI~1\apps\qgis-ltr\python\plugins\fTools\tools', '//fileserver1/Uniform/Mapping/QGIS/Projects']

This error doesn't occur all the time and so far only happened twice out of the twenty times I've tested this. Are you able to advise what the cause of this may be please?


texts: Take first match when encountering a 1-Many relationship when joining with ArcPy Data Access UpdateCursors

texts:
My current code is as follows:
import arcpy

sourceFC = r"C:\Path\SourceFeatureClass"

sourceFieldsList = ["JoinField", "ValueField1", "ValueField2", "ValueField3", "ValueField4", "ValueField5" ]

valueDict = {r[0]:(r[1:]) for r in arcpy.da.SearchCursor(sourceFC, sourceFieldsList)}

updateFC = r"C:\Path\UpdateFeatureClass"

updateFieldsList = ["JoinField", "ValueField1", "ValueField2", "ValueField3", "ValueField4", "ValueField5"]

with arcpy.da.UpdateCursor(updateFC, updateFieldsList) as updateRows:
    for updateRow in updateRows:
        keyValue = updateRow[0]
        if keyValue in valueDict:
            for n in range (1,len(sourceFieldsList)):
                updateRow[n] = valueDict[keyValue][n-1]
            updateRows.updateRow(updateRow)

del valueDict

After the join, it looks like the cursor keeps updating the joined field until another match isn't found.
I would like to make an adjustment to the code, so the cursor stops updating the joined field when the first joined value is found (i.e. when the cursor goes through, matches the join field, it doesn't keep updating until the last value is found, I want it to stop at the first found value).
I did some testing and it looks like the arcpy.JoinField does what I want (i.e. first occurrence), but I am stuck on which part of the code to update to have this parameter.


texts:
QGIS strategy for denoting point feature changes over time?

My project file has a shapefile layer (point symbols) geolocating 50 historic street name mosaics embedded in concrete sidewalks throughout an urban neighborhood. Label for each point displays the historic street name. Contemporary street names and curb outlines are displayed on background shapefile layers from City GIS portal.

As sidewalks are gradually replaced, tiles at some locations will be removed, then reset at same locations in the new sidewalks. Tiles at some locations may never be removed, others may be removed/reset several times. Attribute table for my original shapefile includes separate columns for creation date and subsequent reset date.

How to indicate on map that a tile was present in original sidewalk, and now is reset in new walk, or even was reset several times? Label text needs to show historic street name.

My ideas:

New shapefile layer with larger point symbol, e.g. a circle with no color fill. Circle symbol on new layer would surround but not obscure original symbol on original shapefile layer. New layer table includes only one attribute column = reset date. For tiles that are removed/reset several times, new shapefile layers can be added with bigger concentric circle symbols.  Example, any mosaics that that are removed/reset 2 times would have a symbol on 3 separate layers (orig layer, reset #1 layer, reset #2 layer) and would display 2 concentric circles around the original symbol.
Multiple, nested symbols, or label colors, or other differentiation on original layer. No clue how to place several different symbols on same layer.  If this is a good strategy, please tell me what this process is called and point me to where to learn about it.
Other strategy?


texts: ST_Geometry vs SDE Binary: seeking recommendations, experiences

texts:
I've got some questions I've collected as we try to come to a decision about whether to use ST_Geometry or SDE native geometry for an Oracle/SDE datastore on a new project. The organization has quite a few legacy systems that use the SDE Binary type, and this would be the first to use ST_Geometry if we went that route.

Apart from the ability to access spatial data in ArcSDE with SQL, what are the benefits and drawbacks of ST_Geometry?

What sort of pitfalls have you experienced in implementing with one type or the other?

Can more than one data type live in the same schema?

Note: for this project the data in SDE will be all points and polygons.


texts: Counting nodes/vertices in large polyline GeoJSON?

texts: I have a large .geojson file (~ 1 GB) of made up of thousands of polylines.  My objective is to count the number of vertices (not end points) in the entire file.  Is there a python module that is suited for this type of analysis?


texts:
I've build an image overlay in Google Earth, but I want to add some other layers via QGIS. In QGIS 2.18.15 on Windows 10, the kmz files just open by adding a raster layer. But...they're clearly wrong. The information in the doc.kml file appears to be a bounding box and the display in QGIS looks like the aspect ratio is wrong, as if the rotation parameter was completely ignored.

So what's the magic to make this work? I have no idea where to start on this and my Google-fu trying to find an answer has failed me.

-- Additional info

I'm trying to add a couple of images, but they're not showing up...

Here is the Google Earth Overlay. Note that aspect ratio of the image.

Google Earth Overlay

Here's what happens when I load the kml into QGIS:

QGIS KML Load

The aspect ratio is wrong, and I really expected north to still be "up".

-- Edit

The QGIS CRS is EPSG:4326. I think I tried forcing EPSG:3857 based on something I found elsewhere about Google Earth, but there is no visible different. If I add the north arrow to QGIS, it points straight up, aligned with the image. Not that the original map is an orienteering map so the north lines on it are in fact magnetic north. What I was aiming for was to pull it into QGIS after getting it calibrated via Google Earth so I could layer a UTM grid over the map.

Here's the doc.kml extracted from the KMZ file:

<?xml version="1.0" encoding="UTF-8"?>
<kml xmlns="http://www.opengis.net/kml/2.2"
xmlns:gx="http://www.google.com/kml/ext/2.2"
xmlns:kml="http://www.opengis.net/kml/2.2"
xmlns:atom="http://www.w3.org/2005/Atom">
<GroundOverlay>
    <name>Untitled Image Overlay</name>
    <color>66ffffff</color>
    <Icon>
        <href>files/CCI_000001.jpg</href>
        <viewBoundScale>0.75</viewBoundScale>
    </Icon>
    <LatLonBox>
        <north>41.44234537860543</north>
        <south>41.4037612242825</south>
        <east>-73.82018638989426</east>
        <west>-73.85321641890307</west>
        <rotation>13.67718681063408</rotation>
    </LatLonBox>
</GroundOverlay>
</kml>


texts:
Are there any resources available where I can download mrsid data in mg2, mg3 and mg4 formats? I've been able to find mg2 with no problem but so far nothing for mg3 or mg4 formats.

thanks,
Jeff


texts:
Possible Duplicate:
  Steps to Start Web Mapping

I have taken on setting up a database which stores information of an annual country-wide survey of alien aquatic weeds, linked to a GIS program to map the sites and attributes.

For my database I am using postgres, with qgis for the mapping. The eventual aim is to publish it online for anybody to access, and get basic maps detailing the surveys. I have been going around in circles trying to understand what process to follow.

What tools do I need to get from my qgis map to an online map, and what is the basic outline of the process I need to follow?

I have managed to derive a map and attributes from my database.


texts: Problem with QgsRasterCalculator.processCalculation() in PyQGIS

texts:
I have a really weird problem with a PyQGIS script which use a "for" loop. My script is doing those tasks :

Go in a given directory and list all the files with os.listdir
Enters a for loop to run tasks for each file if a condition is checked
Process at the end a QgsRasterCalculator.processCalculation() function

The problem is that my for loop doesn't always carry out the same number of loops when I run the program. I have 6 files in my actual folder which fulfill the condition, and when I run the program, once 2 files are processed, the next time it's 5, another time 6...

After some little efforts (put some print('test') to locate where the error is), I know that when the for loop doesn't run 6 loops, it's because the QgsRasterCalculator.processCalculation() function doesn't work. The thing really weird is that it works one in two, so I really don't know what's appening here haha.

Here's my code :

import processing
from qgis.core import QgsApplication, QgsProcessingFeedback, QgsRasterLayer
from qgis.analysis import QgsNativeAlgorithms, QgsRasterCalculator,    QgsRasterCalculatorEntry
from processing.core.Processing import Processing
from BaseFunction import *
import os

qgs = QgsApplication([], False)
qgs.initQgis()
feedback_process = QgsProcessingFeedback()
Processing.initialize()
QgsApplication.processingRegistry().addProvider(QgsNativeAlgorithms())

directory = 'Y:\\1_POLD\\TEMPERATURE\\STD_MN_JOUR_LANDSAT\\sources\\'
file_list = os.listdir(directory)

entries = []

for files in file_list:

    if files.startswith('clip') and files.endswith('.tif'):
        # Loading the raster layer for the QgsRasterCalculator instance
        layer = QgsRasterLayer(directory + files, files)

        print(files)

        lyr = QgsRasterCalculatorEntry()
        lyr.ref = layer.name() + '@1'
        lyr.raster = layer
        lyr.bandNumber = 1
        entries.append(lyr)

        # Processing the calculator
        frequ_expr = lyr.ref + ' >= 2'
        frequ_calc = QgsRasterCalculator(
            frequ_expr,
            directory + '1_' + files,
            'GTiff',
            layer.extent(),
            layer.width(),
            layer.height(),
            entries
        )

        frequ_calc.processCalculation(feedback = feedback_process)
        print('test')

QgsProcessingFeedback doesn't return anything by the way.


texts: Why is 'Count Points in Polygon' returning zero values for all NUMPOINTS?

texts:
I'm using QGIS and I have a file of points (Coyote_Clipped).  I'm trying to count the number of points in each grid cell. I ran "Count Points in Polygon" on the grid and it worked without errors, but when I look at the attribute table for the count shapefile the value is zero for all NUMPOINTS, even though there are clearly points that fall inside the grid cells.  What might be happening here?
This is the file I'm working with:https://drive.google.com/file/d/1dNiUPOhNNqFHrJ3dpZ0zHi95Vn3D1iIx/view?usp=sharing


texts: Calculate bounding box coordinates of a selected polygon with QGIS

texts: I just manually (and somewhat inaccurately I'm sure) retrieved the bounding box coordinates for a county polygon, which is a selection from a larger layer of county shapes. I'd like to be able to use a QGIS plugin or other function to calculate this quickly and not-by-hand. This is almost certainly possible but I couldn't find anything with a Google search


texts: Synchronizing views in QGIS print composer?

texts:
I found a question very similar to my problem: How to link dataframes in ArcGIS or QGIS?. But the solution for QGIS doesn't work now (Dual View plug-in is not available anymore).

I want to sync three views:

The idea is to move the three map views at the same time with a master map in QGIS print composer.

Is there a way or plug-in to do this?


texts: I am trying to open a raster image file in Arcmap 10.3 . I have the AUX File, a JGW file , 2 XML Document files ,and a RRD file. It imports into Arcmap 9. However, Arcmap 10.3 does not import the image. Can someone help me?


texts:
I just installed QGIS 3.6 on a new Macbook pro running Mojave, and I'm trying to produce geospatial PDFs of some basic maps of property boundaries overlaid onto aerial imagery base maps to use in Avenza on an iPhone.

I'm new to open source GIS, and what I'm used to with Arc being able to export directly to a geospatial PDF seems to be much more difficult with QGIS. It seems that the best workaround mentioned online is to export a map to a TIF with world data, then to convert that TIF to a PDF. On a PC running QGIS 3.2.3, I can do that just fine--I export a map image which creates a .tif file and a .tfw file, then I can use the "Translate (Convert Format)" tool in the raster drop-down menu and use that TIF as an input, select the projection, and select .pdf as the output file type to generate a geospatial PDF, but on my Mac PDF isn't listed as an output file type option.

I'm struggling to find a way to just get a map formatted so I can open it in Avenza, and I have no idea why this file type is an option on my older PC version of QGIS but not on the newer version on my mac. Does anybody has an idea of how I can take a map I produce in QGIS and easily export it as or convert it to a geospatial PDF? It seems most of the info I've found on either is outdated or involves a writing Python scripts that I just don't understand.


texts: Taking substring of one field and using value to field calculate or update cursor another in ArcPy?

texts:
I have a field that has values such as DIMP 1234567 N. Eucla Ave/ W. Bonita Ave or DIMP Seg #1234567- Spinnaker Dr. Ventur or DIMP SEG 1234567 PYRITES AND ESMERALDA.  I want to extract only the numerical values and update another field with those numerical values of the same row.   The field that has the numerical values is called WORK_ORDER_REF and the field that I want to update with the numerical value is DIMP_SEGMENT.  A variation of the code below updated the DIMP_SEGMENT field only with the first value of my list.

import arcpy

a = []
table = "table"
with arcpy.da.UpdateCursor(table, ["WORK_ORDER_REF_DESCRIPTION", "DIMP_SEGMENT"]) as cursor:
     dimpsegs= [row[0] for row in cursor]
     for seg in dimpsegs:
        line = seg
        for word in line.split():
           try:
              a.append(int(word))
           except ValueError:
             pass
        for row in cursor:
            row[1] = a
            cursor.updateRow(row)


texts: How to use Domain Codes set in a Feature Class when opening the Feature Class in MS Access

texts:
If I set domain codes on a feature class, does this translate into MS Access as well?

For example, I set a domain code and use it on a field in the feature class, I can then go into ArcMap, start an edit session and click on cells in the field and I get the drop down set by the domains. What I want to know is, can you open the feature class in MS access and use that domain drop down? I tried this and the drop down doesn't populate. I'm wondering if there's something I need to enable to get this?

Thanks,
Mike


texts: How to tile and mosaic spatial vector data in R?

texts:
Processing large vector data requires a way to tile the input data and mosaic the output.
The following example only targets the tile creation and final mosaic step. Tiling the data has the reason that some processing steps might be CPU intensive or RAM exhausting if the entire data is processed. Processing steps between making tiles and mosaics are skipped and are not part of this question.
The following code creates polygons, tiles and a try to mosaic:
# create sample data
library(sf)
pts <- st_multipoint(matrix(c(10, 15, 15, 20, 20, 10, 20, 15), nrow = 4, byrow = TRUE), dim = "XY")
pts <- st_cast(st_geometry(pts), "POINT")
pol <- st_buffer(pts, 3)

# create tiles
library(GSIF)
tiles <- getSpatialTiles(as(pol, "Spatial"), block.x=6, overlap.percent=10)
for (i in 1:length(tiles)) {
  tile <- st_crop(pol, st_as_sf(tiles[i]))
  st_write(tile, paste0("path/to/tile_", i, ".shp"))
}

plot(pol, col="green")
plot(tiles, add=TRUE)

# mosaic
file_list <- list.files("path/to/", pattern = "*shp", full.names = TRUE)
shapefile_list <- lapply(file_list, st_read)
pol.union <- do.call(st_union, shapefile_list)
plot(pol.union)

How to mosaic these tiles?
I tryed several approaches (st_union, merge, st_join,...), but my best result yet does not remove the buffer and does not join parts of polygons into single polygons.
pol.union <- shapefile_list%>%
  lapply(st_union) %>%
  do.call(c, .)
plot(pol.union)

btw.: if there are better ways to create tiles (without sp), a hint would be nice.


texts: ogr2ogr SQL Expression Error, From starting with number

texts:
I'm merging shp files with ogr2ogr in a python loop and adding the specific filename as coloumn. A command looks like:

ogr2ogr -f "esri shapefile" -update -append final/merged.shp 3044/2019_0158_Example.shp -sql "SELECT '2019_0158_Example.shp' AS filename, * FROM 2019_0158_Example"

Which fails like:

ERROR 1: SQL Expression Parsing Error: syntax error, unexpected integer number, expecting string or identifier. Occurred around :
2019_0158_Example.shp' AS filename, * FROM 2019_0158_Example

This is clear as FROM should not start with a number but a character.

Is there something I can do to get around this except renaming all my files? For example some dialect or parameter I can update my command with?


texts:
i have ArcGIS Image Server for publishing my raster data to my user. In traditionally way i have been tiling with gdal and it takes a long time besides this it requires more effort for getting true result.

Do you know ArcGIS Image Service support for openlayers like Openlayers ArcGIS Server REST support?


texts: Automated mapping of cost paths from multiple origin points

texts:
I am currently mapping Lahar flows on Montserrat Island (as many have done before!). I have the Cost Distance analysis below, the circle is the origin and the triangle the destination. With a known destination a Cost Path analysis can easily be performed, however I would like to use the Cost Distance map below to automatically find the destination point, based purely on the the direction of flow from the origin/circle.

I was thinking of using a filter to include only edge cells of the raster, and then find the one with the lowest value, is there a way to do this?

Software: Arc GIS Desktop 10.5


texts: PostGIS: Update child geometries based on parent

texts:
I'm using Leaflet and postgis in a case where using Leaflet draw plugin I draw a polygon first like this

Now, I draw a line which cuts the polygon and using line and polygon geometry, based on postgis ST_Split function, I store the geometry of splitted polygons and consider them as children of previous polygon like this.

Coming towards the issue, I want to update the geometries where If I increase/decrease the length or width of black polygon (the parent), red polygons (the children) lengths should also increase/decrease. Is there any combination of techniques and/or PostGIS functions which can help me in this case


texts: georeference and orthorectification of aerial photos

texts:
I am trying to process photos taken with a tetracam ADC micro camera placed on a UAV.

I have GPS coordinates and mavlink logs providing roll, pitch and yaw of the UAV.

I am having a real hard time georeferencing and orthorectifying the pics using tools as ossim or grass.

Has anyone done this before?


texts: How to update one column in file GDB by using FME?

texts:
I have one point shape file (in FME: reader) file which I exported to Geodatabase (in FME: writer). Meanwhile I edited only some features (lines on table) on my point shape file in one column. Now I want only this updated column to be re-written on the GDB. I don't have unique ID but I have coordinates of the points that are unique (X, Y) from both reader and writer.
In other words, I am looking for transformer that compares coordinates of reader and writer and for those features where match is found re-writes one column.
I found "DatabaseUpdater" and I tried it but when I run it an error pops up without any further explanation
Anyone experience with DatabaseUpdater?


texts: Copy ESRI Map by Value (not Reference) in JavaScript

texts:
Any way that you know of in JavaScript to copy an ESRI map object by value, not by reference?

I know that I can do this:

var foo = new esri.Map("foo");
var bar = foo;
foo.id; //"foo"
bar.id; //"foo"
foo.id = "oof";
bar.id; //"oof"

what I want is to make a copy of foo so that

bar.id = "bar";
foo.id; //"foo"
bar.id; //"bar"

Any ideas? So far the best I've come up with to to construct two different maps, and then copy foo's graphics, features, etc over to bar:

var foo = new esri.Map("foo");
//do some stuff to foo, adding layers etc
var bar = new esri.Map("bar");
bar.addLayers(foo.getLayersVisibleAtScale());
//etc


texts: FAIL - Tried to use command /upload via a GET request but POST is required Deploy Geoserver.war

texts:
I'm trying to deploy geoserver.war on tomcat7 but a message error is shown:

  FAIL - Tried to use command /upload via a GET request but POST is required

What should I do?


texts: Using Esri Desktop Advanced 10.2.2 is it possible to create an attribute based on a classification method, i.e. Jenks?  To elaborate - if you use the symbology tab in a layers properties, you can set a classification method that displays the class breaks, but I want to also do some statistical analysis based on those class breaks.  Previously, I have done this by performing the classification using the layer properties and then using the class breaks to do a SQL and manually editing a attribute field with the class number.  This works, but is labor intensive.  Is there a progamatic way to do this?


texts: How to estimate the standard error on the mean difference between two DEMs?

texts:
I have two DEMs of the same area measured on two different dates. Differencing these DEMs gives a map of the change in height between these two dates. I would like to calculate the mean change in height and the standard error on the mean over a region spanning a large number of pixels.

Calculating the mean is straight forward, but the standard error seems more complicated in that the noise is correlated between neighbouring pixels - as a result, just dividing the standard deviation of the height change in my region of interest by the root of the number of pixels gives an unrealistically low standard error.

Is there a way to calculate the standard error that takes into account this spatial autocorrelation?


texts: Mask raster with second, categorical raster in python

texts:
I have two rasters. One is a "normal" Sentinel-2 L1C and the other one is a categorical raster with levels from 0-5. They are both of the same extent. Is there any way to set pixels to NA that are of a specific level in the other raster without doing some kind of polygonization?
No asking for any specific solution, but more for ideas;)


texts: MergeShapes gives AttributeError: 'NoneType' object has no attribute 'isEmpty'

texts:
I try to merge two polygon shapes in QGIS 2.18.4 with the MergeShapes plugin.

The error message shown below is what I get in return.
Does anyone know a solution?

I want to make one table out of two with different fields.
I want QGIS to put all different fields on a row, and combine the fields that are equal in both tables.

AttributeError: 'NoneType' object has no attribute 'isEmpty'

Traceback (most recent call last):
  File "C:/Users/rob/.qgis2/python/plugins\mergeshapes\mergeshapesdialog.py", line 87, in inputFile
    if files.isEmpty() or encoding is None:
AttributeError: 'NoneType' object has no attribute 'isEmpty'

Python version: 2.7.5 (default, May 15 2013, 22:44:16) [MSC v.1500 64 bit (AMD64)]
QGIS version: 2.18.4 Las Palmas, 59c1d21


texts: Matchlines with Data Driven Pages in ArcGIS

texts:
I have a map book created via the Data Driven Pages tool in ArcGIS.  Specifically, I used the Strip Map option when I created the book. And I have employed the 'Angle' attribute of the index layer in order to rotate each page so that the route is (generally) shown west-to-east, with the north arrow turning dynamically with each page.

Now, I understand the labeling options for the Data Driven pages.  But I don't want to merely label the pages - I want to create vertical matchlines that appear on each page.  Can this be done, even as each record of the index is turned 'upright'?

Graphic matchline example:

Credit/source: https://courses.cit.cornell.edu/arch465/arch367f07-workdwgs/notes/04.html


texts:
I am working on an atlas project in qgis 3 that is basically roads, towns, and other features for a quite large area. What I would like to add, is an index of all the towns in the atlas area with at least the page number, if not also a grid coordinate. My towns are all in a single layer, so I'd like to use that for my data.

Is there any way besides doing it manually? I can't seem to find any documentation or examples that would show how to do it.


texts: Merge attribute table with conditional statements in QGIS

texts:
I have two attribute tables with different values in QGIS.
The merge value between both tables is an ID number, each of this ID number have more than one values related with crops and area in the second table (imagine an id number for a polygon and each polygon have for example 4 values by crop and area)
My idea is to choose only one crop by the maximum area for every ID number, I need a conditional expression for that but I don't know how to write it.


texts: How to properly implement EPSG: 27700 using vector tiles in Openlayers 4 and Geoserver

texts:
I am trying to get pbf vector tiles which are served by Geoserver to properly render in Openlayers 4 using the code below:

<!DOCTYPE html>

  <title>Vector tiles</title>

Vector tiles

The problem is that the vector tiles are always offset, see image below:

I'm pretty sure it has something to do with the extent. I have tried both
"Compute from data" and "Compute from SRS bounds" for the Native Bounding Box in Geoserver and inputting those extents into the Openlayers extent with no improvement. Does anyone have any suggestions/solutions?


texts: Undoing Operation while editing in ArcMap with Python Add-ins?

texts:
Using ArcPy Add-ins, I can issue a message, but I don’t understand how to cancel or undo deleting an object. Is it possible?

class ControlUniq(object):
    """Implementation for addin1_addin.ControlUniq (Extension)"""
    def __init__(self):
        # For performance considerations, please remove all unused methods in this class.
        self.enabled = True
    def onDeleteFeature(self):
        chek=pythonaddins.MessageBox(u"Delete object?",u"Deleting...",1)
        print chek
        if chek=="OK":
            pythonaddins.MessageBox(u"The object was deleted",u"Deleting...")
        else:
            ### I don't understand how to undo operation Deleting
            ### Undo()
            pass


texts: Making multiple layers editable in QGIS?

texts:
I am using QGIS and I want to eliminate fields from 300 shapefiles.

All layers are loaded and selected.

Is there a way to enable all 300 for simultaneous editing?


texts:
I'm loading a shapefile in an app through a frontend with OpenLayers and trying to save the geometry to an Oracle Database with Hibernate Spatial. The shapefiles that have valid geometry (validated with JTS in the backend) don't have any problem but those that have invalid geometry (self-intersecting polygons and duplicated nodes) are not saved. Is not possible to save invalid geometry in Oracle?
On the other hand I would like to fix the invalid geometry, I've tried this (according to JTS documentation):
 polygon.buffer(0)

After execute that code, the geometry is still invalid. Is there a way to fix an invalid geometry  with javascript (OpenLayers) or Java (JTS or GeoTools)?


texts: How to select centroids within buffer range?

texts: I want to select the red circled centroids. These centroids must be within the blue buffer.


texts:
I have a linestring on my postgis.
and I want to add a new point linestring. I mean I want to update my linestring?
how can ı do this?

UPDATE roads SET geom = ST_AddPoint(geom, 'POINT(0,0)') WHERE id =9;


texts: Remove Holes Existing Within a Single Polygon for a specific area in Qgis

texts:
I have a shapefile that consists of a single polygon. I have gaps within the polygon that needs to be eliminated, but only the ones that are less than 200m2.

I tried different ways but none of them work, I need to run this process as a batch process.

Any ideas?


texts:
I am looking for a way to convert point cloud collected through a Riegl LiDAR scanner to a perfect 3D model of a town.

The sample data is from vehicle mounted scanner, but final data is going to be from aerial and vehicle mounted scanners.

The final 3D model of the town is supposed to be used as a terrain for a FPS (first person shooting) video game.

I am not sure where to begin. I have some sample data. The data is available in the following formats:

SDW
SDC
LAS/LAZ

Is it even possible? The point cloud has color information. The sample point cloud that I have has RGB values for each point, but as far as I have tested the City Engine doesn't work with point clouds, it needs meshes, am I right?

For classes, can't I use height function for the points, to define all points at a certain height to be the ground, and everything higher buildings? I need to take the cars out anyway, so maybe I can take the trees out too, then I have just two classes, the ground class and the not ground class (the buildings).


texts: Flip or mirror a group of features with QGIS 3.10.4

texts:
I'm trying to duplicate and flip a bunch of features you can see on yellow in order to complete the right part of this map (Abattoirs de Montmartre in Paris, 1850):

I tried with QGIS docs with no success, then with this plugin Shape Tools but it flips each features with their own center as origin
As you can see here:

Is there a way to flip the full layer or full group and to chose the origin?


texts: Multilevel b-spline crs issue in modeler

texts:
I'm struggling with a processing modeler tool. I have a point shapefile (WGS84) which I want to interpolate with the SAGA spline tool. The model is reprojecting the wgs84-points to UTM32N (EPSG:25832). The next step is the spline tool. But I always get this error code:

  Error executing algorithm Multilevel b-spline interpolation 'NoneType' object has no attribute 'crs' See log for more details

When I use the points (reprojected manually to UTM32N) with this model (removing reproject layer tool, because the layer is projected to 32N) it works fine.

Is this a bug?


texts: Which data providers support loading attribute tables without geometry?

texts:
I have a table in Memory from a cx_Oracle query (a list of rows, and a list of field names and their types) in my python plugin and I want to load it as a geometry-less layer (Attribute table) in qgis.  Which data Provider will Support that?

QgsOracleProvider has a constructor that Returns false, disconnects and sets valid to False if the table queried has no geometry column.

QgsMemoryProvider only allows Point, Linestring, Multipoint, MultilineString, and Multipolygon geometries. Any other geometry is interpreted as WKBUnknown, which causes the QgsMemoryProvider.isValid() to return False.  If I cheat and create a QgsVectorLayer("Point", "myLayerName", "memory"), and add features with no geometries, then I get a seg fault when I add the layer to the map layer registry. Also important to note, is that when a layer is imported using QgsMemoryProvider, it is not possible to Join against it (A Join from the layer does work.) Since mostly people want to load non-geometric tables into Qgis in order to join them to an existing vector, this reduces the Memory provider's suitability for non geometric tables.  (Unless I am wrong.  I would be delighted to be wrong.)

QgsDelimitedTextProvider requires an Input file, and I would much rather not write out a csv to the file and then Import it again when I already have the data.  But maybe I have to?

Thanks


texts: How can I mosaic raster data which have different number of bands ( 3 and 4)?

texts:
I would like to mosaic 100 raster data.
several raster data's band are 3, though last data have 4 bands.

How can I mosaic raster data which have different number of bands??

If somebody knows how to solve this issue, let me know.

OS:Windows XP
GIS: ArcInfo 10.1  or QGIS1.8.0


texts: Which value fits the best for Parameter SIZEUNITS for Googles Web Mercator Projection EPSG 3857 to show a radius in meters

texts:
In my mapfile there is a WMS definition, which shows a circle buffer and it's radius in meters. I offer UTM-Systems and Googles Pseudo Mercator. In UTM and SIZEUNIT meters it works fine. Please, could someone give an advice which fits the best for Pseudo Mercator Projection?

SYMBOLSCALEDENOM 24000
SIZEUNITS meters

CLASS
    NAME 'Test'

    STYLE
        OUTLINECOLOR 0 0 0
        PATTERN 40 10 END
        WIDTH 1
    END

    STYLE
        GEOMTRANSFORM "centroid"
            COLOR 0 0 0
            SYMBOL "vertline"
            SIZE [radius]
            WIDTH 3
            ANGLE 60
        END
    END
END

I've tried an Approximation in my POSTGIS SQL-Request:
radius / (2*PI()/10) AS r_google
this would be enough for me.

How can I tell MAPSERVER which radius to use, when the projection system changes? Is there a Switch-Statement for this?


texts: Geoserver WFS-T insert no response. bug?

texts:
I use Geoserver + Openlayers and want to add and save some new Points via WFS-T. I created insert transaction (correct) and sent it, but no response (include errors) was returned. In console was only:

01 dec 17:34:16 INFO [geoserver.wfs] -
Request: getServiceInfo

I tried to execute demo insert transaction but had same result.
what am I doing wrong?
Geoserver 2.5.3, Tomcat 7


texts: How do I read a layer from a .gpkg file whilst selecting on an attribute?

texts:
I have a GeoPackage (.gpkg) database with only one layer:

ogrListLayers("data/LM_data/zipped/ay_riks.gpkg")
[1] "ay_riks"
attr(,"driver")
[1] "GPKG"
attr(,"nlayers")
[1] 1

There are too many objects and I'd like to read only a subset of them based on an attribute:

library(sf)
library(rgdal)
library(tidyverse)
dta<-st_read("data/LM_data/zipped/ay_riks.gpkg",
                 layer = "ay_riks",
                 EWKB = T,
                 geom_column="geom",
                 query= "SELECT * FROM ay_riks WHERE KOMMUNKOD='2480';")

I get an error:

Reading layer `ay_riks' from data source `/Users/chihayakenji/Documents/spatial_r/data/LM_data/zipped/ay_riks.gpkg' using driver `GPKG'
Error in st_sf(x, ..., agr = agr, sf_column_name = sf_column_name) :
  no simple features geometry column present

What is the correct way of doing this?


texts: What coordinate system should be used for projecting a non-geographical scan?

texts:
This is not precisely GIS problem but still spatial.

I would like to measure certain features of a thin wood slice. I have a scan of that slice.

What coordinate system should be used to obtain a view "just as it is" and to be able to measure features (by ruler tool or by calculating geometry in ArcMap) of a object from scan?

Is that even possible in ArcMap?

Example of a scan:


texts: Different databases using different ways to export xml file

texts:
I tried to use different ways to get road osmid database, for same region.
The bounding box of that region with north latitude, south latitude, east longitude, west longitude are equal to 30.7302492849 30.6552968187 104.127065029 104.039651731 (in WGS84 system) respectively.

1st approach: export .xml file (name it as map) from "https://www.openstreetmap.org/export (use Overpass API)", install "https://github.com/pgRouting/osm2pgrouting", create postgresql database using the command

osm2pgrouting --f map --conf /usr/share/osm2pgrouting/mapconfig.xml --dbname road_network --username postgres --password *** --addnodes --clean

In the created database, there are 2 tables (at least from name) that I think might contains all the roads in that region, one table named osm_ways, which is an empty table, the other table named ways, which is a table with 3682 lines.

2nd approach: use osmnx.core.graph_from_bbox function in python package osmnx to get map and corresponding database, the database contains 4665 lines.

Not alone the 2 databases have different lines, some road osm id in table in 1st approach can't be found in database in 2nd approach. I thought since this 2 approaches both use openstreetmap on same region, and also will all updated if openstreetmap updated, they should get exact same database. Can someone tell me why is this the case?