texts: How to find network distances between two sets of origins and destinations?

texts:
I'm learning how to use Network Analyst, and I'm not entirely sure what to do (or search for). My objective is to measure the nearest distance from the red circles to the blue circles along a road network:

I have made a network in ArcCatalog using the network wizard and I'm just interested in the simplest network distance. My network looks like this:

Is this calculation possible to do? My red points are generated as centerpoints of a grid as I am just interested in the average distance for each grid. I realize that these red points are not on the network, so I'm not sure if that is a problem for the calculation.

What should be my next steps? (after reading more of the help documentation!)


texts: Programmatically uncheck layer in AxTOCControl

texts:
I am running my map in dynamic display mode, and I want to programmatically uncheck a layer that is not a IDynamicLayer that is visible in the AxTocControl -- its just a regular layer file.

I thought that layer.Visible = false; with a map refresh would work, but bafflingly it does not.

I'm using ArcEngine 10 VS2010 & C#


texts: How to customize Directions widget of ArcGIS Viewer for Flex 3.4 using Flash Builder 4.6?

texts:
I want to customize the Directions widget in order to be suitable for my application.

However, I do not know how to show a Directions widget without signing in to resources in http://www.arcgis.com.

Is this possible to do?


texts: Selecting data for a region from Google Earth Engine

texts:
I am using Python to work with Google Earth Engine. I would like to extract some data for Paris.
I got the shapefile of Paris and I extracted its boundary as shown below:
import ee
ee.Initialize()
import geopandas as gpd
from shapely.geometry import Polygon
import pylab as plt

gdf = gpd.read_file('drive/MyDrive/Colab Notebooks/shapefiles/France/COMMUNE-shp/')     ## Paris metro area
gdf = gdf.to_crs(epsg=4326)

def returnCityBoundary(gdf):
  ## from shapefile it returns city boundary as a rectangle in different format
  #### polygon_geom: shapely polygon
  #### polygon: geopandas dataframe
  #### geometry: Earth Engine polygon
  bounds = gdf.total_bounds
  lon_point_list = [bounds[0], bounds[0], bounds[2], bounds[2]]
  lat_point_list = [bounds[1], bounds[3], bounds[3], bounds[1]]
  polygon_geom = Polygon(zip(lon_point_list, lat_point_list))
  crs = {'init': 'epsg:4326'}
  polygon = gpd.GeoDataFrame(index=[0], crs=crs, geometry=[polygon_geom])
  geometry = ee.Geometry.Polygon(
          [[[bounds[0], bounds[1]],
            [bounds[0], bounds[3]],
            [bounds[2], bounds[3]],
            [bounds[2], bounds[1]]]])
  return polygon_geom, polygon, geometry

polygon, gdf_polygon, geometry = returnCityBoundary(gdf)
fig, ax = plt.subplots()
gdf_polygon.boundary.plot(ax=ax, color='red')
gdf.boundary.plot(ax=ax)

Now I would like to extract some images in this boundary. This is what I am doing but it does not work.
cover = ee.ImageCollection('MODIS/006/MOD11A2').filter(ee.Filter.date('2018-01-01', '2018-05-01'))\
  .select('LST_Day_1km')\
  .filterBounds(geometry)


texts: Can ArcPy fill forms in Add-in interface?

texts:
I'm trying to build an ArcPy tool, which will ask the user to input some information (ID, Name, Address, Zip, etc.) first. My add-in interface is shown below, and hope to implement that once the user inputs the ID, if all the other related info (Name, Address, etc.) exist in another known table, they could be shown in the following blanks at the same time, instead of letting the user type in everything.

In short, can ArcPy fill forms in the Add-in interface, not in the result window?

The validation works, but extremely slow when I run a search cursor in a .dbf which has more than 160,000 records. How can I improve the code below, or is there a better solution other than using python script tool?
It seems the form will go through the cursor again even after I fill other unrelated blanks.

import arcpy, datetime
import os
import sys
class ToolValidator(object):
  """Class for validating a tool's parameter values and controlling
  the behavior of the tool's dialog."""

  def __init__(self):
    """Setup arcpy and the list of tool parameters."""
    self.params = arcpy.GetParameterInfo()
    fc = "C:\\test\\vectorDBO.dbf"
    field = "PARCEL"
    cursor = arcpy.SearchCursor(fc)
    row = cursor.next()
    n = 0
    while row:
        if row.getValue("PARCEL") == self.params[0].value:
            self.params[1].value = row.getValue("LASTNM")
            self.params[3].value = row.getValue("ADDRESS")
            self.params[4].value = row.getValue("CITY")
            self.params[6].value = row.getValue("ZIPCODE")
            break
        row = cursor.next()

  def initializeParameters(self):
    """Refine the properties of a tool's parameters.  This method is
    called when the tool is opened."""
    self.params[10].value = datetime.datetime.now()
    return

  def updateParameters(self):
    """Modify the values and properties of parameters before internal
    validation is performed.  This method is called whenever a parameter
    has been changed."""
    return

  def updateMessages(self):
    """Modify the messages created by internal validation for each tool
    parameter.  This method is called after internal validation."""
    return


texts: How to make a fixed block of text in cartodb instead of infowindows?

texts:
I started to use CartoDB today for a school project, it is here : https://elenab.cartodb.com/viz/3c08e76a-ec4a-11e5-8d14-0ecfd53eb7d3/public_map

But what I would like to do is to have a block of text below my map, where my textes would appear when I click on country. I don't want a pop up info window, because my texts are long.

I saw this map, did by a french newspaper : http://www.liberation.fr/apps/2015/03/carto-fn/

When you click on the different markers, the text change in a fixed block. I would like to do that with my own map.

Is there any way to do that?


texts: How to assess the visibility of a feature from anywhere on the landscape (DEM)

texts: I'm wanting to assess the visibility of a feature of a given height from anywhere on the surrounding landscape. It seems like the viewshed tools wont work for this since they show the visibility of the landscape from the feature itself or from set observer points. It seems like I need something like the inverse of the viewshed tool.


texts:
I own about 200 acres of forest with about 350' of elevation change. I have a DJI 4 Advanced and my goal is to produce elevation maps to plan the building of one or more ponds.

I've tried and failed, at researching how to do this. It appears to be fairly easy if you don't have a Mac to use the LAS file and the LAS toolkit in QGIS to perform these functions. I've tried to run the LAS tools through WINE with no success. I tried GRASS, but it's LAS toolkit also doesn't work on a Mac.

So my question is "is there a way to filter out trees (it's that time a year again when the leaves are down, so I can get shots through to the floor) with either the LAS file or the GeoTIFF on a Mac?".

If there is no easy way, I'll try again at using WINE, or as a last resort might resort to buying a Windows computer. I've been using MapsMadeEasy to process my pictures.


texts: Cannot run Grass algorithms in QGIS 2.18 on Ubuntu 16.04

texts:
I have QGIS and GRASS installed through the UbuntuGIS repositories. If I try to execute any GRASS 6 or 7 algorithm I get the following error message:

"This algorithm cannot be run :-( It seems that GRASS GIS 7 is not correctly installed and configured in your system. Please install it before running GRASS GIS 7 algorithms."

I checked other questions on this issue however most deal with Setup on Windows systems:

"Missing dependency.This algorithm cannot be run :-(" QGIS + GRASS Windows 7
Grass 6 extension installation in Grass 7 error
How to make GRASS 7 addon available in QGIS 2.12?

On linux I do not have to configure the path of the GRASS folder as the QGIS help on configuring external applications tells me. Grass runs from the shell by typing grass. The Grass-version called by the shell is GRASS 7.2.0 (2016).

Furthermore i checked that GRASS 6 and 7 is activated in the additional providers menu in the toolbox of QGIS.

Any suggestion what could be done to fix this?

Edit

Meanwhile I completly uninstalled and reinstalled qgis python-qgis qgis-plugin-grass from the oficial QGIS repositories...but the error continues.

Edit 2
Checking out the logs for qgis in the Log Messages Panel I found out under the "plugins" tab that the grass plugin does not load properly:

2017-01-23T14:10:22 1   Failed to load /usr/lib/qgis/plugins/libgrassplugin7.so (Reason: Cannot load library /usr/lib/qgis/plugins/libgrassplugin7.so: (libgrass_gis.7.0.3.so: cannot open shared object file: No such file or directory))
2017-01-23T14:10:22 1   Failed to load /usr/lib/qgis/plugins/libgrassprovider7.so (Reason: Cannot load library /usr/lib/qgis/plugins/libgrassprovider7.so: (libgrass_gis.7.0.3.so: cannot open shared object file: No such file or directory))
2017-01-23T14:10:22 1   Failed to load /usr/lib/qgis/plugins/libgrassrasterprovider7.so (Reason: Cannot load library /usr/lib/qgis/plugins/libgrassrasterprovider7.so: (libgrass_gis.7.0.3.so: cannot open shared object file: No such file or directory))


texts: Does OpenGeo Suite no longer exist? What is a good alternative?

texts:
I've been trying to get a copy of open geo suite to work through some tutorials and create a few things. From what I can tell, the most recent version is 4.8. Theres plenty of documentation out there related to the suite and to the tutorials. However, when I try to go to any of boundlesses sites related to this, it seems like they have changed their format entirely and no longer offer any support of this (even previous versions). Anyone have a suggestion for where I can get it now?

"Impossible. Perhaps the archives are incomplete."

"If an item does not appear in our records, it does not exist."

EDIT: I'd also like to add that I personally find it a bit irritating that a company that frames itself as open-source friendly would eliminate such a tool.


texts:
Is there an ArcGIS-like implementation of Flow Accumulation for R?

What I am trying to avoid:
Using ArcGIS; Using Grass; Using RpyGEO (but will if need be)

I am aware of the landscape/hydrology packages of RSAGA, enaR, and vegan, but none of these seem to have a flow accumulation function.  I have seen previous reference to rsaga.flow.accumulation, but that does not appear to be support any longer.

I have tried ??flow, but only the RPyGEO and a flow path from the Raster package show any promise.


texts: Interpolation of 17 data points - temperature data

texts:
I have 17 location points here in Iceland: weatherstation points.

It's a 10-year, month-to-month, composite of temperature values for each station which means 120 values for each station. 2040 values overall. Since I only have 17 station I only have 17 values visible on my map.

What I want to do is interpolate these points to create a temperature map for the area and I want to be able to compare them, f.i. compare june for every year and make a map of that - compare june-sep for each year - compare december for 5 of those 10 years etc. I basically want to be able to compare what I want.

The distance between each point varies but the area is fairly small so only a few kilometers between stations.

When the interpolation is done I hope to be able to do a time animation from the data.

This is a school project. I am attending The University of Iceland.

So, which method should I use and how?

Thanks in advance.

Einar Ingi


texts: Javascript API Scalebar halo CSS Styling

texts:
I am very new to HTML, CSS, and JS.

I have added a scalebar to my map and figured out how to add a halo around the labels so they can be seen better against dark or light backgrounds, but I can't seem to replicate this on the line itself. Below is the sample I found online that i'm using for the labels.

  .esriScalebarLabel, .esriScalebarLineLabel, .esriScalebarSecondNumber, .esriScalebarLabelDiv {
      text-shadow: -1px -1px 0 #fff,
      1px -1px 0 #fff,
      -1px 1px 0 #fff,
      1px 1px 0 #fff; }

I have tried applying this same code (along with other CSS properties like box-shadow, outline, text-shadow) to the .esriScalebarLine selector but the closest I have gotten is a boarder around the lines container, and not the actual line itself.

My bar currently looks like this:

And I need it to look more like this:

Notice the 1px white outline surrounding the line? I need that! Any ideas on how to do that? Am I using the wrong selector?


texts: Desktop ArcGIS Region Group outcome interpretation

texts:
I am new to GIS and I sincerely cannot explain some results.

I have been trying to run the Region Group tool on an unsupervised classification output with 23 classes that looks like this:

The tool used is as follows:

and the outcome is finally as follows:

Could anyone please explain why is the outcome like this?


texts: Using wildcard within ArcGIS ModelBuilder to select raster files?

texts:
I already checked (almost) all other related topics I was not able to figure out how to solve my problem.

I want to process certain rasters in many different subfolders. Within these folders I have also rasters I don't want to process. I want to identify the rasters to process by name. They're all named like 84z93bhjgd_DEM.tif. I thought that putting *DEM into the wildcard field of the iterate rasters menu would solve the problem but it doesn't. Why? Is it not possible to use the * in front of the definition?


texts: Setting Up Snapping Environment Snapping Type using ArcObjects

texts:
I try to setup snapping environment snapping type by referring to this source (http://resources.arcgis.com/en/help/arcobjects-net/conceptualhelp/index.html#//0001000001s1000000) but I can't combine the snapping type like the example. It gives me Operator '+' cannot be applied to operands of type 'ESRI.ArcGIS.Controls.esriSnappingType' and 'ESRI.ArcGIS.Controls.esriSnappingType'

What could be the problem or is there another approach to achieve this?

I try to do this using C# and ArcGIS Engine 10.1.


texts:
I've read many forum posts on the function ST_Transform but unfortunately was unable to resolve my doubts so I posted a new question.

I have a spatial table SRID = 101010 and want to transform it to SRID = 101020, both already on the table my spatial_ref_sys

Initially I did the following SQL and it worked.

SELECT
    ST_Transform (geom, 101020)
FROM
    lotacao

Later I tried to update my geom column

UPDATE lotacao
SET geom = ST_Transform (geom, 101020)

And the error occurs:

[Err] ERRO:  Geometry SRID (101020) does not match column SRID (101010)


texts: Relation aggregate function error in QGIS using Geometry by expression

texts:
I have a point layer, consisting of 12 points with fid 1 to 12. I want to render points with fid > median (fid) - thus in this case, point 7 to 12.
That works perfectly if using geometry generator with this expression - see screenshot, where only 7 to 12 are rendered, whereas 1 to 6 only have labels/callouts:
if (
    "point_id" > median (
        "point_id"
    ),
    $geometry ,
    ''
)

However, when I use the same expression with Menu Processing / Toolbox /Geometry by expression on the same point layer, I get an error: Cannot use relation aggregate function in this context. Why is that  - or asked the other way round: how to correct the expression? It's the same expression as in geometry generator in rather the same context.
I am asking specifically about this relation aggregate function error, as I already encountered it several times in more complex cases. In this simple case, of course there would be easier ways to achieve the goal. However, I wanted to present a simple case to make the error easier to reproduce.


texts: How to get percentage coverage of raster values in an overlapping polygon?

texts:
I have an NDVI image in ArcGIS 10.1 which is from -1 to 1. I have reclassified that image in 0 (range all negative value - non-green) and 1 (all positive value - green).

I have overlayed a 'fishnet' (rectangular grid) of 3 x 3 on image of 5 m resolution. In this case there are green and non green (1 and 0 values) in one fishnet cell.

I want to calculate the percentage coverage of the 'green' pixels within each cell of the fishnet grid (("AREA GREEN" / "AREA of CELL") * 100).


texts: Why don't some layers line up afer converting Mapinfo TAB to Shapefile?

texts:
I have 6 MapInfo layers / TAB files which I need to convert to ESRI ShapeFiles, preferably using QGIS. When a colleague opens these in MapInfo 10, they all look fine and line up correctly. When I open them in QGIS, however, some layers are at completely different X / Y coordinates, and are tiny when compared to other layers.

When my colleague uses the universal translater in MapInfo 10 it creates all the shape files (it actually creates more, there are polyline, rectangle, region and ellipse files), but again when I open up the resulting shape files in QGIS (I click on "user defined projection" during opening) it has the same result.

In MapInfo under Map -> Options -> Projection I can see the projection as "Longitude / Latitude".

Any ideas how I can export all my MapInfo layers to ESRI shape files, and retain their shape / size and location?  I can get access to MapInfo 10 if changes in the source TAB files are required.

Thanks in advance, really stuck on this one!

Dominik

EDIT: I know the process to convert MapInfo to ShapeFiles, and have tried this both in QGIS and MapInfo 10. The problem remains, no matter which application I use to make the conversion, there are problematic layers which are not in the correct location, and are tiny compared to the source layers. If someone has some insight as to how to see the projection for each TAB file, and how to change it that would be the starting point, I guess.  This is NOT a question on how to convert from one format to the other!


texts:
I'm trying to use the Raster Calculator inside ModelBuilder to get rid of null values, but I'm hitting a snag. I'm using the expression:

Con(IsNull("%temp_ras%"),FocalStatistics("%temp_ras%",NbrRectangle(6,6),"MEAN"),"%temp_ras%")

I keep returning Error 000539, which I'm guessing is bad syntax? Is there another way I should be doing this?

This is being used in ModelBuilder as a tool, where temp_ras is a raster layer.

edit: ArcMap 10.0


texts: Why do people pronounce GDAL (Geospatial Data Abstraction Library) like goo-dal?

texts:
Ok, a bit of a silly question, but:

I've heard people pronouncing the library GDAL (Geospatial Data Abstraction Library) like goooo-dal.  My natural inclination was to pronounce it geeee-dal, as in the letter 'g'-dal.

Any explanation on how goo-dal came to be popular?  Was goo-gle involved?


texts: GLCM function for RasterStack using FOR loop in R

texts:
I would like to calculate glcm function in R for all (52) bands in my RasterStack. Using this code:

LS<-stack(img1,img2,img3,img4,img5,img6,img7,img8,img9,img10,img11,img12,img13)

names(LS) <-  paste0("B", c(1:52))

band<-brick()

for (i in 1:52){
     band[i]<-glcm(LS$B1,window = c(3, 3),na_opt = 'center')
    }

I get 1 layer RasterBrick with NA values and an error message:
1: In x@data@values[i] <- value :
  number of items to replace is not a multiple of replacement length

Do you know how to get the result for all the bands in my RasterStack?


texts: QGIS heatmap classification and symbology?

texts:
I want to create a heatmap that doesn't show raster cells that have the value 0. I have 5 classes of which the 0 values have no colour (transparent). Now values between 0 and the first value (124.66...) seem to be shown semi-transparent in the map which looks "cloudy".

How can I create a sharp border between values that are 0 and values that are not 0?

In ArcGIS Desktop you can set the break values for each class but it seems in QGIS you just set one value and the colours are interpolated inbetween?

I want 0 to be 100% transparent and everything else 100% non-transparent.

After adding 0 as no data value and changing the interpolation mode to discrete:


texts:
I want to test performance of ArcGIS Desktop 10.0 based applications. I tried GDBT but it not give me any information about request time etc. I tried also a Workflow Profiler, but this tool usually crash my ArcMap installation. I think about use of ArcPy but I am not sure how to do that properly.

Any Ideas? Maybe someone of you made this kind of tests. For Manual tests I can use AutoIT but I am  not sure what use for loading tests of ArcGIS Desktop. And how to get any queries which are send by ArcMap to Database?


texts: Using Nested Iterators in ArcGIS ModelBuilder?

texts:
I have two polylines (say two paths/roads/etc. for example) that I want to iterate against another set of iterated polylines (ie. geologic fault lines).

I want the first iterator to take Path #1 and iterate against each "Fault Line" (the 2nd iterator).  Then iterate Path #2 against the "Fault Line" polylines.  How do I do this?  I tried but cannot get it to work.


texts: ESRI Table error: does not exist but it does exist in GDB

texts:
I am attempting to iterate through all tables in a Enterprise GDB and get all values in second column but I am getting an IOError error:

processing...
SJLU.FIRE_HYDRANT_PIPE_DIAMR_LK_MV
Traceback (most recent call last):
  File "C:/Users/bailey/Desktop/example_lookups.py", line 14, in <module>
    fieldsInTable = arcpy.ListFields(lu)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.3\ArcPy\arcpy\__init__.py", line 1131, in ListFields
    return gp.listFields(dataset, wild_card, field_type)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.3\ArcPy\arcpy\geoprocessing\_base.py", line 344, in listFields
    self._gp.ListFields(*gp_fixargs(args, True)))
IOError: "SJLU.FIRE_HYDRANT_PIPE_DIAMR_LK_MV" does not exist

Here is my script:

import arcpy
import os

aWS = r"C:\Users\bailey\AppData\Roaming\ESRI\Desktop10.3\ArcCatalog\Connection to GDB.sde"
arcpy.env.workspace = aWS

luList = set(arcpy.ListTables("*LK_MV") + arcpy.ListTables("*_LK"))

print "processing..."

for lu in luList:
    print lu
    fieldsInTable = arcpy.ListFields(lu)
    secondField = fieldsInTable[1].name
    print secondField

Any ideas here?


texts: Starting ArcMap Python Add-in as service?

texts:
Is it possible to start Python add-in as a service which keeps running while in ArcMap?

I am designing something using a add-in to receive socket messages from out of ArcMap and do some operations based on this message.


texts: SDO Oracle Table not recognized as Geometry in QGIS

texts:
I can connect to the Oracle Database with QGIS but spatial table/layers are not recognized as such.

I have read the following in the QGIS documentation:

  Normally, an Oracle Spatial layer is defined by an entry in the
  USER_SDO_METADATA table.
  Source: QGIS Oracle connection documentation

But I have not found this in the database and I cannot edit it.

The layers were created with ArcMap (so with ArcSDE). I know under ArcMap the spatial layers are recognized.
I found that the following tables has informations about spatial tables:

Is it an option to overwrite USER_SDO_METADATA table ?


texts: Comparing rows value of same field in QGIS?

texts:
I have a field in attribute table names "Route", with  rows

I want to comapre row one of Field"Route" with row two of the same field if the value bigger or smaller then insert a value in another field :

WHEN rownumber 1  of " Route" = rownumber 2 of "Route" then insert 1 n " class"

Any idea?


texts: PostGIS operations extremely slow when called in a script

texts:
OK, I am aware that this is rather vague, but unfortunately I was not able not reproduce my problem in simplified code.

I am running a series of functions within another function (main()). All functions do something within a Postgresql 9.4 / PostGIS 2.1.7 database and work well when executed separately.

When running main() basically a couple of large tables and corresponding views with geospatial data will be dropped and recreated. The function in question will then do some polygon analysis using the following functions: st_makeline, st_touches, st_azimuth, st_length. main()uses the same connection which is being committed after each function call.

Now when I run the function from within the main function, it will take take forever (i.e. hours instead of seconds) to execute and I have no clue why that is. If I run the function right after main() has finished, it will also be slower than usual. If I wait 1-2 minutes it will be fast again. So it feels to me like some housekeeping is required. I tired (will be updated):

VACUUM ANALYZE of the involved tables
COMMIT after each step

Do you have ideas what the reason could be or how I can further debug my problem?

Update: Well, I found out how to solve my problem (see my answer below). But if anyone feels like explaining why that works, I will be happy to accept their answer.


texts: Transfering shapefiles from ArcMap to ArcPad makes all form fields required

texts: I've figured out how to create shapefiles in the field and record data.  Sometimes I merge these shapefiles in ArcMap and then transfer them to the GPS device in ArcPad using Get Data for ArcPad.  When I do this, every single field in the data collection screen is now Required, which is just unworkable.  Things like FID, Shape, Latitude, Longitude, y_proj, x_proj, and time (from my most recent example) were carried over from merging shapefiles created on different devices and are now all Required despite being fields created and filled by the devices not the users.  I can delete unwanted fields in the Merge process but there are sonme we may want the data from so I didn't delete them.  How do I stop fields from being Required in this process?  The same thing happens when I create or edit a form in ArcPad Studio and then transfer it.


texts: QGIS polygonize TIFF results in strange vector

texts:
I used an open access data set of GeoTIFFs and used polygonize to transform them into vectors. For the majority of TIFFs this went well, except for two, which create strange longitudinal bars (see images). The TIFF used in the image can be downloaded here (it is the 160_60 TIFF tile/file): https://storage.googleapis.com/uq-intertidal/v1_1/datamask/datamask.zip

Can somebody explain how to solve this problem?


texts: Enabling ArcSDE Archiving with Python and C API using ArcGIS 10.0?

texts:
I have a script that  does some register as versioned and updating privileges.

How would I call into the c api for ArcSDE to enable archiving with python?

http://help.arcgis.com/en/geodatabase/10.0/sdk/arcsde/api/capi/versioning/archcreate.htm

Any direction would be helpful.


texts: Seeking open source Geocoding tool which can be used commercially?

texts:
Is there any Geocoding tool which is open source and can be used in commercial application?

We tried Google geocoding API's but cannot continue with it because of some licensing limitation.

Now we are exploring GisGraphy but haven't finalized it.

Is any other better alternative available?

Currently we want to GeoCode US data only but in future we may extend it to other countries.


texts:
I've tried to make a MouseOver Markers Popup since a couple of days, I'm using Leaflet. I know some solutions I've read in many forums, but I wanna make it work with GeoJSON-Points and I haven't found a solution yet.

Here is a working example, for click popups only:

$('#map').each(function () {
    // map initialization
    map = new L.map('map');

    var gjsonMarker = [{
        "type": "Feature",
        "properties": {
            "name": "start",
            "popupContent": "Portoferraio"
        },
        "geometry": {
            "type": "Point",
            "coordinates": [10.3163552, 42.8072368, 0.0]
        }
    },
    {
        "type": "Feature", "properties": { "name": "stop", "popupContent": "Marciana Marina" },
        "geometry": { "type": "Point", "coordinates": [10.1970935, 42.8069535, 0.0] }
    },
    {
        "type": "Feature", "properties": { "name": "stop", "popupContent": "Bastia" },
        "geometry": { "type": "Point", "coordinates": [9.4520402, 42.6943563, 0.0] }
    },
    {
        "type": "Feature", "properties": { "name": "stop", "popupContent": "Campoloro" },
        "geometry": { "type": "Point", "coordinates": [9.5417333, 42.3393521, 0.0] }
    },
    {
        "type": "Feature", "properties": { "name": "stop", "popupContent": "Cavo" },
        "geometry": { "type": "Point", "coordinates": [10.422377586, 42.859476028, 0.0] }
    }];

    gpsMarker = new L.geoJson(gjsonMarker, {
        onEachFeature: function(feature, layer) {
            if (feature.properties && feature.properties.popupContent) {
                layer.bindPopup(feature.properties.popupContent, {closeButton: false});
            }
        },
        pointToLayer: function (feature, latlng) {
            return L.circleMarker(latlng);
        }
    });

    map.addLayer(gpsMarker);
    map.fitBounds(gpsMarker.getBounds(), {padding: [0, 0]});
});

Any ideas how to make the Popup with MouseOver effect, and when I leave the point the popup must be hidden again.

Greets!


texts: Layer Definition for previous 30 days (JSON)

texts:
I am writing a JSON query into my layer definition of a web map on ArcGIS Online to only show features that were edited in the past 30 days using current timestamp.
What I'm using now (seen below) is not working. No errors are given but it simply hides all features even if features were in fact edited within the 30 day range

"layerDefinition": {
    "definitionExpression": "last_edited_date BETWEEN CURRENT_TIMESTAMP - 30 AND CURRENT_TIMESTAMP"

  },

Any ideas how to correct this code?


texts:
I have a problem trying to import GeoPandas. Since it was impossible to install the package in the normal channels within Anaconda, because of possible incompatibility with other spatial packages, I create a new channel called geo_envas it was sugested in other chats, as follows:

conda create -n geo_env
conda activate geo_env
conda config --env --add channels conda-forge
conda config --env --set channel_priority strict
conda install python=3 geopandas

by this way it was possible to install successfully Geopandas, however, now I have a mixing of environments and by importing geopandas as pdg I get:

import geopandas as gpd
Traceback (most recent call last):

  File "<ipython-input-28-a62d01c1d62e>", line 1, in <module>
    import geopandas as gpd

ModuleNotFoundError: No module named 'geopandas'

and even if I apply

  python -c "import geopandas; print(geopandas.version)"

Now, I feel that I see the end of the path!! I have not idea now how to make this work, and I really need Geopandas installed for some spatial analysis. Important to be mentioned, that I also have installed shapely and fiona. Both of these are sometimes recommended to have installed before Geopandas.


texts: What can POLYCOLOR mean and is there a way to convert to RGB?

texts: It is not for the first time I found in the tables of some data column called POLYCOLOR. Is there any way to figure out what that means and is it somehow convertible to RGB values? If I get it right it should be connected to ArcView 3.0a color system.


texts: Reading file geodatabases without using ArcObjects or any other licenses

texts:
I am interested in writing C# class libraries that can read/modify file geodatabases.

The issue I have with the File Geodatabase API from ESRI is that this is for commercial use, and I don't want to pay licensing fees. So I'd rather write something myself. Has anyone written documentation regarding the makeup of these files, similar to how there's documentation on the different types of shapefiles from ESRI?


texts: What is ST_Geometry, ST_Curve and ST_Surface on PostGIS?

texts:
I dont't have any experience with geographics databases so I hope that you be patiente with me, so I opened this question (PostGis and datatypes WKT) in the past week about the datatypes in the postgis, but later I continued to read the documentation and then I found the following.

ST_Curve = http://postgis.net/docs/manual-2.3/ST_Length.html
ST_Geometry = http://postgis.net/docs/manual-2.3/ST_GeomFromEWKB.html
ST_Surface = http://postgis.net/docs/manual-2.3/ST_Area.html

What are these datatypes? they are instanciable? and why they aren't here https://github.com/postgis/postgis/blob/svn-trunk/liblwgeom/liblwgeom.h.in#L82 (that link is the link that was used to answer my question that I mencioned on the begining of the question)


texts:
I would like to write a script that includes a selection from the attribute table.
But my data is in Hungarian, so it includes special characters like: öüőűáéí

So when there is non of these in the query it executes, otherwise I get this:

Is there a setting or something to solve this?


texts: When using react-leaflet, for Marker it shows

texts:
When using react-leaflet with create-react-app, for Marker it shows <img src=[object Module]> in console after changing the default icon as given in here - https://github.com/PaulLeCam/react-leaflet/issues/453#issuecomment-541142178
. Hence I am unable to view the marker Icon on Map.
delete L.Icon.Default.prototype._getIconUrl;
L.Icon.Default.mergeOptions({
    iconRetinaUrl: require('leaflet/dist/images/marker-icon-2x.png'),
    iconUrl: require('leaflet/dist/images/marker-icon.png'),
    shadowUrl: require('leaflet/dist/images/marker-shadow.png')
});
return (
    <>
        <Marker position={location}> </Marker>
    </>
)


texts:
I am attempting to fill out a VPAT for QGIS so that it can be made more accessible in the federal IT environment. I could use some help and wanted to explore this venue as a means of getting it. If you are familiar with the QGIS codebase, federal IT accreditation, or rely on assistive technology you may have insight which could be useful for greatly expanding the audience for this wonderful software.

Format: I will post the section of the VPAT in question below. Below that I will try and translate the section into a more easily digestable format. I am looking for anybody who can speak to the level of support QGIS provides to meet this item specified. Specific technical information about what features QGIS has that helps it to do this are also useful. If QGIS does not meet this and you know why I would also like to know this and have it documented.

Section 1194.21(c): "A well-defined on-screen indication of the current focus shall be provided that moves among interactive interface elements as the input focus changes. The focus shall be programmatically exposed so that Assistive Technology can track focus and focus changes."

In other words: "The point of action (e.g., cursor) is easy to find and is programmatically exposed to AT(Assistive Technology) so that AT can track focus and focus changes."

Anyone have any insight? If this works out I may be posting similar items from the VPAT that I have struggled answering myself.


texts: How to use sqlite symdifference with ogr via command line

texts:
I am trying to use the sqlite symdifference function in windows command line with ogr but can't seem to get the syntax right. Here is what I have:

ogr2ogr -dialect sqlite -sql "select symdifference(a.geometry, b.geometry) from shp1 a, shp2 b" output.shp E:\temp\shp1.shp E:\temp\shp2.shp

I think there is an issue with how am I referencing the 2 source files (shp1 and shp2) at the end. Is it possible to run this function from command line, and if so how should the syntax be adjusted?


texts:
I've seen varying questions relating to this question. Except mine revolves around C#.  I'm trying to validate that the projection I'm expecting matches what is in all the .shp files I'm reading.

GEOGCS["GCS_North_American_1927",DATUM["D_North_American_1927",SPHEROID["Clarke_1866",6378206.4,294.9786982]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]]

For example I have this prj file, is there anyway to verify that is the SRID that I'm looking for (In the .Net library)?

I'm looking for an answer relating to the .Net Framework:
I don't have ArcObjects or any Esri Libraries. Any chance there is something baked in to the .Net Framework?


texts: Defining point location at EPSG 4326 image

texts:
I am new in the field and I would like to ask which is the correct way to define a point location for an image with EPSG 4326 crs.

I used QGIS software, and created a .csv file with the latitude - longitude information in decimal degrees (I was given the information of those coordinates both in Degrees - Minutes - Seconds, and decimal degrees. I checked these from Google Maps to see if they represent my locations of interest.) I added this .csv file I created as "delimited text layer" in QGIS to define the points in QGIS. And I used these points to extract remote sensing data from the image.

My question is; is defining the points in decimal degrees is correct for EPSG 4326 crs images?


texts: Import CSV file into Postgis DB with Location-field formatting?


texts: r raster:: extract band (integer) values from raster buffer without taking the mean

texts:
I want to extract values from my RGB raster, which are in the buffer(zones) of my 14 station points data.
So far the values are extracted with the raster::extract by taking the mean of some values instead of returning integer values.

For the next step I need the first band values as the original integers from the buffers of my points data. The values vary between 0 and 255.
My progress so far:

# Loading the raster
rasterfile <- raster("rasterfile_example.tif")
rasterfile <- projectRaster(rasterfile, crs = "+init=epsg:3035")

# My station_points:
> head(station_points)
    station_latitude_d station_longitude_d buffer
120           52.54304            13.34933   1000
128           52.48581            13.34877   1000
137           52.39841            13.36810   5000
142           52.47319            13.22514   5000
144           52.48945            13.43084   1000
157           52.64417            13.48308   2500

structure(list(station_latitude_d = c(52.543041, 52.485814, 52.398406,
52.473192, 52.489451, 52.644169, 52.447697, 52.463611, 52.467534,
52.481709, 52.514072, 52.5066, 52.513606, 52.438115), station_longitude_d = c(13.349326,
13.348775, 13.368103, 13.225144, 13.430844, 13.483079, 13.64705,
13.31825, 13.441654, 13.433967, 13.469931, 13.332972, 13.418833,
13.38772), buffer = c(1000, 1000, 5000, 5000, 1000, 2500, 5000,
150, 150, 150, 150, 150, 1000, 150)), row.names = c(120L, 128L,
137L, 142L, 144L, 157L, 162L, 167L, 169L, 170L, 171L, 173L, 174L,
175L), class = "data.frame")
# read coordinates
coordinates(station_points) <- ~station_longitude_d+station_latitude_d

# assign CRS
proj4string(station_points) <- "+init=epsg:4326"

# make sure raster and station points are in the same projection
station_points <- spTransform(station_points,
    crs(rasterfile))

# The extract function I tried:
rasterfile_extract <- raster::extract(
  rasterfile,
  station_points,
  buffer = station_points$buffer,
  method = "simple",
  small = FALSE,
  cellnumbers = TRUE,
  weights = FALSE,
  normalizeWeights = FALSE,
  df = TRUE,
  layer=1,
  nl=3
) %>% as.data.frame()

The results from the extract:

> head(rasterfile_extract)
  ID  cells first_band
1  1 272716   255.0000
2  1 272717   255.0000
3  1 272718   244.3879
4  1 272719   213.4301
5  1 272720   209.0000
6  1 272721   209.0000


texts: Unable to find an inherited method for function ‘writeRaster’ for signature ‘"SpatialGridDataFrame", "character"’

texts:
I have the following error when trying to use the writeraster function in R:

  Error in (function (classes, fdef, mtable)  :
    unable to find an inherited method for function ‘writeRaster’ for signature ‘"SpatialGridDataFrame", "character"’

This is the code:

library(rgdal)
x <- readGDAL("LC81360442014322LGN00_B5.TIF")
z <- reflconv(x,2.0000E-05,-0.100000)
writeRaster(z,filename="test.tif",format="GTiff",overwrite=TRUE) #ERROR LINE

I want to store variable z in a raster file. Can any one help me out from this error?


texts:
I'm looking into performing a spatial correlation test of 2 matrices to see where they overlay. To do so I need to create a grid layer of 2.2km2 around the coordinate points. I've tried the vector - research tools -vector grid but the parameters are in XY, while the UTM locations set to WGS84 - UTM36N, same as the project projection. How can I create a grid layer of 2kmx2km for the extent of my trap layout layer?
Using QGis 2.8.


texts:
I need to export layer features to csv file. I import csv and pyqgis lib.

When I get attribute for each feature I have weird format for each feature.attributes():

[32L, u'value_001', u'value_002', NULL, NULL, NULL, NULL, NULL, 3496L, u'value_003', u'value_004']

Could you explain me why "32" integer value is return as "32L" by feature.attributes()  ?

This is my code :

layers = qgis.utils.iface.legendInterface().layers()
layer = layer[0]
features = layer.getFeatures()
attributes= [el.attributes() for el in features]

QGISv 2.14.2 ESSEN -> use python 2.7


texts:
I am creating maps of various section of the Appalachian Trail.

I have approximately 40 GeoTiffs of STRM data that I would like to merge together.

As I understand these images each gradation is related to an elevation so 0 is 0 meters above sea level and 200 would be 200 meters above sea level. Correct me if I am wrong.

In the individual GeoTiffs, the lowest value that I have is 0 and the highest value is 1324.

When I merge the files the range goes from 0-746.

How would I go about merging them to make sure they all come together fine


texts: Filtering WFS exposed data in GeoServer?

texts:
I want to create web map of the recycling facilities that will display the nearest facilities with the information on what can be recycle as well as enable the user to search for facilities that recycle specific material(s).

I have a shapefile with point data of recycling facilities and materials are represented in table Materials and are comma separated.

I've uploaded the shapefile and configured it on GeoServer and now I want to create the web map. I'm wondering what is there any way that I can filter data according to one field:

For example (this is a link created by GeoServer):

http://localhost:8080/geoserver/Recycling/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=Recycling:RecycleSites_latlong&maxFeatures=1&outputFormat=application%2Fjson&

And I wonder if I can get something like this (filter and get only Materials field)

http://localhost:8080/geoserver/Recycling/ows?service=WFS&version=1.0.0&request=GetFeature&SomeField=SomeValue&typeName=Recycling:RecycleSites_latlong&maxFeatures=1&outputFormat=application%2Fjson&

I'm a new GeoServer user and do not have links to where I could read more about it.


texts: Clip labels to the extent of another data frame

texts:
I'm making a townwide map of a water distribution system with a zoomed in inset (a seperate data frame) for the densest part of town. For the townwide map, I'd like to show the labels for pipe size and diameter everywhere except within the bounds of the inset, as that information is redundant and clutters the townwide map.

The red line on the town map shows the extent of the inset data frame. Is there a way that I can clip or hide labels within the bounds of the extent of the inset data frame?


texts: Creating vector (with lines/polygons) from raster based on paper map in QGIS

texts:
For a task I have several maps created with GIS that look like the following:

I am mainly interested in the colored zones that aren't roads (red) or rivers (blue) and would like to turn these into vectors for further spatial analysis. The extraction of this by a simple Raster to Vectors gives me incomplete lines/polygons because the edges of these zones aren't sufficiently sharp.

I have tried some approaches using reclassify values (with SAGA) to clean up the edges somewhat, edge detection using the Orfeo Toolbox, and raster calculator to make the non-white values stronger, but this doesn't seem to lead to sharper edges. Because these zones have a lot of edges, manual digitizing of the edges would be too time-consuming.

Are there any ways to automatize this process?


texts: Need to join attributes from one point feature class to another based on line connectivity

texts:
I need points representing houses to have the ID number of the electrical pole sourcing the house.
I would use a spatial join "WITHIN_A_DISTANCE", but some houses have sources that are actually closer in proximity but not connected by a line.
I must also account for houses that are connected to power sources by lines made up of multiple polyline features (I think I could handle this with a merge if necessary)
Is there a way to perform a join that joins points from two separate point feature class based on their connectivity via a third line feature class?

The homes need the attributes of the green power source that is connected to them by blue line.


texts: Speeding up ArcGIS arcpy.Merge_management for MANY files?

texts:
I'm currently running a script that is using the arcpy.Merge_management function to merge roughly ~1500 shapefiles together.

When I tried this with much smaller test runs, it ran perfectly.

Now that I'm running all the files it is taking a VERY long time.

I let it run over the weekend, and it is not even half way complete, even though it is still running.

Does anyone have any suggestions on how I could possibly speed up the process? Or would able to suggest another method to merge the files together in a timely manner?

All the files being merged have the exact same schema. None of the files have any more than 150 records. The majority of the files being between 1-20 records.

Here is the code snippet of this part of the script:

# Set the directory and name for the fully merged MESQ_Streets File

outputfile = "Q:\GIS\Admin\Shapefiles\Streets\MESQ_STREETS_SINGLELINE\MESQ_STREETS_SINGLE_LINE.shp"

# This merges the Roads back together

workspace = "Q:\GIS\Admin\Shapefiles\Streets\MESQ_STREETS_SINGLELINE\Rejoined\\"
arcpy.env.workspace = workspace
try:
     arcpy.Merge_management(arcpy.ListFeatureClasses(), outputfile)
except Exception as e:
     print "This did not work", e


texts: Using OpenLayers 3 with feature overlay?

texts:
I left out some very relevant information in my first post. I realized the issue I'm having isn't because of the counties layer itself but of the Feature overlay that is created when the user hovers over a feature from the counties layer.

My issue is this: I have two vector layers (USA states & USA counties) and what I am trying to do is have the map zoom to the extent of the state clicked, and also when a county is hovered over with the cursor it creates a new feature overlay (to highlight it) and pulls county information from the counties layer. So the clicking and zooming action should be limited to just states, and the hover and highlighting/information pulling should be limited to the county. But while the zoom to extent works for the state clicked it will also zoom to the extent of the county clicked (since it has to be hovered over to do so).

I think what is happening is that since the feature overlay is independent of a layer, the layer filter in the foreachfeatureatpixel does not ignore it? So the problem seems to be clicking on the feature layer results in the same zooming to extent action that I want to limit to only clicking on a state. Is there a way to do this such as disabling the on.click event when the feature overlay exists?

I have looked over Using OpenLayers 3 forEachFeatureAtPixel to specify layer? and tried implementing the same attempt at filtering selectable layers with the foreachfeatureatpixel method. But is there away to also filter a feature overlay too?

I am new to both JavaScript and OpenLayers so perhaps I am going about this in the wrong way.

Here are my layers:

var statesLayer = new ol.layer.Vector({
      title: 'usa states',
      source: new ol.source.GeoJSON({
          projection: 'EPSG:3857',
          url: 'simple_states.geojson'
      }),
      style: new ol.style.Style({
          stroke: new ol.style.Stroke({
              color: [44, 127, 184, 1],
              width: 2.5
          })
      })
  });
  var countiesLayer = new ol.layer.Vector({
      title: 'cpo counties',
      source: new ol.source.GeoJSON({
          projection: 'EPSG:3857',
          url: 'simple_counties.geojson'
      }),
      style: new ol.style.Style({
          stroke: new ol.style.Stroke({
              color: 'white',
              width: 1.75
          }),
          fill: new ol.style.Fill({
              color: [127, 205, 187, 0.6]
          })
      }),
      maxResolution: 4000
  });

Here is my attempt to use a singleclick interaction on the map.

      map.on('singleclick', function (evt) {
          var stateFeature = map.forEachFeatureAtPixel(evt.pixel, function (feature, layer) {
              return feature;
          }, null, function (layer) {
              return layer === statesLayer;
          });
          if (stateFeature) {
              var extent = stateFeature.getGeometry().getExtent();
              map.getView().fitExtent(extent, map.getSize());
          } else {
              stateFeature = null;
          }
      });

Here is the overlay:

    var featureOverlay = new ol.FeatureOverlay({
          map: map,
          style: styleFunction
      });

And here is the method for creating the feature overlay for the highlight and getting the county info:

      map.on('pointermove', function (evt) {
          featureOverlay.getFeatures().clear();
          var countyFeature = map.forEachFeatureAtPixel(evt.pixel, function (feature, layer) {
              return feature;
          }, null, function (layer) {
              return layer === countiesLayer;
          });

          if (countyFeature) {
              featureOverlay.addFeature(countyFeature);
              var countyName = countyFeature.get('NAME');
              var pop = countyFeature.get('POP');
          } else {
              countyName = "\u00A0";
              pop = "\u00A0";
          }
          document.getElementById('summaryLabel').innerHTML = '
' + countyName + '
' + pop + '
';
      });


texts: Should I use my DEM or ortho file to create my Hillshade layer?

texts: I've just finished processing my DEM and ortho files, and have been asked to make my hillshade layer. I'm not sure whether to use my DEM or ortho file with this tool. What is the difference between those two files, and which one should be used with the hillshade tool?


texts:
I'm newbie with Linux and QGIS. I'm using Ubuntu 18.04 LTS and QGIS 3.0 Girona and I can't do mosaic from a raster pile because the following log error message returns:

 Processing algorithm… AlgorithmDialog Input parameters: { 'DATA_TYPE'
 : 2, 'INPUT' : '', 'NODATA_INPUT' : 0, 'NODATA_OUTPUT' : 0, 'OPTIONS'
 : '', 'OUTPUT' :
 '/tmp/processing_b9cc3491592e4af584caddae75abe2cf/d6cd346d9c534b2a9bb38c222bbc6cf4/OUTPUT.tif',
 'PCT' : False, 'SEPARATE' : False }

 GDAL command: gdal_merge.py -n 0 -a_nodata 0 -ot UInt16 -of GTiff -o
 /tmp/processing_b9cc3491592e4af584caddae75abe2cf/d6cd346d9c534b2a9bb38c222bbc6cf4/OUTPUT.tif
 --optfile /tmp/processing_b9cc3491592e4af584caddae75abe2cf/mergeInputFiles.txt
 GDAL command output: /bin/sh: 1: gdal_merge.py: not found

 Execution completed in 0.02 seconds Results: {'OUTPUT':
 <QgsProcessingOutputLayerDefinition
 {'sink':/tmp/processing_b9cc3491592e4af584caddae75abe2cf/d6cd346d9c534b2a9bb38c222bbc6cf4/OUTPUT.tif,
 'createOptions': {'fileEncoding': 'System'}}>}

I presume that the gdal_merge.py isn't installed or activated. How can I fix this problem?

OBS.: Maybe other GDAL functions have the same problem, but I have not tried yet.


texts:
I am trying to convert a 13 character UNIX timestamp to datetime and I am getting a year that is in the 44000's. I'm not sure what I am doing wrong, hopefully someone can help.

library(maptools)
GMT <- readShapePoints("~/Documents/Data/Cyber/US/CONTIG_US/UNF01_CONT")
GMT1 <-   as.POSIXct(as.numeric(as.character(GMT$Encoded_Ti)),origin='01-01-1970', tz = 'GMT')
Y <- sample(GMT1, 10, replace = FALSE)
print(Y)

I am grabbing a random sample just to show the results because my file contains 100K+ points. The results I get are as follows:

[1] "44086-01-29 06:33:40 GMT" "44071-01-06 17:21:10 GMT"
[3] "44058-02-21 04:30:10 GMT" "44065-10-06 14:21:40 GMT"
[5] "44061-07-25 17:34:40 GMT" "44038-03-28 21:28:10 GMT"
[7] "44047-07-02 14:27:00 GMT" "44090-04-01 22:37:20 GMT"
[9] "44040-02-16 12:33:40 GMT" "44068-10-08 14:17:20 GMT"

Does anybody have any ideas as to what I am doing wrong?


texts: Google Earth Engine - Reprojecting Vector FeatureCollections

texts:
I am calculating zonal statistics in Google Earth Engine, using an uploaded shapefile(EPSG:4326) and raster dataset (GHSL - EPSG:54009), so I need to reproject one of the layers.

Whilst I have reprojected the raster dataset to the EPSG:4326 and this solves the problem, I still would like to know how to reproject the vector layer to EPSG:54009 for future calculations.

I have tried several solutions, one of which results in an error of:Cannot read property 'transform' of undefined.

Please can someone update my script?

// Define an arbitrary region in which to compute random points.
var region = ee.Geometry.Rectangle(-119.224, 34.669, -99.536, 50.064);

// Create 1000 random points in the region.
var randomPoints = ee.FeatureCollection.randomPoints(region);

// Display the points.
Map.centerObject(randomPoints);
Map.addLayer(randomPoints, {}, 'random points');

// Define mollweide projection as wkt
var wkt = ' \
  PROJCS["World_Mollweide", \
    GEOGCS["GCS_WGS_1984", \
      DATUM["WGS_1984", \
        SPHEROID["WGS_1984",6378137,298.257223563]], \
      PRIMEM["Greenwich",0], \
      UNIT["Degree",0.017453292519943295]], \
    PROJECTION["Mollweide"], \
    PARAMETER["False_Easting",0], \
    PARAMETER["False_Northing",0], \
    PARAMETER["Central_Meridian",0], \
    UNIT["Meter",1], \
    AUTHORITY["EPSG","54009"]]';

// Define wkt as projection
var proj_mollweide = ee.Projection(wkt);

// Apply transformation
var points_moll = ee.randomPoints.transform(proj_mollweide, 0.001);

print(points_moll);


texts: ArcGIS Online - How to show field value on line symbol

texts:
I have a question, similar to this one:

Showing field value on top of line symbol in ArcMap?

How do I reproduce this effect (field value on top of line symbol) in ArcGIS Online?

EDIT:

I looked through the Change Symbol  option but could not find anything for dropping labels on top of the lines, just various ways to change the line width, color, and field to represent.


texts: Strange URL returned when GeoServer WMS layer added to OpenLayers

texts:
I want to add a WMS layer from GeoServer to an OpenLayers website. However, when I try to load it into the website, I get following error in the console:

GET localhost:8080/geoserver/test/wms?%20LAYERS=test.put_25_vlak_1_vulling&SERV…=
&FORMAT=image%2Fjpeg&SRS=EPSG%3A4326&BBOX=-90,0,0,90&WIDTH=256
&HEIGHT=256 net::ERR_UNKNOWN_URL_SCHEME

I have tried to follow this link, but I suspect that I have made a syntax error somewhere.

Details are as follows:

workspace name: test
layer name:put_25_vlak_1_vulling

    var map = new OpenLayers.Map('map');
    var wms = new OpenLayers.Layer.WMS( "World Map", "localhost:8080/geoserver/test/wms",
    { layers:'test.put_25_vlak_1_vulling'
    } );
    map.addLayer(wms);
    map.zoomToMaxExtent();

Can somebody help me with this problem?


texts: Getting QGIS to work with Python2.7 Anaconda Distribution?

texts:
I use Python 2.7 from Anaconda installed in a python environment on my computer (MacOS Sierra 10.12). QGIS is installed from the KyngChaos framework, and I know that there are problems with mixing other python distributions (besides the standard mac installed one) with QGIS, but I had the Anaconda Python2.7 working with QGIS up until today when I randomly received an error that the python console could not be loaded.

I then, 1) Uninstalled QGIS 2.16. 2) Installed QGIS 2.18. 3) Added this to my path: /Users/ron/anaconda/envs/python2/bin/, which is the location of my Python distribution from Anaconda.

These are the errors I am receiving from QGIS.

I am basing my actions off of the accepted answer to Setting QGIS Python version? However I'm having trouble where the accepted answer says "you must modify the symlinks in usr/bin". I am just learning about the file system and I am not really sure how I should go about modifying the symlinks.

I do have VirtualBox with Linux on my computer, so if all else fails I could just install QGIS with Conda on Linux (as right now the Conda QGIS install is only for Linx), but I feel like that would be taking the easy way out and I want to try and figure this out first.


texts: Rescaling LAS file with lidR package in R?

texts:
I am trying to process a large .las file received, which appears to have been scaled in non-standard ways for some reason (likely data transfer portability). The scales between x, y, z are different in and in the range of 10^-5 instead of 0.01 or so. The odd scales in non-standard values prevents me from processing or writing the LAS files with writeLAS or other commands, with an error of:

Error: Invalid header: X scale factors should be factor ten of 0.1 or 0.5 or 0.25 not 1.06228168611949e-05

Is there a good way to actually rescale the point cloud with the lidR package?

I know that the scale factors can be accessed through the header information, as per the code and output from sample LAS file below.

I also seek additional discussion on what it means to rescale the point cloud, and what should be used for physical interpretation and writing to DEMs rather than storage.

Sample reading in and displaying header info in R:

    LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
    las <- readLAS(LASfile)
    las@header
    # las@header@PHB

Output includes:

Scale factor X Y Z:       0.01 0.01 0.01

Currently using lidR version 2.02 and rlas version 1.3.1.


texts: Database based solution for TDSP (Fastest)

texts:
I have tried pgRouting in PostGres for finding shortest path

select d.*, bs.starttime, bs.endtime, bs.trip_id from pgr_dijkstra(
    'select id, source, target, traveltime as cost from my_table',
    19,
    21,
    directed := true
) d JOIN my_table bs on d.edge = bs.id;

The above doesn't consider time as parameter. Before selecting the next edge, I want to check the time of the transport as well.

Eg:

A ---- > B (03.00PM) (travel time 5)

B ---- > C (3.30PM) (travel time 15)

B ---- > C (11:00AM) (travel time 10)

Current algorithm consider the path with travel time 10 (Without considering the time).

Is there a standard solution for this. I have seen lot of implementations like A* algorithm. I'm talking about half million edges.

Or is there any way to get a hook on previous node or next node while computing pgr algorithms ?


texts: Maintaining attributes in feature class in memory using ArcPy?

texts:
i am trying to feature to line conversion  that are not multipart lines and Merge output feature class but its loses the attributes.

Here is my code:

arcpy.env.workspace="in_memory/"
with arcpy.da.SearchCursor(LINEFEATURECLASS,['shape@']) as sip:
    i=0
    for row in sip:
        geom=row[0]
        if not geom.isMultipart:
            op="polyLine"+str(i)
            arcpy.FeatureToLine_management(geom,op,"0.001 Meters","ATTRIBUTES")
            i=i+1
    fclist=arcpy.ListFeatureClasses("","Polyline")
    Merge=arcpy.Merge_management(fclist,"D://Merge_Polyline")

Tool working but it losing attributes.


texts: Using SAGA reclassification with tabled categories in QGIS script

texts:
I am using QGIS, and I would like to write a script to help processing the data.

To do so, I have to reclassify heights values from a raster. I have written the following code, where retab is the path to the file containing the reclassification rules (see later) :

##Layer1=raster

from qgis.core import *
from qgis.gui import *
from qgis.utils import iface
from PyQt4.QtCore import *
from PyQt4.QtGui import *
import os
import processing

retab = "E:\\Divers\\Qgis\\Scripts\\tab_reclass.txt"
#processing.runalg('saga:reclassifygridvalues', input, method, old, new, soperator, min, max, rnew, roperator, retab, toperator, nodataopt, nodata, otheropt, others, result)

Layerout=processing.runalg('saga:reclassifygridvalues', Layer1, 2, 0, 0, 0, 0, 0, 0, retab, 0, 0, 0, 0, 0, 0, 0)

The code gives the following error :

  Error: Wrong parameter value: E:\Divers\Qgis\Scripts\tab_reclass.txt

My rules are the following :

minimum maximum new
0   0.5 1
0.5 2   2
2   1000    3

I could not find a clear example on how to use the table rules, so I guess that it is where the error comes from.

If possible, I would prefer using the rules from the script rather than from the external file.


texts: Is it possible to do a GetFeatureInfo request using WMTS?

texts:
I currently have an interactive map with multiple layers based on the WMS protocol, which has been working fairly well. I use GeoServer to serve the layers and OpenLayers for the front-end. My problem has been that the tiles have been loading very slowly for large data-sets. I decided to use GWC, which is built into GeoServer, to seed the tiles to improve performance. All my layers are in EPSG:27700 and it turns out this CRS does not play nicely with GWC as I found out in my previous question, How to use geowebcache seeded tiles with a custom projection (EPSG: 27700) in openlayers.

I decided to use WMTS, which has solved the problem I was having in the previous question, however I have noticed that there is no GetFeatureInfo call with WMTS as with WMS. I have found some mention of using the WMTS REST API to get feature information, however I am having a lot of trouble finding a good example of how this is done.

Does anyone know of a good resource or have a code example to help with my problem?


texts: Gstat krige() error: Covariance matrix singular at location [917300,3.6109e+06,0]: skipping

texts:
Recently, I am facing a strange problem when using krige() in gstat library. This error was not there a few months ago because I am using the same code with the same data that I successfully used before. But now I am getting singular covariance matrix at grid points at which I want to krige the unknown values. I changed nothing in that old code. Just opened and ran, and now it would not run. Upon searching I found one solution set=list(cn_max = 1e10) from a communication between Edzer and Hengl under title "[R-sig-Geo] 'LDLfactor' error in 'krige' function". However, then krige() gives syntax error about cn_max. I also read "Error in In predict.gstat in R" question on this site, but it also does not work for me.

Why I am getting this error suddenly?


texts: Failure to pick up attributes following use of ArcGIS Pro ModelBuilder: Iterate Feature Classes

texts:
I have a folder containing numerous point shapefiles
(Species occurrence records) and (in a separate folder) a multipolygon
shapefile of small parcels of land. For each species, I need to
extract all the polygons which contain a species occurrence record. As
there are many species, and the process needs to be re-run whenever we
get new occurrence records, I have set up a model that iterates
through an input folder (containing species occurrence record
shapefiles) and runs the process.
The simple model (see attached diagram) is set up as follows: Using
the Select Layer By Location tool I select the polygons which contain
the species occurrence points. I then use the Spatial Join tool to add
all the information contained in the attribute table of the point
shapefiles (e.g. species name, observation date etc.) to this. The
final output is a multipolygon shapefile, named after the input
species containing the detail of the original species occurrence
records.
If I do this manually one at a time this workflow works perfectly;
however, when I run the model using the Iterate Feature Classes
command to iterate through numerous species files in the input folder
the final output does not contain the information from the species
occurrence points. When I open the Field Map of the Spatial Join tool,
the only fields present are those from the original multipolygon
shapefile! I have tried to troubleshoot and for whatever the reason the “Feature Class” output from the
Iterate Feature Classes command does not have an attribute table…


texts: How to set keyboard shortcuts for 'Layers' and 'Browser' panels in QGIS 2.12.1 Lyon?

texts:
I have been using the python solution presented in this answer, which has worked well from 2.6 to 2.10 for me, and still works in the LTR 2.8.4. I can't get it to work with 2.12 (and recent point release 2.12.1).

I have my .qgis2/python/startup.py configured based on the answer linked above.

from functools import partial
from qgis.utils import iface

from PyQt4.QtCore import *
from PyQt4.QtGui import *

mappings = {"Layers": Qt.ALT + Qt.Key_1,
            "Browser": Qt.ALT + Qt.Key_2}
shortcuts = []

def activated(dock):
    dock = iface.mainWindow().findChild(QDockWidget, dock)
    visible = dock.isVisible()
    dock.setVisible(not visible)

def bind():
    for dock, keys in mappings.iteritems():
        short = QShortcut(QKeySequence(keys), iface.mainWindow())
        short.setContext(Qt.ApplicationShortcut)
        short.activated.connect(partial(activated, dock))
        shortcuts.append(short)

bind()

From the python console I am able to toggle the visibility of each panel using .setVisible(True/False). However I can't seem to get the keys to bind.


texts: How to find feature in OpenLayers exsisting layer

texts:
I have a vector layer in my OpenLayers applicaiton.

  var layer = new ol.layer.Vector({
      source: new ol.source.Vector({
          url: "my geojson url"
      })
  })

my GeoJSON data fills the vector layer.

I have an AJAX query that gets single feature from server. But I need to check my layer to control if feature is in layer.

$.ajax(settings.url).then(function(response){
     how to check if response feature in my layer?
})

Response is GeoJSON formatted data.


texts:
I have a shapefile representing the US freight highway network available from Oak Ridge National Lab. I'm trying to prepare this as a .xml network for processing with MATSim. The output structure is very simple, something like

<?xml version="1.0" encoding="UTF-8"?>
<!DOCTYPE network SYSTEM "http://www.matsim.org/files/dtd/network_v1.dtd">
<network name="VISUM export national network 2007-11-28">

<!-- ====================================================================== -->

    <nodes>
        <node id="1000" x="730065.3125" y="220415.9531" type="2" origid="1000" />
        <node id="1001" x="731010.5" y="220146.2969" type="2" origid="1001" />
        .....
    </nodes>

<!-- ====================================================================== -->

    <links capperiod="01:00:00" effectivecellsize="7.5" effectivelanewidth="3.75">
        <link id="100365" from="226" to="227" length="921.0" freespeed="33.3333333333333" capacity="5600.0" permlanes="2.0" oneway="1" modes="car" origid="183" type="10" />
    ...
    </links>
<!-- ====================================================================== -->

</network>

Most users of MATSim pull data from OpenStreetMap using osmosis, but I have to use this network for a variety of reasons.

I have seen lots of suggestions in different places as to what I might try. Some of these include:

v.net in GRASS
shp2pgsql suggested by pgRouting tutorial
manually with pgRouting, though this advice may be outdated?
writing code through the sp and rgdal libraries in R to build a node and link structure.
find someone with a TransCAD license to export the network ORNL provides along with the shapefile, though I don't know TransCAD exports to .xml.

It seems that with options 1 - 3, I'll then have to figure out a way to write the .xml from the database. Option 4 seems difficult, but mainly in terms of my time. I'm happy to learn something new, but I don't want to invest a day towards option 1 only to find that option 3 would have been better.

Or is there a simpler way?

Post-script:

Someone asked the MATSim developers if they had a tool to do the shp -> xml conversion, but the answer was that shapefiles are too wild to build one.


texts: How to calculate area of a contour in QGIS

texts: I have prepared a contour map of an area in a shapefile in QGIS. I want to know how can I calculate area of a particular contour?


texts: How to enable GRASS plugin in QGIS 1.8.0 (Linux Mint 13 - 64 bit)?

texts:
Every time I enable GRASS plugin in QGIS (Plugins -> Manage Plugins -> GRASS), I got a message: Exception - G_getenv(): Variable GISDBASE not set.
I googled and found some information here but don't know how to "unset GISRC".
Anyone please show me how to do it?


texts: How to load WMT_MS_Capabilities using Javascript

texts:
I would like to load and parse an WMT_MS_Capabilities using Javascript. I have unsuccessfully tried severaly ways. I am trying with this file but I always get this error:

uncaught exception: [Exception... "Component returned failure code: 0x80004005 (NS_ERROR_FAILURE) [nsIXMLHttpRequest.send]" nsresult: "0x80004005 (NS_ERROR_FAILURE)" location: "JS frame :: file:///home/frowing/Projects/XMLParsing/testXML.js :: <TOP_LEVEL> :: line 10" data: no]

I hope someone can throw some light on this because I have no idea where to go next.

Thanks!


texts: Can File Geodatabase API be Redistributed?

texts:
Esri's file geodatabase API can now be used with GDAL/OGR - http://www.gdal.org/ogr/drv_filegdb.html however each individual seems to have to login to Esri and download the API files individually. For example this is referenced in a blog post here and in the GDAL build notes.

The license agreement before downloading the API shown by Esri is at: http://resources.arcgis.com/node/agreement/3193

Clause 47 states:

  Licensee may develop and distribute software or web applications that use the Esri File Geodatabase API to Licensee's end users.

The best discussion about the license is at Paul Ramsey's blog post from January 2011. This finishes with a final update that:

  Update 4: Received email from Esri confirming that the final license
  will be reviewed to ensure there are no ambiguities and that it
  reflects their intent that the API be usable by any application in any
  application category and the derived product freely redistributable
  and royalty free. So to the extent that the current license has any
  ambiguity it shouldn't be considered a red flag that the final one
  will.

So has or when will the Esri API be available as part of GDAL downloads or installations, and when can we distribute application which include the API?


texts: Hydro-Enforced/surface constraint DEM Creation using GRASS/QGIS

texts:
I am wondering if there is a way to create a hydro-breakline enforced DEM by means of any open source solution?

I am able to create a DEM from Lidar data using GRASS, but would like to enforce some 3D hydro breaklines into the DEM.  I am looking for something similar to creating a DEM from a Terrain available from ESRI.

Any ideas?

-Jeremy


texts: Change opacity using leaflet without plugin

texts:
I have a geojson and try change the opacity chaging the button, but it dont work anyway.

Where statesData is my geojson.js, style and onEachFeacture are others functions that I have.

Here is my button:

0.5
<input type="range" id="sldOpacity" min="0" max="1" step="0.1" value="0.5" />

And here is my JS

$('#sldOpacity').on('change', function(){
  $('#image-opacity').html(this.value);
  geojson.setOpacity(this.value);
});

var geojson = L.geoJson(statesData, {
        style: style,
        onEachFeature: onEachFeature
    }).addTo(map);

I tried put opacity: opacity and create a function opacity() but dont work too. What I have do?


texts:
I have a polygon file with 90 adjacent polygons. I also have a lat/long point file that overlaps the entire surface of the polygon file (170,000 points). My problem is that working with all of the data at once is far too demanding and slows the project considerably. If I could have the individual polygons separated and have only the points that are intersecting with those polygons on those separated layers, I could work smoother as I wouldn't have to deal with the entire data set at once.

My first idea (which is still in progress) is to do a simple join attributes by location. I say it's still in progress because with 90 polygons and 170,000 points I have been looking at the rainbow loading wheel for a while...

Just wondering if this is the best option or if there's a better idea out there. I'd like to be able to work on one polygon at a time and also be able to interact with the points which are associated with that polygon.

QGIS on a Mac.

Thanks =)


texts: QGIS-3.10 raster calculator crashing on Windows

texts:
Trying to add several rasters, all multiplied by certain factors (like (A@1*.8)+(B@1*.003)+(C@1*.2)...) in raster calculator in QGIS 3.10 the program crashes on Windows (that doesn't happen on Linux, at least on 3.8, which is what I use mostly); exploring further, GDAL's plugin raster calculator outputs something like "missing module" trying to do the same thing; in Osgeo4w gdal_calc was not listed when I checked the available commands.
 You may just add two rasters on QGIS's calculator without crashing.


texts: Getting selected layer name in QGIS Modeler

texts: If I have a geoalgorithm in my model and I want the value from one the parameters, then how can I get it? For example I want the layer name in the parameter layer in below figure so I get the value Pts1.


texts: How to reproject MODIS L2 data with lat/lon

texts:
I downloaded the MODIS level 2 data for precipitable water vapor (MOD05_L2.A2015193.1535.006.2015194020321.hdf) and the corresponding geolocation file (MOD03.A2015193.1535.006.2015193213734.hdf). When viewed in Panopoly the Water_Vapor_Infrared field is of type Geo2D and I can plot it on a map. However, the Water_Vapor_Near_Infrared is labelled as type 2D and can't be projected onto a map. I Downloaded the ModisSwath tool and reprojected the MOD05_L2.A2015193.1535.006.2015194020321.hdf using the geolocation file MOD03.A2015193.1535.006.2015193213734.hdf which seems to have worked. But when I look at the reprojected data in Panoply it is still only 2D!

How do I convert the Water_Vapor_Near_Infrared data so that it can be projected onto a map? Is there a specific setting I need to use in the ModisSwath tool or is there a way to do it straight from Panoply (I've also tried it in a few other programs like SeaDAS with no luck).


texts: VBA not able to fill all FNODE and TNODE because of large number (1,200,000) of records?

texts:
I need to populate FNODE and TNODE fields for my roads dataset.

I have an intersections point feature class that has unique ids, and I need each road segment to have an FNODE and TNODE value that corresponds to the intersection id.

The roads are broken into segments at each intersection so that each segment only touches 2 intersection points.

I thought that maybe Network Analyst generated these, but it doesn't look like it actually writes these values.

Does anyone know how to generate these?

In addition, I have used the VBA code but because of its large records(1,200,000 records), VBA is not able to fill all the FNODE and TNODE.


texts: Installing pl/r pl/python pl/v8 in Ubuntu

texts:
I am trying to find instructions on installing the PostGIS Procedural Languages (PL) for R, Python and V8 as required for the "PostGIS in Action" 2ed. book. I searched on google and couldn't find instructions.

Does anyone know if this can be done using the apt-get system without compiling from source?

I am using Ubuntu 14.04 LTS and a daily built version of the upcoming 16.04 LTS for testing.


texts: Merging ShapeFiles with QGIS for UK counties

texts:
I want to draw a map for all the counties of the UK including Northern Ireland. Since I cannot find a source of counties including both Northern Ireland and Britain, I tried to merge two shapefiles with QGIS from different sources (Northern Ireland and Britain) by myself.

On QGIS, it looks good when I look at the map with two different layers.

But when I try to merge the two shapefiles, Northern Ireland overlays Britain as shown below.

Is there anything I can do to get my ideal result like the first image?


texts:
I enabled labels for a streets layer in my QGIS 3.6.0 project via Properties > Labels > Single labels. The following screenshot shows an example of such labeling:

How can I hide individual labels if the there are too many as in the example?


texts:
Is there a way to calculate a field in ArcMap to be one letter up from another field? I'll never need to go beyond letter 'T', so no worries there.

   Field Ltr| Field NLtr

          a       b
          b       c


texts: Is there a way to create custom fills from a hi res raster format (e.g., png file)?


texts:
Suddenly the button for opening the toolbox is gone. I have searched really long already.
Do you have any idea how to get it back?


texts: Excluding several layers from forEachFeatureAtPixel function using OpenLayers 3?

texts:
I have a map with some vector layers that deliver popups when the user clicks on the markers. I use a forEachFeatureAtPixel function to identify these layers. However, there are also several background layers which should do nothing if you click on them. After looking through this forum for methods, I am able to exclude one layer from the function, but how would I exclude several layers?

Background1 and background2 are the layers that should not deliver any popups. My code is like this:

var background1 = new ol.layer.Vector...

var background2 = new ol.layer.Vector...

var feature = map.forEachFeatureAtPixel(evt.pixel, function(feature, layer) {
if (layer != background1) return feature;
});

if (feature) {

//etc.

That works fine, but I also want to exclude background2. I thought I could write "if (layer != background1 || background2) return feature;" but that doesn't work. What should I do?


texts: Change default date format Dynamic Text or add option to Dynamic text Menu

texts:
I frequently use the Dynamic text for date saved in my map templates but my new boss prefers a different format (February 2nd, 2016).  I have the short code he prefers but would like to change it to the default format or add it as an option to the dynamic text menu under the insert menu.  The code is:

<dyn type="date" format="MMMM"/> <dyn type="date" format="yy"/>, <dyn type="date" format="yyy"/>


texts: I am making a map in QGIS to show different program interventions in a number of districts. There are  programs overlapping in some districts. How can I make solid fill for a give program and diagonal pattern fil to show in a given district that we are running multiple interventions


texts:
I'm projecting some features using JSON format on OpenLayers. I need to get the id from the feature (not the OpenLayers' id but the one in the JSON file). The file looks something like this.

{"type": "FeatureCollection","features": [{ "type": "Feature", "id": 7655199, "properties": {"power":"line","wires":"NA","frequency":"NA","voltage":"380000","operator":"NA","cables":"6"}, "geometry": {"type":"LineString", "coordinates":[[8.5310699,49.8531139],[8.5289255,49.8563936],[8.5268405,49.8594383]]}} etc..

Now in OpenLayers I'm trying to do something like this.

var feature_col = gjson.read(ReadFile(json_path + file[i]) , "FeatureCollection");
//Readfile is just a function for reading the file, what a surprise... :P
var power_res = feature_col[0].attributes.power;
var id_res = feature_col[0].id) //for the first feature

So the value in power_res gives line which is correct. However, the value in id_res is something like OpenLayers.Feature.Vector_362 and I need the id in this case to be 7655199.


texts:
I have never done this before, but how do you find the linear mileage of microseismic lines?  I tried researching it online in different forums/tutorials with no luck. I have ArcView 9.3, no other extensions.  I went into the Toolbox Linear Referencing Box and from there, I tried a few of the tools, but no luck.

This is for PostPlot source points (Microseismic survey). Simply put, as above, I need to find the linear mileage with the postplot points that goes on/through the surface owners.  How do I do this?

Be advised, I am not a programmer nor a techie.  If at all possible, keep any answers as simple as possible unless it is unavoidable.

If you need more information, please let me know as I did my best to provide pertinent information.

I have included a pdf of the data if this helps any.

A cheap and easy way to do this is use the Measure Tool.  But, I figure that HAS TO BE A better way to do it then individually for each parcel.  Is there?


texts: Performing image classification from Nearmap WMS in ArcGIS Desktop?

texts:
I have a subscription to Nearmap that provides me with WMS Server link to their images. I would like to perform image classification using this imagery however I can't work out how to use the WMS feed for classification in ArcGiS (v10.3.1) short of exporting a series of maps as a raster.

I've tried Arc's "make image server layer" tool however I can't seem to select an image service layer. I get to the point in the attached image and can't select anything.

Does anyone know how I do this?


texts:
I am trying to get started with GDAL on a windows machine using Python2.7.
I have tried commands such as:

gdalinfo RasterFileName.TIF and gdal_translate.py -projwin ulx uly lrx lry -of GTiff input.TIF output.TIF

and they run without producing any results. What am I missing that makes these not to work.

For example, I am doing the gdal_translate as follows:

    import os
    infile = path+'InFileName.TIF'
    outfile = path+'OUTFileName.TIF'
    string ='gdal_translate -projwin %.0f %.0f %.0f %.0f -of GTiff %s %s' %(ulx,uly,lrx,lry,infile,outfile)
    os.system(string)

Many thanks


texts: Creating georeferenced and projected data from hand-drawn map using ArcGIS Desktop?

texts:
I am preparing datasets for students in a spatial statistics course for a module on spatial autoregressive models. We will use ArcGIS 10.1. The first dataset will correspond to a mythical pre-industrial village in which houses closest to the tannery will be smaller and more tightly packed and as the distance from the tannery increases, houses will become larger and less tightly packed. Likewise, the property values of houses will increase as distance from tannery increases. NB: Tanneries were historically noxious points that smelled dreadfully and often polluted local waterways.

This example is purposely designed with multicollinearity between distance from tannery and size of house. I want my students to struggle with the question of what is explaining property values so they can understand what multicollinearity means in a concrete example. Later modules will include larger houses near the tannery, as well as smaller houses further from the tannery, so that my students can then use a spatial autoregressive model to identify the effects of positive externalities impacting a small house's property values when it located in a neighborhood with large houses, and the negative externalities impacting a large house located near a tannery.

I currently have a hand-drawn map of this mythical village, and have scanned it as jpg. Each building will have an ID number, a property value, a distance from the tannery and a size in square meters. While I could simply plug these data into Minitab, etc., I want to show my students how to use the spatial autoregressive model function in ArcGIS 10.1 (the proper name from the command escapes me as I type this, but it is under the Spatial Statistics Toolbar). For this function to work, however, I must also project my data. I plan on doing this by selecting a 300m x 300 m polygon somewhere on the globe and simply add control points that arbitrarily line up with my mythical village map (also scaled to capture a 300m x 300m area).

The added benefit of georeference this map is that I can also easily calculate the centers of these structures as well as their areas --- so long as I digitize the buildings.

Does this sound like a reasonable solution in terms of creating georeferenced and projected data to fit an image for an imaginary village?


texts:
I would like to update QGIS from 2.2 to 2.4 on Mac OS X (10.8) but I cannot find how to do so.

Do I have to uninstall and re-install all from KyngChaos? If so, I would therefore lose the extensions I have added manually.

Nothing is really clear about this on the Net.


texts: 'reportError' in subclass of QgsProcessingFeedback never running

texts:
I have two layers, layer_a and layer_b. Both have invalid geometries. When I run a processing tool, it throws an error as expected.
Based on the answer posted here by @ndawson, if you subclass QgsProcessingFeedback, you can implement your own custom logic dealing with the error.

class MyFeedBack(QgsProcessingFeedback):

    def reportError(self, error, fatalError=False):
        print(error)

...

res = processing.run("TOOL_NAME", params, feedback=MyFeedBack())

When I run the following script in QGIS Python Editor
class MyFeedBack(QgsProcessingFeedback):
    def reportError(self, error, fatalError=False):
        print("ERROR! ERROR!")

layer_a = QgsProject.instance().mapLayersByName("a")[0]
layer_b = QgsProject.instance().mapLayersByName('b')[0]

feedback = MyFeedBack()

params = {'INPUT': layer_a, 'OVERLAY': layer_b, 'OUTPUT': 'TEMPORARY_OUTPUT'}
rslt = processing.run("native:intersection", params, feedback=feedback)

print("\nDone!")

I get this error:
Traceback (most recent call last):
  File "/usr/lib/python3.8/code.py", line 90, in runcode
    exec(code, self.locals)
  File "<input>", line 1, in <module>
  File "<string>", line 9, in <module>
  File "/usr/share/qgis/python/plugins/processing/tools/general.py", line 108, in run
    return Processing.runAlgorithm(algOrName, parameters, onFinish, feedback, context)
  File "/usr/share/qgis/python/plugins/processing/core/Processing.py", line 159, in runAlgorithm
    ret, results = execute(alg, parameters, context, feedback, catch_exceptions=False)
  File "/usr/share/qgis/python/plugins/processing/gui/AlgorithmExecutor.py", line 71, in execute
    results, ok = alg.run(parameters, context, feedback, {}, False)
_core.QgsProcessingException: Feature (0) from “b” has invalid geometry. Please fix the geometry or change the Processing setting to the “Ignore invalid input features” option.

According to documentation, reportError method should run in case of any error. But the method never runs although the processing tool throws an error.
What am I doing wrong?
I use
QGIS 3.16.2-Hannover, Ubuntu 20.04.1 LTS
QGIS 3.16.0-Hannover, Windows 10 x64


texts: ValueError error from lyr.replaceDataSource()?

texts:
I'm trying to write a script that allows me to look through map documents and replace each layer's data source in each map document. Everything is working logically I believe, and all my loops are working as they should. The datasets that are referenced exist. Yet, I am getting this error :

  Traceback (most recent call last):   File
  "C:\Python27\ArcGIS10.1\Lib\site-packages\pythonwin\pywin\framework\scriptutils.py",
  line 326, in RunScript
      exec codeObject in main.dict   File "C:\Users\Daimon Nurse\Desktop\DFMPROJECT\Scripts\editmapdocument8.py", line 18, in

      Lyr.replaceDataSource(workspace_path, workspace_type, dataset_name)   File "C:\Program Files
  (x86)\ArcGIS\Desktop10.1\arcpy\arcpy\utils.py", line 181, in fn_
      return fn(*args, **kw)   File "C:\Program Files (x86)\ArcGIS\Desktop10.1\arcpy\arcpy_mapping.py", line 680, in
  replaceDataSource
      return convertArcObjectToPythonObject(self._arc_object.replaceDataSource(*gp_fixargs((workspace_path,
  workspace_type, dataset_name, validate), True))) ValueError: Layer:
  Unexpected error

This is the code I am using:

 import arcpy
import os

PATH2 = r"C:\Users\Daimon Nurse\Desktop\DFMPROJECT"
arcpy.env.workspace = PATH2
arcpy.env.overwriteOutput=True

for file in arcpy.ListFiles("*.mxd"):
    filePath = os.path.join(PATH2,file)
    MapDoc = arcpy.mapping.MapDocument(filePath)
    lyrList = arcpy.mapping.ListLayers(MapDoc)
    i = 1
    for Lyr in lyrList:
        print Lyr
        workspace_path =  r"C:\Users\Daimon Nurse\Desktop\DFMPROJECT\DFM"
        workspace_type = "FILEGDB_WORKSPACE"
        dataset_name = 'Zone{0}'.format(i)
        print dataset_name
        print workspace_type
        print workspace_path
        print os.path.join (workspace_path, dataset_name)
        Lyr.replaceDataSource(workspace_path, workspace_type, dataset_name)
    i = i + 1

These are the printed results on the first iteration for dataset_name, workspace_type, workspace_path and os.path.join respectively:

Zone1
FILEGDB_WORKSPACE
C:\Users\Daimon Nurse\Desktop\DFMPROJECT\DFM
C:\Users\Daimon Nurse\Desktop\DFMPROJECT\DFM\Zone1


texts: Extracting pixel time series from Google Earth Engine

texts:
I have written a short program in Python to extract a time series for any given pixel for MODIS data stored in the Google Earth Engine. The code is working fine and returns a data frame containing the relevant band value and date.

import pandas as pd
import numpy as np
from datetime import datetime as dt
import ee

def extract_time_series(lat, lon, start, end, product_name, band_name, sf):

    # Set up point geometry
    point = ee.Geometry.Point(lon, lat)

    # Obtain image collection for all images within query dates
    coll = ee.ImageCollection(product_name)\
        .filterDate(start, end)

    # Get list of images which correspond with the above
    images = [item.get('id') for item in coll.getInfo().get('features')]

    store = []
    date_store = []

    # Loop over all images and extract pixel value
    for image in images:

        im = ee.Image(image)

        # Obtain date from timestamp in metadata
        date = dt.fromtimestamp(im.get("system:time_start").getInfo() / 1000.)
        date_store.append(np.datetime64(date))

        # Extract pixel value
        data = im.select(band_name)\
            .reduceRegion(ee.Reducer.first(), point, 1)\
            .get(band_name)

        store.append(data.getInfo())

    # Scale the returned data based on scale factor
    store = [x * sf if isinstance(x, int) else np.nan for x in store]

    # Convert output into pandas data frame
    df = pd.DataFrame(index=date_store, data=store, columns=[band_name])

    return df

if __name__ == "__main__":

    ee.Initialize()

    latitude = 9.95
    longitude = -10.09
    start_date = '2018-01-01'
    end_date = '2018-01-31'
    product = 'MODIS/006/MOD11A1'
    band = 'LST_Day_1km'
    scale_factor = 0.02

    # Extract data and obtain pd.DataFrame
    output = extract_time_series(latitude,
                                 longitude,
                                 start_date,
                                 end_date,
                                 product,
                                 band,
                                 scale_factor)

However, I am aware of the way GEE handles reducing regions and how it can sometimes be difficult to know precisely what is happening. I decided to compare my extract MODIS values with the same time series requested from AppEEARS.

I noticed that the values are different and thus my code must be performing some spatial average, probably during the:

data = im.select(band_name)\
        .reduceRegion(ee.Reducer.first(), point, 1)\
        .get(band_name)

line. I set scale == 1 here so that the returned value I receive should be equal to the value at the native resolution (info on this).

Can anyone shed any light on how I can extract the MODIS band value for ONLY my point of interest? Any amendments to the code is also very welcome!


texts:
Currently I am getting only a white dotted line in my map using the code below:

style: new ol.style.Style({
                                    stroke: new ol.style.Stroke({
                                      color:[255, 255, 255, 1.0],
                                      opacity: 1,
                                      width: 4,
                                      lineDash: [4, 4]
                                      }),
                                      fill: new ol.style.Stroke ({
                                          color: [155, 155, 155, 0.4]
                                      })

                      })

But my requirement is to get black and white dotted lines in the form black white black white...


texts: Drawing selection custom selection box using ArcObjects?

texts:
I have some code that I'm using to selecting features in a certain extent. I already rigged the tool properly, so it selected the objects as needed, but I want to draw interactively the selection box the user is making.

I'm using IScreenDisplay.DrawRectangle method, but it pollutes the screen with many small rectangles. Is there a way to clear up the old drawings in the screen?

This Draw method is being called from a mouse move event. Any ideas?

Here are two methods I'm using to do the dirty work:

public override void OnMouseMove(int button, int shift, int x, int y)
{
    if (button != 1)
        return;

    if (!_IsSelectOperation)
        return;

    endPoint = focusScreenDisplay.DisplayTransformation.ToMapPoint(x, y);
    envelope.PutCoords(startPoint.X, startPoint.Y, endPoint.X, endPoint.Y);
    DrawEnvelope(envelope, BuildSelectionSymbol(), focusScreenDisplay);
}

private void DrawEnvelope(IEnvelope envelope, ISymbol symbol, IScreenDisplay screenDisplay)
{
    screenDisplay.StartDrawing(screenDisplay.WindowDC, 0);
    screenDisplay.SetSymbol(symbol);
    screenDisplay.DrawRectangle(envelope);
    screenDisplay.FinishDrawing();
}

Here is the code for BuildSelectionSymbol:

private ISymbol BuildSelectionSymbol()
{
    IColor fillColor = new RgbColorClass();
    fillColor.Transparency = 0;

    IRgbColor outlineColor = new RgbColorClass();
    outlineColor.Red = 0;
    outlineColor.Green = 0;
    outlineColor.Blue = 0;

    ISimpleLineSymbol lineSymbol = new SimpleLineSymbolClass();
    lineSymbol.Color = outlineColor as IColor;
    lineSymbol.Width = 1;
    lineSymbol.Style = esriSimpleLineStyle.esriSLSDash;

    ISimpleFillSymbol fillSymbol = new SimpleFillSymbolClass();
    fillSymbol.Color = fillColor;
    fillSymbol.Outline = lineSymbol;
    fillSymbol.Style = esriSimpleFillStyle.esriSFSHollow;

    ISymbol symbol = fillSymbol as ISymbol;
    symbol.ROP2 = esriRasterOpCode.esriROPNOP;

    return symbol;
}

I already have the answer for this. I was using a normal envelope. The deal is that ESRI already has something to handle that for you and it's contained in INewEnvelopeFeedback class.

http://help.arcgis.com/en/sdk/10.0/arcobjects_net/conceptualhelp/index.html#/Sample_Display_feedback/000100000m9s000000/

Actually, this interface got me where the code snippet in the question was. It stills draws all the rectangles, but I still need a way to clean them.

Any ideas?


texts: Spatially join points to polygons in ArcGIS Desktop?

texts:
I am trying to join a point feature class to a gridded polygon feature class in arcGIS so that all the points are matched to the grid cell that they fall within. I would like to end up with an attributes table that contains the count of points that fall within each grid cell and the average, maximum and stdev value of all the points falling in each grid cell. I have multiple points falling within many of my grid cells.

I thought that this could be achieved by using the 'join' tool in the 'joins & relates' option or the 'spatial join' tool in the 'Analysis' toolbox in arcGIS. Unfortunately, these tools do not give me the output that I am looking for. When I use a point layer with 20024 points as the join_feature, the count column in the output joined feature class exceeds this value (20203). This count column should equal the number of rows in the original point layer.

I had wondered whether the tool was matching points to multiple grid cells if points spanned >1 cell, but apparently points do not have any width associated with them so this cannot be the problem. My grid was created using the 'create fishnet' tool and the boundaries of grid cells do not overlap.


texts: Pivoting Tables and using Frequency in ArcGIS Desktop?

texts:
I have a points layer of specimen samples. They have a field called "RESULT" which is a string field containing the results of a test, "Positive" and "Negative"(you can only see negative in the pic but trust me they're there).

I would like to determine how many positives and negatives are in each county, with each county having a COUNT_POSITIVES and COUNT_NEGATIVES column.

I've tried the Summary Statistics, Table Intersect, Pivot, and Frequency tools in ArcToolBox. The closest I've gotten is the output from the Frequency Tool:

but I'd like to "unstack" the results column so that for each county the table has a column for negatives and positives with their respective counts by county.

You can do this in about 2 seconds in Excel. Is there an easy way to do this in ArcGIS that I'm missing? This is the only response I could find but the case field in summary statistics won't work due my fields being strings (Calculate frequency of a attribute field per polygon)


texts:
I have two polygon shapefiles,each with different data.  One has US counties  It's fields are long, Lat, area_meters, perimeter, state, and county name.

The other shapefile shows wilderness boundaries and has State.  My goal is to figure out every county that these wilderness areas are in.  Some do encompass two or more counties.  I was thinking of doing a table join based on state but that does not seem to work.  It puts each record in only one county.

When I display the counties shapefile and boundary shapefile they show nicely in ArcMap and I could do it manually but their are a lot.

Wanted to know if way to do it with ArcMap or FME?

Any thoughts?


texts: Are there any open source packages to help process network links?

texts:
I'm trying to find an open source package that can process KML-based network links.  I am looking for more than a parser -- more along the lines of setting up timers for interval refreshing, handling <Update> elements, handling the <NetworkLinkControl> element, etc.

I'm am aware that Google Earth API does this for you, but I am in need of an alternate mechanism to mirror that.