texts: How to get point x and y geometry from a shapefile?

texts:
I have a shapefile of US Cities and I want to get the X and Y coordinate of EACH city in that shapefile. I've tried this:

for city in city_cursor:
    geom = city.Shape
    point = geom.getPart()
    citylist_City_Name.append(city.CITY_NAME)
    citylist_Country.append(city.CNTRY_NAME)
    citylist_Admin.append(city.ADMIN_NAME)
    citylist_Population.append(city.Population)
    citylist_X_Coor.append(geom.point.X)
    citylist_Y_Coor.append(geom.point.Y)

but I get the error:

Traceback (most recent call last):
  File "C:/Users/workd.py", line 43, in
    citylist_X_Coor.append(geom.point.X)
AttributeError: 'PointGeometry' object has no attribute 'point'

I don't really understand the error message? How can I fix it? Thanks!


texts: Most common unmet GIS need for businesses

texts:
What would you say is the most common unmet GIS need for businesses, as far as technology goes?

I'm a programmer and I'm trying to think of a business idea. I like GIS, so I think I might want to do something related to GIS. In order to identify opportunities, I'd like to know what the most common unmet GIS needs are for businesses.


texts: how to define user role relationship in geoserver?

texts:
I created the user 'Auser' in geoserver and define the roles Aread, Aadmin and Awrite where Awrite is inheriting the rights of Aread. but when I click Auser to pick the roles from service it gives me only default- how will I configure this? Essentially My roleservice.xml looks like this:

<?xml version="1.0" encoding="UTF-8" standalone="no"?>
<roleRegistry xmlns="http://www.geoserver.org/security/roles" version="1.0">
<roleList>
<role id="admin"/>
<role id="read"/>
<role id="write" parentID="read"/>
</roleList>
<userList/>
<groupList/>
</roleRegistry>

instead of the one given under default

<?xml version="1.0" encoding="UTF-8" standalone="no"?>
<roleRegistry xmlns="http://www.geoserver.org/security/roles" version="1.0">
<roleList>
<role id="ADMIN"/>
<role id="GROUP_ADMIN"/>
</roleList>
<userList>
<userRoles username="admin">
<roleRef roleID="ADMIN"/>
</userRoles>
</userList>
<groupList>
<groupRoles groupname="Role_administrator"/>
</groupList>
</roleRegistry>


texts: How to draw polygons in leaflet using x and y coordinates instead of lat lng?

texts:
I am currently working with leaflet. I was able to draw the polygons with geojson using lat,lng with below code snippet:

var geolayer=L.geoJson(geodata);

But above code uses the latlng for drawing the polygons. Now is there any way possible to draw the polygon using x and y coordinate(eg:1,1)?


texts: Selecting polygons which contain majority of area from intersected layer by polygons using QGIS?

texts:
I have two polygons layers (a squares grid (blue in image) and a rectangles grid (green in image)). Because the different polygon sizes, different squares overlap with a single rectangle. I need to find the square with the highest overlapping area for each rectangle, so at the end, I could assign each square to one rectangle. Somehow, I would like to filter just the rectangles with squares assigned after applying this rule. Notice that two squares could be assigned to the same rectangle.

Just to have in mind: the squares don´t have the same area, as in the borders it may be an incomplete square (see image), but still need to know in which rectangle fall the majority of its area.

I am using QGIS.

Any ideas?


texts:
I'm having a problem with the python support of my QGIS version. When I start QGIS, I get the error message:

Traceback (most recent call last):
   File "", line 1, in
ImportError: Module use of python34.dll conflicts with this version of Python.

Then the python version and the QGIS version are displayed, followed by a bulk of path descriptions.

Since the python version displayed in the error message (2.7.5) doesn't match the python version I've installed (3.4), I tried to remove the python34.dll in the system32 folder, but that leads to a similar error, something like:

Traceback (most recent call last):
   File "", line 1, in
ImportError: Couldn't find module

Do you have any ideas how to fix that? I don't want to downgrade my python version. Is there a way to tell QGIS to use a a newer version of python? I'm using windows 7, 64bit.


texts: QGIS Georeferencer returning blank output

texts:
I cannot get QGIS's (2.8.2 Wein/Win 7) georeferencer to work for the life of me. I select points, set the CRS, run the function,

http://i.imgur.com/lsQUOU3.jpg

and the output looks like this.

http://i.imgur.com/pGgoBlF.png

Blank on the map, set to 0 to 0 in the red band and blank in the others in Layer properties (hitting Load changes these values, but apply does not change the screen and the changes are not maintained after I close the properties window), but a garbled thumbnail of the original image still appears. Zoom to layer does not change the display location. Original image in .jpg.

I've seen bug reports for this (or similar) issue stretching back to 1.8 but never any resolutions. Where am I going wrong?

edit:
Andre, using your method produced the same result, unfortunately. I cannot seem get a useable image out of the georeferencer, even with adjusting the stretch values, unloading/reloading, et c.

http://i.imgur.com/VyZCu5T.png

http://i.imgur.com/NC2v239.jpgg


texts: Converting XYZ files (with regularly spaced xy) to ESRI GRID format?

texts:
Normally consultants (LiDAR survey) provide me with data to my specifications; already in ESRI GRID or ASCII GRID format. This way LiDAR data can be easily used with ArcGIS. This time I only have xyz files (i.e., with regular spaced xy coordinates).

Previously, I would use a utility by Min-Lang Huang named GridBatch that would convert XYZ to ASCII GRID then to ESRI GRID files. Unfortunately, the utility is crashing since I upgraded my machine to W7 64bit and ArcGIS to version 10.

Is there another way to batch convert LiDAR surface XYZ data to ESRI GRID files?

(Convert to points -> Interpolate is an option I would like to avoid considering there are 300+ tiles to convert).


texts:
I'm using this code which I got from U-PENN geog 585:

https://www.e-education.psu.edu/geog585/node/718

                  var response = JSON.parse(event.text);
                  if(response.features.length !== 0){
                    var returnedFeature = response.features[0];
                    // Configure the popup
                    map.addPopup(new OpenLayers.Popup.FramedCloud(
                      "marketInfo",
                      map.getLonLatFromPixel(event.xy),
                      null,
                      "" + returnedFeature.properties.NAME + "
" + returnedFeature.properties.ADDRESS,
                      null,
                      true
                    ));

I have the webmap up and all of the layers and basemap are rendering ok but when I click on a feature the cursor just spins and spins with no end. I can still navigate on the map and pan around, zoom in and out, but as soon as I stop navigating the spinning cursor resumes and a pop up never generates. I'm a complete nube to web mapping. Tried looking at a few other posts with pop up issues but couldn't see the same cause.


texts: Using computed min/max image values to set min/max visualization parameters in GEE

texts:
I would like to use the calculated min/max values of an image region to set the min/max visualization parameters of that image.  I can retrieve the actual numbers, but can't use them as variables in the visualization parameters block.  I keep getting the error:

  Shaded: Layer error: Can't encode object: abs() Computes the
  absolute value of the input. Args: this:input (Number): The
  input value.

Here's the code I'm using:

var fc = ee.FeatureCollection('TIGER/2018/States')
    .filter(ee.Filter.and(ee.Filter.eq('NAME', 'Utah')));

var DEM = ee.Image('USGS/NED').clip(fc);
var Terrain = ee.Terrain.products(DEM).select('slope');
print(Terrain);

var visPct = Terrain.reduceRegion({reducer: ee.Reducer.percentile([5,95]).setOutputs(['min','max']),
    geometry: fc,
    scale: 10,
    bestEffort: true
  });

var Min = visPct.getNumber('slope_min');
var Max = visPct.getNumber('slope_max');
print(Min,Max);

var vizParams = {
  min: Min,
  max: Max,
  palette: ['blue','green','Yellow','red']
};

Map.centerObject(fc);
Map.addLayer(Terrain, vizParams, 'Shaded');

The only thing that I can figure out is that GEE is trying to display the image before the min/max values are calculated.


texts: The same labels show up on different points across a layer

texts:
I am making a map of a bicycle tour I took on ArcGIS.

I have created a layer called "Sleep Stops".

Furthermore, I have created 3 label classes called "Name", "Favors", and "NotableFacts"

I have quickly realized that trying to enter a text string for the label class expression results in every single point on the layer being labeled the exact same thing. For instance, every single point on "Sleep Stops" gets labeled "Friendly Crossways Hostel" as its name.

One possible solution is to create 3 new fields for this layer: "name", "favors", and "NotableFacts" and to edit these for each point. Then, I can simply label each point with its value for each of these fields.

How do I create my own text values for fields as I create points?


texts: Getting Image (Error) User memory limit exceeded in Google Earth Engine

texts:
I am coding a simple RF classification algorithm with 2 categories, burned and unburned. My training data are at .geotiff format. 3 purely burned and 6 purely unburned areas, at sentinel-2 10 meters resolution.
I still get this error, and it keeps me back from having the output I want to keep progressing with my thesis .
Does GEE supports ImageCollections as for training data and not just FeatureCollections?
This is my code:
var image = image;
var image2 = image2;
var image3 = image3;

Map.addLayer(image, {bands:['b1', 'b2', 'b3'], min:0, max:1000}, "burned");

Map.addLayer(image2, {bands:['b1', 'b2', 'b3'], min:0, max:1000}, "unburned1");

Map.addLayer(image3, {bands:['b1', 'b2', 'b3'], min:0, max:1000}, "unburned2");

var S2 = ee.ImageCollection('COPERNICUS/S2_SR')
                  .filterDate('2019-11-01', '2020-02-27')
                  .filterBounds(geometry);

var S2 = S2.filterMetadata('CLOUD_COVERAGE_ASSESSMENT', 'less_than', 10);

var S2_image = S2.first();

  var S2_clipped_collection = S2.map(function (S2_image){
  return S2_image.clip(geometry);
});

Map.addLayer(S2_clipped_collection, {bands:['B4', 'B3', 'B2'], min:0, max:3000}, "Sentinel 2 Collection");

//create the training data
var burned = image;
var unburned1 = image2;
var unburned2 = image3;

var image_bands = ee.Image(image.select('b1', 'b2', 'b3', 'b4', 'b5', 'b6','b7','b8','b9'));

var S2_clipped_collection = S2_clipped_collection
  .reduce(ee.Reducer.mean())
  .reproject(image_bands.projection());

var S2_clipped_collection = ee.Image(burned.select('b1', 'b2', 'b3', 'b4', 'b5', 'b6','b7','b8','b9'));

var bands = ['b1', 'b2', 'b3', 'b4', 'b5', 'b6','b7','b8','b9'];

// Make a FeatureCollection from the hand-made geometries.
var polygons = ee.ImageCollection([
  ee.Feature(burned, {'class': 0}),
  ee.Feature(unburned1, {'class': 1}),
  ee.Feature(unburned2, {'class': 1}),
]);

// Get the values for all pixels in each polygon in the training.
var training = S2_clipped_collection.sampleRegions({
  // Get the sample from the polygons FeatureCollection.
  collection: polygons,
  // Keep this list of properties from the polygons.
  properties: ['class'],
  // Set the scale to get Landsat pixels in the polygons.
  scale: 10
});

var classifier= ee.Classifier.smileRandomForest(30).train({
  features: training,
  classProperty: 'class',
  inputProperties: bands
});

// Train the classifier.
var trained = classifier.train(training, 'class', bands);

// Classify the image.
var classified = S2_clipped_collection.classify(trained);
print(classified);

var palette = [
  'D3D3D3', // unburned (0)  // red
  '0000FF', // burned (1)  // green
];

// Display the classification result and the input image.
Map.addLayer(classified, {min: 0, max: 2, palette: palette}, 'classified');


texts: Producing a print layout using 2nd map view in QGIS

texts:
I am trying to produce a map like this completely in QGIS:

I can generate each half of the map in the main program by using View->New Map View but I can't then find a way to make a map in the print composer refer to my 2nd map view (Map 1).

I thought maybe I could just change the map projection in the 2nd map on the composer but then it remains blank (may just need a method for setting the extent?).

I've also seen reference to using map themes and follow map theme but I don't seem to be able to create more than one map theme at a time so that doesn't seem to help either.


texts:
I'm struggling to understand how to manage print layouts using python for QGIS. For example, I am starting with Zero exisiting print layouts in my layout manager. I create a new print layout and add a label to it with this code:

project = QgsProject.instance()
manager = project.layoutManager()
layout = QgsPrintLayout(project)
layout.initializeDefaults()
layoutName = "PrintLayout"
layout.setName(layoutName)
manager.addLayout(layout)

layout.setName(layoutName)
manager.addLayout(layout)

title = QgsLayoutItemLabel(layout)
title.setText("Title Here")
title.setFont(QFont("Arial", 28))
title.adjustSizeToText()
layout.addLayoutItem(title)
title.attemptMove(QgsLayoutPoint(10, 4, QgsUnitTypes.LayoutMillimeters))

This is simple enough, it creates a new print layout called "PrintLayout" and it looks like this:

But if I run this code another time, I start getting problems. If I run exactly the same code, except I change the following line:

title.setText("Title There")

I get the following result. As you can see, "Title There" and "Title Here" are now overlaid on top of each other. This isn't good.

So I started thinking if an existing print layout exists with the same name, maybe I can delete that layout and create a new one. So I try the following code:

project = QgsProject.instance()
manager = project.layoutManager()
layout = QgsPrintLayout(project)
layout.initializeDefaults()
layoutName = "PrintLayout"
layout.setName(layoutName)
manager.addLayout(layout)

layouts_list = manager.printLayouts()

for layout in layouts_list:
    print(layout.name())

for layout in layouts_list:
    if layout.name() == layoutName:
        layouts_list.remove(layout)

for layout in layouts_list:
    print(layout.name())

title = QgsLayoutItemLabel(layout)
title.setText("Title There")
title.setFont(QFont("Arial", 28))
title.adjustSizeToText()
layout.addLayoutItem(title)
title.attemptMove(QgsLayoutPoint(10, 4, QgsUnitTypes.LayoutMillimeters))

Notice, I created a variable "layouts_list" which stores my layouts and for loop to check the names of existing layouts, if the name equals the new layout I want to create, I remove before creating a new one. Also, I print the layout names from the layout manager before and after checking the layout names.

However, it seems the .remove() method doesn't actually remove a layout from the layout manager.

J.Monticolo suggested I try the .removeLayout() method. So I modified my for loop like this:

for layout in layouts_list:
    if layout.name() == layoutName:
        manager.removeLayout(layout)

This does actually remove the layout from the layout manager. But I'm met with this error:

This makes sense. I did delete the layout. But a new layout is not created. They should be as my remaining code makes a new layout and adds a title to it.

Why does this error message stop the code, before creating a new layout?

Here is the final code all together:

project = QgsProject.instance()
manager = project.layoutManager()
layout = QgsPrintLayout(project)
layout.initializeDefaults()
layoutName = "PrintLayout"
layout.setName(layoutName)
manager.addLayout(layout)

layouts_list = manager.printLayouts()

for layout in layouts_list:
    if layout.name() == layoutName:
        manager.removeLayout(layout)

layout.setName(layoutName)
manager.addLayout(layout)

title = QgsLayoutItemLabel(layout)
title.setText("Title There")
title.setFont(QFont("Arial", 28))
title.adjustSizeToText()
layout.addLayoutItem(title)
title.attemptMove(QgsLayoutPoint(10, 4, QgsUnitTypes.LayoutMillimeters))


texts: Raster Layer Statistics Producing NAN Values From Input Raster Layer in QGIS?

texts:
I've just tried to determine the maximum value for a raster layer using Raster Layer Statistics.

However when I use the tool it produces nan values for the statistics, even though the layer clearly has values.

What would be the likely cause of this?


texts: Arcgis API for Javascript: how to change color of simple point marker using colorpicker

texts:
as title states im having a hard time changing the color of my simplemarkersymbol using javascript.

I have read that the format of color needs to be [r, g, b, a] which I have a color picker giving me the correct format but whenever i try to set the value with a variable, or the .val(); it will not work! i have spent far too much time trying to tweak this Ill post my code.

    function toRGBA(color) {

        var colorToChop = color;
        var result;

        var findme = "a";
        if ( colorToChop.indexOf(findme) > -1 ) {

        colorToChop = colorToChop.substring(5, colorToChop.length -1)
        result = "[" + colorToChop + "]"

        console.log(result);
        var tempOpacity = result.substring(result.lastIndexOf(",")+1,result.lastIndexOf("]"));
        var jAlpha = parseInt(tempOpacity * 255,10);

        var color = result.substr(0, result.lastIndexOf(",")+ 1) + jAlpha + result.substr(result.lastIndexOf("]"));
        console.log(color);

        return color;

        } else {
          colorToChop = colorToChop.substring(4, colorToChop.length -1)
          result = "[" + colorToChop + ", 1" + "]"

          var tempOpacity = result.substring(result.lastIndexOf(",")+1,result.lastIndexOf("]"));
          console.log(tempOpacity);
          var jAlpha = parseInt(tempOpacity * 255,10);
          var color = result.substr(0, result.lastIndexOf(",")+ 1) + jAlpha + result.substr(result.lastIndexOf("]"));

          return color;
        }
    }

the above function takes any rgb or rgba color value and will convert it to the 4 digit array like i need... now the following function is where the simplemarkersymbol is defined and set up.

            colorToUse = toRGBA($("#colorPicker").val());
            console.log(colorToUse);
            var outlineColor = toRGBA($("#outlinePicker").val());

            markerSymbol = new SimpleMarkerSymbol({
              "color": colorToUse,
              "size": $("#sizeSpinner").val(),
              "angle": 0,
              "xoffset": 0,
              "yoffset": 0,
              "type": "esriSMS",
              "style": $("#pointIcon").val(),
              "outline": {
                "color": outlineColor,
                "width": $("#widthSpinner").val(),
                "type": "esriSLS",
                "style": "esriSLSSolid"
              }
            });

Now i am using the spectrum color picker to get the colors in the first place.
If i manually enter "color": [255, 0, 0, 255] I will get the color red obviously to display on the map.  My problem is when I pass it as a variable which I need to as it is a dynamically changing marker color.


texts: how to add grid coordinates in QGIS 2.12.0-Lyon in the print composer canvas?

texts:
The attached image is the municipal basmap created using QGIS.

I was not able to add the grid coordinates in the print composer canvas.

How do I add those coordinates?


texts: Pop-up hyperlink to local network in Portal for ArcGIS

texts:
I am working in ArcMap 10.5 and published a feature service to the Portal for ArcGIS. One of the features I have published has a field which contains a hyperlink to a folder on our local network.

Is there a way to click the hyperlink in the pop-up window that will access the folder on our local network?


texts: Relationship between WMS styles and layers parameters for GetMap request in WMS 1.3.0?

texts:
Here is the description which I coped from geoserver

layers  Yes Layers to display on map. Value is a comma-separated list of layer names.
styles  Yes Styles in which layers are to be rendered. Value is a comma-separated list of style names, or empty if default styling is required. Style names may be empty in the list, to use default layer styling.

However, I have questions:

Can 1 layer have multiple styles applied on it for just 1 GetMap request? e.g: layers=world_map&styles=style1,style2,style3 (with style1, style2 and style3 belong to world_map layer)
Can multiple layers have only 1 style (i.e: this style is applied to all input layers or just the layer1)? e.g: layers=layer1,layer2,layer3,...&styles=style1 (with style1 belong to layer1)
Can the number of layers and the number of styles be different? e.g: layers=layer1,layer2&style=style1,style2,style3,style4,style5? (with style1 and style2 belong to layer1, style3, style4 and style 5 belong to layer2).

I think these 3 cases are valid?


texts: Why increasing decimal points not working for XY point features with sufficient precision length?

texts:
I am working with ArcMap 10.3.1 on a Lenovo ThinkPad (Windows 10). I want to display and copy x-y field values with more decimal positions (from 5 to 8; i.e. 34.3456 to 34.34561739), but after instructing and applying the numeric attribute changes under the field Properties there is no evident change in the table.

I am using a 'double' format.


texts: Google Maps (street, physical, satellite) layers stay fixed on zoom/pan in QGIS

texts:
I'm using QGIS 1.8 and Openlayers plugin 1.1.0 to work with Google Maps as basemap (I need to use this old versions for both QGIS and plugin). All was working OK, but from some days Google Maps (street, physical, satellite) loaded with OL plugin stay fixed when I zoom or pan, so my other layers move but map stays fixed. Other maps (Bing, Apple, ...) work perfectly. Any help on this? Is it a googlemaps issue? Is it due to the plugin? I do not how to tackle this issue.

I have been searching for information about it and this seems to be related:
Google Satellite layer does not allow zooming in and out in QGIS
 but in any case I cannot change to QuickMapServices plugin instead as proposed, as it needs QGIS 2.x. I am restricted to QGIs 1.8 and OL plugin 1.1.0.

Note: I have a quasy-migrated version of the same project on QGIS 2.18 and OpenLayers 1.4.1, and it still works, so I do not know if the issue is only with the old version of QGIS, the old version of Openlayers plugin or is a GoogleMaps issue. Old project stop working in August.


texts: how to compare geometries (contains, within) in arcpy?

texts:
I am trying to use 2 geometry objects to check if a line is within a (rather complex) polygon. The geometry check always fails, even if the line is clearly within the polygon. I use the following code:

#get the polygon
sC = arcpy.SearchCursor(polygon)
for row in sC:
    polygonshape = row.Shape
del sC

#creation of the polyline(s):

feature_info = [[[p1.x, p1.y], [x, y]], [[p2.x, p2.y], [x, y]], [[p3.x, p3.y], [x, y]]]
features = []
for feature in feature_info:
    array = arcpy.Array()
    for point in feature:
        P = arcpy.Point()
        P.X = point[0]
        P.Y = point[1]
        array.add(P)
    features.append(arcpy.Polyline(array))

#the actual check:
for item in features:
    in_poly = polygonshape.contains(item) #this is always false

Any ideas why this fails?
I considered, that one of the layers (polygon) has a coordinate system, while the other doesn't but according to another answer here (which I can't find right now) coordinate systems are supposedly ignored for this kind of checks.
Could the size/complexity of the Polygon be an issue?

If I create a Feature Class of the line and do a spatial selection (within) afterwards, it works as expected.
Also: this behaviour returns unexpected results for contains as well as within.


texts: Why is my GeoJSON polygon only displayed as a line?

texts:
I'm trying to convert a simple shapefile drawn in GIS to GeoJSON. However, it's not working and I don't know why.

The file converts with QGIS (Save As), but the http://geojsonlint.com/ doesn't plot the shape, nor does my Mapbox application - it plots it as a line.

I've set the CRS to EPSG 3857.

Here's the GeoJSON, http://pastebin.com/93e0wWzA

Do you need the Shapefile to figure this out?


texts: Using supervised vs unsupervised classification in identification of region of interest?

texts:
I was introduced to machine learning and remote sensing recently.
My task was to classify the satellite images into vegetation and non vegetation.
We were introduced to two approaches.

Supervised learning - where we had wkt or geojson files made from ground truth. These files had polygons which were used to train the model.
satellite images from WorldView-3 Satellite Sensor
Unsupervised classification where the pixels were classified based on NDVI values using clustering models such as K-means, Fuzzy C-means clustering.
satellite images from landsat 8

While all of these things were virtually spoon fed and I took the code samples from here and there. I still fail to understand which method is used where, specifically with context of crop forecasting.

What is the advantage of collecting the ground truth, when we can use the unsupervised learning to classify the images?

If it is about accuracy, then are there any specific examples as to how ground truth helps in accuracy in crop forecasting?


texts: QGIS project with PostGIS layers freezes irregularly for minutes

texts:
I'm using about 30 more or less complex PostGIS layer in my QGIS project. Some of them are PlanetOSM layer generated with osm2pgsql.

The Project loads fast (ca 5 sec). The canvas update by moving or zooming is also fast (2 sec per update). My computer have got enough resources.

Overall, everything seems okay, but sometimes after doing or saving changes on a layer, switching between edit mode or by adding or removing a layer, QGIS freezes for a moment (up to 2 minutes). While hanging, one of my CPU cores is at 100%. It seems, that it does not matter, how big the layer is. Sometimes also adding a empty Postgis layer causes that freezing.

EDIT: Sometimes changes are not saved or the DB does not update the changes.

After restarting the Project, everything works fine... for the first edits.

What could be the problem?
How can I monitor whats going on in my database?

System configuration: Xubuntu 14.04, QGIS 2.6, PostgreSQL 9.3.5, PostGIS 2.1.2
Intel i5-4570 CPU @ 3.20GHz, 8GB RAM


texts: Uploading lots of raster files/georeferenced maps to ArcGIS Online?

texts:
The end goal is to build a web application where people can explore the historical maps I've georeferenced - among other things. ArcGIS online is the ideal venue for everything else I want, but for whatever reason it seems I can't easily upload raster files as in ArcMap. And I have approx. 12 rasters I want to upload (about the size of Bermuda).

There seems to be some confusing work arounds, but my concern is that a lot of credits will be used. Also, even more importantly, it seems I have to publish my rasters as layer files to ArcGIS online as a work around, but I feel this makes my work public to use by other users. I want only those in my organization to be able to access and use my georeferenced work in their maps. The work I've performed isn't mine to share, it belongs to my educational institution.

What's the best route for me to go? Is there a way to publish these layers to ArcGIS Online in a way that makes the layers accessible only to my organization?


texts:
I'm trying to install ArcGIS WEB ADF. How can I install this tool?

My project contains ESRI.ArcGIS.ADF.Web.UI.WebControls, Version=10.0.0.0 but I cannot run project because I must install ArcGIS WEB ADF.

I found the following page, but I don't understand how can get administrative privileges.

http://help.arcgis.com/en/arcgisserver/10.0/install_guide/arcgis_server_net_install_guide/index.html#//008900000010000000.htm


texts:
I've been trying to convert a map from ArcMap to KML without exit. The issue is that I want to export the complete map, not only the layers in the data frame but also legend, symbols, grid ticks and so on but the Map to KML tool does not have an option (as far as I know) to keep them.

I must do it in ArcMap so does anyone know how to do it?


texts: ArcGIS 10.3 Dissolve tool splitting up lines, not dissolving them

texts:
I have noticed strange behaviour with the Dissolve tool so tried a test on three line segments, snapped to the same central point. The line segments each consist of exactly two vertices i.e. the end vertices. There are no vertices in between.

When I ran the dissolve tool on this (multipart deselected, for testing purposes), I got a new feature class consisting of five records: two of the original line segments had been split in two! I was expecting to get a new feature class with a single record, as all the lines shared a vertex. At the very least, I didn't expect to get more records in the output feature class than were in the input feature class.

Anybody have any ideas as to why it might be behaving this way?


texts: ModelBuilder layer name output is not unique

texts: After using inline variable substitution in ArcGIS Pro ModelBuilder it generated the outputs that I wanted, but their layer names in the output were not unique. In the metadata tab, I can see that the feature class is named as I wanted it, but it's not displayed as the layer name. Is it possible to make corrections to the model so it would rename them correctly or create a script that would give the layer name the feature class name?


texts:
I am trying to write a series of raster files into a new folder (after doing some image processing). However, when using my code, the files are but are not in the folder that I have created for that purpose.

I say they are but they are not because when I run my code and take a look inside the folder that I have created, a file pops-up, disappears, pops-up, disappears and this goes until I hit stop the R console.

Here is my code

Img.Raster <- list()
for (k in 1:length(list.pics)){
  Img.Raster[[k]] <- raster(list.pics[k])
  extent(Img.Raster[[k]])=c(coordinate.one.Image$CenterLeftx[[k]], coordinate.one.Image$CenterRightx[[k]],
                            coordinate.one.Image$CenterBottomy[[k]], coordinate.one.Image$CenterTopy[[k]])
  proj4string(Img.Raster[[k]]) <- "+proj=utm +zone=15 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
  writeRaster(Img.Raster[[k]],filename =paste0("/Users/Juan/Desktop/Images_Folder/"
                                                , "Img_",names(list.pics[[k]])), format="GTiff",overwrite=TRUE)
}

list.Image is a variable that contains 1745 images

Those anyone know what could be wrong?

Suggestion by @Jul

for (k in 1:length(list.pics)){
  Img.Raster[[k]] <- raster(list.pics[k])
  extent(Img.Raster[[k]])=c(coordinate.one.Image$CenterLeftx[[k]], coordinate.one.Image$CenterRightx[[k]],
                            coordinate.one.Image$CenterBottomy[[k]], coordinate.one.Image$CenterTopy[[k]])
  proj4string(Img.Raster[[k]]) <- "+proj=utm +zone=15 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
  writeRaster(Img.Raster[[k]],filename =paste0("/Users/Juan/Desktop/Images_Folder/"
                                                , "Img_",stringr::str_pad(string = k,width = 4,side = "left",pad=0),
                                               format="GTiff",overwrite=T)
}


texts: Error Using PostGIS Address Standardizer

texts:
I am trying to use the PostgreSQL 9.5 and PostGIS address standardizer to parse an address field. The particular field in the table I am working with only has the street name and the street type (ex. 'MAIN ST'). Ultimately what I would like to do is create an insert statement from the existing table into a new, empty table, but for now I will just settle for getting this to work.

I have created the 'address_standardizer_data_us' extension. But when I try to use the information directly from the Postgre documentation:

SELECT house_num, name, suftype, city, country, state, unit  FROM standardize_address('us_lex',
           'us_gaz', 'us_rules', 'One Devonshire Place, PH 301, Boston, MA 02109');

I get this error:

ERROR:  function standardize_address(unknown, unknown, unknown, unknown) does not exist
LINE 1: ..., name, suftype, city, country, state, unit  FROM standardiz...
                                                             ^
HINT:  No function matches the given name and argument types. You might need to add explicit type casts.
********** Error **********

ERROR: function standardize_address(unknown, unknown, unknown, unknown) does not exist
SQL state: 42883
Hint: No function matches the given name and argument types. You might need to add explicit type casts.
Character: 67

What am I missing here? Can I substitute the string 'One Devonshire Place, PH 301, Boston, MA 02109' with a field name in an insert statment?


texts: Warping raster to new domain/resolution/projection causes data loss

texts:
I am attempting to warp a GTiff from a lat/lon projection to a higher (~2x) resolution Polar Stereographic grid.

In the output GTiff any cells north of 64 deg lat become NaN instead of their actual values.

I am using gdalwarp thusly:

  gdalwarp -t_srs 'wkt string here' -tr 5000 5000 -te -4252339.5,-5302590.0,4597660.5,2747410.25 in.tif out.tif

What could be causing the output raster to have the blank area?

Input extents are [-84.84, -50.19, 55.07, 79.96]. [minLon,maxLon,minLat,maxLat]

Output Coordinate system:

proj4 (+proj=stere +lat_0=90 +lat_ts=60 +lon_0=-100 +k=1 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs) .

It's standard Polar Stereographic but with Central Meridian at -100 deg and latitude of origin at 60 deg.

The issue was caused by the large difference in coordinate systems. (The domain was self-intersecting in one of the coord systems). It was largely solved by explicitly entering the output extent instead of having GDAL attempt to guess it.

i.e -te -180 0 180 90 was added to the gdalwarp line.


texts:
I am looking for the land use data of London in vector format, preferably in polygons. I am interested in whether a building / an area is residential, office or retail use, like this one (https://www.researchgate.net/figure/The-spatial-distribution-of-three-key-land-use-categories-in-Inner-London_fig5_301482105).

But I tried to download several datasets including Strategi® which includes developed land use areas. But unfortunately, no land use attributes are available in the shapefile. I could only find the data from UMMAP product from Verisk: http://www.geoinformationgroup.co.uk/ukmap But the data is not an open data.

Where can I find the open free land use map for London and the UK?


texts: How to change position of CartoDB legend in Editor?

texts:
I have a map here that I'm building in Editor with some custom CSS.
I'm having trouble positioning my legend in the place I want it on the screen.
I referred to this post, and tried following the instructions from iriberri, but my solution isn't working.

My goal is to get the legend on the lower left of the screen.

Here's my code:

and here's the result. The CSS I wrote in Editor is not overwriting correctly:


texts:
I have a raster and I'd like to underlie it with a basemap in tmap. In the view-mode this kind of happens automatically. But when I set it to plot-mode I really struggle to get a basemap. Anyone maybe has a pointer on where to look?

I thought about something like this:

tm_shape(raster1) + tm_raster("col1.VALUE", palette = "viridis") +
tm_basemap("Stamen.Watercolor") + tm_tiles("Stamen.TonerLabels")


texts: Converting 4 Band RGBA to 2 Band LA (Gray + Alpha) mode while keeping georeferenced data in GeoTiff

texts:
I created some GOES-16 images where the dark/night side of the visible image is transparent (alpha channel). I used satpy's DayNight Compositor to make the night side of the image transparent, while keeping the day time in gray scale. Unfortunately, using this method has converted the original L (gray) image into an RGBA image.
I tried using PIL to convert the image from RGBA mode to LA mode but it removes the georeferenced data in the GeoTiff.
img = Image.open(filename_of_RGBA).convert('LA')
img.save(new_filename)

How can I keep the georeferenced data while converting from RGBA to LA?


texts:
In PostGIS I have a large database of linestrings which I have obtained using ST_makeline function.

I would like to cluster similar lines based on where the lines start and end.

Here is an example of what I would like to obtain:

From this post Grouping lines together to form one line after a distance threshold is met I know that in ArcGIS Desktop there is a tool that identifies the mean direction, length, and geographic center for a set of lines.

Do you know if there is an equivalent tool in QGIS?

Do you think it would be easier to cluster start and end coordinates instead of lines?

Basing on this post How to merge connected lines with same direction (PostGis), I wrote the following query that seems to work although it is very slow:

select a.routeid as route1, b.routeid as route2,a.the_geom as geom1, b.the_geom as geom2 from prova a, prova b
WHERE (ST_DWithin(st_startpoint(a.the_geom),st_startpoint(b.the_geom),0.1) and ST_DWithin(st_endpoint(a.the_geom),st_endpoint(b.the_geom),0.1) and (a.routeid <> b.routeid))

To prevent a line from being included in more than one separate cluster, I used distinct on as follows where line1 and line2 are the same tables as prova:

 SELECT * FROM(
 SELECT DISTINCT ON (line2.routeid) line1.routeid as route1, line2.routeid as route2,line1.the_geom, line2.the_geom
 FROM line2 INNER JOIN line1 ON
 ST_DWithin(st_startpoint(line1.the_geom),st_startpoint(line2.the_geom),0.1)
 AND ST_DWithin(st_endpoint(line1.the_geom),  st_endpoint(line2.the_geom),  0.1)
 ORDER BY line2.routeid) as clusters
 ORDER BY route1;


texts: How to make line size scale dependent in OpenLayers webmap made using qgis2web?

texts:
I've made a map in QGIS v3.2 with a vector layer named "Roads." In QGIS, the line widths representing roads are set to "meters to scale." When viewed in QGIS, this behaves exactly as I want. When I zoom in, the roads appear wider, and when I zoom out, they appear more narrow.

However, after exporting this as an OpenLayers webmap, the opposite is true. The roads appear widest when zoomed out, and become more narrow as I zoom in. I've played around with the line width in the JS files and have made them more narrow, but still not scale dependent. What I really want is to make it scale dependent, just like when viewed in QGIS.

Here's a sample of the code I'm working with. On my last attempt at fixing this, I set the width to 'auto.' This just made the lines consistently narrow, same as if I set it to '1'.

 var style = [ new ol.style.Style({
    stroke: new ol.style.Stroke({
        color: 'rgba(237,206,174,1.0)',
        lineDash: null,
        lineCap: 'round',
        lineJoin: 'bevel',
        width: 'auto',
    }),
    text: createTextStyle(feature, resolution, labelText, labelFont,
                          labelFill, placement)
})];

return style;

So, can somebody help me correct this or direct me to a possible solution? I haven't been able to find something that works yet. FWIW, I'm new to Javascript and am learning as I go here...


texts: Intersecting great circle arcs and polygons

texts:
I am trying to find the point of intersection between a great cirle arc and a polygon.  I am working in long/lat, where I am given the endpoints of the great circle arc, A and B, and the vertices, C, D, E and F of the polygon in long/lat.  The point A is inside the polygon and point B is outside of the polygon.  It is obvious by construction that the great circle arc intersects the polygon on the CD side of the polygon.  When I compute the intersection of the great circle arc AB, represented by a LINESTRING, with the POLYGON CDEF, I get a LINESTRING from a point on CD to A;  when I compute the intersection of the great circle arc AB with the LINESTRING CD, I get a different point on CD.  Here are the details:

A:  -122.3 47.4
B:  -84.4 33.6

C: -107 49
D: -107 31
E: -130 31
F: -130 49

(1) intersect LINESTRING and POLYNOMIAL

select st_astext(st_intersection(st_geographyfromtext('linestring(-84.4 33.6, -122.3 47.4)'), st_geographyfromtext('polygon((-107 49, -107 31, -130 31, 130 49, -107 49))')));

====>
LINESTRING(-107 42.1770201452269, -122.3 47.3999999999286

(2) intersect LINESTRING and LINESTRING

select st_astext(st_intersection(st_geographyfromtext('linestring(-84.4 33.6, -122.3 47.4)'), st_geographyfromtext('LINESTRING(-107 49, -107 31)')));

====>
POINT(-107.040211660479 43.9998296744469)

As you can see, this point is not the same as the endpoing in the output of (1).

What is going on?

I am running PostgreSQL9.3.4 on x86_64-unknown-linux-gnu, compiled by gcc (Ubuntu 4.8.2-16ubuntu6) 4.8.2, 64-bit POSTGIS="2.1.2 r12389" GEOS="3.4.2-CAPI-1.8.2 r3921" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.10.1, released 2013/08/26" LIBXML="2.9.1" LIBJSON="UNKNOWN" RASTER

Finally, a colleague of mine computed the intersection of the two great cirlce arcs directly calling the Java version of GDAL/OGR and got -107 41.8290237496.


texts: Filling shape with labels or attributes in QGIS?

texts:
In QGIS, I would like to set the fill of a shape to be a pattern repetition of label values (i.e. an attribute of the item).

Is this possible, looking a little like the sketch below (for shapes with ID = 1,2,3)?


texts: WMS in Google Chrome (ArcGIS API JavaScript) not working

texts:
I am doing a small web with WMS, but I have a problem with Google Chrome (and IE... but with this one I don't mind). Down of this lines there is two website:

Galicia: This one contains WMS from Galicia (north of Spain) and it loads correctly on Firefox and Chrome. Right

http://signo-geo.com/alejandrowms/mapacorrecto2.html

Andalucia (south of Spain): In this second App, the WMS on Chrome is not loading but in Firefox it works fine.

http://signo-geo.com/alejandrowms/mapaprueba2.html

I was looking for information about that and only find this thread.

http:// forums.arcgis.com/threads/98149-Chrome-WMS-and-content-type-application-vnd.ogc.wms_xml-fails..

Here is an Image which shows the problem:

Why is the WMS not appearing in the second App, even though the code is similar to the first one? And why does this problem occur only in Google Chrome and not Mozilla Firefox?


texts:
I'm running a GeoServer 2.17.1 that is sitting behind a proxy that handles the TLS termination and send the HTTP traffic to GeoServer.
The issue is that when I'm accessing trough HTTPS to the web admin interface, I can't make any change that contains a form as I get an HTTP 403 response. But doing the same trough plain HTTP the requests return OK
My first thought was that it had to be an issue with the Origin beeing now https://mydomain.com instead of http://mydomain.com and also found this link in the documentation about CSRF Protection and issues with proxies.
Adding the environment variable GEOSERVER_CSRF_WHITELIST that it mentions, solved this issue.
But now I wanted to add a CORS policy, so I added
<filter>
    <filter-name>CorsFilter</filter-name>
    <filter-class>org.apache.catalina.filters.CorsFilter</filter-class>
    <init-param>
        <param-name>cors.allowed.origins</param-name>
        <param-value>*</param-value>
    </init-param>
    <init-param>
        <param-name>cors.allowed.methods</param-name>
        <param-value>GET,OPTIONS</param-value>
    </init-param>
    <init-param>
        <param-name>cors.allowed.headers</param-name>
        <param-value>Origin,Access-Control-Request-Method,Content-Type,X-Requested-With,Accept,Access-Control-Request-Headers</param-value>
    </init-param>
</filter>

<filter-mapping>
    <filter-name>CorsFilter</filter-name>
    <url-pattern>/*</url-pattern>
</filter-mapping>

To /usr/local/tomcat/conf/web.xml. The cors filter worked, but once added the issue with the forms and receiving 403 returned.
I've been testing this but it seems to that these two configurations are conflicting with each other. AS they work separately but not in tandem
I've also tried with other GeoServer version (2.16.4) but the issue remains
Anyone can help me understanding why this can be happening? I would like to enable HTTP -> HTTPS redirection but this issue is preventing me from doing that

Edit:
Adding logs:
logs for configuration with CORS in /usr/local/tomcat/conf/web.xml
and the env variable
logs for config with only the env variable
And the tomcat access log that show the 403 when trying to a add a new workspace
[...]
172.20.0.5 - - [06/Aug/2020:08:12:19 +0000] "GET /geoserver/web/wicket/bookmarkable/org.geoserver.web.data.workspace.WorkspaceNewPage HTTP/1.0" 302 -
172.20.0.5 - - [06/Aug/2020:08:12:19 +0000] "GET /geoserver/web/wicket/bookmarkable/org.geoserver.web.data.workspace.WorkspaceNewPage?6 HTTP/1.0" 200 4970
172.20.0.5 - - [06/Aug/2020:08:12:20 +0000] "POST /geoserver/web/wicket/bookmarkable/org.geoserver.web.data.workspace.WorkspaceNewPage?6-1.IFormSubmitListener-form HTTP/1.0" 403 -
172.20.0.5 - - [06/Aug/2020:08:12:22 +0000] "GET /geoserver/web/wicket/bookmarkable/org.geoserver.web.data.workspace.WorkspaceNewPage?6 HTTP/1.0" 200 4970

To summarize, I want to enable CORS and to use the web admin interface trough HTTPS, for this I tried different configs:

Only setting the env variable GEOSERVER_CSRF_WHITELIST with my domain (the web admin works)
Only adding the Cors config in tomcat root conf /usr/local/tomcat/conf/web.xml (Works and send the Allow-Origins headers)
Set the env variable and the cors config (Only the cors works, but the web admin don't and receive 403)
Cors config in the webapp geoserver dir /usr/local/tomcat/webapps/geoserver/WEB-INF/web.xml (doesn't work and can't log in as I get a 403)


texts:
Is it possible to make something like this using python ?

How should i best structure my data in order to make this map and what should i use to draw it. Code is not necessarily required. I just need some advice at this point.


texts: Calculating damage from a flood at the village level in percentage

texts:
I am using ArcGIS for my analysis. I want to calculate the damage from a flood at the village level in terms of percentage area damaged.

I have done the following steps.

My village file is in point shape ad I have overlayed it with a district-level polygon in order to identify it (by using "Join Relates" option, the resulting output is in point shape).
I have flood polygon, I have used the output at step 1 for intersection with flood polygon, I get the results (the resulting output is point shape).

Now I want the guidance for area calculation given output file is point shape. How can I do that?


texts: Getting coordinates on click and pass it to dialog using PyQGIS?

texts:
I am creating a plugin in which i need to get coordinates by clicking on QGIS canvas so i tried the solution mentionned at Getting coordinates by clicking on QGIS Canvas with PyQGIS? but it didn't work for me , but this one Programatically check for mouse click in PyQGIS? seems working well now i want to pass the coordinates to GUI LineEdit, using QGIS 2.8.6, Python 2.7 and PyQT4 but it's not working so here's my code :

from PyQt4.QtGui import *
from PyQt4.QtCore import *
from qgis.core import *
from qgis import *
from PyQt4.QtCore import QSettings, QTranslator, qVersion, QCoreApplication,QObject
from PyQt4.QtGui import QAction, QIcon

import os.path
from qgis.core import QgsApplication,QgsMapLayerRegistry
from qgis.gui import QgsMapCanvas,QgsMapToolEmitPoint,QgsMapTool

# Initialize Qt resources from file resources.py
import resources
# Import the code for the dialog
from Proejt_dialog import ProjetDialog
from Proejt_dialog import ProjetDialog as ProjetDialog1
import os.path
class Projet:
   def __init__(self, iface):
        """Constructor.

        :param iface: An interface instance that will be passed to this class
            which provides the hook by which you can manipulate the QGIS
            application at run time.
        :type iface: QgisInterface
        """
        # Save reference to the QGIS interface
        self.h=0
        self.iface = iface

        # initialize plugin directory
        self.plugin_dir = os.path.dirname(__file__)
        # initialize locale
        locale = QSettings().value('locale/userLocale')[0:2]
        locale_path = os.path.join(
            self.plugin_dir,
            'i18n',
            'Projet_{}.qm'.format(locale))

        if os.path.exists(locale_path):
            self.translator = QTranslator()
            self.translator.load(locale_path)

            if qVersion() > '4.3.3':
                QCoreApplication.installTranslator(self.translator)

        # Declare instance attributes
        self.actions = []
        self.menu = self.tr(u'&Projet')
        # TODO: We are going to let the user set this up in a future iteration
        self.toolbar = self.iface.addToolBar(u'Projet')
        self.toolbar.setObjectName(u'Projet')
   def functions(self):
        do stuffs
   def run(self):
        """Run method that performs all the real work"""
        tool = PointTool(self.iface.mapCanvas())
        self.iface.mapCanvas().setMapTool(tool)
        # show the dialog
        self.dlg.show()
        # Run the dialog event loop
        result = self.dlg.exec_()
        # See if OK was pressed
        if result:
            # Do something useful here - delete the line containing pass and
            # substitute with your code.
            pass

class PointTool(QgsMapTool):
    def __init__(self, canvas):
        QgsMapTool.__init__(self, canvas)
        self.canvas = canvas
        self.ui = ProjetDialog1()

    def canvasPressEvent(self, event):
        pass

    def canvasMoveEvent(self, event):
        x = event.pos().x()
        y = event.pos().y()

        point = self.canvas.getCoordinateTransform().toMapCoordinates(x, y)
        self.ui.X.setText(str(point[0]))
        self.ui.Y.setText(str(point[1]))
    def canvasReleaseEvent(self, event):
        #Get the click
        x = event.pos().x()
        y = event.pos().y()

        point = self.canvas.getCoordinateTransform().toMapCoordinates(x, y)

        P=Projet(self) #error

        P.dlg.X.setText(str(point[0]))
        P.dlg.Y.setText(str(point[1]))

    def activate(self):
        pass

    def deactivate(self):
        pass

    def isZoomTool(self):
        return False

    def isTransient(self):
        return False

    def isEditTool(self):
        return True

any help please, i'm still new to PyQGIS

Update :
i couldn't link the class 'PointTool' to the dialog initialised in the 'Projet' class, so i added a parent parameter which contain the dialog to the 'PointTool' then call it in run function from Projet Class


texts: qgis2web incorrect leaflet popup position

texts: I have encountered some problems with positions of popup containing image in leaflet format (version) of this great tool. Sometimes has the pointer of popup (in a picture marked with red circle) incorrect (to the right moved) position. Is it easy to fix it? In an extreme case could help remove pointer from popup. Is it possible? qgis2web 2.34.0, QGIS 2.14.2-2


texts:
I am trying to add a WebTileLayer to my ESRI map and the code runs without any problems. I can see the protobuf files being received by looking at the network inspector, however they do not appear on the map.

<script>
require([
  "esri/Map",
  "esri/views/MapView",
  "esri/layers/WebTileLayer",
  "esri/layers/support/TileInfo",
  "dojo/domReady!"
], function(Map, MapView, WebTileLayer, TileInfo) {
  var map = new Map({
    basemap: "streets"
  });

  var view = new MapView({
    container: "viewDiv",
    map: map,
    zoom: 18,
    center: [-89.37318187952042, 45.61155180826724]
  });
  var vtlayer = new WebTileLayer("http://localhost:3000/featurelayers/{level}/{col}/{row}.pbf");

  map.add(vtlayer);
  vtlayer.visible = true;
});
</script>

I have tried styling the layer hoping that they might just be hidden but that didnt seem to help. I am using ArcGIS API for JavaScript 4.


texts: LS7 filling the gaps image with google earth engine

texts:
I would like to work with Landsat 7 images, but most of them have gaps... I would like to fill the Gaps by aggregate images bands by month. So I'm looking for a way to deal with this problem.
I have tried to create a composite image with

var collectionLS7 = ee.ImageCollection('LANDSAT/LE7_L1T_TOA') //10m
    // Select the bands of interest to avoid taking up memory.
    .select(bandNames)
    // Filter to get only imagery at a point of interest.
    .filterBounds(roi);
var composite = collectionLS7.median();

But it doesn't work

I have also tried

var customComposite = ee.Algorithms.Landsat.simpleComposite({
  collection: collectionLS7,
  percentile: 75,
  cloudScoreRange: 5,
  maxDepth: 4
});

But with an error Layer 1: Layer error: ImageCollection.reduce: Error in map(ID=LE72040492015361ASN01):
Landsat.TOA: Band 'B1' is Type<Float>, expected Type<Short<0, 255>>.

Any ideas are welcome !


texts:
We would like to use sDNA to calculate accessibility measures, in addition to the other measures sDNA outputs.

I've been going over the documentation a fair bit, but cannot find whether one of the outputs equals an accessibility measure like

Where the accessibility A of location i is the sum of the destinations j with destination weight d, accessible from i multiplied by the costs of reaching destination j from origin i. The cost function is often referred to as the impedance function or distance decay function.

In its simplest form, the cost function

Would equate to a rectangular function like

Of course, other cost functions are imaginable. For instance, a commonly used distance decay function is the exponential function.

Using the example network below, where every link is 100 meters and the link weight is depicted as a label, the accessibility within 100 meters from A is 30 (6 links accessible with weight 5), the accessibility from B is 50 (1 link accessible with weight 50) and from C is 40 (2 links accessible with weight 20) .

Which sDNA measure would equate such an accessibility measure?


texts: ogr2ogr from file to PostgreSQL 12 table fails

texts:
Using CentOS 8.0 with PostgreSQL 12, PostGIS 3, and GDAL/OGR 3.0.1.
I want to append a file to a database table.

Previously (before upgrading to the new system) I used

ogr2ogr -f "PostgreSQL" PG:"dbname=dev user=admin password=123 host=0.0.0.0 port=5432" -lco "SCHEMA=myschema" -lco "GEOMETRY_NAME=geog" -append -sql "select ID as gid from input" input.gen

However, that fails now with

  ERROR 1: ERROR:  column s.consrc does not exist
  LINE 1: ...nrelid = c.oid AND a.attnum = ANY (s.conkey) AND (s.consrc L...
                                                               ^

  HINT:  Perhaps you meant to reference the column "s.conkey" or the column "s.c                                                                                onbin".

  ERROR 1: ERROR:  column s.consrc does not exist
  LINE 1: ...nrelid = c.oid AND a.attnum = ANY (s.conkey) AND (s.consrc L...
                                                               ^

  HINT:  Perhaps you meant to reference the column "s.conkey" or the column "s.c                                                                                onbin".

  ERROR 1: Layer myschema.input already exists, CreateLayer failed.
  Use the layer creation option OVERWRITE=YES to replace it.

  ERROR 1: Terminating translation prematurely after failed
  translation from sql statement.

If I make a PGDUMP with

ogr2ogr -f "PGDUMP" -lco "SCHEMA=myschema" -lco "GEOMETRY_NAME=geog" -append -sql "select ID as gid from input" test.dump input.gen

I can run that file in psql and import the data, but I would like to avoid that extra step of writing a large pgdump file.

Notes

It also fails without -append
It fails with other input data types
and files as well (Shapefile, ArcGEN, and GeoJSON)
Columns mentioned in the Error message are not part of my input data.
I suspect that it has something to do with PostgreSQL 12 dropping OID


texts: ECW Exception when exporting a map in QGIS

texts:
Since monday Qgis has been unstable, everytime I try to export a map (or worse an atlas) it has crashed.

I get the folowwing message :

---------------------------
ECW Exception
---------------------------
(thread) ECW Version : 5,3,0,862
 Unknown Error:
An error has occurred in this application.  If you receive this message,
please report it to Hexagon Geospatial Support to assist with fixing this error.

An Access Violation Exception occurred at code address 0x1805f79df attempting to
WRITE to memory address 0x18.

Dump: 488942f8488bc24883c4205bc3cccccccc40534883ec204533c033d233c9ff15

Stack Trace:

---------------------------
OK
---------------------------

I tried to contact Hexagon Geospatial without much succes. I guessed the problem was  caused by an ecw file, but most the maps didn't have one.


texts: Using vegetation cover and DEM elevation data together

texts:
I'm using Terrain data (DEM) from USGS and I would like to build a map in which I can keep the information on the elevation in 3D and on the vegetation cover. I have obtained the 3D raster for the elevation and I have the shapefile for the vegetation cover but I do not know how to mix them together.

Do you have any suggestion?


texts: Add file name as attribute field in QGIS

texts:
I have a folder with several shapefiles (around 75). I would like to create a new field into each .shp file which has to be filename. i.e: 20201003.shp --> 20200103.shp (New field 'Fecha' = 20200103)
I have tried this but nothing is done and no errors neither:
import os, datetime

`shapefolder = r'C:\Users\juavafer\Proyectos\Costas_CoastSat\Output\CMillor\Comparacion\matlab\Qgis'
newfieldname = 'test3'

for root, folder, files in os.walk(shapefolder):
    #Find the shapefile
    for file in files:
        fullpath = os.path.join(root,file)
        if os.path.isfile(fullpath) and fullpath.endswith('.shp'):
            print(fullpath)
            #Extract numbers from filename
            datename = ''.join([num for num in file if num.isdigit()])
            #Create a datetime object and convert to string in format 'YYYY-MM-DD', which is what a shapefile wants
            QgsField(name='test3', type=QVariant.String, len=20)
            newdatename = dt.strftime('%Y%m%d')

        print(datename)
        # dt = datetime.datetime.strptime(datename, '%Y%m%d')
        # newdatename = dt.strftime('%Y%m%d')
        print(datename)
        #Add the field https://opensourceoptions.com/blog/pyqgis-adding-and-deleting-vector-layer-fields/
        layer = QgsVectorLayer(fullpath)
        layer.dataProvider().addAttributes([QgsField(newfieldname, QVariant.Date)])
        layer.updateFields()

        #Update each feature with date https://gis.stackexchange.com/questions/381174/setting-all-the-values-in-a-attribute-field-to-the-same-value-with-pyqgis)
        attr_map = {}
        field_index = layer.fields().indexFromName(newfieldname)
        for f in layer.getFeatures():
            attr_map[f.id()] = {field_index: newdatename}
        layer.dataProvider().changeAttributeValues(attr_map)`


texts: What server and API to use for displaying complex polygons with dynamic symbology?

texts:
Background

I need to interface a map into an existing website which allows the user to query a database and which currently results in a table listing population statistics per county.

For example, display the population under 15, for 2003 - the result is a table listing each county's population values.

I need to display a map showing county polygons which are symbolised with those values. If the user modifies the query, the map needs to update (eg, if the year is changed to 2004, or the statistic is changed to population over 65).

Performance will be a factor, so the map should update rapidly when the query is run/modified. The map should support mouse-over events on the counties, so that statistics can be displayed in a popup as the user's mouse moves.

My questions

What are the options for the GIS server (to serve the polygons) and the mapping API (to embed the map on the website)?

Options

I like the look of Leaflet, but it seems that the county polygons would need to be supplied as coordinates/rings. Would this be too slow when zoomed to national scale? Can Leaflet handle on-the-fly simplification (as with feature layers in ArcGIS Server)?

Using OpenLayers with GeoServer could work too, but the performance is quite slow for hover events, and it appears that the symbology is hard-coded in the WFS service(?) I'd need the ability to update the symbology when the query was modified.

ArcGIS Server 10.1 has advanced on-the-fly symbology options, but hasn't been released yet. The client will be wary of deploying beta software on a production site.

Thanks for any advice.


texts: I have some point data in PostgreSQL / PostGIS. How do I export the data to a GPX file? I looked for commands similar to ST_AsGML and ST_AsKML but could not find anything to convert to GPX (like a ST_AsGPX).


texts: Create new table with a subset of fields from the original table in arcgis desktop?

texts:
I am using ArcGIS Desktop/ArcMap. I have a table with over 500 columns. I want to create a table off of that table that contains only 3 out of the 500 columns. I see a "Make Table View" tool in the toolbox, but I have to manually uncheck all of the fields that I don't want visible. There appears to be no option to unselect/select all.

Is there a better way to do this?


texts: Heat Map to find 10 restaurants, 10 retail, and 10 stores open after 6pm within 3 linear blocks

texts:
I would like to find where in the City there is a combination of 10 restaurants, 10 retail, and 10 stores open after 6pm within 3 linear blocks. What is the best way to do this? What data would I need and how would I analyze this? I would assume I need hours of operation, not sure how to get this (business license data, business associations, scraped from websites)? I will need land use classification, zoning data, block parcels, …

Any ideas would be helpful. I was also thinking of do the same analysis using zoning data to compare where this combo is possible vs where it exists. I would like to use QGIS or GeoMedia.


texts: Plotting a CSV with coordinates one row at a time via GeoPandas

texts:
I have a CSV of all earthquakes in my country for the last 10 years arranged by date and time. My goal is to plot each line chronologically one by one, so that I can compile these maps and produce a "time-lapse" of the earthquakes. The CSV is already sorted but I'll code the "sorting" function when I get there.

So if I have 10 earthquakes, the end products are 10 maps with each a single plot of where the earthquake happened.

My code:

from shapely.geometry import Point
import geopandas as gpd
from geopandas import GeoDataFrame
import pandas as pd
import os

os.chdir(r'path')

def plotPoint():
    df = pd.read_csv('earthquakes.csv')
    geometry = [Point(xy) for xy in zip(df['Longitude'], df['Latitude'])]
    crs = "epsg:32651"
    gdf = GeoDataFrame(df, crs=crs, geometry=geometry)
    basemap = gpd.read_file('basemap.shp')
    for i in range(len(gdf)): # Counts how many rows there are then loops through each line of "gdf"
        if i == 0: # Skips over first line since it is the column header
            continue
        pt = gdf.plot(ax=basemap.plot(figsize=(10,6)), marker='o', color='red', markersize=15)
        print(pt[i]) # Executes "pt" one line at a time
plotPoint()

When I run this, I get

File "C:\path\file.py", line 26, in plotPoint
    print(pt[i])
TypeError: 'AxesSubplot' object is not subscriptable

I'm still new to Python. Any leads?


texts: Processing algorithm fails from Python script

texts:
I'm trying to write a Python script performing a geoalgorithm.
What is surprising me is the following:

I first test the algorithm by QGIS (2.8) interface. In my case, the GRASS interpolator v.surf.idw
I see that the solution is sufficiently good using a certain setting of parameters.
Then, I run the same algorithm with the same parameters, froma Python script. In my case:

out_ras = processing.runalg("grass:v.surf.idw", vl,12,2,"field_3" ,False, "%f , %f, %f, %f "% (xmin , xmax , ymin , ymax), 0.5, -1, 0.001, fileoutput)

where:

vl is the point vector layer
field_3 is the filed where vlaues to be interpolated
fileoutput is the raster file in output
(xmin, , xmax , ymin , ymax) are the Extent of my layer

This setting (perfectly working when launched from QGis interface) produices a Nodata value Raster (only 1 cell).
It seems that the algorithm does not recognize the vector in input.
I've also checked the CRS of the layer (with vl.crs().authid() ) and everything sounds good.

Any help?
Any experience in detecting different behaviour of the SAME algorithm run by Python through processing instead of from QGIS UI ?


texts: How to calculate monthly maximum values from daily data in a 365 band grd?

texts:
I have a MINT2005.grd file which has 365 bands representing days of year from 2005. Each band is a 32x35 matrix representing the temperature at different latitude and longitude pair.

How do I calculate maximum value at each pixel position,so that the output contain 12 array of dimension 32x35 representing the month wise maximum value ?


texts: Replacing two characters of field in attribute table using QGIS Field Calculator?

texts:
Working in QGIS I have an attribute table containing details of polygons I have digitised, including the two fields, COLUMN_1 and COLUMN_2, each of which is a 4 character Integer that I entered.

However, having entered the last two characters of some entries in COLUMN_2 incorrectly, I now need to correct hundreds of records. Specifically the last two characters of the entry in COLUMN_2, for all COLUMN_1 records = 7777.

I could do this in a Spreadsheet but it's cumbersome and requires too many operations.

Opening the attribute table and starting the Field Calculator, I think the solution lies in the expression...

Case When 'COLUMN_1' = 7777 Then regexp_replace(input_string, regex, replacement)
End

or in my 'pseudocode' …

Case When 'COLUMN_1' = 7777 Then regexp_replace(“COLUMN_2”, right(contents of the field ,2), 99)
End

Am I working the right direction and if so how do I write this expression correctly?

The example below illustrates the problem.


texts: Python GDAL/OGR: check lines geometries connectivity

texts:
I'm working with shapefiles representing roads, so basically it's stored as (poly)lines.
I need to perform a road connectivity checking. In other words, I want to ensure that roads meant to communicate are really in contact, because I don't wanna see a little gap between my 2 segments of road.
I don't really see how to proceed.
Should I compute the distance from the current line to the closest ones and then look if this distance is under a threshold? Or is there an easier solution?


texts:
I get a Python error when I run any plugins as batch processing in QGIS 3.6.3. The problem does not solve if I change layers. Any idea how to solve it?
The error I get is:
WARNING    Traceback (most recent call last):
              File "C:/PROGRA~1/QGIS3~1.6/apps/qgis/./python/plugins\processing\gui\BatchAlgorithmDialog.py", line 88, in runAlgorithm
              parameters[param.name()] = wrapper.parameterValue()
             AttributeError: 'NoneType' object has no attribute 'parameterValue'


texts: Create table from results of ST_Intersection queries in PostGIS

texts:
I create a query in PostGIS using 'where' clause, and with that a table. now I want to use those attributes to create a intersection with other table... but I don't know how to call the attributes, because is always asking for the original table and like that it will skip the query

First query to get Rail features:

CREATE TABLE RAIL AS
SELECT rw.gid, rw.name, rw.type, rw.geom
FROM railways rw
WHERE "type" = 'rail';

Here I have my new table (layer in qgis) with only the type "rails". NOw when I want to do the intersection with rivers table (layer):

CREATE TABLE RW_W_RESULT AS
SELECT rw.gid as rwgid, rw.name as rwname, rw.type as rwtype, w.gid as wgid, w.name as wname, w.type as wtype,
ST_COLLECTIONEXTRACT(ST_INTERSECTION(w.geom, rw.geom),1) as rw_w_geom,
COUNT(DISTINCT rw.gid)
FROM waterways w, rail
WHERE ST_INTERSECTS(w.geom, rw.geom) AND w.gid < rw.gid
GROUP BY rw.gid,w.gid;

I got this error: ERROR:  missing FROM-clause entry for table "rw", that is clear to me, but how im suppose to call the attributes of the new created table(layer) to work with this query? If I eliminate rw. from the attributes its says is ambiguous.


texts:
I have a polygon shapefile (buildings) with an 'appart' attribute (see file), how may select the rows whose sum will satisfy a condition?
in other words, how may select all rows whose 'appart' sum is 'x'?
I tried to use QGIS's aggregates functions, but haven't been able to impose a condition on 'group-by'.
The solution is not univocal since many rows may have the same regrouping; but once the first group is found, I assume it should be removed from the following selection.
geojson: https://gist.github.com/pigreco/f439d9f3f15b1941466f9791ff393b6c
EDIT 1:

each building contains 'n' apartments (field 'appart'), how to select the buildings for which the sum of the apartments is >= 15?
EDIT 2:
here is a possible output, create for each select a different layer
https://pigrecoinfinito.files.wordpress.com/2016/10/output_u4o6tu.gif


texts: What GPS unit can I use to track myself all day for a 10 day trip?

texts:
Possible Duplicate:
  Logging GPS and Accelerometer Data - Device Suggestion?

I will be visiting China this October as part of a field studies course and I want to be able to take the opportunity to try out something unusual (maybe usual?). I want to be able to track my movement throughout China using the GPS unit attached to my bag/clothing. I want something that will be reliable/accurate and where the battery lasts a decent amount of time (might go without being able to charge 8 hours at a time). Bonus, this GPS will stay with me for the years to come so I will think of it as an investment, so my budget is fairly loose at the moment.

I've drawn some inspiration from these maps:

http://www.youtube.com/watch?v=iC56P-Ro4PY
http://www.flickr.com/photos/walkingsf/sets/72157629014750905/
http://www.flickr.com/photos/walkingsf/sets/72157628738161697/

Anyways, just thought I'd poke around in gis.se because there should be a fair amount of geographers in here with GPS usage experience!


texts:
I have created a DEM, DSM, and nDSM from a LiDAR Point Cloud using SAGA-GIS.  I extracted the appropriate returns from the point cloud, then converted to grid, and then closed gaps to fill in any nodata regions.

I'm finding that my resultant rasters are still pretty noisy. For example, a region of pixels in the DSM that is obviously a tree or cluster of trees has several pixels with a very similar value to ground instead of the surrounding higher value pixels.

I've tried doing a majority filter on the DSM, but it didn't help.

Any suggestions for other filters that might help create more contiguous regions of pixels?


texts: Removing frame/plot borders from plot.im output

texts: I am plotting point density using the plot.im command in spatstat, and my output always has a frame around the plotted image that becomes thicker when I increase resolution for export. I tried frame=F, axes=F, plot.frame=FALSE, bty='n' but none of them seem to fix the issue. Anyone know a solution to this?


texts:
I am working in ArcMap 10.2.2. I have a shapefile with line features which I labeled with some vb code. The problem is when I zoom out there are labels that disappear, I have around 4000 features but I don't know how to keep all the labels. I would like to create some DDP (data driven pages) and it is important that all labels are in the created maps.
Any suggestions how to do this?


texts: How to export list of feature classes in multi gdb's to a csv file?

texts:
Using the codes mentioned here Listing feature classes in multiple geodatabases in folder using ArcPy? i was able to list all the feature_classes in multiple geodatabases within a folder. Problem is i don't know how to save the print results to a csv file. I knew, in theory - i have to create a csv file -> open it -> Append the "Print_results" to that file. Since i don't have any experience in python i don't know how to execute it.

Currently i'm trying with this code

import arcpy

dir = r'D:\Test\gdb'
arcpy.env.workspace = dir

gdbList = arcpy.ListWorkspaces('*','FileGDB')

filename = r'D:\Test\FC_List.csv

with open(filename, 'wb') as myfile:
    wr = csv.writer(myfile,quoting=csv.QUOTE_ALL)

for gdb in gdbList:
    arcpy.env.workspace = gdb
    datasetList = arcpy.ListDatasets('*','Feature')
    fcList = arcpy.ListFeatureClasses()
    for fc in fcList:
        print arcpy.env.workspace,fc
    for dataset in datasetList:
        arcpy.env.workspace = dataset
        fcList = arcpy.ListFeatureClasses()
        for fc in fcList:
            print arcpy.env.workspace,fc
        arcpy.env.workspace = gdb
        for each_fs in fcList:
            wr.writerow(each_fs.split(","))
            print (each_fs)

I've refered to the code mentioned here https://www.snip2code.com/Snippet/163571/Write-a-list-of-personal-geodatabase-fea. Actually, this code is creating the Fc_list.csv file. But only as an empty csv file


texts:
I used python add in to make a tool bar in Arcmap that edits features on a map. It works on a file geodatabase, but not on an ArcSDE version one. I've been getting an FDO error: 0. The solutions that I've found all involve opening a workspace, but I'm getting an error doing that as well. Here is the relevant portion of my code:

        outFolderPath = r"U:\tempStuff"
        outName = "Stuff.sde"
        databasePlatform = "SQL_Server"
        instance = "SQL Server"
        arcpy.CreateDatabaseConnection_management(outFolderPath, outName,
                                                  databasePlatform,
                                                  instance)
        print "Centroid of buidling:(%r,%r) Side length: %r" % (xCoordinate,
                                                                yCoordinate,
                                                                bob)
        print "workspace is %r" % outName #outName is workspace

        #Start edit session with an undo/redo stack then start
        #edit operation
        print "Starting editing"
        edit = arcpy.da.Editor(outName) #outName is workspace Line 86
        edit.startEditing(True, True)
        edit.startOperation()

        #stuff happens

        cursor = arcpy.da.InsertCursor( fc, ["SHAPE@"])
        cursor.insertRow([stuff])
            #stop the edit operation and then save changes
        edit.stopOperation()
        edit.stopEditing(True)

        del cursor
        print "Building should be drawn"

This is the exact error I've been getting:

Traceback (most recent call last):
  File "C:\Users\153289\AppData\Local\ESRI\Desktop10.3\AssemblyCache\{0BEF4E61-6332-4BCA-A908-959CA3C0B4E7}\leetScripts_addin.py", line 86, in onClick
    edit = arcpy.da.Editor(outName)#outName is workspace
RuntimeError: cannot open workspace

What's am I doing wrong in my code?


texts: Calculating attribute for network distance between multiple points in ArcGIS Desktop?

texts: There is a road network and two point feature classes. For one feature class, I need to calculate a field to hold road (not Euclidean!) distance from every point to a certain point in the other feature class. Is there a way of doing it in ArcGIS 10? If I understand correctly, what I need to do is broadly similar to generating a Near table, but for road distances.


texts:
I am using Boundless Geo-server & Ckan on Cent-OS & I have configured Geo-server auth-plugin successfully. We have nginix server as well as reverse-proxy server.
When I am trying to access map URL outside our environment then using 'auth-key' then I am able to access it but interesting part is to access all layers within our environment, In 'Ckan' using Get-Map request I am not able to display map.
For Ckan : GetMap request : https://domain-name/geoserver/workspacename/layername/ows?service=WMS&version=1.1.0&request=GetMap
I checked this on my side and looks like there is some setting is missing in 'Filter-Chain' part.If I remove this then internal layer works fine but auth-key layers are not working.
Anyone has configured this before ?


texts: How to create a raster "grading" several different layers (polygons)

texts:
So I have a bunch of different shapefiles loaded into qgis:

i.e.:
 - shape of poor soil quality
 - shape of cheap land areas
 - shape of poor weather areas

now I would like to create a map which "grades" the whole map, considering each one of these layers. For example: if the area is under poor soils it gets 0 points for that criteria, but if it is under cheap areas, it gets 10 points.. etc, etc..

any ideas?


texts: Geoserver CQL filter with many “IN” conditions

texts: I am using Geoserver to generate colored regions on a map a wms link, CQL condition can contain a long string of "IN". Max number of INs could be up to 1-2 thousands. I tried to divide it into 3 requests, but such filters are loading the processors very much. Is there any way to get the result and not use long CQL filters?


texts: Clusters disappearing when using minpoints > 2 in ST_ClusterDBSCAN

texts:
I'm using ST_ClusterDBSCAN to cluster polygons. I'm familiar with the algorithm (or so I thought) but the results of the PostGIS implementation are confusing me. When running it with minpoints larger than 2, many of the clusters simply vanish, despite these very same clusters being present in a previous minpoints=2 run while having more than the specified minpoints number of features.
Here is an example using a subset of my complete data. The top image is a clustering done with eps=5 and minpoints=2. The bottom image is the same clustering, but with minpoints=5:
minpoints=2

minpoints=5

Notice how the dark blue cluster in the first image disappears when minpoints is raised to 5, despite there being 12 features in the first case. eps was not changed.
Am I missing a crucial detail in the algorithm or its implementation? Why are these clusters vanishing? I know that DBSCAN is not technically 100% deterministic but I would think that obvious clusters like this just straight-up vanishing is wrong.
So far as a workaround I just run it with minpoints=2 and then filter the results afterwards so exclude clusters with less than my actually desired minpoints value. But this is cumbersome and inefficient when scaling up.
I am using:

Postgres 11.2
PostGIS 2.5.2
GEOS 3.5.1

Here is some sample code to try it out on this data sample:
Create table
Too large for Stack Exchange, here is a pastebin link.
ST_ClusterDBSCAN with different minpoints
drop table if exists buildings_cluster_test_result_mp5;
create table buildings_cluster_test_result_mp5 as (
    SELECT
        geom,
        ST_ClusterDBSCAN(b.geom, eps := 5, minpoints := 5) over () AS cid
    FROM
        buildings_cluster_test b
);
drop table if exists buildings_cluster_test_result_mp2;
create table buildings_cluster_test_result_mp2 as (
    SELECT
        geom,
        ST_ClusterDBSCAN(b.geom, eps := 5, minpoints := 2) over () AS cid
    FROM
        buildings_cluster_test b
);


texts: How to use GeoServer securely in production with Tomcat Policies?

texts:
I want to use GeoServer (first testing, then production) but we have to use Tomcat policies for security reasons.

I have already tried to deploy GeoServer with the security manager enabled, but GeoServer will only work if granting explicit access to the root (/) of the server.

Is this so? I didn't find anything in the docs.

Does anyone know how to get GeoServer working with Tomcat policies?


texts: OpenLayers Map Click won't fire on Android Chrome (Samsung S6 device)

texts:
I'm having an issue getting a map click event to fire properly on a Samsung Galaxy S6 device using Chrome for Android. I have not seen this issue on other devices (yet). I'm using v3.19.1 of the OpenLayers API.

Here's how I create my map and register my click event

        CreateMap : function() {
            this.bLayer = new ol.layer.Tile({ source: new ol.source.OSM() });
            this.vLayer = new ol.layer.Vector({ style:styleFn, source:new ol.source.Vector() });

            var controls = ol.control.defaults({ attributionOptions: { collapsible: true }});
            var view = new ol.View(this.options.view);

            this.map = new ol.Map({
                layers  : [this.bLayer, this.vLayer ],
                controls: controls,
                view    : view
            });

            this.map.on("singleclick", this.onMap_Click.bind(this));

            function styleFn(f) {
                return Styles[f.getGeometry().getType()];
            }
        },

I tried debugging my device, the code never goes into the onMap_Click function. I've also tried a regular click event.

One thing I noticed is that most of the time, when I'm trying to tap on the map, the map pans a little. If I start tapping furiously on the map, at some point one click event will register properly.

Would moving to the most recent version of OpenLayers help at all?


texts:
Following the example Draw Feature

in openlayers 4 i was able to addInteraction and successfully drew a line.

My questions here are

How do i capture all the lines drawn and populate them on the map again?
How do i remove the latest drawn line. (Something like an undo button that keeps removing the latest line)?


texts: Populating a tool's field in Qgis Graphic Modeler

texts:
I am trying to pouplate a tool's X and Y coordinates from a point file and I want to do this in the graphic modeler.
-But how do i connect the coordinates with the appropriate field in the Target X and Y coordinates of the tool, in this case when calculating the Upslope Area with a SAGA tool?


texts: Setting the time step interval in ArcGIS Online

texts:
I'm trying to publish a time aware layer with a time span of about 5 years to ArcGIS Portal, and it works, but instead of using the time step interval I have set in ArcGIS Desktop (1 week), it seems to use it's own (1 year).

Is there any way I can set the time step interval to 1 week in ArcGIS Portal?


texts: How to edit imported shapefile description data in google earth desktop

texts:
I have imported a line shapefile into google earth pro desktop which is displaying fine and in the right geographic location.

The issue I have is that when I click on the layer in google earth the description balloon that pops up has data in it which I cannot edit (see image below).

I would like to find out how to edit that data in the balloon so that data I put into the description shows up instead.


texts: Routing engine that supports right sided pickups

texts:
I am using open street maps and I need to incorporate routing in to my application.  I know there are various routing engines that I can use, but I have a very specific requirement for the direction of travel along the road when the vehicle approaches a stop location.

Specifically this:

I am routing school buses.  Some roads are considered hazardous, which means that students are not allowed to cross the street to get on the bus.  In this situation, I would need the bus to be traveling in the proper direction to prevent the child from crossing the street.  On other (non-hazardous) roads, students are allowed to cross the street to load the bus.

I am currently using asp.net, c#, and javascript.  I would prefer an open source solution but would be willing to pay for something if I know it will fulfill my requirements.


texts:
This may be a question to post on a GIS forum, but I was wondering if among the spatial data experts here, someone could help me calculate the extension of Moran's I to 3D? Our data is about 3D habitat complexity, specificlaly, the 3D complexity of the habitat provided by algae to other marine organisms. The idea is to estimate Moran's I using 3D data we have in which each voxel (3d pixel) is either water or algae. So, the neighbours for each point are distributed in 3d space.

Norbert Marwan used this to calculate bone complexity (medical journal) but there is something about how poitns are accounted for in the index calculation avoid repetition of neighbour pairs that I don't understand

Any help would be most appreciated.

All the best,

Ana


texts: Iterate through date columns in CARTO (frames, VL, etc.)

texts:
I currently have a data file that has columns for each day over an 8-month period.
I am seeking of a way to create an interactive CARTO map that displays the change over time. I want to know if there is a way(or even possible) in Cartoframes or CARTO-VL that the animation over time can change by iterating through the date columns.
I know there is a way to display this change over time by making all date columns unique to each geography (polygon). The issue is I have nearly 850 geographies. Even limiting my dataset to 3 months, we are approaching 1 GB of data (due to the geometries of each field needing to be stored). The size won't generate a time series map on any library I've found or import to CARTO.
With the wide format (each date has its own field) the dataset is only 23MB for 8 months of data. I have not seen any examples iterating through columns but I feel there has to be a way (maybe I'm crazy).
Image of column sample:


texts: CartoCSS & change colours on road-colors.yaml

texts:
I am looking at just having one shading of colours for OSM roads like red.  When I update the hue values to something like 10-45.  I try the generate the generate_road_colours.mss file.  But this usually fails with the one of the values being outside the sRGB.

I have tried changing the lightness and chroma, but I am having no joy, can anyone help?

This is the input for scripts/generate_road_colours.py

# All road classes colours will be generated for, in order of importance
# (biggest first).
roads:
  - motorway
  - trunk
  - primary
  - secondary

# The starting and ending hue. The range goes from 0 to 360, with 0 and 360
# representing red.
hue: [10, 40]

# The lightness and chroma (also known as saturation) for each type of colour.
# Lightness ranges from 0 to 100; dark to bright.
# Chroma ranges from 0 to 100 too; unsaturated to fully saturated.
classes:
  # Colours for output into the MSS file
  mss:
    fill:
      lightness: [80, 97]
      chroma: [35, 29]
    casing:
      lightness: [50, 50]
      chroma: [70, 55]
    low-zoom:
      lightness: [62, 92]
      chroma: [50, 40]
    low-zoom-casing:
      lightness: [50, 70]
      chroma: [50, 65]
    shield:
      lightness: [20, 25]
      chroma: [40, 42]
  # Colours used by the road shields script
  shield:
    fill:
      lightness: [85, 95]
      chroma: [12, 14]
    stroke_fill:
      lightness: [70, 80]
      chroma: [22, 24]


texts: Line through polygons spanning 180-degree line

texts:
In ArcMap 10.2, I am trying to render a shapefile with polygons in the Pacific Ocean, which span the 180-degree line. The resulting image shows a line going through the polygons at the meridian, and I would like to know how to avoid this. Note: the polygons are not split into separate features in the attribute table, they are simply rendered with a line through them.

My basemap is in WGS_1984_PDC_Mercator coordinate system, with Central_Meridian: 150.0, so as to depict the Pacific Ocean in one piece (rather than split on opposite sides of the map). I've tried projecting the polygon shapefile into WGS_1984_PDC_Mercator, as I understand this is a Pacific-centered projection. However, this does not get rid of the line.


texts: How to connnect two lines with different width so it looks good?

texts:
Simple question, how do you create a symbol with or without representation to connect two lines with different line width so it looks good?

I need to have a smooth transition between the lines, otherwise I need to create a polygon (hard job).

https://docs.google.com/open?id=0B8sJ0k30rFr2N1ZRQTRVZldRVGs

https://docs.google.com/open?id=0B8sJ0k30rFr2MHNPemRDVTlRbnM

ok so I tried the solution and it works, but I want it to be moore gradient like next link

https://docs.google.com/open?id=0B8sJ0k30rFr2NVdTTGpGbTBMT2c just for an example

we dont use outline and the round cap makes it looks like a half ball, maybe i will need to use a polygon to mask somehow

the problem is when one long straight line change symbology to much bigger width, and I want the transition to be like the third link, and when one line ends, i need it to be flat so i think i have to use polygon to mask and fill the transitions


texts: Constructing spatial queries that include data from non-spatial tables in PostGIS

texts:
I've got a database with the following tables
CREATE TABLE survey (
id_s integer primary key,
eventdate timestamp without time zone,
locality character varying(100)
geom geometry(Point,25830))

CREATE TABLE provinces (
id_p integer primary key,
name character varying(100)
geom geometry(Polygon,25830))

A simple query would be
SELECT locality, eventdate
FROM provinces p, survey s
WHERE st_contains (p.geom,s.geom) and name = 'myprovince'

But, in order to save space (and, hopefully, speed queries), I re-structured the database as follows:

The problem now is that I can't figure out how to construct that simple query. What woud be the syntax?


texts: Why is ArcGIS API for JavaScript Map.setExtent snapping my extent?

texts:
I'm having some problems setting the Map extent using ArcGIS API for JavaScript 3.18.

Whenever I call setExtent with myExtent, the Map ends up at some extent that is close to, but significantly different than myExtent. For example I call setExtent with:

myExtent = Object {type: "extent", xmin: -13486723.42416174, ymin: 6046817, xmax: -11828314.857440159, ymax: 8563452}

Then I handle the extent-change event to log the actual extent:

actualExtent = Object {type: "extent", xmin: -13334853.828803658, ymin: 6277279.110955888, xmax: -11980184.45279824, ymax: 8332989.889044112}

In this example the actualExtent is significantly more zoomed in than myExtent. But sometimes it's more zoomed out, sometimes it's more north, etc. Basically the Map seems to be snapping myExtent to certain zoom levels and positions.

Our server is ArcGIS Server 10.4.1. I have an ArcGIS API for Silverlight client connected to the same server doing the same work but it doesn't have this snapping problem. Can anyone suggest why this is happening and how to stop it?

UPDATE: To clarify how I build and set myExtent, it goes a little like this:

let myExtent = angular.copy(this.map.extent);
myExtent.xmin = -13486723.42416174;  // calculated
myExtent.ymin = 6046817;  // calculated
myExtent.xmax = -11828314.857440159;  // calculated
myExtent.ymax = 8563452;  // calculated
this.map.setExtent(myExtent);


texts: Layer intersection using ArcGIS Server with Flex/REST APIs?

texts:
I'm trying to implement an intersection with ArcGIS 9.3.1, flex and rest api.

I have two layers. I query the first one that returns a Graphic object. With this Geometry I want to query a second layer with an intersection.

The problem is that the second query returns all the features even outside the intersection geometry.

Is there an example somewhere on the internet?


texts: What should a scheme file contain in ManageTilecache_management?

texts:
I am trying to create a big google kml tiles by using several small mosaic. While doing that in a python command I am having problem with a "Import_scheme". Can anybody help me about this import_scheme file?

For following python command,

arcpy.ManageTileCache_management("H:\FV","RECREATE_ALL_TILES",”"Visible", "H:\FV\MosaicData.gdb\FVImages", "IMPORT_SCHEME", " " , "#", "H:FV\Bound.shp", "#", "4224470", "1128.497176")

I am getting this error,

  File "", line 1, in  File "c:\program files
  (x86)\arcgis\desktop10.3\arcpy\arcpy\management.py", line 15656, in
  ManageTileCache raise e ExecuteError: ERROR 000735: Import Tiling
  Scheme: Value is required


texts: QGIS: Using split lines with lines makes unwanted cut on self-overlapping lines

texts:
I have some long lines that get through some points. I need to cut those lines each time they pass on one of these points.
I managed to do it using the graphical modeler but here comes my issue.
In my model I use split lines with lines tool. There, I cut my long lines with smaller ones I have created on each point. But I discovered that when a line I want to cut overlaps itself, it gets cut on this point and not only when it gets through a point.
So here is my point, is there a way to make "split lines with lines" not splitting my lines that overlaps themselves. Or have you got any idea of how I can manage to deal with this.


texts: Passing attributes from join features to target features in 'SpatialJoin_analysis'?

texts:
I'm doing a spatial join between two polygon data sets (A and B). B polygons are always smaller than A polygons and there's only a small number of them. They have 2 attributes that I'd like to pass on to A polygons ('height' and 'area') if A contains B. One A polygon may contain more than one B polygon.

Let's imagine there's a case where 1 A polygon contains 5 B polygons. The 'height' and the 'area' from B that I want to pass on to A are those ones from the B polygon with the maximum 'height'. In ArcGIS Pro -> Geoprocessing I know how to pass the 'height' [Merge Rule -> Maximum]. However, I don't know how to pass also the 'area' from that very same B polygon with max 'height'.

Steps in ArcGIS Pro -> Geoprocessing (to explain what I've done)

Choose Target Features
Choose Join Features
Choose Output Feature Class
Join Operation: Join one to one
Untick 'Keep all target features'
Setting up 'Field Map of Join features'

Output field 1:  'Out_height'
Source 1:  'Height' (from B)    &    Merge Rule: 'Maximum'
Output field 2:  'Out_area'
Source 2:  'Shape_Area' (from B)    &    Merge Rule: 'Join' (what should I use otherwise if not 'join'??)

Match Option: Contains

What I have done in 'Output field 2' will populate 'Out_area' with the area of the A polygon (which will therefore be duplicated) instead of the area of the B polygon with max 'height (which is what I want).


texts:
I am trying to geocode a dataset using the given latitude and longitude. I am then going to use these points and join polygon features to them. However, the projection I am using for the points is incorrect. How do I go about this in a better more efficient way?

I begin by uploading the files and setting the projection. The polygon file I have already has a projection so I used that projection for the points.

fema <- readOGR(paste(fema.path, "/2010_DFIRM.shp", sep=""))
fema@proj4string
#CRS arguments:
 #+proj=lcc +lat_1=30.28333333333333 +lat_2=28.38333333333333
 #+lat_0=27.83333333333333
 #+lon_0=-99 +x_0=600000 +y_0=3999999.9998984 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0
 #+units=us-ft +no_defs

data <- read.table(paste(inf.path, "/samp_HT_ANALYSIS", sep = ""), sep=",", fill=T)
coords <- data[,c(22:23)]
data.sp <- SpatialPointsDataFrame(data = data, coords = coords, proj4string =
            CRS("+proj=lcc +lat_1=30.28333333333333 +lat_2=28.38333333333333 +lat_0=27.83333333333333 +lon_0=-99 +x_0=600000 +y_0=3999999.9998984 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs"))

However, when I do the spatial join using the following code...

  test <- over(data.sp, fema, returnList = T)

No attributes overlap with my points.

I went on to map the points in ArcMap with a base layer map and they aren't mapped at all on the plane. Can someone please explain how I can solve this? Also, general information on projections in R would be great. I thought I could just set the CRS equal of one file equal to the other but I may have to reproject all of my files.

UPDATE:
I have rerun my code using suggestions below. However, they are still not plotting correctly. I checked some of the indivudal lat long coordinates in my data and they are within my study region. Here is the updated code.

# Data temp for map
data.use <- read.csv(paste(inf.path, "/Houston_sample_subset.csv", sep = ""))
coords <- as.data.frame(matrix(0,nrow=nrow(data.use),2))
coords[,1] <- data.use$PROP_LATITUDE
coords[,2] <- data.use$PROP_LONGITUDE
data.sp <- SpatialPointsDataFrame(data = data.use, coords = coords, proj4string =
                                    CRS("+proj=longlat +datum=WGS84 +no_defs"))

# Study area
study.area <- readOGR("./FEMA_FLOOD/tl_2017_us_county/tl_2017_us_county.shp")
study.area@proj4string
study.area <- spTransform(study.area, crs(data.sp))
study.area <- subset(study.area, study.area$CBSAFP == 26420)

tiff("./PLOTS/test.tiff", width = 5, height = 5, units = "in", res = 800)
plot(study.area)
plot(data.sp, add =T)
dev.off()

The resulting plot is

Only the study region is plotted. I switched to plotting study region vs. fema maps because the fema maps take too long to generate.


texts: How to avoid geopy time out problem with nominatim?

texts:
I am using geopy with nominatim to get cities names from geographic coordinates.
I am using the below code :

for index,row in data.iterrows():
    georeverse = geolocator.reverse
    locationdep = georeverse(row[3],language='fr-FR',exactly_one=True)
    address = location.raw['address']
    state=address.get('state','')
    town=address.get('town','')

My problem is that it throws a time out error after treating 256 rows and my data frame contain over 100 000 rows.
So how would I handle this to avoid the timeout error?