texts: How to print legend using GeoServer Print plugin?

texts:
I am using Geoserver 2.4.2 on windows. The print plugin works properly without the legend block. However if I include the legend block in my spec file, it throws an error.

Following is my spec file:

{
    "units" : "degrees",
    "srs" : "EPSG:4326",
    "layout" : "A4 Landscape",
    "dpi" : 150,
    "mapTitle" : "Map Title",
    "comment" : "Map Comment",
    "resourcesUrl" : "http://localhost/portal/img",
    "layers" : [{
            "baseURL" : "http://localhost:8080/geoserver/wms/",
            "opacity" : 1,
            "singleTile" : false,
            "type" : "WMS",
            "layer" : "portal:geo",
            "tileSize" : [256, 256],
            "format" : "png"
        }
    ],
    "pages" : [{
            "center" : [78.800618943104, 23.871591352951],
            "scale" : 100000,
            "rotation" : 0
        }
    ],
    "legends" : [{
            "name" : "Geo",
            "classes" : [{
                    "name" : "",
                    "icon" : "http://localhost:8080/geoserver/wms?version=1.3.0&TRANSPARENT=TRUE&SERVICE=WMS&REQUEST=GetLegendGraphic&EXCEPTIONS=application/vnd.ogc.se_xml&LAYER=portal:geo&FORMAT=image/png"]}]}]"
                }
            ]
        }
    ]
}

I feel the configuration for config.yaml is correct, since it works without the legend block definition in my spec file. This is my legend block definition on the yaml file:

lastPage:
      pageSize: A4
      landscape: true
      rotation: true
      items:
        - !columns
          absoluteX: 60
          absoluteY: 570
          width: 800
          items:
            - !text
              align: left
              font: Arial
              fontSize: 16
              text: 'Map Legend'
        - !columns
          absoluteX: 60
          absoluteY: 515
          width: 720
          items:
            - !legends
                backgroundColor: #FFFFFF
                borders: true
                horizontalAlignment: left
                inline: true
                maxHeight: 700
                maxWidth: 550
                defaultScale: 0.5
                iconMaxWidth: 80
                iconMaxHeight: 80
                iconPadding: 3 0 0 0
                columnMargin: 10
                classIndentation: 3
                classSpace: 5
                backgroundColor: white
                layerFontSize: 9

If I exclude the legend part in my spec file, it is producing the desired pdf file, with a blank second page. Otherwise it is giving the following error:

  Error while generating PDF:

  java.lang.RuntimeException: Cannot parse the spec file


texts: location-allocation with minimum cost optimization

texts:
I am wondering if it is possible in ArcGIS 10.3 to customize the location-allocation network analyst or use other tools to find the location of a facility which minimize the costs and not the distances?

In other words, instead of minimizing the sum of the weighted transport distances from candidate facilities to the demands to find the best location(s), I would like to:

associate a cost to each demand (different between the demands)
associate a cost to the facilities (different between the facilities)
associate a cost which is a function of the distance between them
find the location which give me the minimum of the sum of the cost.


texts: How to make openlayers features externalGraphic reappear after hiding

texts:
I'm displaying a vector layer with pois, generated from a text file:

lat lon icon    iconSize    iconOffset  title   description
56.030325   14.4774116666667    http://www.medinsbiologi.se/hano/Foto.png       -10.5,-25   102_1   UV-foto<br>More text and images

Code to draw layer:

    pois = new OpenLayers.Layer.Vector("Provpunkter",
    {
    projection: map.displayProjection,
    strategies: [new OpenLayers.Strategy.Fixed()],
    protocol:   new OpenLayers.Protocol.HTTP(
        {
        url:        "../poi/DynPOIHano.txt",
        format:     new OpenLayers.Format.Text(
            {
            extractStyles: true,
            extractAttributes: true
            })
        })
    });

The layer shows up fine, and I have now added some checkboxes to hide icons, based on which icon type is used in my text file as so:

        var feature = pois.features;
        for (var i = 0; i < feature.length; i++) {
           if (feature[i].style.externalGraphic == "url to icon") {
           feature[i].style = { display: 'block' }; // or 'none'
           }
        }

Hiding works, but I cannot get the externalGraphic back when I uncheck the box, because I don't know how to restore the externalGraphic style. I've tried setting features.stayle to null and redraw layer but then only the default orange circles appear. If I use fillOpacity=0 to hide, the icons can be brought back, but when hidden they are still accessible by users clicking in the map, which I don't want.

How can I reset the original style with the externalGraphic for some of the features in my layer?


texts: Creating a new feature class with python and then adding attibutes

texts:
I am having trouble adding rows to a shapefile after it has been created. I'm not sure if I need to create the shapefile in memory and add the data or creating then re-open the file.  While I've attempted to re-open the file, it will exit code 1 that file doesn't yet exist.  Any and all suggestion would be useful. (except setting it on fire). I'm getting a "AttributeError: Row: Error in parsing arguments for SetValue" ....  Please help.

for row in arcpy.da.SearchCursor(featureClass, ['FID','SHAPE@JSON','STREAM_ID','REACH_ID','SR_ID','STATION','INTERP','XS_ID','WSP001','WSP002','WSP003','WSP004','WSP005','WSP006','WSP007']):
featureClassData.append([row[0],row[1],row[2],row[3],row[4],row[5],row[6],row[7],row[8],row[9],row[10],row[11],row[12],row[13],row[14]])

outPutPath = "C:\\projectTemp\\FloodAnalysis\\output"
for i in range(0, len(featureClassData) -1):
FirstRowData = featureClassData[i]
SecondRowData = featureClassData[i+1]
FirstRowData[1] = eval(FirstRowData[1])

for cord in FirstRowData[1]['paths']:
    for c in cord:
        point = arcpy.Point(c[0],c[1])
    dataSetCreate = arcpy.CreateFeatureclass_management(outPutPath, 'test_section_' + `i` ,"", featureClass)
    dataSet = arcpy.InsertCursor(dataSetCreate)
    row1 = dataSet.newRow()
    row1.Shape = arcpy.geometries.Polyline(array1,spatialRef)
    row1.setValue(FirstRowData[2], FirstRowData[3], FirstRowData[4],
                  FirstRowData[5], FirstRowData[6], FirstRowData[7],
                  float(FirstRowData[8]), float(FirstRowData[9]), float(FirstRowData[10]),
                  float(FirstRowData[11]), float(FirstRowData[12]),
                  float(FirstRowData[13]), float(FirstRowData[14]))
    dataSet.insertRow(row1)


texts: Export spatial data from SQL Server to ESRI shapefile

texts:
I want to export spatial data from SQL Server to ESRI shapefile and I have an issue :
ogr2ogr -f "ESRI Shapefile" " C:\Users\sqlexport.shp "
"MSSQL:server=PORT_7FMW8H2;database=testdbspatial;trusted_connection=yes;
"-sql "select * from DRShape" -overwrite

ERROR 1: Failed to create directory C:\Users\sqlexport.shp for
shapefile datastore.
ERROR 1: ESRI Shapefile driver failed to create C:\Users\sqlexport.shp


texts: geoserver world image publishing error (georeferenced by ArcMap)

texts:
I want to add store - WorldImage in geoserver 2.5, but I still get this error:

"Could not list layers for this store, an error occurred retrieving them: Failed to create reader from file:data/picture.jpg and hints null"

According to manual at http://docs.geoserver.org/stable/en/user/data/raster/worldimage.html I only have to provide .jpg + .jgw (or .tif + .tfw).

Here they write ( http://en.wikipedia.org/wiki/World_file ) that I have to provide .prj file together with jgw (which is only transformation info, not projection) and jpg.

World files do not specify a coordinate system; this information is generally stored somewhere else in the raster file itself or in another companion file, e.g. Esri's .prj file.

I create my "picture.jpg" in ArcMap with georeferencing panel (update/rectify functions). It is weird that sometimes I get a result as a four files and sometimes only three (but still different) :

Case 1:
picture1.jgw
picture1.jpg
picture1.jpg.aux.xml

Case 2:
picture2.jgwx (different file ending!)
picture2.jpg
picture2.jpg.aux.xml

Case 3:
picture3.jgw
picture3.jpg
picture3.jpg.aux.xml
picture3.jpg.ovr (this file is here new!)

Is problem here (ArcMap 10.2) or is it in geoserver? As a next step I would like to publish my rasters as WMTS using Openlayers. Geotiffs are much larger than jpegs, so I am trying to use world images.


texts:
Running Linux Mint 19, Spyder 4, and QGIS 3.12. I have set my PYTHONPATH to point to where qgis.core is saved. When I try to import qgis.core in Spyder I get the following error:

ImportError: libqgis_core.so.3.12.1: cannot open shared object file: No such file or directory

I have located the file and set the PYTHONPATH to point to the location where it is saved, but that didn't work. PYTHONPATH = /var/lib/flatpak/app/org.qgis.qgis/x86_64/stable/203aec307e8e8084b84920ccf7588cb9182ebbd6d1d134d27a087b5d01d8a3df/files/share/qgis/python:/var/lib/flatpak/app/org.qgis.qgis/x86_64/stable/203aec307e8e8084b84920ccf7588cb9182ebbd6d1d134d27a087b5d01d8a3df/files/lib

Nothing that I have found anywhere on Stack Exchange, Stack Overflow, Reddit, or GitHub has worked to fix this issue. The steps in the PyQGIS Cookbook also failed. Any help on this issue would be greatly appreciated, as I need this to graduate on time.


texts: Writing a shapefile using OGR, from Shapely Geometry - No feature added Error

texts:
I am trying to write a shapefile from a Shapely geometry using PyQGIS in the QGIS 2.6.0 console.

Lines of code before the OGR method to write shapefile:

extnt = poly.extent() # gives a QtRectangle output -> extnt
x_min = extnt.xMinimum()
x_max = extnt.xMaximum()
y_min = extnt.yMinimum()
y_max = extnt.yMaximum()

A_bbox = QgsPoint(x_min, y_min)
B_bbox = QgsPoint(x_max, y_min)
C_bbox = QgsPoint(x_max, y_max)
D_bbox = QgsPoint(x_min, y_max)

where poly is a polygon shapefile.

Source for the Code below :- How to write Shapely geometries to shapefiles? and GitHub post

linelyr = ogr.Geometry(ogr.wkbLineString)
linelyr.AddPoint(x_min, y_min)
linelyr.AddPoint(x_min, y_max)
out_line = ogr.Geometry(ogr.wkbLineString)
out_line.AddGeometry(linelyr)
out_line.ExportToWkb

driver = ogr.GetDriverByName("Esri Shapefile")
ds = driver.CreateDataSource("outputlocation.shp")
layr1 = ds.CreateLayer('', None, ogr.wkbLineString)

layr1.CreateField(ogr.FieldDefn('id', ogr.OFTInteger))
defn = layr1.GetLayerDefn()

feat = ogr.Feature(defn)
feat.SetField('id', 1)

geom = ogr.CreateGeometryFromWkb(out_line.wkb) # or maybe just 'outline'
feat.SetGeometry(geom)

layr1.CreateFeature(feat)
ds = layr1 = feat = geom = None

The shapefile is created, along with all the necessary files, meaning no errors. But the shapefile, when opened in QGIS has no features. Just a constructed attribute table.

Pls assist...

P. S. Any Other solution besides OGR usage also works, as long as it can be coded in the Python console of QGIS


texts: Using Point as the geom type in pgrouting?

texts:
in my table i have the field geom(POINT, srid) and its populated, this table also have the source an target field witch is blank(just like pgrouting guide).

I tried the function pgr_createTopology() and it returns ok, created the vertex_pgr table but it is empty.

Must the geom field be a geom(LINESTRING)?

My problem is: i have some points and need to trace routes between them.

Do i have to make line string between my points?


texts: How get latitude, longitude and altitude from SRTM SWBD and .hgt files in one table

texts:
I know this question may seem repetitive but I am not able find a solution hence I am asking this question .
 1)I have downloaded OSDM data, SRTM SWBD data and also the corresponding .HGT file containing altitude for Sachsen Germany and successfully imported all the data both vector or raster into a postgis database.
 2)I want to combine the Planet_osm_line data along with the swbd data and .hgt data to form a single database which consists of latitude, longitude and altitude and then visualise this using a geoserver.
3) With the help of SQL queries I have found the geometry of all these three data types as planet_osm_line is linestring(xy,xy) , SWBD data is multipolygon(x,y,z,m) and hgt data is polygon(xy,xy)
4) I am not able find the altitude value which is given as raster data for each longitude and latitude points, the altitude file .hgt originally has only 1 row and 3 columns like rid, rast and file name and when transformed using ASTEXT it gives the geometry as polygon(xy,xy) the file name is N50E013.hgt.

What is the correlation between lat/lon and altitude? How to find the altitude for each lat/lon value?

The data of the first rows of each of the datatypes are as follows
   SWBD data :- "MULTIPOLYGON ZM (((12.5719444444444 50.9923611111111 191 0,12.5719444444444 50.9926388888889 191 0,12.5741666666667 50.9926388888889 191 0).

planet_osm_line :- "LINESTRING(13.7424219061054 51.0133491016638,13.7429582886724 51.0134227613141,13.7439974414817 51.0133856445798,13.7448337999365 51.0134121637281,13.7453285287494)

.HGT file "POLYGON((12.9995833333333 51.0004166666667,14.0004166666667 51.0004166666667,14.0004166666667 49.9995833333333,12.9995833333333 49.9995833333333,12.9995833333333 51.0004166666667))"

Can anybody guide me?


texts:
I would like to create a raster where each cell equals the distance to the nearest cell >= the value of the source cell plus a constant. In the example below, the input raster has a grid size of 1m x 1m. Assume we are adding constant of 2 to the cells. For the top left cell, the distance to the nearest cell that contains a value >= 5+2 is 1.0m. I don't think any of the distance tools in ArcGIS will do this.

Input raster:

Contant = 2

Output raster:


texts: Connection Error: Add PostGIS Database on Geoserver through Docker on Geonode:

texts:
I am using GeoNode with connection to Geoserver through Docker.(from here: https://github.com/GeoNode/geonode).
I am now trying to connect to a PostGIS database. I have indicated the host ip address  172.17.0.1 on docker network, using the "ip addr show" in Terminal. The port stays same: 5432.
When trying to connect I get the error message: "Error creating data store, check the parameters. Error message: Unable to obtain connection: Cannot create PoolableConnectionFactory (Connection to 172.17.0.1:5432 refused. Check that the hostname and port are correct and that the postmaster is accepting TCP/IP connections.)".
What could be the solution?


texts: How to use/interpret/visualise a distance matrix with unknown coordinates?

texts: How would you go about using a distance matrix where points have no coordinates, i.e. what useful statistics can be made and how to visualise it in some useful way?


texts: Problem in loading map to Georeferencer in QGIS

texts:
I'm trying to load a map to the Georeferencer in QGIS using the approach from http://www.qgistutorials.com/en/docs/georeferencing_basics.html. I try to get the locality coordinates(lat/lon) of the 'Stetson Pit' and 'Ponzer' on the map.

However, when I set all the points on the map and start georeferencing, it shows "transform failed: failed to compute GCP transform" transform not solvable"

May I know if it is possible to get coordinates of those sites from this map. Am I setting the CRS wrong (WQS84)??


texts: QGIS equivalent of Calculate Adjacent Fields for Atlas

texts:
Is there a QGIS plugin available or a processing model to calculate and populate fields that can be used to label the adjacent pages in a map book?

Although there are tools to create vector grids, I cannot find one to calculate the adjacent polygons. i.e. which adjoining polygon is N,S,E and W of the existing polygon, or even better in all 8 directions surrounding a regular polygon.

This is how it works in ArcGIS;

http://resources.arcgis.com/en/help/main/10.1/index.html#//007000000009000000

I am looking to modify a python script by Ujaval Gandhi on http://www.qgistutorials.com/en/docs/find_neighbor_polygons.html and will post further if I have anything to add. (It currently outputs the surrounding polygons in a single field.)


texts: What is the highest spatial resolution for coronavirus cases?

texts:
I understand that the JHU CSSE Dashboard became the watering hole for public health professionals and data science enthusiasts alike, and it visualizes the spread of novel coronavirus (SARS-CoV-2) infections over countries.

In this phase of the disease, however, a finer gradation of the spread could be highly informative, and I cannot locate any good feed or repository of data for smaller units and locales than countries (or Chinese provinces and US states). Not surprisingly, many countries report data at the regional level, but I am not aware of anyone collecting/collating/scraping it.

Does it call for a crowdsourced entry of national (regional?) reports? Has it not been crowdsourced yet?


texts: Extending a line (e.g. road) using ArcGIS for Desktop?

texts:
I am using ArcGIS 10.3 for Desktop.

I followed the instructions provided in the Help menu i.e. edit, select feature, edit vertices, and then used the continue feature tool to extend the line.

Why is the continue feature tool greyed out?


texts: Adjusting Transparency of Labels on ArcMap?

texts:
I am pretty new to GIS, and I am trying to figure out how to adjust the transparency of my labels.

I know you have to convert the labels to annotations, but I can't seem to do that.

I right click on the layer and go to "Convert Labels to Annotations". According to my many Google searches, I am supposed save the annotations in a database, not in the map. (I tried in the map and it wouldn't let me adjust transparency anyways.)

My biggest problem is that "Destination" is set to "Unknown", and I don't know how to change this or where my Destination should even be. I attached a screenshot so you can visualize what I am saying.


texts:
What is a GIS Professional? Is it the same as GIS Developer or Analyst?

It's just that the term has been used a lot and I'm not really sure what it means.

Update: Voted to close the question as it does not conform to the Q&A format. I'm still a bit confused by how people just use the term as if everybody knows what it means when people can't agree on what it is though.


texts: Clever way of calculating and storing cumulative rainfall for each polygon of a MULTIPOLYGON sf-object

texts:
What I have is this sf-object:
Simple feature collection with 50 features and 8 fields
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 611398.7 ymin: 5122699 xmax: 745144.8 ymax: 5207956
CRS:            +proj=utm +zone=32 +ellps=WGS84 +units=m +no_defs
First 10 features:
   20171225 20171226 20171227 20171228 20171229 20171230 20171231 20180101                       geometry
1         0      0.0    0.600        0        0        0        0 0.600000 POLYGON ((675867.9 5138308,...
2         0      0.0    0.900        0        0        0        0 1.600000 POLYGON ((671435.1 5122705,...
3         0      0.2    0.000        0        0        0        0 1.600000 POLYGON ((688357.4 5152557,...
4         0      0.7    0.000        0        0        0        0 0.900000 POLYGON ((698123.5 5169814,...
5         0      0.0    0.000        0        0        0        0 0.800000 POLYGON ((672429.8 5197674,...
6         0      0.7    0.000        0        0        0        0 0.900000 POLYGON ((698977.7 5168984,...
7         0      0.0    0.000        0        0        0        0 0.500000 POLYGON ((717004.7 5196926,...
8         0      0.0    0.000        0        0        0        0 3.116667 POLYGON ((622728.9 5152046,...
9         0      0.8    0.025        0        0        0        0 0.900000 POLYGON ((695677.5 5172000,...
10        0      0.9    0.300        0        0        0        0 0.600000 POLYGON ((704247.9 5184909,...

Where each column (minus the geometry-column) is a date, and each row is a polygon. The cells are mean rainfall values for that day and polygon. Now I'd like to calculate the cumulative rainfall for each day and polygon and store that in a clean way in some R-object. Up to now I just worked with spatial objects that either were just one point or one polygon. There I could simply do something like:
# get the geometry column
    geom = df %>% st_geometry()

    df_long = df %>%
      st_drop_geometry() %>%
      pivot_longer(cols = everything(), names_to="dates", values_to="precip") %>%
      mutate(dates = as.Date(dates, "%Y%m%d")) %>%
      # get the cumulative count for the days
      mutate(accumulated = cumsum(precip)) %>%
      mutate(geom = geom) %>%
      st_as_sf()

to get something simple like this:
  dates      precip accumulated
* <date>      <dbl>       <dbl>
1 2017-12-25  0           0
2 2017-12-26  0           0
3 2017-12-27  0.600       0.600
4 2017-12-28  0           0.600
5 2017-12-29  0           0.600
6 2017-12-30  0           0.600
7 2017-12-31  0           0.600
8 2018-01-01  0.600       1.20

However, when I have a spatial object of multiple points or polygons, things get quite messy in my head. I could maybe split the sf-object into nrow(df)-sf-objects and apply the same operation and store them in a list. But I don't feel like this is the best way to go. Maybe someone has some kind of idea on how to approach this:)


texts:
I have my map running over GeoServer and OpenLayers, but I want to add a search box to find places in it. Well not specifically places but directions.

The map reads the layers directly over GeoServer (without PostGIS) but I'm a bit lost with the search thing.

How can I implement this feature and what tools are appropriate for it?

Now I will explain how the searching process must be:

1.- I have layers with streets, blocks and house numbers.

2.- the user will search by street name or house direction (house number and street name) and the map will show a list of concordances, if the user select one result on the list it will be zoomed to the destination. Very similar to an ArcGIS search tools.

Any help will be appreciated.


texts: I am using ArcGIS 9.3 without any extensions.  How do you get your map to rotate?  I've done it before but it has been awhile and now the "rotate" does not seem to function now.  I've tried using while in edit and non-edit.  Any suggestions?


texts: ArcGIS remove negative values from DEM, DEM max increases

texts: I have a D.E.M. that has values between -1.13 and 23.33. In the Symbology, I changed the Stretch Type to Minimum-Maximum to get the true values and those previous values are what the raster draws. I wanted to eliminate negative values and found other questions that dealt with getting rid of negative values. One answer said use SetNull(raster <0, raster). I also tried Con. The output from these operations create new files with increased max values. When I used Con, the max went to 64.01 and when I used SetNull, the max increased to 43.4. I don't know why this is happening. The raster is a 32 Bit floating point Tiff.


texts: NA values in dbf when exporting shapefile from R with writeOGR?

texts:
I´m writing an esri shapefile from R (2.14.1) on windows, using writeOGR in the rgdal package, in the resulting dbf file, NA values in the spatialpointdataframe (in numeric columns) are not translated to null values but to 0, and sometime the whole column results in 0 values.

I´ve seen this thread, but nothing more recent,  can this be resolved in other ways, than converting NAs to a known value before and then back again after writing the .shp?

http://osgeo-org.1560.n6.nabble.com/GRASS-stats-Re-R-sig-Geo-writing-shapefiles-DBF-files-when-input-data-contains-NA-td4049981.html

# Example
require(sp)
require(rgdal)

# example dataframe
x <- data.frame(ID = 1:10, Lon = seq (-85,-76,1), Lat = seq(-9, 0, 1), Data=c(1:5, NA, NA, 7, 9, 10))

#create spatialpointsdataframe
coordinates(x) <- c("Lon", "Lat")
proj4string(x) <- CRS("+proj=longlat +datum=WGS84")

#write to esri shapefile
writeOGR(x, dsn=getwd(), layer="example", driver="ESRI Shapefile")

resulting dbf table:


texts: Define polygon colour by a field with CYMK data separated by a comma using QGIS 2.8.2

texts:
I have a similar problem to the following post, but the answers to that question did not assist with my problem.

Data-defined Styles in QGIS

I would like to colour the polygons in a shapefile using the BASECOLOUR field of the shapefile's attribute table.  The BASECOLOUR field is a string that specifies the CYMK colours delimited by commas.


texts: Plotting a map of New Zealand with regional boundaries in R

texts:
I have been asked to make a map of New Zealand using R with regional boundaries. I need something exactly like this example of a map for Switzerland:
http://blog.revolutionanalytics.com/2009/10/geographic-maps-in-r.html

I tried copying the same commands but after I downloaded the data from gadm and tried plotting it with spplot, I ended up with ten tiny plots with New Zealand on the extreme right hand edge. I can't paste the output into this question because R will not let me save it as a file. But what I did was simply download the "level 1" data from gadm in R format and then spplot(gadm).

Please could someone help me to make a map of New Zealand showing regional boundaries, so that I can make the regions different colours? It looks like the main problems at the moment are

I would like the North and South Islands to be in the middle.
I don't want a super-detailed map as I would like to plot it quickly.
I just want a single map, not ten of them.

Alternatively, any pointers to tutorials would be very welcome. I don't necessarily want to use the gadm data if there is a better alternative.


texts: What's the reason for this Invalid Layer "Cannot open GDAL dataset ..." error

texts:
When adding a rasterfile with QGIS python console there appears an error when the first letter of a rasterfile to be added is not written in capital letters. Why this behaviour?

Example case1:

>>> qgis.utils.iface.addRasterLayer("C:\geodaten\natural_earth.tif", "raster")

Error:

  Invalid Layer "Cannot open GDAL dataset C:\geodaten atural_earth.tif:
  `C:\geodaten atural_earth.tif' does not exist in the file system, and
  is not recognised as a supported dataset name."

Example case2:

>>> qgis.utils.iface.addRasterLayer("C:\geodaten\Natural_earth.tif", "raster")

No error!

Rasterlayer is added properly.

With vectorlayers the above mentioned problem doesn't exist.
The QGIS version in use is 2.2.0


texts:
Does anyone know how one can change the Fill style option in QGIS?

I need to use a diagonal pattern that does exist, but I need to adjust it. I have a cartography specification that defines the diagonal pattern to have a specific angle, line width and distance between diagonal lines.

Anyone has experience with this?

Best regards


texts: Creating points with co-ordinates in QGIS to upload to GPS

texts:
I have generated a random point within a set area. I have done this using the Google Maps plugin layer. The CRS is set to WGS 84.

I want to get the co-ordinates for this point so I went into the attribute table and used the $X and $Y functions which gave me some numbers. I read these were in UTM so used an online converter to put them into lat long but they always end up being in the ocean somewhere.

All I want to do is to create a random point in QGIS, get the co-ordinates, save it as a GPX file which I have tried to do also, and then upload it into my Garmin GPS using the GPS Toolbox. When I save it as a GPX file it still never comes up in the upload to GPS dropdown box but I am not sure this would even matter as the co-ordinates never seem to be accurate.


texts: ST_Intersects one polygon with layer - always returns 0 results

texts:
I'm trying to write a query to find if any features in a layer intersect with a certain polygon, but it is always returning 0 rows. I've tried messing with SRIDs, different geom types, sizes, shapes, etc. and have never had it return a single result, making me think there is just a basic misunderstanding I'm having.

The geometry data in the database is stored under the "geometry" column in the WKT format with the 3857 SRID. All I want to is to get returned all the rows whose geometry touches the polygon at all.

The query I wrote looks like this:

SELECT *
FROM tableOfFeatures
WHERE ST_Intersects(ST_GeomFromText('SRID=3857;POLYGON((-118.09636448889168 83.80774305261292, -103.722647355961 -83.72028886793554, -13.87670908426287 -83.7764092545634, -38.85866039455652 84.01409674357751, -118.09636448889168 83.80774305261292))'), geometry);

Am I missing some kind of conversion or cast?


texts: PyProj incorrect transformation from EPSG:25832 to EPSG:4326

texts:
I have the following code to transform coordinates:

from pyproj import Proj, transform

def transform_coordinates(coordinates):
    (x, y) = coordinates
    inProj  = Proj('EPSG:25832')
    outProj = Proj('EPSG:4326')
    lat, lon = transform(inProj,outProj,x,y)
    print(x,y, " -> ", lat, lon)
    return (lon, lat)

transform_coordinates((6559455.75409500021487474, 278682.56987731397384778))

I get the following:

6559455.754095 278682.569877314  ->  1.6901601579676362 56.71184635815745

When according to epsg.io, it should be (5.1333339, 59.1166659)

I could swear it was working yesterday.

EDIT:

inProj: Proj('+proj=utm +zone=32 +ellps=GRS80 +units=m +no_defs', preserve_units=True)

outProj: Proj('+proj=longlat +datum=WGS84 +no_defs', preserve_units=True)


texts: GeoServer & QGIS - Using GeoJSON in a WFS request

texts:
I'm trying to load a WFS Layer in QGIS. I requst the data from my GeoServer and would like to import the data as GeoJSON.

I can export my data in the layer preview section as GeoJSON. Everything works fine.

But i am unable to display my data in QGIS. I tried to add the "&outputFormat=application%2Fjson" at the GetFeature request.

I think my issue belongs to this similar question but the author did not answer the question of @nmtoken if the WFS GetCapabilities request for the service is listed GeoJSON as a GetFeature ResultFormat.
In my case, the format is listed.

Here is the output of the GetCapabilities request:

<ows:Operation name="GetFeature">
   .
   .
   <ows:Parameter name="outputFormat">
     <ows:AllowedValues>
      <ows:Value>application/gml+xml; version=3.2</ows:Value>
      <ows:Value>GML2</ows:Value>
      <ows:Value>KML</ows:Value>
      <ows:Value>SHAPE-ZIP</ows:Value>
      <ows:Value>application/json</ows:Value>
      <ows:Value>application/vnd.google-earth.kml xml</ows:Value>
      <ows:Value>application/vnd.google-earth.kml+xml</ows:Value>
      <ows:Value>csv</ows:Value>
      <ows:Value>gml3</ows:Value>
      <ows:Value>gml32</ows:Value>
      <ows:Value>json</ows:Value>
      <ows:Value>text/xml; subtype=gml/2.1.2</ows:Value>
      <ows:Value>text/xml; subtype=gml/3.1.1</ows:Value>
      <ows:Value>text/xml; subtype=gml/3.2</ows:Value>
     </ows:AllowedValues>
   </ows:Parameter>
   .
   .

Please give me a helping hand here, maybe it's not possible to load my data as GeoJSON through a WFS because it will always provide the data as XML? I tested it with WFS 1.0.0, 1.1.0 and 2.0.0

Versions:

QGIS 3.2.
GeoServer 2.12.1


texts: Getting distances between mutiple points and one "central" point using QGIS?

texts:
I have a GPS data including static measurements and would like to look more in detail on those static ones.

Is it possible to take the point cluster from one measurement (I can extract it by time attribute), manually specify a "central" point or "real position" (as I know where the GPS was placed) and relatively easily to calculate distances between this central point and the measured ones?

Or at least basic distance statistics (max, min, average). I would like to know this without need to calculate distance for each point pair.

Is there any easy way to do it in QGIS?

I prefer GUI if possible.


texts:
I am performing a Near analysis by group (https://blogs.esri.com/esri/arcgis/2010/09/16/nearbygroup/) on points:

Points following each others linearly represent structures. They have same attribute. I would like to do a Near analysis by group, but instead of analyzing groups of same attribute, I want to analyze points with different attributes.

I cannot find any way to do this except maybe changing the core python code of the toolbox?

My goal is to statistically get the averaged distance at which structures meet closely each other the most often. This distance will be used as bandwidth (search radius) for a kernel density analysis on ArcGIS. I try to optimize this bandwidth.


texts: What are Esri's new ArcGIS for Professionals and ArcGIS Pro?

texts:
After looking through all the stuff that was released or talked about at the Esri UC 2013  I am still genuinely confused by what ArcGIS for Professionals is. I will admit that it looks pretty cool, and has some sweet features, but I just do not know how it is different from ArcMap or ArcGIS Online.

Does anybody have any clearer answers beyond the flashy video and homepage?


texts: Openlayers error after removing selectFeature control

texts:
I am using a selectFeature control for a point layer in OpenLayers 2.12. I have a problem after removing the selectFeature control from the map.

The reason why I want to remove the control is because the target layer of the control is destroyed and a new one is added every time the map is moved.

In this process I also remove the selectFeature control from the map and add a new one for the updated layer. But that does not work properly. After the first move/zoom event of the map get the following error:

Error: getFeatureFromEvent called on layer with no renderer. This usually means you destroyed a layer, but not some handler which is associated with it.

Here is a version of the function refresh, which is called on any pan/zoom event of the map:

var mylayer;
var mycontrol;
refresh = function(){

     // Remove previous layers and controls
    if(mylayer != null){
        mylayer.destroy();
        map.removeControl(mycontrol);
    }

    // Get mylayer
    mylayer = new OpenLayers.Layer.Vector("GeoJSON",
        ...
    );

    // Add layer and selectFeature control
    map.addLayer(mylayer);
    mycontrol = new OpenLayers.Control.SelectFeature(mylayer,{
        onSelect: function(e) {
        $("#tmpOutput1").html(e.attributes.name + "\n" + e.attributes.id);
        }
    });
    map.addControl(mycontrol);
    mycontrol.activate();

}

UPDATE

To solve the problem I had to change the layer/control removal to this:

if(mylayer != null){
    map.removeControl(mycontrol);
    mycontrol.destroy();

    map.removeLayer(mylayer);
    mylayer.destroy();
}

Note that here the order matters (at least in my example): first remove the control, then remove the layer.


texts:
How can I ensure that the SQL query will return a field containing integers and not text?

I am using a Spatialite DB for my QGIS project (QGIS 2.18). One of the task is to count the occurrences of values in an attribute (group), to join the resulting view table ("countvalues") to a spatial layer, and display the values based on the joined values.

Example of the grouping query in QGIS DB Browser:

CREATE VIEW countvalues AS
SELECT "items"."group", COUNT("items"."group") AS "cnt"
FROM "items"
GROUP BY "items"."group"

However, the problem is that the type of the type of the newly created column ("cnt") is unknown. So, the symbology in QGIS will allow to symbolize the values as "Categorized" but not as "Graduated". (I assume this requires numbers)

I tried the follwing to cast the field values to integer numbers:

CREATE VIEW countvalues AS
SELECT "items"."group", CAST(COUNT("items"."group") AS INTEGER) AS "cnt"
FROM "items"
GROUP BY "items"."group"

This did not work. If using REAL instead of INTEGER I can see that the values in the table contain a period. However, even then the data type of the "group" and "cnt" column remain empty in the field descripion of the DB Browser .


texts:
We want to periodically poll various WMS or WFS services and do something if the source data has changed. Is there any way to do any of the following:

find the date of the most recent update to a data layer
find whether data has changed since a previous request
get some kind of hash indicating a summary of the data at present

I don't see anything in the specs, but perhaps there is some other mechanism I'm not aware of?


texts: st_transform on planet_osm_line to get lat and lon coordinates - extra data returned


texts: Removing lines from inside polygon in QGIS

texts:
Through a "dissolve" process in QGIS, I achieved some defective polygons generated after dissolve process that have some fails, like these in the screenshot:

Is there any automatic process I can make in QGIS to remove those lines inside the polygon and leave only the outer borders? I can manually delete using the "Node Tool", but I was looking for an automated solution.

The polygon seems to be closed, as it can be seen:


texts:
I've been looking for a way to convert txt file columns to shapefile's. With the code below I am able to divide txt columns, what I don't know is how to store all those columns in a point shapefile, any idea? I am still looking for the answer.

f = open(txtRuta)
contents = f.read()
f.close()
new_contents = contents.replace('\n', ' ')
new_contents = new_contents.strip(' ')

f = open(txtRuta, 'w')
f.write(new_contents)
f.close()

new_contents = re.sub(r" + ", " ",  new_contents)

itemsInString = new_contents.split(" ")
lon= ((len(itemsInString)/9))

idx = []
aa = []
gpstime = []
x = []
y = []
z = []
bb = []
cc = []
dd = []

numbeer = 0

for index, item in enumerate(itemsInString):
    if numbeer == 9:
        numbeer = 0

    if numbeer == 0:
        idx.append(int(item))
    elif numbeer == 1:
        aa.append(float(item))
    elif numbeer == 2:
        gpstime.append(int(item))
    elif numbeer == 3:
        x.append(float(item))
    elif numbeer == 4:
        y.append(float(item))
    elif numbeer == 5:
        z.append(float(item))
    elif numbeer == 6:
        bb.append(float(item))
    elif numbeer == 7:
        cc.append(float(item))
    elif numbeer == 8:
        dd.append(float(item))
    numbeer = numbeer + 1
for i in range(0, lon):

    line = ogr.Geometry(ogr.wkbLineString)
    line.AddPoint(x[i], y[i], z[i])
    print line.ExportToWkt()

EDIT:

From GDAL/OGR cookbook, code below works for .txt to .shp, configure it as you need it.

f = open(txtRuta)
contents = f.read()
f.close()
new_contents = contents.replace('\n', ' ')
new_contents = new_contents.strip(' ')

f = open(txtRuta, 'w')
f.write(new_contents)
f.close()

new_contents = re.sub(r" + ", " ",  new_contents)

itemsInString = new_contents.split(" ")
lon= ((len(itemsInString)/9))

idx = []
aa = []
gpstime = []
x = []
y = []
z = []
bb = []
cc = []
dd = []
larg = len(idx)

numbeer = 0

for index, item in enumerate(itemsInString):
    if numbeer == 9:
        numbeer = 0

    if numbeer == 0:
        idx.append(int(item))
    elif numbeer == 1:
        aa.append(float(item))
    elif numbeer == 2:
        gpstime.append(int(item))
    elif numbeer == 3:
        x.append(float(item))
    elif numbeer == 4:
        y.append(float(item))
    elif numbeer == 5:
        z.append(float(item))
    elif numbeer == 6:
        bb.append(float(item))
    elif numbeer == 7:
        cc.append(float(item))
    elif numbeer == 8:
        dd.append(float(item))
    numbeer = numbeer + 1

shpout = shpRuta + "\\" + pNombre + ".shp"
driver = ogr.GetDriverByName("ESRI Shapefile")
data_source = driver.CreateDataSource(shpout)
srs = osr.SpatialReference()
srs.ImportFromEPSG(int(epsgUse))

layer = data_source.CreateLayer(pNombre, srs, ogr.wkbPoint)

layer.CreateField(ogr.FieldDefn("id", ogr.OFTInteger))
layer.CreateField(ogr.FieldDefn("aa", ogr.OFTReal))
layer.CreateField(ogr.FieldDefn("gpstime", ogr.OFTInteger))
layer.CreateField(ogr.FieldDefn("x", ogr.OFTReal))
layer.CreateField(ogr.FieldDefn("y", ogr.OFTReal))
layer.CreateField(ogr.FieldDefn("z", ogr.OFTReal))
layer.CreateField(ogr.FieldDefn("bb", ogr.OFTReal))
layer.CreateField(ogr.FieldDefn("cc", ogr.OFTReal))
layer.CreateField(ogr.FieldDefn("dd", ogr.OFTReal))

for i in range(0, lon):

    feature = ogr.Feature(layer.GetLayerDefn())
    feature.SetField("id", idx[i])
    feature.SetField("aa", aa[i])
    feature.SetField("gpstime", gpstime[i])
    feature.SetField("x", x[i])
    feature.SetField("y", y[i])
    feature.SetField("z", z[i])
    feature.SetField("bb", bb[i])
    feature.SetField("cc", cc[i])
    feature.SetField("dd", dd[i])

    wkt = "POINT(%f %f)" % (float(x[i]) , float(y[i]))

    point = ogr.CreateGeometryFromWkt(wkt)

    feature.SetGeometry(point)
    layer.CreateFeature(feature)
    feature.Destroy()

data_source.Destroy()


texts: Modify features not updating on vector layer

texts:
I have point vector layer in my application.Now I am trying to modify its coordinates. I am able to do it but its not reflecting on layer.
Here is code I am using :

var selectInt = new ol.interaction.Select({
    wrapX: false
});

var modify = new ol.interaction.Modify({
    features: selectInt.getFeatures()
});

modify.on('modifyend', function (e) {

    var features1 = rdprvectorlayer.getFeatures();
    var featureToUpdate;
    var features = e.features.getArray();
   var id=e.features.getArray()[0].get('id');
   for (var i = 0; i < features1.length; i++) {

       if (features1[i].get('id') == id) {

           featureToUpdate = features1[i];
           break;
       }
   }

   var coordinates = e.features.getArray()[0].getGeometry().getCoordinates();

   var lon =coordinates[0];;
   var lat = coordinates[1];

   featureToUpdate.set('geometry', new ol.geom.Point(getPointFromLongLat(lon, lat)));

   var fea = new ol.Feature({
       geometry: new ol.geom.Point(getPointFromLongLat(lon, lat))
   })

   rdprvectorlayer.addFeatures([fea]);

})

function getPointFromLongLat(long, lat) {
    return ol.proj.transform([long, lat], 'EPSG:4326', 'EPSG:3857')
}

but after adding to vector laying nothing is displaying on map.
can any one help me to solve it.


texts: Importing Shapefile to File GDB with domains and coded values

texts: If I have a shapefile with attributes and values that I need. What happens when I import it into a GDB with domains for that attribute? Say in my shapefile I have a feature with attribute 'TYPE' with a value of 'HOME' and import it into the GDB the value is 'HOUSE', does it convert my value 'HOME' to 'HOUSE' to maintain consistency?


texts: Is there a faster way in python of finding the smallest number in a field?

texts:
Using arcgis desktop 10.3.1 I have a script which uses a search cursor to append values to a list and then use min() to find the smallest integer. The variable is then used in a script. The Feature class has 200,000 rows and the script takes a very long time to complete. Is there a way to do this quicker? At the moment I think I would just do it by hand rather than write a script due to the length of time it takes.

import arcpy
fc = arcpy.env.workspace = arcpy.GetParameterAsText(0)
Xfield = "XKoordInt"
cursor = arcpy.SearchCursor(fc)
ListVal = []
for row in cursor:
    ListVal.append(row.getValue(Xfield))
value = min(ListVal)-20
print value
expression = "(!XKoordInt!-{0})/20".format(value)
arcpy.CalculateField_management (fc, "Matrix_Z" ,expression, "PYTHON")


texts: spatial queries by webgis and virtualizing

texts:
i have a question regarding this tutorial:

http://www.gistutor.com/postgresqlpostgis/10-intermediate-postgresqlpostgis-tutorials/75-understanding-postgis-spatial-queries.html

As u see, it demonstrates how spatial sql statements are working. But my question is, how this tutorial with this webgis was built up.

I know that openlayers,json,php,wfs-t,ajax and postgis comes into the play but how all these componentes are working together. Is there any good tutorial out there or do somebody can gives workshops? thank you

regards
ulan


texts: Get layer name by click on features OpenLayers3

texts:
I want to add selec interactions to my OL3 map.
Adding two different select interactions like this is not working. Only the second one works.
Like this:

        var select = new ol.interaction.Select({
        layers: [clusters]
    });
select.getFeatures().on('add', function (e) {alert("orange")});
    map.addInteraction(select);
    var select2 = new ol.interaction.Select({
        layers: [clustersSarv]
    });
select2.getFeatures().on('add', function (e) {alert("blue")});
    map.addInteraction(select2);

Now I am trying to do something like this, but I am unable to get layer on which is clicked, to define different actions for different layers.

var select = new ol.interaction.Select({
        layers: [clusters,clusters2]
    });
select.getFeatures().on('add', function (e) {
   console.log(e.element);

// if layer is clusters  alert("orange")});
    map.addInteraction(select);

Any suggestions? Main question is how to get a layer name by clicking on layer's features.


texts: Changing nodata pixels of -3.40282347e+38 to a different number in QGIS

texts:
I have a raster file with two different nodata values: -32768 and -3.40282347e+38.  This happened because I am using -32768 as my nodata value, and the georeferencer in QGIS adds nodata pixels as -3.40282347e+38.  I would like to have all my nodata pixels as a single value, which I would like to be -32768.

If the no data value is -32768 (or any other reasonable number), it is easy to change the value of the no data pixels into any other number using the raster calculator in QGIS.

NOTE: For this to work, you first have to uncheck "No data value" in the Transparency window of Layer Properties, otherwise the raster calculator won't recognize the no data pixels.

To change nodata pixels with a value of -32768 to -999 using the raster calculator, you can use the following formula:

("RasterFileName@1"=-32768)*(-999)+("RasterFileName@1"!=-32768)*("RasterFileName@1")

However, this will not work for nodata pixels with a value of -3.40282347e+38.  The raster calculator does not recognize these pixels, and they do not get changed.

I thought that this might because the raster calculator did not like the format of -3.40282347e+38, so I tried using the raster calculator to change any pixel that was less than -10000, using the following formula:

("RasterFileName@1"<-10000)*(-999)+("RasterFileName@1">=-10000)*("RasterFileName@1")

This correctly finds all the pixels that are -32768, but still ignores the pixels with a value of -3.40282347e+38.

Is there any way to get the raster calculator to recognize pixels with a value of -3.40282347e+38?  Or some other strategy of getting rid of this noxious number?


texts: Geoserver CQL + SQL Server Plugin - .Filter Statement

texts:
Does anyone know why GeoServer's SQL Server Plugin when making a WFS request with a spatial cql/ecql predicate ends up adding a SQL .Filter to the generated query in addition to the spatial predicate included in the request? And/or how to make it not do so?

Example:

cql_filter=WITHIN(the_geom, POLYGON((-117.86681324243546 33.63613597564614,-117.86551505327225 33.63613597564614,-117.86551505327225 33.635347676345916,-117.86681324243546 33.635347676345916,-117.86681324243546 33.63613597564614)))

Translates to the following SQL Query:

SELECT count( * ) FROM (SELECT * FROM dbo.Whatever) as "vtable"
WHERE  "the_geom".Filter(geometry::STGeomFromText('POLYGON ((-117.86681324243546 33.63613597564614, -117.86551505327225 33.63613597564614, -117.86551505327225 33.635347676345916, -117.86681324243546 33.635347676345916, -117.86681324243546 33.63613597564614))', 4269)) = 1
AND "the_geom".STWithin(geometry::STGeomFromText('POLYGON ((-117.86681324243546 33.63613597564614, -117.86551505327225 33.63613597564614, -117.86551505327225 33.635347676345916, -117.86681324243546 33.635347676345916, -117.86681324243546 33.63613597564614))', 4269)) = 1

It seems that .Filter is superfluous and in many cases drastically increases the query execution time (compared to the same query with just one or the other). I can't seem to figure out why this would be necessary nor how to tell it to stop.

From my testing I'm finding that every request that has a spatial cql/ecql predicate that executes a layer from a SQL Server DataStore does this. Meaning that all of our spatial WFS queries have an additional WHERE clause that potentially increases their execution time. We've found this to be true on several different versions of GeoServer (2.10.0, 2.12.1), all versions of WFS (1.0, 1.1.1, 2.0) and changing the spatial predicate (CONTAINS, WITHIN, INTERSECTS, etc.) doesn't seem to change anything. When the spatial predicate is removed, so are BOTH the spatial SQL statements.

Full Example Request:

http://localhost:8080/geoserver/wfs?service=wfs&REQUEST=GetFeature&OUTPUTFORMAT=application/json&version=2.0&typeName=Whatever&cql_filter=WITHIN(the_geom, POLYGON((-117.86681324243546 33.63613597564614,-117.86551505327225 33.63613597564614,-117.86551505327225 33.635347676345916,-117.86681324243546 33.635347676345916,-117.86681324243546 33.63613597564614)))

UPDATE
Based on Andrea's information below, perhaps the "superfluous" categorization is going a bit far. I personally don't know if SQL 2016 still requires a .Filter to force spatial index usage...so it's entirely possible they're right on the money and it does need to be there for that purpose. In any case, this explains WHY this happens and hopefully will help someone else who is scratching their head on this one in the future.


texts: QGIS Atlas generation - 2 panels with different sections both Atlas generated

texts:
I have a polyline that I split in 5000 m sections and I made it to generate an Atlas for the first step, but since the sections that I want to show are long and not wide I would like to add sections, e.g. Atlas Map page1 shows section 1 and 2 (instead of the same green section1, the second panel should show orange section2.
I tested filtering the Atlas, but I got stuck and wasnt sucessfull. I created 2 virtual layers containing every second section, but after I added the second panel I could not find out if it is possible to let the panels be controlled by different layers for  Atlas generation. Any idea how to get 2 panels with continuing  sections?


texts: Area behind polygon starting from a point (radial distancem like viewshed analysis)

texts:
I'm looking for a way to create a polygon (grey area in picture) behind another polygon (green), just like a shadow from a point / another polygon (blue).

To make it clearer here's what I mean:

It seems obvious to create a 100 m buffer around the green polygon, but what would be the next step? If it's only one polygon, I'd cut it manually, but it's +100.

I'm using ArcMap Basic licence without any extensions. And no coding skills.


texts: MongoDB to store vehicle tracking information

texts:
I am building an app similar to Uber for tracking vehicles. Since the update frequency is so high (accounting for several users), I want to know the general practices involved for making writes faster to mongodb collection.

I am maintaining a database to store historical location information from all vehicles but it is bound to grow very fast once we go into production. I need to get the list of vehicles closest to a point. For this should I implement a separate table (with one row per vehicle) which gets updated after every update or there is better/faster way to do this using the existing table?

Also I would like to know if there is a more efficient & faster database out there which can server my purpose.

Right now both the tables mentioned above have same structure

id
vehicle_id
loc
timestamp

Links to articles with best practices are welcome

Is there a limit to number of inserts / updates mongodb can handle per second?
In this video Curtis from Uber mentions that they discarded mongodb as they reached its limits soon after they went into production. I couldn't find any information on what they replaced it with.

I am expecting around 10000+ updates every second. Will this design still work for 100000 updates?


texts:
I was trying to include the v.in.wfs tool in my processing modeler however whenever I run the model this tool gives me an error.

While trying to connect the same link through layer--> add layer --> add WFS layer it works

WFS link I'm using: http://202.159.223.250/erdas-apollo/vector/INDIASTATES


texts: Write layer style .qml as predefined within a gpkg using R

texts:
I have this code to export a layer (Pbm) to a gpkg called "Atlas":

writeOGR(Pbm, paste0(normalizePath(OutputFolder),'\\', "Atlas.gpkg"),"Atlas",
             driver="GPKG", encoding = "UTF-8", overwrite_layer = TRUE)

In other hand I have a .qml file with the QGIS style I like and which I want to save within the .gpkg. Does it exist the possibility to write this style within the .gpkg directly from R ?

PD: I already know that copying the .qml file with the same name as the .gpkg layer will add that style directly when loading into QGIS. The issue here is that I don't want that extra file out of the .gpkg.


texts: PyQGIS setSubsetOfFields function in QGIS 2.4

texts:
I found the function setSubsetOfFields to filter fields in a layer in the PyGIS cookbook:

http://docs.qgis.org/testing/en/docs/pyqgis_developer_cookbook/vector.html#iterating-over-a-subset-of-features

But it seems like that these function doesn't exist in QGIS 2.4. Neither has the layer a fields() function. There is only a setSubsetOfAttributes function.

I tried to use it like this:

layer = qgis.utils.iface.activeLayer()
request = QgsFeatureRequest()
request.setSubsetOfFields(['id', 'name'], layer.fields())

features = layer.getFeatures(request)
for feature in features:
    print feature.attributes()

Am I using this wrong or is there another possibillitie to set a filter on the layer fields?


texts: Make QGIS use a conda virtual environment for python libraries

texts:
I wanted to install a plugin, but it requires scipy. Thus I realized that QGIS use the default python3 installation. I tried to launch it from a conda virtual environment, as advised in this related stackexchange question, but I end up with the "Couldn't load SIP module." Issue#37783 (which isn't resolved imho).
My QGIS installation is system-wide (/usr/bin/qgis).
How can I make QGIS use a conda virtual environment for python libraries?
QGIS and OS versions

QGIS version
3.16.1-Hannover
QGIS code revision
b381a90dca

Compiled against Qt
5.12.8
Running against Qt
5.12.8

Compiled against GDAL/OGR
3.0.4
Running against GDAL/OGR
3.0.4

Compiled against GEOS
3.8.0-CAPI-1.13.1
Running against GEOS
3.8.0-CAPI-1.13.1

Compiled against SQLite
3.31.1
Running against SQLite
3.31.1

PostgreSQL Client Version
12.5 (Ubuntu 12.5-0ubuntu0.20.04.1)
SpatiaLite Version
4.3.0a

QWT Version
6.1.4
QScintilla2 Version
2.11.2

Compiled against PROJ
6.3.1
Running against PROJ
Rel. 6.3.1, February 10th, 2020

OS Version
Ubuntu 20.04.1 LTS

Active python plugins
openlayers_plugin; MetaSearch; processing


texts: QGIS 3x rotated label returns angle value instead of initial one

texts:
I would like to rotate the text value assigned to my point in QGIS.
For this purpose I used the following threads:
Rotating labels in QGIS 2.8.2
QGIS Labels Rotation
from which the best appears to be this one:
https://gis.stackexchange.com/questions/232727/is-there-a-way-to-interactively-rotate-labels-in-qgis-using-any-of-the-edit-tool#:~:text=Yes%2C%20you%20can%20rotate%20labels,rotation%20using%20the%20label%20tools.
Anyhow, I think I got the mechanism of rotation. Unfortunately, instead of presenting my initial values allocated to the given point, I am getting the angle value.

I would like to admit, that maybe my situation isn't as straightforward as in the given example, because I used the rule-based leveling, where I set the text character based on the "Visible" column in my data attribute table. The string still comes from the "Text" column, but the rule from the "Visible" column. Is that something, which might affect the rotation?
I've picked up the "Text" for the rotation rule. Shall I use the "Visible" variable then?
UPDATE:
Manual changing the label properties doesn't help at all.
When I put the right value, the angle immediately equals this value and other way round.

Is there a way to make the label value independent from the angle value?


texts:
I'm applying a Carto named map layer to an existing Leaflet map similar to:

cartodb.createLayer(map,{
    type: 'namedmap',
    user_name: **username**,
    options: {
        named_map: {
            name: **name**
        }
    }
},{
    https: true
}).on('done',function(layer){
    layer.setAuthToken(**token**);
    layer.addTo(map);
}).on('error', function(err){
    console.log('CartoDB Layer Error: '+err);
});

Our named map templates each cover a city and our UI allows users to change the city without leaving the page, which means I need to update the source being used by this layer.

With a normal Leaflet tileLayer this was as easy as:

layer.setUrl(**URL**)

But is there an equivalent with a Carto layer?

I'd really prefer not completely destroy/recreate the layer from the ground up due to the challenges of keeping its state (hidden/shown), ensuring it stays in the same order in the Leaflet controls, etc.


texts: Slow Elevation Batch Queries in PostGIS?

texts:
I've been exploring PostGIS and have found it to be quite slow and I thought I would ask the community if they notice any blatant flaws in my queries. First, let me say that I need to query elevation points as fast as possible and I am very close at this point to switching to an in-memory data store like Redis to perform this operation if I cannot speed my queries up. My first naive query is the following:

SELECT ST_Value(rast, 'SRID=4326;POINT(-111.5000000 38.5000000)'::geometry)
FROM srtm
WHERE ST_Intersects(rast, 'SRID=4326;POINT(-111.5000000 38.5000000)'::geometry);

Initially, I thought that this would be a bad query because I need to get 10s of thousands of points in a single query, this means that I would have call this query 10s of thousands of times in my software. The first test that I ran showed that I could get a single elevation in just under 100ms (I'm running postgis in a Docker container too so that might slow things down a bit, but I'm looking at relative benchmarks at the moment). So I then posted Querying Thousands of Points with ST_Value()? asking the community if there was a better way to do a batch call to postgis (in hopes of speedings things up since I would only have to make a single query). The result is the following:

SELECT gv.val
FROM srtm , LATERAL ST_Intersection(rast, 'SRID=4326;MULTIPOINT(-111.5000000 38.5000000)'::geometry)  AS gv
WHERE ST_Intersects(rast, 'SRID=4326;MULTIPOINT(-111.5000000 38.5000000)'::geometry);

This query took nearly a second in my current configuration, YIKES! So then I thought, well maybe if I have more points it will scale better, say like 10,000 multipoints. So right off the bat I tried 10,000 points. My system sat there for a very long time and never finished the query before I gave up and ctrl-c'd the process. This was very bad news. So then, I created a multipoint with only 4 entries as shown in the following query:

SELECT (ST_Dump(gv.geom)).geom, gv.val
FROM srtm , LATERAL ST_Intersection(rast, 'SRID=4326;MULTIPOINT(-111.305048568 38.0601633931,-111.822991286 38.6025320796,-111.136977796 38.3631992596,-111.206470006 38.971228396)'::geometry)  AS gv
WHERE ST_Intersects(rast, 'SRID=4326;MULTIPOINT(-111.305048568 38.0601633931,-111.822991286 38.6025320796,-111.136977796 38.3631992596,-111.206470006 38.971228396)'::geometry);

This query took nearly 3 seconds to complete. As a result, I found that my initial query was the fastest possible way to get the data, which is still rather slow.

Can someone explain to me why the multipoint query did not speed things up?

Is there a better batch type call I can make that could get me closer to thousands of queries per second?

If not, I think I may have to switch to using something like Redis and come up with some kind of clever scheme to load billions of elevation points into memory.


texts: How to create a single raster based on the values of three rasters (Conditional)?

texts:
I have three rasters of soil granulometric fractions (clay, sand and silt), all in %. I want to make a single raster of the soil texture classes. For that, I need to make some conditionals (i. e.: when sand > 90 and clay < 10, so the class in the new raster must have the value 1).

Any suggestion? I'm using Arcgis 10.6.


texts: How to fix 'Could not find the Qt platform plugin “Windows” in “”' when implementing qgis in custom Python3 application?

texts:
I am implementing qgis 3.4.3 in a custom application, but when I instantiate the QgsApplication() class, I receive the error "Could not find the Qt platform plugin "Windows" in "".

I've attempted using various installs of qgis 3.x+ including the standalone installer and OSGeo4W web installer. I currently am sticking with the OSGeo4W web installer installation. I am using a Python 3.7 installation separate from what comes with OSGeo4W and attempting to integrate qgis functionality.

I have followed the instructions in the following QGIS help doc under the section: "Using PyQGIS in Custom Applications"
https://docs.qgis.org/2.18/en/docs/pyqgis_developer_cookbook/intro.html#run-python-code-when-qgis-starts

After attempting to run the script I realized a dll plugin was missing. After some research, I found its the qwindows.dll that qt uses. The qwindows.dll is included in the OSGeo4W installation under: C:\OSGeo4W\apps\Qt5\plugins\platforms

I changed the QT_PLUGIN_PATH variable on the local command prompt to include the directory above, but the same error prevailed. I also changed the QT_DEBUG_PLUGINS variable to 1 which printed out the locations Qt is looking for plugins. Interestingly, it wasn't looking for plugins in the path I specified in the QT_PLUGIN_PATH variable.

Python Code:

import sys
#add qgis python libraries to instance python path in case not added
#at the environment variable %PYTHONPATH% level
sys.path.extend(['C:\\OSGeo4W\\apps\\qgis\\python',
'C:\\OSGeo4W\\apps\\Python37\\Lib\\site-packages'])

from qgis.core import *

# supply path to qgis install location
QgsApplication.setPrefixPath('C:\\OSGeo4W\\apps\\qgis', True)

# create a reference to the QgsApplication
# setting the second argument to True enables the GUI, which we need to do
# since this is a custom application
qgs = QgsApplication([], True)

# load providers
qgs.initQgis()

# Write your code here to load some layers, use processing algorithms, etc.

# When your script is complete, call exitQgis() to remove the provider and
# layer registries from memory
qgs.exitQgis()

Batch File Code to Start Cmd Prompt with Proper Variables:

@echo off
SET OSGEO4W_ROOT=C:\OSGeo4W
::Include environment variables set by qgis on startup
call %OSGEO4W_ROOT%\bin\o4w_env.bat
call %OSGEO4W_ROOT%\bin\qt5_env.bat
call %OSGEO4W_ROOT%\bin\py3_env.bat
@echo off
path %PATH%;%OSGEO4W_ROOT%\apps\qgis\bin
path %PATH%;C:\OSGeo4W\apps\Qt5\bin
path %PATH%;C:\OSGeo4W\apps\Python37\Scripts
path %PATH%;C:\OSGeo4W\apps\Python37

set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\qgis\python
set PYTHONHOME=%OSGEO4W_ROOT%\apps\Python37
set QT_PLUGIN_PATH=%OSGEO4W_ROOT%\apps\Qt5\plugins;%OSGEO4W_ROOT%\apps\qgis\qtplugins
set QT_DEBUG_PLUGINS=1

cmd.exe

I expect to run the python script and receive no errors since I'm pointing to the directory where the qwindows.dll missing plugin is stored. However, I am still receiving the missing windows plugin error.

Here is the actual message with QT_DEBUG_PLUGINS set to 1:
QFactoryLoader::QFactoryLoader() checking directory path "C:/OSGeo4W/apps/qgis/qtplugins/platforms" ...
QFactoryLoader::QFactoryLoader() checking directory path "C:/OSGeo4W/apps/Python37/platforms" ...
qt.qpa.plugin: Could not find the Qt platform plugin "windows" in ""
This application failed to start because no Qt platform plugin could be initialized. Reinstalling the application may fix this problem.


texts: Why isn't View>View Settings working in ArcGIS 10 ArcScene?

texts: I am running ArcGIS 10 desktop verison with the 3d Analyst extension to view DEMs in ArcScene. I want to view the DEM from a specific xyz location, which should be achievable using the View Settings dialogue box in the View menu. When I click on View>View Settings nothing happens except the menu closes. I can't find any hint of why this isn't working in the ArcGIS help, forums or blogs. Does anyone know where I'm going wrong?


texts: Reduce.region odd result in GEE, however the map seems to be correct

texts:
I'm using Google Earth Engine to calculate the 5th and 95th percentiles of the intraannual coefficient of variation (CV) of the NDVI within a polygon (my area of interest). My code actually runs and the results seem to be OK when I look the CV map in GEE, but percentiles  within my aoi are exactly the same number when I use a reduce.Region function.

What could be wrong with my code?

//study area
    var geometry = ee.Geometry.Polygon(
              [[-60, -30],
               [-60, -28],
               [-62, -28],
               [-62, -30]],null, false);
    Map.addLayer(geometry, null, 'geometry');
    Map.setCenter(-61, -29, 6);

//temporal subset of MODIS NDVI
    var NDVI = ee.ImageCollection('MODIS/006/MOD13Q1')
    .filterDate('2001-01-01', '2001-12-31')
    .select(['NDVI']);

// NDVI functional attributes
    var NDVImean = NDVI.mean();
    var NDVIsd = NDVI.reduce(ee.Reducer.stdDev());
    var NDVIcv = NDVIsd.divide(NDVImean); // Intraannual Coeficient of variation of NDVI
    Map.addLayer(NDVIcv, null, 'NDVIcv');

//Percentiles of NDVImean within geometry
    var quartM = NDVImean.reduceRegion({
              reducer: ee.Reducer.percentile([5,95]),
              geometry: geometry,
              crs:'EPSG:4326',
              scale: 231.656358,
              bestEffort: true,
              maxPixels: 100000
              });
    print ('percentiles of NDVImean', quartM)

//Percentiles of NDVIcv within geometry
    var quartCV = NDVIcv.reduceRegion({
              reducer: ee.Reducer.percentile([5,95]),
              geometry: geometry,
              crs:'EPSG:4326',
              scale: 231.656358,
              bestEffort: true,
              maxPixels: 100000
              });
    print ('percentiles of NDVIcv', quartCV)

code link: https://code.earthengine.google.com/de6d677c879e4a9b5240e849ff620bdd


texts:
I am trying to export a DEM from QGIS into a SDTS format for somebody to open it in SketchUp. Anyone has any ideas? My understanding is you can do it in Arc, but I don't have an Arc license. Sketchup will take a .dem or a .ddf file. Here is a link for the help section for ways to bring in a DEM into sketchup.
https://help.sketchup.com/en/sketchup/importing-dem-files-terrain


texts: QGIS how to access to constrained clustering plugin? I have downloaded the plugin but I'm wondering where to find it?


texts: How to project feature class located into 3 UTM Zones?

texts: i have created a silverlight API web map application,for Saudi Arabia country, i have created a 3 search button, once search for "X,y", the problem with me is to how to define a probability projection for this country, KSA located into 3 UTM Zones, i i have already added UTM Zone 38 only for all feature classes but when i entered number for another area located into 37 zone the search result not correct it showing me into zone 38 only.


texts: Deleting fields if name contains underscore ("_") with ArcPy?

texts:
I have several duplicate fields with under score characters in a feature class that I would like to delete.  How can this be accomplished using ArcPy?  My code below lists the fields, what steps should be executed next to delete anything in field_names?

fc= "featureclass"
field_names = [f.name for f in arcpy.ListFields(fc)]
print field_names

Output:

[u'OBJECTID', u'JOBID', u'MOBILEID', u'COMPANY', u'ADDRESS', u'CITY', u'ZIPCODE', u'WORKSTATUS', u'WTBDTDATE', u'WTBDTTIME', u'WORKTYPE', u'LATTITUDE', u'LONGITUDE', u'DAYS_OVERDUE', u'OVERDUE_STATUS', u'SHAPE', u'STATUS', u'ZONEID', u'PRIORITY', u'ORIGPRIORITY', u'NEXTJOBID', u'MOBILE', u'DESCRIPTION', u'REGION', u'DISTRICT', u'COMPANY_1', u'SHAPE_LEN_1', u'MOBILE_1', u'DESCRIPTION_1', u'REGION_1', u'DISTRICT_1', u'COMPANY_12', u'SHAPE_LEN_12', u'WORSTOFFENDER', u'SHAPE.AREA', u'SHAPE.LEN']


texts: Geocortex essentials Identity server credentials?

texts:
I am working with Gecortex Essentials 4.4 and HTML5 viewer 2.5. Feature layers are getting loaded to the site from ArcGIS Online account (created for my organization). I have set up geocortex identity server with few users and roles (admin and general)  as GE site should work based on user authentication.

But issue is if i log in with GE Identity server credentials, feature layers never load (as arcgis online log-in is never done).

Is there a way to store arcgis online credentials somewhere in config file and call that setting on site startup to load feature layer?

This way i can use both identity server and can also load feature layers from arcgis online account.


texts: How to find parent-child relationships between geometric features? Or is there another way of dealing with GIS data?

texts:
I have bunch of shapefiles from OpenStreetsMap. A simple understanding of human is that we know which feature is inside which feature. And it's mathematically constrained. In other words, it's impossible for a given tree to be in two cities. Thus if we might define inclusion attribute for geometric and geographical features.

I need to make my data more meaningful by creating relationships between features. Like for example I need to store that this city belongs to this state, or this college is in this city, and so on and so forth.

This seems to be a very natural requirement, thus I searched and to my surprise I found not much about how it works and how to achieve it. For example, Managing MySQL spatial polygons in parent child style? is unanswered 2 years later with no upvotes.

Maybe I'm not thinking GISly and I need to change my understanding.

Is it OK to think hierarchically in GIS?

If so, how to related unrelated data?

If not, what other approaches do we use?


texts: Geoprocessing Service works in ArcMap but not via Server Manager or JavaScript API

texts:
I've published a simple geoprocessing service that takes an uploaded zip file and unzips it. The geoprocessing service was published with "uploads" enabled. I have also enabled "Info" level messaging, but nothing gets logged so that's been useless. When I run the GP service in ArcMap, it runs perfectly fine and does what I expect it to do. When I try running from the rest API, I get a failure message. When I try running from the JavaScript API, I get a failure message also that's a little more detailed, it says I have an Invalid URL but I'm not sure what URL it's evaluating...

Rest API input (I'm not entirely sure I'm inputting it right so it could be that):

After clicking "Submit Job (POST)", I get this:

In the JavaScript API, this is my code:

require([
"esri/map",
"esri/request",
"esri/arcgis/utils",
"esri/layers/FeatureLayer",
"esri/tasks/Geoprocessor",
"dojo/dom",
"dojo/on",
"dojo/domReady!"
],
function(
Map,
esriRequest,
arcgisUtils,
FeatureLayer,
Geoprocessor,
dom,
on
) {

var map = new Map("map", {
  basemap: "hybrid",
  center: [-77.069, 36.950],
  zoom: 14
});

var featureLayer = new FeatureLayer("https://www.myurl.org/xxx/rest/services/HbMonitoringTest/HabitatData/MapServer/1");

map.addLayer(featureLayer);

gp = new Geoprocessor("https://www.myurl.org/xx/rest/services/Web_Map_GP_Services/ZipTest/GPServer");

on(dom.byId("upload"), "click", upload);
function upload(){
  //upload the zip file and get back the itemID (via uploadSucceeded)
  var upload = esri.request({
    url: "https://www.myurl.org/xxx/rest/services/Web_Map_GP_Services/ZipTest/GPServer/uploads/upload",
    form: dojo.byId("uploadForm"),
    content: {f: 'json'},
    handleAs: 'json',
  }).then(uploadSucceeded, uploadFailed);
}

function uploadSucceeded(response) {
  var itemID = response["item"].itemID;
  console.log("File upload successful, item ID: ", itemID);
  var params = {"Input_Zip_File": "{'itemID':" +itemID+ "}" };
  gp.submitJob(params, gpJobComplete, gpJobStatus, function(error){
    console.log(error.jobStatus + '(' + error.jobId + '):' + dojo.toJson(error));
  });
}

function gpJobComplete(result){
  console.log(result.jobStatus + '(' + result.jobId + ')')
}

function gpJobStatus(result){
  console.log(result.jobStatus + '(' + result.jobId + '):' + dojo.toJson(message));
}

function uploadFailed(response){
  console.log("Failed: ", response);
}

});

For the error I get a HTTP 400 error, "Invalid URL":


texts: How to force one layer to stay attached to another layer when moving feature in QGIS

texts:
QGIS 2.18.13
Layer one) Sewer Main  [Layer Type: "MultiLineStringZ"]
Layer two) Sewer Lateral [Layer Type: "LineString"]

When moving the Sewer Main, the sewer laterals do not move with the main.  Instead, the user has to go back to each lateral and manually snap the sewer lateral node to the sewer main.

How can I force the lateral layer to stay attached to the main layer when moving the main layer?

Below I have tried to show a simple example of our project.  When the main is moved, it leaves behind the lateral (as shown in the Not Corrected example).  We are then forced to manually snap every lateral to the main (as shown in the Corrected example).


texts: Is there a built in way to reduce the amount of points in an IElement

texts: I have a polyline element and I want to reduce the number of points (vertices?) that make up the line. I know it will make the curved sections jagged and that is fine. I just want to control (perhaps by a percentage of the total?) how many points are decimated. I looked at the ITopologicalOperator.Simplify method but this does not seem to be what I need.


texts: GeoTIFF to ENVI GDAL translation causes rotation

texts:
I have a GeoTIFF(you can download it here ~16MB) which I want to convert to ENVI raw. I'm using this GDAL command:
gdal_translate -ot UInt16 -scale -of ENVI -outsize 1025 1025 dgm1_32372_5680_2_nw.tif dgm1_32372_5680_2_nw.raw
The command runs fine, but gives a warning:
Warning 1: Geotransform matrix has non rotational terms
When I compare the *.tif and the *.raw file in QGIS (3.14), the *.raw file got somehow rotated.
So I did some googling and found this question.
Ok, let's have look:
grep "map info" dgm1_32372_5680_2_nw.hdr gives:
map info = {Transverse Mercator, 1, 1, 32371999.5, 5679999.5, 1.95121951219512, 1.95121951219512, rotation=90}
But according to this blog post the bug have been fixed with GDAL 2.20. I'm using GDAL 3.0.4.
So does anybody know what the problem might be?


texts:
How can I generate turn-by-turn directions using PgRouting?

At the moment, I have a web-page where a user can enter and start/end location and a route is generated on a map (using PgRouting and OpenStreeMap). I know how to list all of the street names, but I don't know what the best approach for stating which direction the user should take at intersections.

I have read through previous questions on this site:
Does pgRouting support turn-by-turn routing? and How to emulate Google Maps driving directions using pgRouting? but both of these questions are over two years old, and several of the answers link to sites that no longer exist.


texts: Count Lat Lon Data Points Lying Between 2 Azimuths

texts:
I need to count the number of points (LAT, LON) from a cvs file that exist inside an arc defined by 2 azimuths.  Ultimately it will need to be iterative where the arc defined by the 2 azimuths rotates a full 360 degrees counting the # of points in each arc.
I have installed geopandas as a tool box.


texts: Problems with ST_ClusterDBSCAN cluster sizes

texts:
I want to cluster all villages (and towns, ...) around the world. So instead of having millions of them, I'd like to reduce them by combining villages that are close enough (like 10km or so) to each other.

So I was looking into ST_ClusterDBSCAN and it's doing quite a good job.
This is how I made my table:

INSERT INTO villages_clustered
  SELECT
    name,
    way,
    ST_ClusterDBSCAN(way, eps := 2000, minpoints := 1) over () AS cluster_id
  FROM villages;

Now when I look at my data, I get lots of small clusters which is what I wanted. But there are also some, which are just suuuuper large that should split into smaller ones.

Would really love to know how to improve my query to get better results.
What I basically want is clusters of like 10km or so.


texts:
I am new to all of this and I would like to know if it is possible to convert .xyz to .dxf in QGIS? I heard that it is possible in ArcGIS but I am not sure about it in QGIS. (I hope I am right asuming that an Ascii file is .xyz?)
My problem is that I have to convert a digital terrain model (DGM in German) to .dxf but  I have actually no idea how to do it and how it should look like.
When I just open it, I get the following image (and my computer gets pretty slow...):

I saw that I can save it as .tiff, but there are no more options to choose from.


texts: Halo mask in symbol property editor not working in ArcMap

texts:
I am trying to add a halo around symbols i have mapped in ArcMap 10.4 to make them more visible on the map. The symbols are plotted as points. I went to the shapefiles properties -> symbology tab -> clicked on symbol -> edit symbol -> mask tab -> changed "none" to "halo" and set halo to be 5 in size. When i apply this and press okay, it does nothing to the symbol. This is happening with my other shapefiles as well.
Am i doing something wrong?


texts: Is there a way to digitize an ROI in ENVI 5.2 that is larger than the area shown in the viewer?  I'm setting up labs for an intro to RS class, and the next step is to subset the image to get rid of the edge pixels.  As far as I can tell, the ROI tool in the new ENVI doesn't allow you to pan your view and keep digitizing, which means that you can't get an ROI that is larger than whatever subset of the image you are looking at...I hope I'm describing this problem well enough, let me know if it's unclear.


texts: Using Python to set text for a label in a QGIS3 print layout label

texts:
I am able to assign a QGIS 3.2 print layout label item to a Python variable using:

MyMapLabel = MyPrintLayoutObject.itemById('MyLabelName')

Which is all well and good but I need to be able to change the text of this label. In QGIS 2.18 this was possible using:

MyMapLabel = MyComposerObject.getComposerItemById('MyLabelName')
MyMapLabel.setText('A text string')

When trying the setText method in QGIS 3.2 the following error is returned:

  AttributeError: 'QgsLayoutItem' object has no attribute 'setText'

I notice that QGIS 3 includes a new class QgsLayoutItemLabel but this is a constructor and doesn't seem to include a method for identifying existing label items.

Has anyone had more luck than I in interpreting the QGIS 3 Python API?


texts: I'm loading four .ASC files into QGIS, each containing topographical data about a tile of 10 x 10 km. On the map, this data is displayed as black/white and different shades of gray. My problem is that the edges of each tile are very clear where the four tiles meet. Is there any way to get QGIS to polish this up, so as to achieve a smooth transition?


texts: How to convert a vector map to raster in Qgis?

texts: I am creating a thematic map (vector data) with different classes. How to convert this vector map to a raster map, keeping the classification?


texts: I am trying to import data from an existing ArcMap document to ArcGIS Online but can't seem to get it to work. I have tried zipping the GDB from my map document but this will not import to ArcGIS Online. Any suggestions on what I should do or what I am doing wrong?


texts:
I'm trying to claculate the NDWI values for some plots that I research in Africa using LANDSAT 7. it seems like the code work, but when I try to visualize the polygons, I get only parts of them (see image)

the dark area is the RGB median image from Landsat 7, the orange part is the NDWI calcualtion for the same geometry as the geometry of the RGB polygon.

i'm worried because my end goal is to get statistics for this area for many years, but I found out  that when I try to visualize specific images from the Image Collection, some of them can't be visualized, even though they have values.

this is the code I have written:

var geometry=plot2;
/**
 * Function to mask clouds based on the pixel_qa band of Landsat SR data.
 * @param {ee.Image} image Input Landsat SR image
 * @return {ee.Image} Cloudmasked Landsat image
 */
var cloudMaskL457 = function(image) {
  var qa = image.select('pixel_qa');
  // If the cloud bit (5) is set and the cloud confidence (7) is high
  // or the cloud shadow bit is set (3), then it's a bad pixel.
  var cloud = qa.bitwiseAnd(1 << 5)
                  .and(qa.bitwiseAnd(1 << 7))
                  .or(qa.bitwiseAnd(1 << 3));
  // Remove edge pixels that don't occur in all bands
  var mask2 = image.mask().reduce(ee.Reducer.min());
  return image.updateMask(cloud.not()).updateMask(mask2).divide(10000)
  .copyProperties(image, ['system:time_start']);
};

//Create LANDSAT7 dataset

var dataset = ee.ImageCollection('LANDSAT/LE07/C01/T1_SR')
                  .filterDate('1999-01-01', '2019-07-01')
                  .select('B1','B2','B3','B4','B5','pixel_qa')
                  .filterBounds(geometry)
                  .map(cloudMaskL457);

//RGB visualization

var visParams = {
  bands: ['B3', 'B2', 'B1'],
  min: 0,
  max: 1,
  gamma: 1.4,
};
// Map.centerObject(geometry);
// Map.addLayer(dataset.median(), visParams);

//clip the dataset according to the geometry
var clippedCol=dataset.map(function(im){
  return im.clip(geometry);
});

//test if clipping the image collection worked
Map.centerObject(geometry);
Map.addLayer(clippedCol.median(),visParams, 'RGB TEST');

// Get the number of images.
var count = dataset.size();
print('Count: ',count);
//print(clippedCol);
print(dataset,'dataset');

//function to calculate NDWI in LANDSAT7
var addNDWI = function(image) {
  var NDWI = image.normalizedDifference(['B4', 'B5'])
  .rename('NDWI')
  .copyProperties(image,['system:time_start']);
  return image.addBands(NDWI);

};

//NDWI to the clipped image collection
var withNDWI = clippedCol.map(addNDWI).select('NDWI');

// Test the addNDVI function on a single image.
var ndwi1 = withNDWI.select('NDWI').mean();

var NDWIcolor = {
  min: -1,
  max:1,
  palette: ['0000ff', '00ffff', 'ffff00', 'ff0000', 'ffffff'],
};

// Map.addLayer(ndwi1,NDWIcolor,'test');

// Test the addNDWI function on a single image.
var listOfImages =(withNDWI.toList(withNDWI.size()));
// var firstImage = ee.Image(listOfImages.get(0));
// var secondImage =ee.Image(listOfImages.get(1));
var randomImage=ee.Image(listOfImages.get(7));

Map.centerObject(geometry);
// print('NDWI Images',withNDWI.select('NDWI'));
// Map.addLayer(firstImage,NDWIcolor,'test1');
// Map.addLayer(secondImage,NDWIcolor,'test2');
Map.addLayer(randomImage,NDWIcolor,'Random_Image_test');

//Long-Term Time Series

print(ui.Chart.image.series(withNDWI, geometry, ee.Reducer.mean(), 30));

I also get graph which looks nice but now I'm not sure that the data there is correct.

the end goal is to create time series chart graph with the NDWI mean value  for each image.


texts: Setup Mapcache (ms4w package) with FastCGI on IIS

texts:
When I try starting MapCache using the IIS FastCGI module, I get an Internal Server Error

  .. mapcache.fcgi.exe - The FastCGI process exited unexpectedly

I've already got MapServer running like this without any problems.

I'm using the MS4W v3.2.1 CGI binaries (includes compiled executables for FastCGI) and setting MAPCACHE_CONFIG_FILE in environment variables for FastCGI to the path of the xml configuration.

I wasn't able to find a good explanation of how to get this running. (Only things like "I got it running" - not very helpful)

I'm running Win10 with IIS10 and VC++2012 redistributable installed.


texts: Cannot get the correct OS grid references on the outside of the map

texts: I am using both 50K and 25K Rasters from OS to try to produce smaller scale maps, eg Exmoor, but try as I might I cannot get the OS grid references on the outside of the map grid, the best I have got is on the centimetre or millimetre button to produce figures on the outside of the map, but both X and Y start at 0 not the OS grid reference, I am using CRS ESPG:27700, can anyone help?


texts:
I've got a WMTS service that work well, but I would like to display a grid in qgis over it with X, Y, Z information how may I do that ?

If QGis can't generate it for us, how may I generate this kind of grid ?

My problem is that I would need very high zoom (17/ 18) and the generation of geopackage is impossible (30 000 000 000 000 polygon is way to much with today PC...).

Is there a way to have some kind of "stream" of grid for this level f zoom ?

Like this: http://openlayers.org/en/latest/examples/canvas-tiles.html


texts: Wrong area using Albers Equal Area Projection in GRASS for Java, Indonesia

texts:
I am just getting to know GRASS 7.0 and I am struggling finding a good projection to work on a spatial analysis of Java, Indonesia. I am new to the forum, so if information's are missing for a better understanding of the problem please let me know. I built my own projection using the Albers Equal Area Projection (which is similar to http://spatialreference.org/ref/esri/asia-south-albers-equal-area-conic/ only the central meridian and the standard parallels are changed):

["Java Equal Area Conic",

GEOGCS["GCS_WGS_1984",  DATUM["D_WGS_1984",
SPHEROID["WGS_1984",6378137.0,298.257223563]],
PRIMEM["Greenwich",0.0],     UNIT["Degree",0.0174532925199433]],
PROJECTION["Albers"], PARAMETER["False_Easting",0.0],
PARAMETER["False_Northing",0.0], PARAMETER["Central_Meridian",110.0],
PARAMETER["Standard_Parallel_1",-6.38],
PARAMETER["Standard_Parallel_2",-8.3],
PARAMETER["Latitude_Of_Origin",-0.0], UNIT["Meter",1.0]]

I checked if the projection works by letting GRASS calculate the area of the Java Island. The "Query results" show Sq_Meters=10.371, but it should be around 128.297 km², though.

Is some conversion not working? Is the projection file wrong? How can I let GRASS calculate the actual area of the Java Island in meters? Or is this simply not possible with Albers?


texts: Identifying overlapping polygons in single layer using QGIS?

texts:
I have a shapefile of buffered points in QGIS. I need to display all buffer polygons in this layer which overlap. I have experimented with the intersect tool, but this only seems to work if I am looking at the intersection of 2 separate layers otherwise it just says that all polygons intersect.

Does anyone have a QGIS solution?


texts: Can I project a single set of coordinates that are not in a table or feature class with arcpy?

texts:
I have been searching for the easiest way to convert a single set of State Plane Coordinates (in units of "Feet") to Decimal Degrees using arcpy. All the methods I have seen so far involve one of the following ...
 1. Installing a Library like pyproj
 2. Creating (or starting with?) a feature class or table to project with "Project_management" or "ConvertCoordinateNotation_management".
 3. Using more math than I care to ... and some sort of Haversine or Vincenty algorithm?

I really only want to convert one set of coordinates that I have already arrived at in my python script (using arcpy.mapping). If I have to I suppose I can create a table or feature class out of them in order to then read in the table or feature class and project that but that seems ridiculous. Having recently done exercises using the "utm" module where this kind of thing was really easy, I don't understand why I should have to jump through hoops to convert a single set of coordinates. Shouldn't arcpy provide a way to do this? If one of the above is the only way for me to do it, then that is what I will do. Installing the library seems easiest but I'm not sure how to install a module in the version of Python that I use in ArcMap (I have done it with PIP in the python I have in the OSGeo4W shell).

Reference Links:
Search Cursor Solution w/ pyproj link also
Adding XY Data
UTM module


texts:
I successfully find coding for calling data from postgresql into QGIS layer panel.

The coding as below;
import qgis
root = qgis.core.QgsProject.instance().layerTreeRoot()
groupName = root.addGroup('Cerun')

uri = QgsDataSourceURI()
uri.setConnection("localhost", "5432", "amsdb", "postgres", "postgres")

uri.setDataSource('public', 'slope', 'geom', '', 'name')
layer = qgis.core.QgsVectorLayer(uri.uri(), 'slope', 'postgres')
print layer.isValid()
QgsMapLayerRegistry.instance().addMapLayer(layer,False)
groupName.addLayer(layer)

uri.setDataSource("public", "slopeinfo", None)
layer = QgsVectorLayer(uri.uri(), "slopeinfo", "postgres")
print layer.isValid()
QgsMapLayerRegistry.instance().addMapLayer(layer,False)
groupName.addLayer(layer)

uri.setDataSource("public", "services_facilities", None)
layer = QgsVectorLayer(uri.uri(), "services_facilities", "postgres")
print layer.isValid()
QgsMapLayerRegistry.instance().addMapLayer(layer,False)
groupName.addLayer(layer)

uri.setDataSource("public", "sign_of_distress", None)
layer = QgsVectorLayer(uri.uri(), "sign_of_distress", "postgres")
print layer.isValid()
QgsMapLayerRegistry.instance().addMapLayer(layer,False)
groupName.addLayer(layer)

uri.setDataSource("public", "surface_subsurface", None)
layer = QgsVectorLayer(uri.uri(), "surface_subsurface", "postgres")
print layer.isValid()
QgsMapLayerRegistry.instance().addMapLayer(layer,False)
groupName.addLayer(layer)

uri.setDataSource("public", "engineering_judgement", None)
layer = QgsVectorLayer(uri.uri(), "engineering_judgement", "postgres")
print layer.isValid()
QgsMapLayerRegistry.instance().addMapLayer(layer,False)
groupName.addLayer(layer)

uri.setDataSource("public", "structure", None)
layer = QgsVectorLayer(uri.uri(), "structure", "postgres")
print layer.isValid()
QgsMapLayerRegistry.instance().addMapLayer(layer,False)
groupName.addLayer(layer)

enter code here

The problem is I want to adding all the table layer under Group>SubGroup>Layer just like picture below. All the Group>SubGroup>Layer already in Layer Panel. Need to find coding to adding the table under Group>SubGroup>Layer.

Can anyone help?


texts: Distance from one point to many polygons

texts:
I'm using ArcMap 10.1 and have a shape file of polygons of municipalities. I also have the coordinates of the centroids of the polygons (if I need to use them).
I want a column in my attributes table which gives me, for each polygon, the distance to one specific point (lat 37N long 16E). I've been fiddling around with the Proximity Toolbox (under Analysis Tools), but don't know where to enter the coordinates of the specific point.


texts: PostGIS - incorrect result distance between two points

texts:
I am quite confused the way PostGIS is returning the result for this basic query.

So I am trying to find the distance between two points
45.2714,71.2087 and 42.3739,66.39

and the distance should be ~505 km

but when I ran the query (considering 1 degree = 111km)

select st_distance (st_setsrid(st_makepoint(71.2087,45.2714),4326),
      st_setsrid(st_makepoint(66.39,42.3739),4326))*111 as d

I get the result as 624 km.

Can anyone please tell what I am doing wrong here ?

Thanks.


texts: Differences between ArcSDE 9.3 and 10.0?

texts:
I am curious that what is differences between ArcSDE 9.3 and ArcSDE 10?

Is there ArcSDE 10 to use it?

Or ArcSDE 10 is coming with ArcGIS for Servers?

Can anyone explain me?

Edit : Another helpful answer to this question is here ; ArcSDE performance and table changes on upgrade from 9.3 to 10.1?


texts: How to symbolize line jumps/crosses in ArcMap?

texts:
Where two line features cross but are not physically connected (or related), I would like to symbolize the crossing using a line jump symbol. Visio is capable is using such symbology:

In 2008, Esri said "we might get this on the list of enhancements for future releases." I've seen no update to indicate that the enhancement was ever made.

Are line jumps possible in ArcMap today, either natively or via a workaround? If you've employed a workaround, what is it?

EDIT 3/7/2016: I would like to do this without modifying the underlying geometry, if possible.


texts: What's is the difference between slicing and stacking in a time-series analysis?

texts:
My adviser told me that there are two categories of time-series in remote sensing: slicing and stacking.

What is the difference between time-series stacking and slicing pertaining to landcover change in the field of remote sensing


texts: Arrange overlapping point labels so they are all aligned on one side (ArcMap 10.0)

texts:
I know by changing the conflict detection settings in placement properties it's usually possible to force overlapping labels to different corners or sides of a point, but would it be possible to have all of these labels arranged in a stacked list on one side of the point without converting them to an annotation?

It would look something like this:

If it was simply 3 labels each time I could set one to top right, one to center right, and one to bottom right but there could be different number of labels for each location where an overlap occurs.


texts: Converting polygon geometries to polygonz geometries using GeoTools?

texts:
I have a shapefile which contains a lot of polygon records. These polygon geometries contains points which has x and y coordinate pairs. I want to add z values (z=0) for every polygon points thus I can convert polygon geometry type to polygoz. But I don't know how do that using geotools. Here which I tried so far;

public List<Geometry> getGeometries(String filePath) throws IOException {
    List<Geometry> geometriListesi = new ArrayList<Geometry>();

    File file = new File(filePath);
    Map<String, Object> map = new HashMap<>();
    map.put("url", file.toURI().toURL());

    DataStore dataStore = DataStoreFinder.getDataStore(map);
    String typeName = dataStore.getTypeNames()[0];

    SimpleFeatureSource source = dataStore.getFeatureSource(typeName);
    SimpleFeatureCollection collection = source.getFeatures();
    System.out.println(collection.size());

    SimpleFeatureIterator itr = collection.features();

    while (itr.hasNext()) {

        SimpleFeature f = itr.next();
        Geometry g = (Geometry) f.getDefaultGeometry();
        geometrilistesi.add (g.boundary());
     }
     return geometriListesi;
   }

Now I am stuck and can't decide how to go on.


texts: I have installed geonetwork 3.10.4. When I print a raster WMS, a VERSION = 1.1.1 is put in the getmap, while the service is 1.3.0. This throws an error in the EXCEPTION field. Can you help me?


texts:
Has anyone already used IGeometryServer2 or IGeometryServer within C# code? I'm looking for a sample since I have not succeeded in using it.

I have already selected:
One Polyline (IPolyline) : for cutting
One Polyline (IGeometry) : the shape to be cut


texts: QGIS 3.10 drop-drown list (value relation)

texts:
I am using QGIS 3.10. I created a field with a drop-drown list of 10 country names (using an imported csv file)  . Is there an option to add a country that does not inculded in the list manually without updating the list itself? The field is a widget type "value relation".
@sn1ks I want to add a new value when I open the feature form.


texts: select special rows from feature using condition in ArcGIS ModelBuilder?

texts:
I am using ModelBuilder and would like to delete specific rows of a feature where the value of two columns are equal.

For instance: if there are columns a,b,c,d and rows 1,2,3,4. I would like to delete the rows where a=c. For that, it is first needed to select those features. When I enter the syntax in the query builder a = c I receive error (000229) and it does't work. I have tried another selection method (Make feature Layer, Selection by Attribute) I am using ArcGIS Desktop 10.2.


texts: Find the fastest route depending on slope

texts:
I have two ways drawn on a map and I'm supposed to calculate which one is the fastest from start to finish. The assumption is that you run with a speed of 4 min/km on flat ground, each percent increase of average slope increase the time with 15 sec/km.

As indata I have an elevation raster and an orthophoto in  jpg-format with corresponding jgw-file.

Does anyone know how to solve this?


texts:
I often switch between graphical modeler and the main QGIS application and wonder why I am not able to operate QGIS when modeler window is open. Saving my models and closing / reopening modeler annoys me, so...

Is the a deeper reason for modeler being modal?
Is there a simple way to turn off this behaviour (python console)?


texts:
I was working on an example provided in the  pyqgis developer cookbook  and performed the code in the QGIS version 2.8.9 but i am getting the following error:
ImportError: No module named MapComposer
I had tried the link provided in the book but that is not working.
My code is as follows:
from PyQt4.QtCore import *
from PyQt4.QtGui import *
from qgis.core import *
from qgis.gui import *
import MapComposer

lyr = QgsVectorLayer("/qgis_data/ms/mmississippi.shp",
"Mississippi", "ogr")
reg = QgsMapLayerRegistry.instance()
reg.addMapLayer(lyr)
mr = iface.mapCanvas().mapRenderer()
qc = MapComposer.MapComposer(qmlr=reg, qmr=mr)

setGridAnnoPos = qc.composerMap.setGridAnnotationPosition
setGridAnnoDir = qc.composerMap.setGridAnnotationDirection
qcm = QgsComposerMap

qc.composerMap.setGridEnabled(True)
qc.composerMap.setGridIntervalX(.75)
qc.composerMap.setGridIntervalY(.75)
qc.composerMap.setGridStyle(qcm.Solid)

qc.composerMap.setShowGridAnnotation(True)
qc.composerMap.setGridAnnotationPrecision(0)

setGridAnnoPos(qcm.OutsideMapFrame, qcm.Top)
setGridAnnoDir(qcm.Horizontal, qcm.Top)
setGridAnnoPos(qcm.OutsideMapFrame, qcm.Bottom)
setGridAnnoDir(qcm.Horizontal, qcm.Bottom)
setGridAnnoPos(qcm.OutsideMapFrame, qcm.Left)
setGridAnnoDir(qcm.Vertical, qcm.Left)
setGridAnnoPos(qcm.OutsideMapFrame, qcm.Right)
setGridAnnoDir(qcm.Vertical, qcm.Right)

qc.composerMap.setAnnotationFrameDistance(1)
qc.composerMap.setGridPenWidth(.2)
qc.composerMap.setGridPenColor(QColor(0, 0, 0))
qc.composerMap.setAnnotationFontColor(QColor(0, 0, 0))
qc.c.addComposerMap(qc.composerMap)

qc.output("/qgis_data/map.jpg", "jpg")