texts:
At the moment I have a WFS with the following cql_filter:

'cql_filter: 'dwithin(bag:geometrie, point(lon lat), 100, meters)'

The problem is that this filter does not sort the features from closest to farthest.

I would like to know if anybody is either able to sort this list appropriately, or if anybody knows how to get the closest feature in another way.


texts:
I use QGIS and I want to work with custom forms using python logic and I try to follow this manual.  In my custom qu I add one new pushButton and I want for example for educations rules to open some image.  But that dont work dont show me anything,any idea where I have wrong in my python code ?

Here the python code :

from PyQt4.QtCore import *
from PyQt4.QtGui import *
import os
from PyQt4.QtGui import QWidget

nameField = None
myDialog = None
def my_form_open(dialog, layer, feature):
    geom = feature.geometry()
    control = dialog.findChild(QWidget, "MyLineEdit")
    global myDialog
    myDialog = dialog
    global nameField
    nameField = dialog.findChild(QLineEdit,"Name")
    buttonBox = dialog.findChild(QDialogButtonBox,"buttonBox")
    pushButton = dialog.findChild(QPushButton,"pushButton")

    # Disconnect the signal that QGIS has wired up for the dialog to the button box.
    pushButton.accepted.disconnect(myDialog.accept)

    # Wire up our own signals.
    pushButton.accepted.connect(validate)

def validate():
    filename='C:/Users/username/Desktop/image.jpg'
    os.startfile(filename)


texts: Why does means not appear in the Feature Collections after applying a ee.Reducer.mean() through the image.reduceRegions

texts:
I am trying to aggregate the MODIS ET data into 1km by 1km grid squares. I have tried applying the ee.Reducer.means through the image.reduceRegions to get the mean value of MODIS-ET pixels within the grid square. However, the resultant Feature somehow did not include a resultant mean value for each feature (grid square).
var AWD = table // 1km by 1km grid feature collections across Ayerwady region in Myanmar;

var ETcollection = ee.ImageCollection('MODIS/006/MOD16A2');

//print(ETcollection);

var ET_first = ee.Image('MODIS/006/MOD16A2/2001_01_01');

print(ET_first);
print('Scale:', ET_first.projection().nominalScale());

var ET = ET_first.select('ET');

print(ET);

var evapotranspirationVis = {
  min: 0.0,
  max: 300.0,
  palette: [
    'ffffff', 'fcd163', '99b718', '66a000', '3e8601', '207401', '056201',
    '004c00', '011301'
  ],
};

Map.addLayer(AWD, {}, 'Ayerwady');
Map.centerObject(AWD,8);

Map.addLayer(ET.clipToBoundsAndScale({geometry:AWD,scale:463.3127165275}), evapotranspirationVis, 'ET');

var awdgrid = ee.FeatureCollection(table2);

Map.addLayer(awdgrid,{}, 'Grid');

var stats = ET.reduceRegions({
  collection: awdgrid,
  reducer: ee.Reducer.mean(),
  scale: 463.3127165275
});

print('1km mean', ee.Feature(stats.first()));

print('mean', ee.Feature(stats.limit({max:10, property: 'mean' }))); // there are no mean values in the mean column

function getCols(tableMetadata){
  print(tableMetadata.columns);
}

//stats.limit(0).evaluate(getCols);

Map.addLayer(stats.select('mean'),evapotranspirationVis,'1km');

var empty = ee.Image().byte();
var outline = empty.paint({
  featureCollection: stats,
  color:'mean'
});
//var palette = ['#000000', '#F6FD0C', '#FF9A00', '#FF0000'];

Map.addLayer(outline, evapotranspirationVis,'ET_1km'); // I could still visualise the 1km by 1km means
//Map.addLayer(grid);

I can't seem to find where I went wrong. Is there a way to retrieve the mean value of each grid square?
UPDATE
'1km mean'
Feature 000000000000000033e8 (Polygon, 5 properties)
type: Feature
id: 000000000000000033e8
geometry: Polygon, 10 vertices
properties: Object (5 properties)
id: 492075
bottom: 16.924260574
left: 94.205277649
right: 94.214277649
top: 16.933260574

'mean'
FeatureCollection (10 elements, 7 columns)
type: FeatureCollection
columns: Object (7 properties)
features: List (10 elements)
0: Feature 000000000000000034e2 (Polygon, 6 properties)
type: Feature
id: 000000000000000034e2
geometry: Polygon, 5 vertices
properties: Object (6 properties)
id: 503135
bottom: 15.979260574
left: 94.250277649
mean: 828.770462633452
right: 94.259277649
top: 15.988260574


texts: How to convert SHP/POSTGIS data into x, y, z, "cost" data?

texts: Is there a way to convert a line layer (shp or PostGIS) data into (x, y , z, cost) i.e. (From_node, To_node, optional_Elevation, Cost/Weight) ?


texts: Creating shapefile of S2 cells for given level

texts:
I want to create a shapefile that has polygons matching with the outlines of the S2 cells (http://s2geometry.io/) of a given level, within a certain bounding box.

Roughly what this website does: Sidewalk Labs Region Coverer, but giving me the data as an actual shapefile.

I'll accept any working option:

generated through a website
though QGIS (and/or GRASS) if properly explained how to (as in if you're giving me a Python script, please tell me HOW to use it within QGIS, rather then just listing me all the options for the parameters)

Whatever works and is understandable and doable without having to purchase the software.

EDIT: if it requires hard-coding the S2 cell level, I'm particularly interested in levels 14 and 17, if the code requires manual altering the level please make sure it's clearly visible what value needs editing.


texts: Vector maps comparison (looking for a Kappa-like index for vectors) in R?

texts:
I am currently trying to compare two shapefiles (representing buildings) in order to validate simulation results. One shapefile contains real data and the other is the product of an urban simulation model.

The thing is that in most of the papers I found, researchers are using raster data (i.e. http://onlinelibrary.wiley.com/doi/10.1111/tgis.12047/full ) and apply grid-oriented indexes like the Kappa statistics. The only paper I found that compares two vector maps is the following one : http://link.springer.com/chapter/10.1007/978-3-319-46840-2_10
In this paper, "the error is computed as the rate of new real buildings that are not intersected (with a given buﬀer) by buildings obtained by simulation, weighted by their area" (Taillandier et al., 2016). This is a good start but I'd rather use more than one statistic.

Thus, I was wondering if anyone had any idea of how to compare two vector maps, both containing one feature class (buildings), in order to express/measure the agreement between the two maps.

I mainly work with R.


texts: Removing grid from showing up in shapefile

texts: I am pretty new to QGIS. I have the following result generated (map shown). The issue I am having is I am not sure how to get rid of gridlines inside the map in a manner that it just shows the outside boundary of the country?


texts:
With python-addin-wizard I created a combobox which enables the user to select a specific ID of features in multiple featureclasses of the Table Of Content and then it zooms to the selected feature. So far so easy. It's not really fast, but it works okay.

My problem is, if a featureclass is selected in the TOC e.g. ,

items in the TOC appear to flicker and the processing is delayed considereably. I searched a lot in the arcpy functionalty, but can't find a bit of code to solve this Problem (quasi unselect the feature in the TOC). The only thing related to my Problem, that I found is pythonaddins.GetSelectedTOCLayerOrDataFrame() > But I have no idea how to unselect the layer.

Here is the part of the code which makes the selection and zoom:

def onSelChange(self, selection):
    self.value = selection
    fcsites2 = [fc_site_pol, fc_site_punkt, fc_site_linie, fc_neg_pol, fc_neg_punkt, fc_neg_linie, fc_unter_pol, fc_unter_punkt, fc_unter_linie]
    for x in fcsites2:
        arcpy.SelectLayerByAttribute_management(x, "NEW_SELECTION", fields_site + "='" + selection + "'")
        d = arcpy.Describe(x)
        n = len(d.FIDset)
        if n > 0:
            break
    self.df.zoomToSelectedFeatures()
    self.df.scale = 5000

@artwork21 helped me start to solve my problem using the comtypes module. I have no much experience with integrating ArcObjects in Python. I reached to reference the TOC and the selected layer once. Then I can unselect the layer with the code below.But when I try a second unselect it doesn't work. It seems like the first selected layer is still known as selected like in a cache.

mU = GetModule(r"C:\Program Files (x86)\ArcGIS\Desktop10.2\com\esriArcMapUI.olb")
mF = GetModule(r"C:\Program Files (x86)\ArcGIS\Desktop10.2\com\esriFramework.olb")

app = CreateObject(mF.AppROT, interface=mF.IAppROT)
pDoc = CType(app.Item(0).Document,mU.IMxDocument)
pMap = pDoc.FocusMap
ptoc = CType(pDoc.ContentsView(0),mU.IContentsView)
pselLayer = pDoc.SelectedLayer

if pselLayer != 0:
    print pselLayer.Name
    ptoc.RemoveFromSelectedItems(pselLayer)
    ptoc.Refresh


texts: Timeout to load ArcGIS Server hosted layers

texts:
I'm sharing a map and a geoprocessing service in ArcGIS Server. I configured the CMV(Configurable Map viewer) with the rest URLs to consume them. I deployed the CMV in a IIS server, then call the CMV in a iframe in my asp.net MVC 5 application.
When I'm running the CMV in local, the layers and the geoprocessing services are well loaded so usable. However when I'm running from my application it never loads and it gives me a timeout error.
I googled a bit and I found that I have to add this line in the web.config file, so the timeout would be extended :
<httpRuntime executionTimeout="21600" />

Yet nothing changes and the layers are not added to the map.
Do you know what may cause this  issue right here?


texts:
I have QGIS 2.18.16, in Arch Linux; I wanted to use the qgis2web plugin and I realized it is not working anymore (it used to work last time I used it a few months back)
I get the following error

  cannot import name QtWebKit

I have QtWebkit installed - not sure why it cannot find it;
as a side note, with the developer version qgis2web works fine (but not all needed plugins are available, so I cannot use it)


texts: Is it possible to use esri.dijit.Popup for editing

texts:
In the new arcgis js api ESRI added new nice infowindow - esri.dijit.Popup. We can use it this way:

 //create a popup to replace the map's info window
 var popup = new esri.dijit.Popup(popupOptions, dojo.create("div"));

 map = new esri.Map("map", {
      extent: initExtent,
      infoWindow: popup
 });

The code from samples.
Is it possible to use it for editing? I like this new infoWindow more than standart editing infoWindow


texts: Moving Pie Chart in QGIS and creating leader line

texts:
I'm using QGIS 2.18.3 (Las Palmas).
I have created pie charts for a layer in QGIS and in order to be able to move the pie charts around i processed the centroids of my polygons. By moving the (Centroid) points I am now able to move the pie charts around.

I still have some overlapping in my charts, so i'm looking into creating leading lines, so I can allocate the pie chart to its polygon. I have looked into previous questions here (for example How to create dynamic leader lines? or Creating leader lines for moved labels) but I dont seem to be able to use it for the charts.

How can I create a leading line for a pie chart?


texts: Method for Spatial Analysis of Raptor Predation

texts:
I currently have presence data on two birds, a falcon species and the grouse that it prefers to hunt. These data are in two point shapefiles, one for each species. The grouse shapefile represents all of the individuals that I've seen during observational transects, while the falcon shapefile represents all falcon nests in my study region.

I have a hypothesis that the grouse species are located in areas where they are near the fewest possible number of falcon nests (there is almost nowhere on the study site completely without falcon nests, but there are locations where nests are more concentrated). I'd like to examine this with some sort of statistical GIS analysis, but I'm not sure what steps I should take. I've been looking into kernel analysis and the near distance table.

How do I put things together into a cohesive method?

My grouse data is biased toward 3 roadways that run for miles through my study site. If possible, I'd like to also entertain method that may take this bias into account -- but only if possible!


texts: Change HTML markup used it PopupTemplate in ArcGIS (javascript api)

texts:
everyone. I'm sorry if it is common question but I'd googled it for a couple days and still haven't find the answer. Everything I see in documentation is style an exist HTML layout, which looks like:

            Zoom to

But I need to setup totally different popup, put style on it and, surely, link it to data. I see there's the PopupTemplate class that allow to put data in 'content' block or 'title' block, but how to define own sections, change markup and do stuff like this?

What I tried:

delete sections in popup.domNode and then add my ones - data stoped
synchronize.
changed markup directly in popup.domNode.innerHTML - same thing
thought that it is crucial to update ._title and ._content properties to show where's now these blocks placed, but have same result.


texts: Adding multiple addresses - ArcGIS Javascript

texts:
I am trying to use ArcGIS Javascript API to add all of the service orders that are open to a map. I can get this to work with a single address, but it is not working when I try to loop through a list containing addresses. I am not too familiar with using javascript, but I do not see anything wrong with the actual code.

I have already looked at "Geocoding Multiple Addresses using ArcGIS JS API," but it did not work for me.

ESRIMapView (cshtml)- This contains the code to map all of the service orders:

 window.onload = function locateAll() {
               map.graphics.clear();

               var addresses = Viewbag.Addresses;
               var addressCount = Viewbag.Addresses.length;

               for (i = 0; i < addressCount; i++)
               {
                   locator.outSpatialReference = map.spatialReference;
                   var options = {
                       address: addresses[i],
                       outFields: ["Loc_name"]
                   }
                   locator.addressToLocations(options);
               }
           }

HomeController.cs:

 public ActionResult ESRIMapView(string address)
        {
            ViewBag.Message = "Your contact page.";
            MFSRepository _repository = new MFSRepository();

            List<ServiceOrder> orders = new List<ServiceOrder>();
            orders = _repository.GetServiceOrderAddresses();

            List<string> addresses = new List<string>();

            foreach (ServiceOrder order in orders)
            {
                addresses.Add(order.StreetNumber + " " + order.StreetName + " " + order.City + " " + order.State + " " + order.Zip);
            }

            ViewBag.Addresses = addresses;
            ViewBag.Address = address;

            return View();
        }

As I step through the code in the controller, everything is being set correctly. I also know that it is getting the list of addresses from the controller because I was able to put the first one in the list into a textbox through that ESRIMapView. Nothing is being plotted on the map, and I am not getting any errors. My guess is that my syntax for the loop is probably wrong.


texts: Converting 3d profiles to points in ArcGIS Desktop?

texts:
I have created some 3d profiles using the Interpolate Line button (in ArcGIS 10.1). Which looks like:

I would like to convert these lines into points and in order to do that I am under the impression that I need to have a shapefile of these lines. I am not really sure how to do that and was wondering whether anyone could give me some pointers.

Also, do I have to create a new feature class and re-draw these lines?

If so is there a way to avoid that?


texts: Change the resolution and size of a raster with GDAL in Python

texts:
I am quite new to GDAL and Python but want to use it to create groundwater models in FloPy.

I have a raster with rivers (roughly drawn) and a DEM raster. The resolution of the river-raster is higher than that of the DEM. Thereby the river-raster only covers a part of the total area of the DEM-map:

What I want is the following:
1) Apply the resolution of the DEM to the river-raster
2) Extend the river-raster so it has the same extent as the DEM by filling it with no-data values within the boundaries of the DEM.

In QGIS it works easily by saving the river raster with a different extent and resolution but since I want to do this for whole batches it would be nice if I learn how to do this with GDAL in Python.


texts: Using Value from Iterate Feature Selection to Calculate Field in ArcMap 10.1 ModelBuilder

texts:
I have two polygon feature classes, an alignment layer and a parcels layer. I am trying to use ModelBuilder to create a tool that will do the follow:

iterate through groups in a single field of a polygon alignment layer (ie select and iterate through certain attributes in a field with same names such as "tunnel", "portal", etc)
Clip a parcel layer by each group
Populate a new field in the parcel layer called "EsmtType" with the name of each group (group is called "Value" in model)
Merge the clipped groups into one feature class

I have gotten steps 1 and 2 to function but I am stuck on step 3. I am not sure how I need to configure the "Calculate Field" tool in order to recognize the "Value" variant of the "Iterate Feature Selection" tool. I tried using the advice from this post Using a value from "iterate feature selection' to calculate a field but I couldn't figure out how to implement it for my situation.

Here is the model I have thus far:
http://imgur.com/ji2nAJJ

I would also like advice on how to implement step 4, I am not sure if it would require integrating a model within a model.

Thank you so much!!!
-kkaszas


texts: Merge nearby items into a single representation

texts:
I'm investigating OpenLayers and I'm trying to see if it's possible to filter down the number of displayed feature items.

Is it possible to display many nearby items as a single point item on a map when zoomed out.  So instead of hundreds of nearly overlapping points I can filter them down a few dozen consolidated points?  Or is this something that has to be in the Spatial server?


texts: Gdal_calc outputs image with no projection

texts:
I am trying to convert a pixel value above 0.02 to be nodata. However, the output file comes out with no projection. When I assign projection, QGIS loads it in the wrong place.

gdal_calc.py -A "000-IMG_0037_4_gcp.tif" --outfile="000-IMG_0037_4_gcp_noglint8.tif" --calc="((A<=0.02)*(A>0))*A+((A>0.02)*0)" --NoDataValue=0

The gdalinfo of the original file:

Files: /drive/000-IMG_0037_4_gcp.tif
Size is 1048, 1422
Coordinate System is `'
GCP Projection =
PROJCS["WGS 84 / UTM zone 16N",
    GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9122"]],
    AUTHORITY["EPSG","4326"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",-87],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
    AUTHORITY["EPSG","9001"]],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH],
    AUTHORITY["EPSG","32616"]]
GCP[  0]: Id=1, Info=
      (0,0) -> (280191.450981,3352631.04297,0)
GCP[  1]: Id=2, Info=
      (1048.9984,0) -> (280298.687984,3352631.04297,0)
GCP[  2]: Id=3, Info=
      (1048.9984,1422.4381) -> (280298.687984,3352484.8363,0)
GCP[  3]: Id=4, Info=
      (0,0) -> (280191.450981,3352631.04297,0)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 1422.0)
Upper Right ( 1048.0,    0.0)
Lower Right ( 1048.0, 1422.0)
Center      (  524.0,  711.0)
Band 1 Block=1048x1 Type=Float32, ColorInterp=Gray

The output file gdalinfo:

Driver: GTiff/GeoTIFF
Files: /drive/000-IMG_0037_4_gcp_noglint8.tif
Size is 1048, 1422
Coordinate System is `'
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 1422.0)
Upper Right ( 1048.0,    0.0)
Lower Right ( 1048.0, 1422.0)
Center      (  524.0,  711.0)
Band 1 Block=1048x1 Type=Float32, ColorInterp=Gray
  NoData Value=0

Is there a way to make sure that the output image comes with the same srs of the input raster?


texts: How to use SpatialOperator in filter expression in QGIS

texts:
In a python script I want to specify a request filter expression using a spatial operator. I am following the general syntax "field name" operator condition. For example:

setFilterExpression('$geometry soOverlaps QgsRectangle(...)')

This gives no results and no error.

Please can somebody provide an example on how to build a valid expression using these spatial operators.


texts: Earth Engine: How to input band values of image to a tensorflow model trained on array of values

texts:
I have a TensorFlow model that is trained on an array of reflectance values (red, green, blue, nir) to predict sediment (SSC).

Now, I want to use Landsat 8 image over a region to get the SSC for each pixel in the image.One way that I'm aware of but not sure of implementation is to export each band of image to Google Drive as a GeoTIFF and then read in the file to Colab using the rasterio package.

Here is my code:

L8imgs = (ee.ImageCollection('LANDSAT/LC08/C01/T1_SR').filterDate(startDate, endDate).filterBounds(ROI))
L8imgs = L8imgs.select(['B2', 'B3','B4','B5'], ['blue', 'green', 'red', 'nir'])
L8imgs = L8imgs.map(runModel)

def ann_SSCL8(img):
  I = ee.Image(img)
  ## HOW TO INPUT BAND VALUES AS ARRAY TO tfModel ??
  input=[I.select('blue'),I.select('green'), I.select('nir'), I.select('red')]
  ssc = tfModel.predict(input)    #tfModel stored on drive


texts: Performing Conditional Calculate Field with Many to One using ArcMap?

texts:
Anyone have experience calculating a field based on a many to one scenario?

I have a feature class in a geodatabase with a field of unique names.

Feature Class Table  (let's call this "MAINTABLE") looks like this:

OID   ;  Name ; B
1    ;   Breakfast ; NULL
2    ;   Lunch ; NULL
3    ;   Dinner ; NULL
4    ;   Dessert ; NULL

I also have a stand alone table in a different location that has several entries for each unique name in the feature class.

Stand alone table (Let's call this "SIDETABLE") looks like this:

UID ; Name ; A
1 ; Breakfast ; 0
2 ; Breakfast ; -1
3 ; Lunch ; -1
4 ; Dinner ; 0
5 ; Dessert ; NULL
6 ; Breakfast ; -1
7 ; Breakfast ; NULL
8 ; Breakfast ; -1
9 ; Breakfast ; -1
10 ; Breakfast ; NULL

I need help to search all values ("0" or "-1" or NULL) in SIDETABLE field "A",  and if any entry in SIDETABLE field "A" has a "0", calculate field "B" in MAINTABLE as value "0".

BUT if a unique name entry in SIDETABLE field "A" contains "0" and "-1", calculate field "B" in MAINTABLE as value "0".

If a unique name from MAINTABLE does not in field "A" of SIDETABLE have value "0" but has value "-1", calculate field "B" in MAINTABLE as "-1".

If a unique name from MAINTABLE does not in field "A" of SIDETABLE have value "0" but has NULL value, calculate/keep field "B" in MAINTABLE as NULL.

The resulting feature class table ("MAINTABLE") should look like this:

OID   ;  Name ; B
1    ;   Breakfast ; 0
2    ;   Lunch ; -1
3    ;   Dinner ; 0
4    ;    Dessert ; NULL

***I would assume python is the best option. I am new scripting/ using python or whichever is best practice. I am assuming I would need some sort of conditional statement using Calculate Field on field "B" in the feature class. And I am also assuming I need to relate the feature class to table based on NAME.

***I am using ArcGIS 10.2 for Desktop Advanced

The following is the sort of logic I'm looking for

if A.SIDETABLE = 0
return into B.MAINTABLE the value "0"

*[if there are 10 "breakfast" entries in SIDETABLE,
  look at all values in field "A" and if one of the values is 0,
  then return "0" into B.MAINTABLE]

else if A.SIDETABLE = -1
return into B.MAINTABLE the value "-1"
else if A.SIDETABLE is NULL
return into B.MAINTABLE NULL


texts:
I am developing a plugin which reads MODIS images, and I know that MODIS images are in sinusoidal projection and I need to set that CRS and project the image to WGS84 (epsg 4326) in python. I search and find QgsRasterProjector class but I don't understand how to use this class. Is there another way?

I also tried the processing commands but there isn't Assign Projection command in processing, just the warp.


texts: Dorling Cartogram Map from Shapefile using QGIS

texts:
I would like to create a simple Dorling Cartogram map in QGIS.  My goal is for it to be a similar style to these maps created by the New York Times:

http://www.nytimes.com/interactive/2013/11/20/us/which-states-will-allow-old-health-policies-to-be-renewed.html?smid=pl-share&_r=0

I have a shapefile of jurisdictions, and I would like to size them based on their population.  The shapefile already has equal-sized square polygons for each jurisdiction, with each polygon roughly placed in its correct geographical location relative to the other polygons.

The Cartogram plug-in in QGIS does a good job of creating a Contiguous Cartogram, but I would like to keep the square-shape of each polygon; just size them larger or smaller based on their population.

After they have been sized, I would also like to be able to do some basic styling on them (i.e. have different colored polygons based on the type of jurisdiction, like City, Village, etc.).


texts: Spatial (nearest neighbour) Join Between Points (GPS coordinates) and Line in R

texts:
I want to do spatial join to an OSM line data and GPS point data both in geographic coordinate system, the result which I want is a road layer(has a unique ID) with the attributes of of road and the point file. I don't want to use buffer: sp::over (I get NA values). I tried using  gDistance() it only gives me the distance but does not join the attribute tables. I want to replicate NNJoin of QGIS in R since I can't use batch processing for NNJoin.

I don't think I am satisfied by the answer here https://stackoverflow.com/questions/47675571/r-spatial-join-between-spatialpoints-gps-coordinates-and-spatiallinesdatafra.

longitude <- c(10.86361, 10.96062, 10.93032, 10.93103, 10.93212)
latitude <- c(44.53355, 44.63234, 44.63470, 44.63634, 44.64559)
longlat <- data.frame(longitude, latitude)
my_spdf<- SpatialPoints(longlat, proj4string=crs("+init=epsg:4326"))
shapefile_roads<- spLines(cbind(longitude+.1, latitude), cbind(longitude-.1, rev(latitude)), cbind(longitude-.1, latitude+1), crs=crs("+init=epsg:4326"))
shapefile_roads <- spTransform(shapefile_roads, CRS("+proj=utm +zone=43 ellps=WGS84"))
my_spdf <- spTransform(my_spdf, CRS("+proj=utm +zone=43 ellps=WGS84"))
result <- over(shapefile_roads, my_spdf)


texts: Vector to raster by maximum combined area

texts:
I need to convert a polygon to raster resulting in a raster with cell size of 250 meters and each cell value takes the value of the objects, with equal values, that together sum the largest surface inside a cell. In ArcGIS, there is an option to use the maximum combined area algorithm within the feature to raster tool but now, I do not have ArcGIS installed.

Can it be done using QGIS?

I just want to get something like that:


texts: Add extra metadata to each message of a GRIB file?

texts:
I wonder, is there a tool to allow one add an extra metadata to each message of a GRIB file?
For example: from grib_dump of a GRIB file, I have:
#==============   MESSAGE 6 ( length=976330 )              ==============
GRIB {
  editionNumber = 1;
  table2Version = 2;
  # Rome  (RSMC)  (common/c-1.table)
  centre = 80;
  generatingProcessIdentifier = 131;
  # None Wind speed m s**-1 (grib1/2.0.2.table)
  indicatorOfParameter = 32;
  # Fixed height above ground height in meters  (2 octets)  (grib1/3.table)
  indicatorOfTypeOfLevel = 105;
  level = 10;
  # Average  (reference time + P1 to reference time + P2)  (grib1/5.table)
  timeRangeIndicator = 3;
  # Unknown code table entry ()
  subCentre = 98;
  paramId = 500025;
  #-READ ONLY- cfNameECMF = unknown;
  #-READ ONLY- cfName = unknown;
  #-READ ONLY- cfVarNameECMF = ws;
  #-READ ONLY- cfVarName = ws;
  #-READ ONLY- units = m s**-1;
  #-READ ONLY- nameECMF = Wind speed;
  #-READ ONLY- name = Wind speed (SP_10M);
  decimalScaleFactor = 0;
  dataDate = 20201026;
  dataTime = 0;
  # Hour (stepUnits.table)
  stepUnits = 1;
  stepRange = 0;
  startStep = 0;
  endStep = 0;
  #-READ ONLY- marsParam = 32.2;
  shortName = sp_10m;
  GDSPresent = 1;
  bitmapPresent = 1;
  numberOfVerticalCoordinateValues = 0;
  Ni = 928;
  Nj = 448;
  latitudeOfFirstGridPointInDegrees = 60;
  longitudeOfFirstGridPointInDegrees = -18;
  earthIsOblate = 0;
  uvRelativeToGrid = 0;
  latitudeOfLastGridPointInDegrees = 32;
  longitudeOfLastGridPointInDegrees = 40;
  iScansNegatively = 0;
  jScansPositively = 0;
  jPointsAreConsecutive = 0;
  #-READ ONLY- alternativeRowScanning = 0;
  jDirectionIncrementInDegrees = 0.063;
  iDirectionIncrementInDegrees = 0.063;
  #-READ ONLY- numberOfDataPoints = 415744;
  #-READ ONLY- numberOfValues = 308091;
  missingValue = 9999;
  tableReference = 0;
  #-READ ONLY- binaryScaleFactor = -19;
  #-READ ONLY- referenceValue = 0.529586;
  sphericalHarmonics = 0;
  complexPacking = 0;
  integerPointValues = 0;
  additionalFlagPresent = 0;
  packingType = grid_simple;
  bitsPerValue = 24;
  values(415744) =  {
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999,
  9999, 9999, 9999, 9999, 9999
  ... 415644 more values
  }
  #-READ ONLY- numberOfCodedValues = 308091;
  #-READ ONLY- maximum = 17.9256;
  #-READ ONLY- minimum = 0.529586;
  #-READ ONLY- average = 5.81554;
  #-READ ONLY- numberOfMissing = 107653;
  #-READ ONLY- standardDeviation = 3.87567;
  #-READ ONLY- skewness = 0.861793;
  #-READ ONLY- kurtosis = -0.237433;
  #-READ ONLY- isConstant = 0;
  gridType = regular_ll;
  #-READ ONLY- getNumberOfValues = 415744;
}

There are 100 messages from this GRIB file. I want to add an extra metadata processingLevel=x with values for x from 1:100 (each message has an unique number).


texts: v.kernel is missing under Processing Toolbox GRASS GIS 7 commands

texts:
I installed the QGIS Standalone Installer Version 2.14 (64bit), that automatically contains GRASS 7.0.4. I created a map with several layers and rasters in the program called QGIS Desktop 2.14.3 with GRASS 7.0.4.

I need one tool from Processing -> Toolbox -> GRASS GIS 7 commands -> Vector -> v.kernel => but is missing there. All the other vectors (v.XY) work well. Because I found v.kernel under GRASS commands, I activated it in the Processing options.

Under GRASS commands -> Vector -> v.kernel appears but is not working, none of the vectors (v.XY) is :

When I click on Enable additional providers at the bottom of the Processing Toolbox I see the Processing options :

The path for the GRASS7 folder is right. I put the same for the GRASS folder. But I cannot found the Msys folder in whole PC.

When i check Plugins -> Manage and Install Plugins -> GRASS7 is installed.

When I go to XY\QGIS Essen\apps\grass\grass-7.0.4\bin, I can see v.kernel.

Did anyone face the same problem?

Basically I have two questions:

1) Why v.kernel is not under GRASS GIS 7 commands? And how can I have it?

2) Why v.kernel and all the other vectors (v.XY) under GRASS commands are not working?


texts:
Could someone help me install QGIS 1.7.3 for a Mac and possibly give step by step instructions?

I have already installed QGIS 1.8 for a Mac, however I have heard that the 1.7.3 will run smoother without any problems.


texts: Exporting to PDF: scale changes print area

texts:
I have read similar questions, but I haven´t been able to do what I want in a quick way.

I want to export the area that I am seeing. For that, for example the scale is 1:8000, but I want the resolution of scale 1:2500.

In both options to export (Export map as PDF or with Composer) when I select the scale of 1:2500, the layout is reduced.

The only solution I have found is this: Export as PDF, copy the four coordinates in a notepad, change scale and, finally, copy the coordinates from notepad to QGIS. But I think there must be a more simple way.


texts: Problem adding layers with REST Services in GeoServer

texts:
This is my first question in this forum and I'm not sure if it is going to be properly structured. I've asked same question in other forum but I did not receive any response. I'll see if I have good luck here!

My problem happens when I want to upload a GeoTiff using the Import Data plugin through REST Services. I can't override same layer, method create a new one. But I need to keep same name for each update.

I'm using cURL commands to call the REST services.

1. Creating a new import

curl -k -u admin:pwd -XPOST -H 'Content-type: application/json' -d @import.json 'http://localhost:8080/geoserver/rest/imports'

where @import.json is:

{
   "import": {
      "targetWorkspace": {
         "workspace": {
            "name": "SIHA"
         }
      }
   }
}

2. Creating a new task (in this case import_id = 0)

curl -k -u admin:pwd -F name=test -F filedata=@file201510.tif 'http://localhost:8080/geoserver/rest/imports/0/tasks'

3. Running the task

curl -k -u admin:pwd -XPOST 'http://localhost:8080/geoserver/rest/imports/0'

At first time, it works perfectly. But when I've tried to call these commands again the name change ending with a sequence. (file2015100, file2015101, file2015102)

I've tried to delete the layer previously, but without success:

curl -k -u admin:pwd -XDELETE http://localhost:8080/geoserver/rest/layers/file201510

Layer is deleted properly but the new layer is created following the same strange pattern (adding a sequence at the end of the name)

Also, I've tried to force the updateMode, but without success as well.

curl -k -u admin:pwd -XPUT -H 'Content-type: application/json' -d @replace.json 'http://localhost:8080/geoserver/rest/imports/0/tasks/0'

where @replace.json is:

{
  "task": {
     "updateMode": "REPLACE"
  }
}

I hope that you understand my problem and someone could help to solve it.


texts: ESRI CreateGeometryFromWkbVariant creates CCW polygon

texts:
I am using ESRI GeometryEnvironmentClass that implements IGeometryFactory3 and I use the method CreateGeometryFromWkbVariant.
For the same wkb, in ESRI 9.3.1 the returned polygon is oriented clockwise, but in ESRI 10.0.0 it seems that it is counter-clockwise.

(I have found this out using the method IsSimpleEx in ITopologicalOperator3 interface: in 9.3.1 the polygon is simple, but in 10 it is not, with reason esriNonSimpleRingOrientation - "A polygon is topologically simple, but its rings may not be oriented correctly (outer rings - cw, inner rings - ccw)".)

Any idea what changed? Why was it working in 9.3.1 and no longer works in 10?

Thank you.


texts: QGIS Autofields expression function editor

texts:
I would like to have an Autofield that returns the length of a line segment plus 20. If that value is less than 150, I would like the Autofield to return 150. I am having trouble with syntax in the Function Editor for the Expression Builder for AutoFields. Below is my code:

from qgis.gui import *

@qgsfunction(args='auto', group='Custom')
def BoreLength(geometry, feature, parent):
    layer=qgis.utils.iface.activeLayer()
    length=geometry$length+20
    if length<150:
        length=150
        return length
    return length

And my Expression is currently this:


texts: Polygon fill colour specified by column R spplot

texts:
I am trying to plot a geological map using the spplot function in R. For every unit I already assigned a specific color in the col-column of the SpatialPolygonsDataFrame.
Now I would like to plot the polygons from the shapefile filling the polygons with the colours according to col-column.

Example Data based on How do I change the polygon fill color and border color for SpatialPolygons objects?

# create data.frame
df <- data.frame(GEOLOGY = c("Quarternary","Tertiary"), col = c("yellow","green"))
# create polygons
Srs1 = Polygons(list(Polygon(cbind(c(2,4,4,1,2),c(2,3,5,4,2)))), "1")
Srs2 = Polygons(list(Polygon(cbind(c(5,4,2,5),c(2,3,2,2)))), "2")
# create spdf
SpDF <- SpatialPolygonsDataFrame( SpatialPolygons(list(Srs1,Srs2)),  df )
# plot spdf
spplot(SpDF, zcol = "GEOLOGY", col.regions = "col")

which results in

Error in .local(obj, ...) : length of col.regions should match number of factor levels

I also tried:

spplot(SpDF, zcol = "GEOLOGY", col.regions = SpDF$col)

which plots the spdf but with wrong colours - red and black?

Any ideas how to plot the correct colours in the correct polygons?

@cengel's suggestion works fine.

After that I added 1 more polygon to my map

df <- data.frame(GEOLOGY = c("Quarternary","Tertiary","Quarternary"), col =  c("yellow","green","yellow"))

Srs1 = Polygons(list(Polygon(cbind(c(2,4,4,1,2),c(2,3,5,4,2)))), "1")
Srs2 = Polygons(list(Polygon(cbind(c(5,4,2,5),c(2,3,2,2)))), "2")
Srs3 = Polygons( list (Polygon(  cbind( c(5,5,4.5,5), c(2,3,2.5,2) ))), "3")

SpDF <- SpatialPolygonsDataFrame( SpatialPolygons(list(Srs1,Srs2, Srs3)), df )
## plot with labels for polygons
spplot(SpDF, zcol = "GEOLOGY", col.regions = as.vector(SpDF$col),
   sp.layout = list(
    list("sp.text", loc = coordinates(SpDF)[1,], txt = SpDF$GEOLOGY[1],cex=0.6),
    list("sp.text", loc = coordinates(SpDF)[2,], txt = SpDF$GEOLOGY[2],cex=0.6),
    list("sp.text", loc = coordinates(SpDF)[3,], txt = SpDF$GEOLOGY[3],cex=0.6)))

removed the typo. Now the Quarterary polgons get a green fill and the Tertiary polygon and yellow fill colour which is not according to the data.frame. It seems to pure coincidence if the colours match?


texts: I want to upload my local web GIS created by GeoServer to the internet

texts: I used GeoServer to create a web GIS in my localhost 8080. I really want to know that can I upload my local web to the the internet.


texts:
In the latest version of Qgis (2.12.3) when trying to create a filter the feature of right clicking and selecting a value does not work. What am I missing ?

I have attached a picture :

At this stage I want to right click and make a selection. However nothing happens. The data refers to a joined table and the attribute table does show values.


texts: Visualize overlapping line vertices in QGIS

texts:
I am trying to inspect the vertices of some linestrings in QGIS 2.4. My goal is to visualize every vertex of a line and easily find those that are duplicated. To clarify this, for an example line composed of vertices a, b, c, d in following order:

  a, b, b, c, d, d

I want to be able to clearly distinguish vertices b and d from a and c.

Right now I have two symbol layers for the line, a simple line and a marker line with markers placed on every vertex. I tried to lower the simple marker's fill's alpha channel so that overlapping markers are less transparent but this doesn't seem to work. (Setting the marker's transparency via the transparency slider seems to be bugged).

Anyone know of any practical solutions? (e.g. inspecting the actual coordinate values is not ideal)

I am after something like this (without labels of course):

To clarify, only points b and d occur twice. Points a and c (not labelled here) occur only once.

EDIT

Although lowering the alpha channel doesn't appear to work initially, the visual effect interestingly appears when you go in editing mode. No idea why this happens.


texts: Generating covariance matrix from ENVI binary file of SAR polarimetric image?

texts:
I have ENVI binary files of SAR polarimetric image.

How can the covariance matrix [c2] be generated from the image?


texts: Does CartoDB automatically fill in Polygon data for city/state?

texts:
I'd like to upload a spreadsheet that has custom city/state information and have CartoDB automatically fill in the polygon data.

If this is possible, then how?


texts: Removing unmatched points from geocoded FeatureClass using ArcObjects?

texts:
I have geocoded some locations but I need to remove the unmatched locations that are in the resulting FeatureClass. How can I do this programmatically with C#? I am able to get the unmatched features but I can't seem to figure out how to actually remove them from the FeatureClass.

public static void RemoveUnmatchedFeatureFromFeatureClass(string workSpace, string featureClass)
    {
        IWorkspaceFactory workspaceFactory = System.Activator.CreateInstance(System.Type.GetTypeFromProgID("esriDataSourcesGDB.FileGDBWorkspaceFactory")) as IWorkspaceFactory;
        IFeatureWorkspace featureWorkspace = workspaceFactory.OpenFromFile(workSpace, 0) as IFeatureWorkspace;
        IFeatureClass unMatchedFeatureClass = featureWorkspace.OpenFeatureClass(featureClass);

        IQueryFilter queryFilter = new QueryFilterClass();
        queryFilter.WhereClause = "Status = 'U'";
        using (ComReleaser comReleaser = new ComReleaser())
        {
            IFeatureCursor newCursor = unMatchedFeatureClass.Search(queryFilter, false);
            IFeature unmatchedfeature = newCursor.NextFeature();

            while (unmatchedfeature != null)
            {
                //Remove feature

                unmatchedfeature = newCursor.NextFeature();
            }
        }
    }


texts: Any issues going back-and-forth between QGIS and ArcGIS Desktop software?

texts:
I have very simple work tasks at the office (digitizing) that I would like to be able to work on at home because we are exceptionally busy.

I don't want to spend thousands of dollars purchasing ArcView, so I've been considering the possibility of downloading QGIS on my home computer.

Bringing the work home is highly preferable to staying late at the office.

Two questions:

Will I experience any issues with my data moving it between the two packages?
Is QGIS difficult to learn for someone who has only had exposure to ESRI software?

Any advice is welcome (ie - interface differences, known QGIS bugs, etc, etc).


texts:
Sorry if this question is too simple/easy for you guys, but I'm still learning QGIS and pretty much teaching myself where I can.

Is there a way to add data from 2 columns into one new column inside  QGIS ?

I have cattle population for each province in a certain country, and buffalo populations for each province in another column (I entered all these manually) - what I'm looking to do is to make a new column (within the same shapefile/database) that combines both cattle and buffalo populations in each province.

I've tried using the field calculator, get error messages all the time - not sure if I have my functions wrong. With the field calculator, I've been trying to enter a simple function of

Catt Pop + Buff Pop

with a new output field as Real number data, but I get error messages such as:

  Invalid expression syntax. The error message of the parser is: 'syntax
  error, unexpected COLUMN_REF, expecting $end'

Secondly - if I can't do this in QGIS, is the only way to open the relevant dbf file in open office and then edit it and then save and load again? From the manual that I read to learn GIS, this can be done but is not recommended (by the guys who wrote the manual) - does anyone have any links to any foolproof guides on how to do this anyway since I've never done it?


texts:
I have two datasets of the same addresses, one geocoded with Google, the other with Nominatim.

Now I want to compare address1 of dataset1 with address1 of dataset2, then address2 of dataset2 with address2 of dataset2 and so on.

Using QGIS, I tried the distance matrix and the Distance to nearest hub-modules. The problem is, that they only give out the nearest point and not the appropriate point/address.

So far, I didn't try ArcGIS, but am not aware of an easy solution there anyway.

Is there a way to write a python script for my issue or a module I am not aware of?


texts: Building height calculation from stereo Satellite images (Satellite Azimuth and Satellite Elevation)

texts:
Can anybody help me in how to calculate building height from stereo satellite images sensor metadata (Mean satellite Azimuth and Mean Satellite Elevation)

Stereo Satellite  Images Azimuth's  (335.6616, 209.5149)
Stereo Satellite Images elevation   (83.8, 83.8)

I have applied some procedure to get multiplication factor with relative displacement but it is not working

First I calculated 1/tan(satele)* Sin(sataz) & 1/tan(satele)* Cos(sataz)

then vector addition to find multiplication factor


texts:
I'm getting the following error - only on one of the vectory layers in my project. Can anyone help? Have tried a restart and reloading the layer and neither worked!

An error has occured while executing Python code:

Traceback (most recent call last):
  File "C:/Users/Claire/.qgis2/python/plugins\EasyCustomLabeling\EasyCustomLabeling.py", line 506, in runLabel
    WKTLine = 'LINESTRING('+ str(geom.centroid().asPoint().x() +0.0001) +' '+  str(geom.centroid().asPoint().y() +0.0001 ) + ' , '+ str(geom.centroid().asPoint().x())+ ' ' +str(geom.centroid().asPoint().y())+ ')'
AttributeError: 'NoneType' object has no attribute 'asPoint'

Python version:
2.7.5 (default, May 15 2013, 22:44:16) [MSC v.1500 64 bit (AMD64)]

QGIS version:
2.10.1-Pisa Pisa, d20c5b7

Python path: ['C:/PROGRA~1/QGISPI~1/apps/qgis/./python/plugins\\processing', 'C:/PROGRA~1/QGISPI~1/apps/qgis/./python', u'C:/Users/Claire/.qgis2/python', u'C:/Users/Claire/.qgis2/python/plugins', 'C:/PROGRA~1/QGISPI~1/apps/qgis/./python/plugins', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\matplotlib-1.3.1-py2.7-win-amd64.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\nose-1.3.3-py2.7.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\tornado-4.0.1-py2.7-win-amd64.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\backports.ssl_match_hostname-3.4.0.2-py2.7.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\certifi-14.05.14-py2.7.egg', 'C:\\PROGRA~1\\QGISPI~1\\bin\\python27.zip', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\DLLs', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\plat-win', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\lib-tk', 'C:\\PROGRA~1\\QGISPI~1\\bin', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\PIL', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\jinja2-2.7.2-py2.7.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\markupsafe-0.23-py2.7-win-amd64.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\pytz-2012j-py2.7.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\win32', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\win32\\lib', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\Pythonwin', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\Shapely-1.2.18-py2.7-win-amd64.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\six-1.3.0-py2.7.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\wx-2.8-msw-unicode', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\xlrd-0.9.2-py2.7.egg', 'C:\\PROGRA~1\\QGISPI~1\\apps\\Python27\\lib\\site-packages\\xlwt-0.7.5-py2.7.egg', u'C:/Users/Claire/.qgis2//python', 'C:\\PROGRA~1\\QGISPI~1\\apps\\qgis\\python\\plugins\\fTools\\tools', 'L:/AC/00018A Northampton Cathedral/GIS']


texts: Writing Python script to give multiple buffers of different sizes in ArcGIS for Desktop?

texts:
I want to write a loop that will cause the output to be multiple buffer shapefiles.

For example, I want a circle buffer of 50ft, 100ft, 150ft, etc.

Can anyone help?

I'm new to this and am trying to figure it out myself.


texts:
I have an application which access same map from two geoserver using loadbalancing concept
the map is accessed using WMS service.
where the two geoserver connected to a single same POSTGIS database.
even the both tomcat is alloted with 3GB of RAM but the map rendering is still too slow

Is this have anything to do with the performance of POSTGreSQL??

is there any way to increase map availability ?

Please help me


texts: Installing third party python library in QGIS 2.18.x Windows

texts:
We are running into a strange behaviour while trying to use pyoidc with QGIS 2.18.x windows standalone version.
Due to pip/pipy no longer working with QGIS 2.18 (see this post), we installed the pyoidc package and its dependencies manually by getting the corresponding wheel packages and running pip from the QGIS python environment though the OSGeo4W.bat console.

Installation runs fine and the pyoidc lib can be used in the python interpreter of the OSGeo4W console.

from oic.oic import Client
from oic.utils.authn.client import CLIENT_AUTHN_METHOD

client = Client(client_authn_method=CLIENT_AUTHN_METHOD)

However, when trying to use the same instructions in a qgis plugin or inside the QGIS python embeded console, we get an import error :

  OSError: Cannot load native module 'Cryptodome.Cipher._raw_ecb': Trying '_raw_ecb.pyd': cannot load library '.......\Cipher_raw_ecb.pyd': error 0x7e

Does anyone have already faced this kind of problem or have any idea what is the reason why the import works well using the python console and does not work inside QGIS ?

We have already checked file permissions and nothing seems different form other libs.
We have also tried with the OSGeo4W network installer, but with the same results...

We are using QGIS 2.18.17, python 2.7.5 and Windows 10 professional (1709 - 16299.371)


texts: Unexpected behaviour when summarizing a field in ArcMap (Un-Equal strings)

texts:
I couldn't figure out why when I use summarize the Field " LandUse_Arabic " it the summary shows two identical values in the output table.

But when I use select by attribute, despite the fact that ArcMAP 10.2.2 displays the 2 values It selects all values which are 47

What might be the issue here ?


texts:
I am using the proximity algorithm to generate a raster proximity map indicating the distance from the coast (points shapefile) of my study area and the center of each planning units (I have divided my study area in almost 5000 square planning units).
I would like to know the difference between the following values I have obtained:

-distance count;
-distance sum;
-distance mean;
-distance minimum;
-distance maximum.

and which one I should consider in order to understand the importance of the distance from the coast to each planning units.


texts: Strategies for managing (and ultimately merging) duplicate business locations in a spatial database

texts:
I am working on a project right now that involves using multiple (very large) "location directories" to identify with the highest possible inclusiveness, a diverse set of business locations.

These locations are often geotagged close to one another but typically not directly on top of one another.

For example, if a walmart was in one database at the NW corner of its parking lot, in another database it might be in the SE corner. In this example the two points might be as far as 2000ft apart. I'm using the extreme examples to create some basic buffer rules (after consolidating brandnames to identify same-brand store locations) for allowing my program to identify duplicate locations.

But there are also some business locations that have closed. Some are geocoded directly on top of the "center of town" (i.e. they had an address of "Cityname, State" and got geocoded incorrectly in the middle of town).

There are a lot of other edge cases and I guess what I'm asking is, is there some python code package or qgis plugin or ruby gem (haha, just kidding) or some other larger abstract object oriented design pattern (way of thinking) that somebody learned as a Computer Science major to handle (entity validation & consolidation) problems like this?

I'm perfectly will and able to cook up my own solutions but don't want to miss out on sophisticated techniques and/or waste my time reinventing potentially less-effective, more time-consuming solutions for managing entity databases.

Am perfectly willing to modify, close this question or migrate it to another forum if a mod or mods thinks this something is not appropriate here. Please comment before taking any administration or negative voting action.

The technology I'm using currently is a ruby on rails database running a postgis database. I extensively utilize activerecord to manage my data. I also use qgis plugins a lot and am building a growing capability set with python though still very much a beginner.


texts: Export background tile map to raster file

texts:
In QGIS I managed to display the ship traffic density from marinetraffic.com as a background tile map using the TileLayer plugin. Here is my TileLayer.tsv file:

#title  credit  serviceUrl  yOriginTop  zmin    zmax    xmin    ymin    xmax    ymax
MarineTraffic   Mtra    http://tiles.marinetraffic.com/ais/density_tiles/{z}/{x}/tile_{z}_{x}_{y}.png   1

for example one of the tiles has the url:
http://tiles.marinetraffic.com/ais/density_tiles/8/164/tile_8_164_109.png

For a selected area, Persian Gulf, I need to convert the map to a GeoTiff raster file in the UTM zone 39N projection. When I enable on-fly projection in QGIS the background map is correctly reprojected:

But how do I export it to a raster GeoTiff file? I couldn't find any export of "save as" option for the background tile layer map in the QGIS menu.


texts: Cropping and transforming geoPDF with GDALwarp

texts:
I got some great advice on how to transform a geoPDF to a geoTIF at Historical U.S.G.S. Georeferencing, which used the available downloads from this site. I would like to take it a little further and not only convert, but remove the border of historical USGS maps in order to be able to tile them in QGIS. Being a basic GIS user, I get lost in the "coding" involved to implement these practices.

How can I convert and crop historical USGS maps so they can be tiled in QGIS?

I am attempting to do this again, it has been a while. I made it through the steps and I get a cropped quad that has back areas along the edges. I don't believe this happened when I did this in the past.


texts: How to sum values of overlapped pixels in thousands of rasters layers and clip it using a polygon in QGIS?

texts:
I have a lot of rasters files (more than 6000), each one represents the presence and absence of a species  (1 for presence and NoData for absence). The set of rasters covers the whole Latin America and they can overlap eventually. I also have a vector polygon layer which bounds an area in which I am interested in.  As a final goal, I would like to sum all the species in each grid in my area of study.
Any ideas on how to achieve this?
I have tried to sum all raster using Grid Sum in Saga and in order to clip it afterwards. Problem is that the output has only NaN and zeros values - which does not make sense... I used a sample of 100 rasters, in order to test it. And it lasts more than two hours! Any ideas on how to do it faster?!


texts: Generating random but believable digital elevation model?

texts:
Is there a way to generate a DEM, programatically or otherwise, that will feed into ArcGIS Desktop for further spatial analysis?

Perhaps this needs to be broken down into smaller incremental steps:

Generate a grid
Fill grid with values where: 0 > value < maxElevation
Neighbouring cells: (x1-x2) < maxSlope


texts: How to get rid of lines within merged poligons?

texts:
Sometimes, when I merge polygons there are still some lines showing within them (see screenshot).

How can I get rid of them?


texts:
I am trying to install Python plugins (to develop QT/Python plugins for QGIS).
I see the following error when I try to run the command:

Compile.bat:

@echo off
call "C:\OSGeo4W64\bin\o4w_env.bat"
call "C:\OSGeo4W64\bin\qt5_env.bat"
call "C:\OSGeo4W64\bin\py3_env.bat"

@echo on
pyrcc5 -o resources.py resources.qrc

Commands:

C:\OSGeo4W64\bin>Compile.bat

C:\OSGeo4W64\bin>pyrcc5 -o resources.py resources.qrc
C:\OSGEO4~1\...PyQt5\pyrcc_main.py: File does not exist 'resources.qrc'

How do I solve this error?


texts: Displaying multiple items at same address using Carto?

texts:
I'm creating a map with specific information (name, DOB, etc) to be displayed for a street address, which is mapped to an exact lat and long. Problem I'm realizing is that for items such as apartment buildings, this isn't going to work within CartoDB using the Simple Visualization Wizard.

Any suggestions in CartoDB to display multiple fields of the same type in a click/hover window?


texts:
I have a series of elevation points from transect surveys completed on a river. I would like to interpolate this data to reveal the bathymetry of the river. The survey points are shown below:

I have tried the default options on topo to raster (cell size 109) and then used extract by mask to focus just on the river channel. It formed this raster below. The transect locations stand out clearly however I am hoping for a more continual surface. [

I then tried IDW with a barrier polyline and this reduced the cell size down to 3. However there appears to be "seams" at the mid points between the transects.

Can anyone recommend a method to interpolate these points into a smoother surface?


texts: GeoServer OpenLayers and WFS Performance

texts:
Currently I'm using GeoServer to dispatch WFS GeoJSON data to OpenLayers using Bbox strategy, the map is made of a WMTS base layer via loaded via MapBox API, although the WFS layers are not so big features are highly concentrated in certain areas and the OpenLayers performance around this areas is a bit slow and laggy, not so bad but it's certainly not the best...

The question is, how could I improve this performance?


texts:
I have a problem that must be solved in automatic manner.
I am using ArcGIS 10.2.

I have a raster files that I have to pair and then make the calculation. I have tried to use ModelBuilder but there is not possible to use two iterations. I have made submodel but then it did not pair the layers correctly.

THe calculation on the pairs is simple. Just minus them than the results must be divided by the maximum value from the resulted layer.

I have :

raster1_yearA, raster1_yearB, raster1_yearC, raster1_yearD
raster2_yearA, raster2_yearB, raster2_yearC, raster2_yearD

But I have to make pairs lik this:

raster1_yearA - raster1_yearB
raster1_yearA - raster1_yearC
raster1_yearA - raster1_yearD


texts:
In a Leaflet app I want to overlay raster data (technically speaking, a field) on a web map.
The data consists of triples < lat,lon,value >, where value is in range [0,1], for example:

23.33,10.55,0.54
23.331,10.55,0.5
23.331,10.56,0.4
...

I tried with a heatmap plugin, but it doesn't really support field visualization (it's a visualization of points, regardless of their values).
This is an example of what I mean from a GIS.

How can I achieve that in Leaflet, without switching to a full-fledged GIS?


texts: Need help understanding why gdal projections aren't giving me what I expect?

texts:
When using gdaltranslate and gdalwarp I specify the following transfomation (in for GCP's and translate, out for warp)

+proj=lcc +lat_1=40.03333333333333 +lat_2=38.73333333333333 +lat_0=38 +lon_0=-82.5 +x_0=1968500 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=us-ft +no_defs

This is intended to be the proj4 equivalent of Ohio South state plane US feet. If I open it in GlobalMapper it is placed exactly correctly (about -83 lon) though it isn't identified as Ohio South.  However when I run listgeo the coordinates in feet are listed, but the lat long are not shown, and if I run gdalinfo the coordinates come back as

Upper Left  ( 1747866.344,  712518.960) ( 99d 6'46.00"W, 38d45'47.64"N)

with the longitude being completely off base. I also notice that the false_easting in gdalinfo come up as

"false_easting",6458320.416666666

which I also don't understand since I set it directly in the projection as 1968500'

In the end what I really want to do is digitize some coordinates in state plane, US survey foot, run them through translate to create a VRT file, then warp that and output Geotiff that has a recognized projection and will read directly into Global Mapper.  I tried just using a prj or wkt file, but once I edit it to change the units to feet then the resulting Geotiff has the following tags either deleted or set wrong.

GTModelTypeGeoKey (Short,1): ModelTypeProjected (ends up user-defined)
ProjectedCSTypeGeoKey (Short,1): PCS_NAD83_Ohio_South (ends up missing)

Any thoughts on what I am doing wrong would be greatly appreciated.


texts: ArcMap Export PDF vs arcpy Export PDF scale differences

texts:
I'm experiencing differences from exporting features to PDF from within ArcMap and using a stand alone python program. Basically, I'm selecting a feature from a layer, zooming to the selected extent, changing the scale of the dataframe, and exporting the image to a pdf from PAGE_LAYOUT view- something pretty simple.

The issue that I'm having is that from within ArcMap, everything works fine. In the page layout view, everything looks as expected, and once exported, it looks great.

I'm trying to automate this so after building a stand alone python program, it exports like it's supposed to, however the scale is off by about 200.0.

Here's a snippet of my code:

mxd = arcpy.mapping.MapDocument("C:\Path\To\MxdFile.mxd")
dataFrames = arcpy.mapping.ListDataFrames(mxd)
df = dataFrames[0]

layers = arcpy.mapping.ListLayers(mxd, "*", df)
for layer in layers:
    if layer.name == "LayerIWantToUse":
        break

if layer.name != "LayerIWantToUse":
    raise Exception("Cannot find layer.")

arcpy.SelectLayerByAttribute_management(layer, "CLEAR_SELECTION")

query = "OBJECTID = 1"
layer.definitionQuery = query
arcpy.SelectLayerByAttribute_management(layer,"NEW_SELECTION", query)

df.zoomToSelectedFeatures();
df.scale = 2400.0
fileName = r"Python_2400.pdf")
arcpy.mapping.ExportToPDF(mxd, fileName, "PAGE_LAYOUT", 1056, 1632)

The scale that it actually exports, appears to be 1:2200, instead of 1:2400 like it's supposed to. I've attached both images to show the differences.

I have noticed that if I save the MXD with a different scale (1:100000 for example) and with ArcMap in "Data View", I have the scale issues. If I make sure "Layout View" is selected, change the scale to 1:2400, and save the MXD, I do not have the scale issues within python/arcpy.

I've seen a couple forum posts regarding something similar, some saying it's a bug that should have been fixed in 10.1, but a) I'm on 10.2.1, and b) no explanations on what could be causing it. Another post recommended that if you change the extent/scale twice, and then export the image/pdf twice, it might fix itself, but that solution didn't work for me either.

Has anyone else encountered something similar, and if so do you have any recommendations?


texts: Creating QGIS Plugin starting from Graphical Modeler

texts:
I'd like to know how can I create a QGIS Plugin starting from a model built with the graphical modeler. Is it possible or i have to use necessarily scripts python made in the console of QGIS?
Any suggestions? I have to use "Plugin Builder"?


texts: Leaflet Buffer value by GeoJSON feature.properties

texts:
Is it possible to add a buffer to a list of points, taking the value of the radius from the properties of the GeoJSON file?

This is my code to take points from JSON files and represent them as cluster markers. I would like to add a buffer equal to the value of the buffer property in the JSON file to the points. I searched the turf documentation but didn't find something like this.

var Ar_clusters = L.markerClusterGroup();
$.getJSON('a.json', function(data) {
  var Arc = L.geoJson(data, {
    pointToLayer: function(feature, latlng) {
      return L.marker(latlng, {
      icon: L.icon({
        iconUrl: "img/ar.png",
        iconSize: [40, 40],
        iconAnchor: [12, 28],
        popupAnchor: [0, -25]
      }),
      title: feature.properties.name,
      riseOnHover: true
    });
  },
    onEachFeature: function (feature, layer) {
            if (feature.properties) {
      var content = "
" + "
Name	" + feature.properties.name + "
" + "
Info	" + feature.properties.popupinfo + "
" + "
";
      layer.on({
        click: function (e) {
          $("#feature-title").html(feature.properties.name);
          $("#feature-info").html(content);
          $("#featureModal").modal("show");
          highlight.clearLayers().addLayer(L.circleMarker([feature.geometry.coordinates[1], feature.geometry.coordinates[0]], highlightStyle));
        }
      });
      $("#feature-list tbody").append('
	' + layer.feature.properties.name + '
');
      ArSearch.push({
        name: layer.feature.properties.name,
        address: layer.feature.properties.ADDRESS1,
        source: "Ar",
        id: L.stamp(layer),
        lat: layer.feature.geometry.coordinates[1],
        lng: layer.feature.geometry.coordinates[0]
      });
      }
      layer.addTo(Ar_clusters);
    }
  });
});


texts: How to choose a latitude-longitude among conflicting sources

texts:
Suppose you use various geocoding facilities (like Mapquest's API, Texas A&M's, etc.) to get a latitude-longitude pair for a given address in the United States. They give you differing answers, sometimes miles apart from one another. How should you choose among them?

Well, I suppose you can inspect the various answers on Google Maps or similar and see which makes the most sense. But now suppose you're doing this for a large dataset, where it is very impractical to inspect all the answers; you want an automated way to do it. What should you do?

I will gladly accept answers of any of the following forms:

a way to choose among the answers obtained

which might simply be a single authoritative data source to believe even if others disagree (in which case there's no need to check multiple sources in the first place, unless the good one has missing data)

a way to find some sort of mean value among the answers obtained
other


texts: How to apply the four colors theorem in a polygon map in ArcGIS/ArcToolBox automatcally?

texts: I need to apply the four colors theorem in a polygonal shape in a way that I do not need to choose manually each color to put in each region. I wish to knows if there is any extension, plug-in, script or data base that may be used with ArcGIS and ArcToolBox to do it mathematically or programatically, so I could use it for now on with every map I come to create.


texts: NoData values in classification accuracy assessment

texts:
I want to assess the accuracy of a settlement classification I made (classification map) by comparing it to another product's settlement classification (reference map). Both maps have two "normal" classes, but also have large areas where the land use is not classified, so they rasters are full of noData values.

For the error matrix, do I convert noData into a new class for which sample points will be created? And do you then also calculate overall, producer, user accuracy etc. for this noData class.

I didn't find any real info how to deal with noData in these assessments, but got the impression that people actually create these sample points but do not include to calculate the performance indices...


texts: OpenLayers 3 how to get only the modified feature?

texts:
I have a problem with OL3 and it's ol.interaction.Modify Event. Whenever I modify a feature, just one, it returns all features inside the layer. Currently I have 77 features (polygons) and I can modify them with the modify interaction. Inside the callback function it returns the event with the feature member array and this seems to be buggy:

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

    console.log('=== Modify end ===');
    var modifiedFeature = e.features.item(0);

    //Why do I get 77 modified features if I have only touched one?
    // e.features size is = 77
    console.log(e.features);
});

Do you know a working way to get the one modified feature?


texts: Choosing bandwidth and scaling factor for Kernel Density Estimate in ArcGIS Desktop?

texts:
I have some points that represent animal locations over the course of a few months at the scale of ~20 sqkm. I am emulating the work of previous researchers in this field, but I'm having trouble understanding why they chose the parameters they did. Specifically, I want to know how they chose the bandwidth and the scaling factor. Furthermore, what is the difference between the scaling factor and the bandwidth? It is my understanding that the scaling factor helps to transform bandwidths between different types of kernels, but I don't understand how that would apply with just one kernel in this case.

Since they didn't provide an explanation of how they arrived at the parameters, I was hoping someone would be able to estimate, back of the envelope, if these parameters look appropriate for the locations and kde provided. If not, how would I determine an appropriate bandwidth and scaling factor myself? The irregular rectangle depicted is roughly 15km x 15km. The citation is Management of pinyon-juniper woodlands at Kirtland Air Force Base: Pinyon Jay summer and winter home ranges and habitat use 2009 final report and below.

Referenced literature:

  For the kernel analysis we used the Hawth’s Tool, Fixed Kernel Density
  Estimator. Kernel analysis is a nonparametric statistical method for
  estimating probability densities from a set of points. When used to
  analyze home range data, kernel density methods describe the
  probability of finding an animal in any one place. The Fixed Kernel
  Density Estimator calculates a fixed kernel density estimate and
  produces contour lines representing the boundary of the area that
  contains a specified percent of the volume of a probability density
  distribution. A 95% volume contour, for example, typically contains
  95% of the points used to generate the kernel density estimate. Kernel
  parameters were set as follows: scaling factor- 1000000, kernel-
  bivariate normal, single parameter smoothing factor (h)-1000, raster
  cell size-100, and percent volume contours- 50, 90, and 95.


texts: Storing GPS location in postgis efficiently?

texts:
I'm getting info from GPS And I have 4 data :

latitude
longitude
altitude
precision

Can I store the precision in M coordinate in PostGIS?

Ex :

 CREATE TABLE gps_pos ( id uuid primary key, position Geometry(Point, 4326))
 INSERT INTO gps_pos (id, position) VALUES (<uuid>, ST_MakePoint(0.0, 0.0; 400.0, 20)

I never used this 4th coordinate and don't know the best practice about it.


texts:
I am using OpenLayers version 3 (answers with ver 2 are also OK) and want to measure how long it takes to render the vector layer. So, I have a vector layer with features:

var vectorLayer = new ol.layer.Vector({
    source: vectorSource
});

And my map object:

var map = new ol.Map({
        renderer: 'canvas',
        target: document.getElementById('map'),
        view: new ol.View({
            center: [550000, 6520000],
            resolution: 500,
            minZoom: 1,
            maxZoom: 10,
            zoom: 1,
            resolutions: [500.0, 250.0, 125.0, 62.5, 31.25, 15.625, 7.8125, 3.90625, 1.953125, 0.9765625, 0.48828125],
            projection: projection
        })
    });

I want to measure how long it takes to render all features from the layer. If I do like this

 var start=performance.now();
 map.addLayer(vectorLayer);
 var time=performance.now()-start;
 console.log(time)

then I see time in console much much earlier than points on map. I want to see the time right after all points are rendered. How can I do that?


texts:
I am working on python plugins for QGIS.I developed the plugins same as RT SQL LAYER plugins.

When we draw the  layers using query builder of "Add Postgis Layer" option,toggle editing option is available for the same.

But when i try to draw layers using RT SQL LAYER plugin ,it is not available for editing.
Same happens with my plugins too.

How do I make toggle editing option available for my plugins.??
My plugins are in python.


texts: Creating Feature Shape Templates in ArcMap?

texts:
Working in oil & gas, I find myself drawing the same polyline shapes over and over again in my company's pipeline database.

For example, a 12" valve set might always have the same dimensions, a 24" pig launcher might have a couple different designs, a 12" drip might have a couple different designs, and so on.

Most of these items require multiple separate polyline features each with its own attributes (such as a 12" mainline pipe and a 4" bypass pipe) and also the possibility of having multiple separate valve and fitting features that would be snapped to the polylines.

Here is an example of a Pig Launcher which I would like to be able to replicate multiple times in the future without having to draw each individual item every single time.  As you could imagine, this process is time consuming, but it is critical to properly locating unique valves and other appurtenances.

Here is an example of a fairly simple bypass valve set that would also be a good candidate for replication.

Is there a way to have pre-made shape templates that I can just select my template, and rotate and snap into place?


texts: Exporting Geodatabase with Subtypes and Domains to Microsoft Access Without Losing Information?

texts: I want to export a geodatabase file with all subtypes and domains values to open in Microsoft Access without losing any information.


texts: How can I get a shapefile of coverage from the Planet Labs scenes API?

texts: I'm using the Planet Labs Scenes API and need to get the results of a (large) search I'm making as a shapefile. Any hints on how to go about it? I'm comfortable with command line tools like the command line interface of the Planet Python Client and OGR, but not quite up for custom python coding.


texts:
I have 1.000 point layers (shapefiles).

How can I find points intersected to a line without defining the points layer?

Because defining thousand layer takes time. And I dont wanna merge those point layers because I need it for visual purpose.

I make a custome line which has branch.
I already have a solution to 'sum' each value in points that intersected in line from previous question but not finished. This is the function :

  aggregate(layer:='point1',aggregate:='sum',expression:="number",filter:=intersects($geometry, geometry(@parent))) +
  aggregate(layer:='point2',aggregate:='sum',expression:="number",filter:=intersects($geometry, geometry(@parent)))

1st problem : i got wrong sum result after the first line i draw ( i want duplicate points only count once)

2nd probelm : imagine i have 1000 layers of points, it is so hard to type them all. So i need generalization for N points layer or any solution which doesnt require to define points layer.

UPDATE

I was thinking of another option to get correct sum number.

The picture above is the attribute table for drawn line to get intersecting point and get the point attribute.

Maybe we can merge these two column, so the yellow-highlighted line one group into one.


texts:
We currently updated to SQL Server 2012 and now we are encountering a slower performance when connecting to the geodatabase on SQL Server 2012.

Here are the things that we tried to fix the issue:
- Re-indexed all the tables.
- Installed the support patch from ESRI.(http://support.esri.com/en/downloads/patches-servicepacks/view/productid/160/metaid/1954).

But none of those fixed the performance issue. The machine hosting the SQL Server 2012 is a VM but it is running with 24 gb of ram. So it shouldnt be under resourced.

I found a similar problem and one of the post was saying to use the SDE Binary Storage type. I didn't do that yet because he did mention that its not recommended.


texts: Create within-country grids for multiple countries at once in Google Earth Engine

texts:
I have a dataset on countries in sub-Saharan Africa (which is set so that anybody can use it):

var SSA = ee.FeatureCollection('users/salem043/Africa_Districts')
.select("ADM0")

I would like to produce grids that cover this entire region but are within-country. That's to say, I do not want any grids to overlap a country's border.

For country-level grids, I usually use the country boundary, create a bounding box, then split that box up into grids. Here, for example, is my code for creating 20,000 grids for Iraq:

//Bounding box over Iraq
var lat_start = 29.0990251735;
var lat_end = 37.3852635768;
var lon_start = 38.7923405291;
var lon_end = 48.5679712258;

//Set number of grid cells
var num_cells = 20000;
var lon_edge = (lon_end-lon_start)/Math.sqrt(num_cells);
var lat_edge = (lat_end-lat_start)/Math.sqrt(num_cells);
// 3) Create the grid
var polys = [];
var polys_line = [];
var cell_id = 0;
for (var lon = lon_start; lon < lon_end; lon += lon_edge) {
  var x1 = lon;
  var x2 = lon + lon_edge;
  for (var lat = lat_start; lat < lat_end; lat += lat_edge) {
    cell_id = cell_id + 1;
    var y1 = lat;
    var y2 = lat + lat_edge;

    polys.push(ee.Feature(ee.Geometry.Rectangle(x1, y1, x2, y2), {label: cell_id}));
  }
}

I imagine there is a function that can call on each ADM0 (country boundary) in my feature collection, but I want the grids to be equal size, and it would be very tedious testing each one to figure out what to set the grid count to in order to ensure they have the same dimensions. What else can I do?


texts: Is there a web mapping server works on a linux shared hosting account?

texts:
Question: Is there a web mapping server software out there that works on a typical linux shared hosting account?

Background: I have been trying find a web mapping server that would allow me to serve WMS services, WFS would be a bonus. I do not have root access and this is a shared linux account. It is a typical shared hosting account that can use php, ruby, python, cgi etc. I don't have root access.

Additional Information:

As with most shared hosts java isn't an option.
I can install binaries locally, I just found out that I can request SSH access which I would assume would open up a few more options. Most shared hosting accounts wouldn't have this option.


texts: Openlayers 2 aand Openlayers 3 Code convertion

texts:
I have an old Opanlayers 2 application. And I need to convert an Opanlayers3 code to Opnelayers2, but could not.

Opanlaysrs 3 Code

  var epsg = "EPSG:900913";
  var projection = ol.proj.get(epsg);

  var projectionExtent = projection.getExtent();
  var size = ol.extent.getWidth(projectionExtent) / 256;
  var resolutions = [];
  var matrixIds = [];
  for (var z = 0; z < 28; z++) {
    resolutions[z] = size / Math.pow(2, z);
    matrixIds[z] = epsg + ':'+ z;
  }

I will convet it to Opanlaysers 2 But could not:

  var epsg = "EPSG:900913";
  var projection = new Opanlayers.Projest(epsg);

  var projectionExtent = ???? // projection.getExtent();
  var size = ???? //ol.extent.getWidth(projectionExtent);
  var resolutions = [];
  var matrixIds = [];
  for (var z = 0; z < 28; z++) {
    resolutions[z] = size / Math.pow(2, z);
    matrixIds[z] = epsg + ':'+ z;
  }

projection.getExtent(); is not a function


texts: Clipping polygons to existing plotting area

texts:
I have an existing plotting area onto which I want to add some polygons.

In order to be sure the polygons' labels are displayed within the plotting area, I'd like to clip the polygons to this area so that the coordinates of the resulting object are guaranteed to be within the existing area:

set.seed(30482)
x = runif(20)
y = runif(20)

base_poly = cbind(c(0, 0, 1, 1),
                  c(0, 1, 1, 0))
shift = .5*expand.grid(c(-1, 1), c(-1, 1))
SP = SpatialPolygons(
  lapply(1:4, function(ii)
    Polygons(
      list(Polygon(
        sweep(base_poly, 2L, unlist(shift[ii, ]), '-')
      )),
      ID = sprintf('Polygon %d', ii)
    ))
)

plot(x, y)
plot(SP, add = TRUE)
text(coordinates(SP), labels = sapply(SP@polygons, slot, 'ID'))

As you can see, the centroids of the SP polygons lie outside the plotting area, so they're not really showing up in the plot (in this example, I could use xpd = TRUE to some success but this is not true in general).

My fix is to clip SP to the plotting area like so:

plot(x, y)
pu = par('usr')
SP = gIntersection(
  SP,
  SpatialPolygons(list(Polygons(
    list(Polygon(cbind(pu[rep(1:2, each = 2L)], pu[c(3:4, 4:3)]))),
    ID = '0'
  )), proj4string = SP@proj4string),
  byid = TRUE,
  id = sapply(SP@polygons, slot, 'ID')
)

plot(SP, add = TRUE)
text(coordinates(SP), labels = sapply(SP@polygons, slot, 'ID'))

Is this the canonical approach to this? The construction of the plotting area as a SpatialPolygons object seems pretty hack-y, am I missing some convenience functions for accomplishing this without as much manual input?


texts: Calling One Plugin from another Plugin in QGIS

texts: I am working on various plugins in QGIS which I have made, so I was wondering if it is possible to make one plugin which would call these plugins in a predetermined order, so that I would not have to call them individually.


texts: Prevent symbols from overlapping in qgis

texts:
I have mutiple of points layers that I want to represent in a map. Because of the reality of the field and some  errors of mapping,  some layers are overlaying so that the printed map is not good and not readable. I'm wondering if there is any possibitilty to manage them so that they would more clear and visible.

you can clearly see the mess on the map. I would be pleased to get some help from you.


texts: Grouping carpool participant addresses for network

texts:
I have a list of addresses of people participating in carpools to a workplace and a list of cars given by the workplace. Both lists change daily.

I need to group people by their addresses and assign them cars by the group size (and then calculate their route to work, but that's the easy part).

How can I group points (worker addresses) by their geographic location, or better yet, by their closeness on the road network?

I'm using ArcInfo 10.1. I tried using the Grouping Analysis (which groups by geographic location and distance from the workplace), but it's giving me weird results.


texts: which map projection should I use to calculate distance in Poland

texts:
I'm using/learning geodjango + postgis, and I'm trying to calculate distance between 2 cities in Poland. First using WGS84/4326 I define 2 points, and then after transforming this points into spherical mercator 900913 i get distance which is about 40km too long. Which map projection should I use then? Where can I find proper srid's for different countries?

sorry for the delay, 2 cities (lat; long): Slupsk(54.465249;17.026062) and Koszalin(54.193378;16.179428)


texts: Calculate frequency of a attribute field per polygon

texts: I have spatially joined a polygon file with census tracts of areas in Houston with a theft occurrence point shapefile. The attribute table contains a field that represents different classes of theft (ie. Class A: $0-$1000, Class B: $1001 -$2000 etc). I want to calculate the frequency of classes in each of my census tracts. So in essence how many different theft classes occur in each of the census tracts. I have figured out how to get the count, but this is not what I want. How would I go about obtaining and representing the frequency?


texts: How do I create a shapefile of circles of varying radii around lat/long points and map the overlap onto another shapefile?

texts:
I have a shapefile of US counties downloaded from the Census, a set of points, and a calculated radius around each point. I'm trying to draw a circle of that radius around each point, then plot the overlap and do various calculations on that overlap.

The calculations come later, but right now, when I run the following code, I get an error. Here is a reproducible example:

library(sf)
library(tidyverse)
library(USAboundaries)

albers_epsg <- 42303 # 102003

# download county boundaries
counties <- us_counties(map_date = "2000-01-01", resolution = "high",
                        states = c("Utah", "New Mexico", "Colorado", "Arizona"))
st_crs(counties)
plot(st_geometry(counties))
counties <- counties %>% st_transform(albers_epsg)

# initialize point and project latitude and longitude to a coordinate system
# that prioritizes accurate distances
lon <- -111.101099
lat <- 40.477140
elev <- 150
tower <- st_point(x = c(lon, lat, elev), dim = "XYZ")
tower <- tower %>% st_sfc(crs = st_crs(counties)$epsg) %>% st_transform(albers_epsg)

# draw the radio horizon around the tower
scalar <- 3.57
radius <- scalar * sqrt(tower[[1]][3]) * 1000
tower_buffer <- st_buffer(tower, dist = radius)

plot(st_geometry(st_intersection(counties, tower_buffer)))

I understand that I have to project the latitude/longitude points from EPSG 4326 (which is how the shapefiles come from the Census) into a form that can measure distance. I care about accurate distance most of all, so I chose the Albers, which is EPSG 102003 or EPSG 42303.

When I run the code above, however, I get this error:

Error in CPL_transform(x, crs$proj4string) : OGR error
In addition: Warning message:
In CPL_crs_from_epsg(as.integer(x)) :
  GDAL Error 6: EPSG PCS/GCS code 42303 not found in EPSG support files.  Is this a valid
EPSG coordinate system?

The error is the same if I use EPSG 102003 instead. Does R not support these projections? I'm using R 3.5.1 that I built from source on a Debian 9 machine.


texts: Spatial join with a calculated WHERE comparison for >1 result

texts:
I am attempting to identify households which are downstream and within 5km of another set of point features (water sources). We are defining downstream as being at a lower altitude than the water source, so water_source_altitude - household_altitude > 0.

In several instances I have water sources with overlapping 5km buffers, and there are households that fall within this overlap. If household A falls within a 5km buffer of water source 1 and water source 2, is closest to and upstream of water source 1, and is downstream of water source 2, then a spatial join will not identify the correct water source for this household. Is there a way to include criteria which would calculate (water source 1 altitude - household altitude) and (water source 2 altitude - household altitude), and choose the water source where result > 0? I've looked at a few other questions about spatial joins and criteria, but they aren't quite what I'm looking for.

By the way, I am running a Basic license of ArcMap 10.2.2 and my knowledge of python is very limited. I feel like the answer could involve PostGres or PostGIS, but I am not quite sure where to start with those.


texts: Why can't I multiply two string fields in the Field Calculator by converting them to int() and float()

texts:
I'm using ArcGIS 10 and new at python scripting.  I have a string field that I need to calculate by multiplying two other string fields.  I can do it in python using two arbitrary values...

CINPT1 = "1.00"
CNETVAL = "8500"

print int(float(CINPT1) * int(CNETVAL))

but...I need this as an arcpy fieldcalculation.  What I have now is below. I'm basically trying to what is above: the CNETVAL field should equal CINPT1 * CUNITV

import arcpy, sys, string, os

MAPCHAR = "MAPCHAR"
CINPT1 = "CINPT1"
CNETVAL = "CNETVAL"
CUNITV = "CUNITV"

expression = "multiStringFunc(!CNETVAL!)"
codeblock = """def multiStringFunc(netvalue):
    value1 = float(CINPT1)
    unitvalue = int(CUNITV)
    netvalue = int(value1 * unitvalue)

    return netvalue"""

# Process: CUNITV = CINPT1 * CUNITV
arcpy.CalculateField_management(MAPCHAR, CINPT1, expression, "Python", codeblock)

What i get in return is "global name 'CINPT1' is not defined".  I am pretty new at using python, so i don't know if I'm just using the value wrong within the codeblock or I'm just totally off.  Anything you can contribute would be appreciated.  Thanks!


texts: Styling a polygon layer using external graphics

texts:
I have a layer which is comprised of a number of different polygons symbolised according to their type.

I have set up a style which successfully filters by type then applies a corresponding fill (in some cases overlays this with a hatch style) and this works great.  Now I want to apply the same filters but in some instances apply an external graphic with a transparent background instead of a hatch style.

I have used the sld code found in the Geoserver SLD Cookbook and this works fine on its own but when I try and use it in conjunction with a filter it fails to work.

Has anyone come across this same issue and if so can you help me please?

Thanks


texts:
Using CARTO VL.

Having a dataset, is there a way to fetch all columns this dataset has before I create the viz object so I can programmatically decide which variables to add to the viz?

Something like:

const source = new carto.source.Dataset('my-dataset');
const columns = source.getColumns();


texts: Count points from two tables in polygons with PostgreSQL/PostGIS

texts:
I've got

CREATE TABLE points1 (
    point1id varchar PK,
    geom1 geom,
    date1 date);

CREATE TABLE points2 (
    point2id varchar PK,
    geom2 geom,
    date2 date);

Both with indexes. I need to count points (also zeros) within a UTM grid and this seems to work fine for one table

SELECT utm.utmid, count(points1.geom1)
FROM utm LEFT JOIN points1
ON st_contains(utm.geom,points1.geom1)
GROUP BY utm.utmid;

But I need to obtain only one table with something like

  utmid / Date / # points1 / # points2


texts: Query from plugin that shows result on Map and in Attribute Table

texts:
I'd like to create a plugin that performs a static (no user input required) query that, based on a condition, selects all the features on the map and displays these in the attribute table.

I have the creation of the plugin itself under control.

However.

1) Querying the PostgreSQL DB with PyQGIS from the plugin and,

2) based on the result, select the features on the map (I'm not out to create a new layer),

3) and showing the selected features in the Attribute Table, is where I could use some hints.

Anyone that has experience with this or knows of one or more plugins where I can get some really good inspiration?


texts: Delete Table Rows Based On Two Field Values

texts:
How can I delete a row in a table if value A from Field 1 does not match value B from Field 4 ?

I attempted to use code from a previous post as I figured it would be a start but didn't work.

import arcpy
shp = r"C:\path\to\your\shapefile.shp"
rows = arcpy.UpdateCursor(shp, "", "", "some_field")
for row in rows:
    if row.some_field == 2:
        rows.deleteRow(row)
del row
del rows


texts:
I need to add a header to my GeoExplorer application in the same way they do here http://demo.geonode.org/maps/new. Where should I modify the code to achieve this?

I try adding Ext.Panel within composer.js, but not the solution I think:

var northPanel = new Ext.Panel({
        region: "north",
        id: "north",
        height: 100,
        border: false,
        split: true,
        collapsible: false,
        collapseMode: "mini",
        collapsed: false,
        hideCollapseTool: true,
        header: true,
        layout: "border",
    renderto: function (){console.log("em composer ");console.log(document.body)},
        items: [{
            region: "center",
            id: "header",
            title: "header",
            layout: "fit"
        }]
    });

Also detected in the GeoExplorer/app/viewer.html there html code file, the one shown at the top of this link http://demo.geonode.org/maps/new, but I can not make this code display my GeoExplorer application.


texts: Mapping polygons into lines and projecting lines on road network

texts:
I have a layer with polygons that describe road damages. I have then another layer with road network lines. The layers are shown in the picture.

I need to project the polygons as lines to the road network. How can I iterate over each polygon, find its two extreme points in relation to the nearest road line and project the line connecting the two extreme vertices onto the road line? Are there any tools for achieving the result in QGIS or is necessary to make an own Python script?
The final output is supposed to be something like the picture below. The colored lines represent the polygons and cover the road network.

Following the steps from the answer, I end up with the points being in a way randomly placed on different roads as shown on picture below (the projected_vertices layer is marked with yellow)


texts: How to delete a feature from attribute table in ArcGIS Online?


texts: How to remove a feature style in openlayers?

texts:
I need to remove a feature style to set up another one. How can I do it? Do I just have to set up the new one and the last one will be removed by default?

I've looked at the documentation and did not found any answer to this.


texts:
I'm playing with QDockWidgets. Here is the piece of code i'm trying :

def dock(self):
    widget = uic.loadUi("C:\dock.ui")
    self.iface.addDockWidget(Qt.RightDockWidgetArea, widget)
    widget.setFloating(True)
    widget.setWindowTitle("Dock ...")
    #QMessageBox.information(None,"Testing dockWidgets ...","Working or not ... ?")

When the QmessageBox line is left commented nothing works (the widget doesn't appear on the screen nor in the panel list of QGIS)

When the Qmessagebox line is uncommented the widget appears, title is set but widget disappear once the QMessageBox is closed ...

I fear i might be something linked with the indentation of my code but i feel i tried everything ...

Am I forgetting something ?

What should i do to have the widget remain after this piece of code is executed ?

Or could it be linked with the dynamic loading of the widget's ui ?


texts:
I have a Raster map with East,West, North and south boundaries as follows
(-118.875,-118.75,37.125,
37) and a Vector data. Now I want to plot a vector lines on this raster map. the vector line is in JSON format  : "geometries": [
    {
      "type": "LineString",
      "coordinates": [
        [
          -124.21855199991967,
          40.725894000326946
        ],
        [
          -124.21783999992073,
          40.72704999928351
        ],
        [
          -124.21660999992264,
          40.72901999928047
        ]
      ]
    },

To plot this line on the map, I was confused how to consider the E,W,N,S values. If I take line vector as (x1,y1), (x2,y2)..... and if I can conver E,W,N,S values to (x,Y) format, with this reference I can start plotting. I am coding in C# to achieve the above. How to consider the E,W,N,S values as (X,Y) ?

Is my approach is right if their are other efficient way please do share.


texts:
I am using FME to fill in some small holes within polygons.  I have used DonutHoleExtractor to split these out, however I have found the output via the HOLE port also brings existing 'filled' holes (valid polys that have an outershell).

How can I differentiate these??


texts: Standalone Web Application Development putting GIS cordinates over layers and Visualising it

texts:
I am trying to understand what architecture I can use for my GIS Application. I am from a programming background, however quite inexperienced with GIS and Web development, yet the task is to create a combination of them.

I am trying to create a simple standalone application in the web (Started using PythonAnywhere as it seemed very easy to understand). A set of GIS coordinates are uploaded in a txt format and they are then plotted using MatPlotLib and saved in PDF format. However, I would like to make the UI more intuitive, by having an option to see these coordinates in the web browser itself and to move around or put layers behind it that shows the actual world map and then move around it a bit like how we manipulate them Google Maps.

Could you provide suggestions to this implementation?


texts: How do I import specific fields from a shapefile?

texts: Using arcpy, how do I import only select fields from a shapefile? Or, if I must import the entire shapefile, how do I delete certain fields after import?


texts: Use enterprise architect to export models with metadata

texts: I am using Enterprise Architect to build a data model, and I need to know if there is a way to write somewhere in EA "in notes or glossary....." the description of each feature and class that I have in my diagram (like metadata) so when I export the model as XML to arc  I would find the earlier description that I wrote using EA in ArcCatalog


texts: QGIS 2.01 Fedora 19 vs Windows SQL-Server Access

texts:
i just installed QGIS 2.01 on Ubuntu / Fedora 19 / windows 8, configured ODBC-support for Fedora 19 and created ODBC-DSNs on the windows machine and the Fedora machine to the same SQL-Server 2012 instance where i stored some layers.

When i load these layers in the windows version of QGIS 2.01 everything works as expected, but in Fedora i get a strange behavior:

My point-Layers are loaded and displayed but every polygon layer i open, i get nothing or just the smallest polygons in that layer. No error message at all.

The connection works perfectly with iODBC, unixodbc and freetds.

Since only some small polygons appear: Is this some kind of resource problem?


texts: Converting MultiPatch to Point Shapefile?

texts:
I have a File GeoDatabase that has a Feature Dataset with MultiPatches.  The MultiPatches are tree models.  I would like to convert the MultiPatch into a point layer.

Is there a tool, plugin, or script that will convert the MultiPatch layer into a point feature class?


texts: I want to select multiple values within a field and I am trying the following method: "xyz" in (1,2,3). I am including several values that may not exist in the attributes table, but I would like to cut and paste the formula into multiple shapes to select the values that may or may not exist in the attributes table for each shape. Since I am including values that don't exist the statement is invalid. Is there a way to enter this "blanket" formula for the shapes that would be acceptable?