texts: Can OpenGeo Suite co-exist with Homebrew-installed PostGIS/PostgreSQL in Yosemite?

texts:
I installed PostgreSQL (9.4.1) and PostGIS (2.1.7) using Homebrew in Yosemite (10.10.3). I'm just starting my way through "PostGIS In Action."  I've been starting and stopping the server in Terminal with pg_ctl, and interfacing with psql, since pgAdmin III wasn't playing nicely with the rest of my system.

I wanted to give the Boundless workshops/tutorials a whirl (as per the recommendations given in response to How do I get started with PostGIS?), but can't quite tell from the other Q&As here if OpenGeo Suite can comfortably coexist with my current PostGIS environment, and whether (as suggested in this response re command line functions in pgAdmin) I'd lose functionality if I operate through OpenGeo rather than the command line.


texts: What is the quickest way to find in which plots a set of points lie?

texts:
I know this sounds easy and lazy to ask, but I'm new to (spatial)database, and I have the feeling my solution is sub-optimal. But for my situation it's important to have good performance.

What do I have now:

Python 3.4 @ Win7-64,
SQLite v3.8.11.1 with Spatialite v4.3.0 via APSW
Two tables; Plots (up to 1m rows) and Points (100 to 200 rows). In real life the points lie on a (almost) straight line.
The position of the points in update constantly, and after each update I need to find out in which plot the points are.
A View containing all plots inside the MBR of the points-table. The view is made by use of a SpatialIndex.

CREATE VIEW 'View' AS SELECT * FROM Plots, idx_Points_geom
WHERE idx_Plots_geom.xmax > (SELECT MIN(xmin) FROM idx_Points_geom)
AND   idx_Plots_geom.xmin < (SELECT MAX(xmax) FROM idx_Points_geom)
AND   idx_Plots_geom.ymax > (SELECT MIN(ymin) FROM idx_Points_geom)
AND   idx_Plots_geom.ymin < (SELECT MAX(ymax) FROM idx_Points_geom)
AND   idx_Plots_geom.pkid == Plots.id

The following query, giving me for each point the corresponding plot.

SELECT * FROM "Points", "View"
WHERE Within("Points"."geom", "view"."geom")

What did I do already:

I reduced the size of the view, by making it from the MBR of subgroups of points (i.e. 1-50, 51-100, 101-150). This reduced the size of the view dramatically, and speeds up the whole process.

What do I want:

A general question I have is: Am I using the right approach to get plots in which the points are? Especially because I know the points are on a line and have a known distance between them.

And more specific: I believer my query is inefficient. Because I know each point is in only one plot, I would like the 'Within-function' to continue to evaluating the next point as soon a it finds the first point-plot match. This could reduce the number of test by 50%.
How do I enforce this? Putting a 'LIMIT 1' at the end of the query will only result in the answer for one point.


texts: How to insert a superscript label in QGIS 2.14 print composer?

texts:
I'm trying to get superscript in my legend label in my print composer and am struggling. I've tried a few of the methods here but don't really understand how to make character map idea work in osx. Specifically how to find the appropriate symbol. I've tried to copy and paste other symbols into the legend title but that hasn't worked

I tried creating what I need in word and copying across but the superscript characters don't reflect in the legend title. So finally I went with a render as HTML text box as a work around which is a little big even though size = 0.

mg/m3

So, could somebody please explain either how to get mapped characters from osx El Capitan into the legend title or how to shrink the html text to useable format?

Apparently you can use basic CSS for the html features. This has a decent workaround for what I need. Anyone know if you can do anything like this in the Legend title?

mg/m3


texts: Placing image between two marker using css in Leaflet?

texts:
In my application where i use leaflet api, i am able to put two gate markers and polyline drawn between them.

Here, i just wanted to find out the way to put the image from css on the polyline which should fit between two markers.I know the image has width and height.

I believe there must be some css trick to do that, although i can't use fill pattern method due to limited understanding about leaflet.

Script(putting two markers and line between them)

var polylineGate2 = new L.Polyline([], { color: 'yellow' });
            function gatePlace2(e) {
                // New marker on coordinate, add it to the map
                new L.Marker(e.latlng, { icon: gateIcon, draggable: false }).addTo(currentFeatureGroup);
                // Add coordinate to the polyline
                polylineGate2.addLatLng(e.latlng).addTo(currentFeatureGroup);
                var decorator = L.polylineDecorator(polylineGate2, {
                    patterns: [
                        // define a pattern of 10px-wide dashes, repeated every 20px on the line
                        { offset: 5, repeat: '20px', symbol: new L.Symbol.dash({ pixelSize: 4 }) }
                    ]
                }).addTo(currentFeatureGroup);
                map.off('click', gatePlace2);
                polylineGate2 = new L.Polyline([], { color: 'pink' });
            }


texts:
I have a stream network and outlets as shown below.

It is overlapped with outlet cells raster which I can use to get sub-watersheds but when I add an outlet (red color cell in current example) then I have to divide stream segment in order to have correct "From Node" to "To Node" relationships with use of Stream to Feature (Spatial Analyst) tool in ArcGIS 10.1

Any suggestions on how to correctly divide stream segment to have correct topological relationships?


texts:
Openlayers doesn't seem to support the URL that Carto prefers for getting UTF Grids.

import TileLayer from 'ol/layer/Tile';
import UTFGrid from 'ol/source/UTFGrid';

const newUTFLayer = new TileLayer({
  source: new UTFGrid({
    url: 'https://PATH.TO.CARTO/LAYERGROUP_ID/0/{z}/{x}/{y}.grid.json';
  })
});

I get the error

"errors":["Invalid coordinates: NaN/NaN/NaN"],"errors_with_context":[{"type":"unknown","message":"Invalid coordinates: NaN/NaN/NaN"}]}

Clearly Openlayers isn't replacing {z}/{x}/{y} with coordinates.

In the official example of using a UTF grid, they seem to be requesting some kind of config from Mapbox to get a UTF grid, but Carto seems to use a different mechanism?

-

Edit: I've found some additional information from people having similar problems to me:

Openlayers GitHub issue
Someone who made a work-around

It seems the issue is that OpenLayers' UTFGrid source is trying to request a TileJSON file (which I believe is a metadata file?) rather than the individual UTF grid tiles themselves.

Frustrating that there needs to be a server-side component to solve it using the aforementioned work-around.

Edit 2: The OpenLayers UTFGrid spec seems to allow you to set a tileJson object, can I use this to get around having to request a TileJson file?


texts: Creating individual points from count field in point feature class using ArcGIS for Desktop?

texts:
I have a point feature class with a count field associated with each of the points.

I would like to create a new feature class where each of the points are duplicated by the value of the count field, so that rather than being an aggregate count, they represent individual events (I need to run a Monte Carlo on them so individual events are required).

Does anyone have any idea how to do this in ArcGIS for Desktop?


texts: Dynamic GeoServer SLD styling template to use on Multiple Polygon shapefiles

texts:
I'm trying to create 1 SLD template for the Geoserver application on any point shapefile I upload.
For example, I want a styling that will automatically color the attributes for the selected layer. If layer 1 has 5 unique records the styling should pull through this.

And if I upload a different point shapefile(another layer) and select the same SLD template with more unique attributes it should automatically pull it through.

Please is this possible??


texts: How to count pixels of a raster in polygons using QGIS?

texts: I need to compute some statistics of the pixels of a given raster that are coincident to the polygons of a vector file. Is there a simple and effective way to do this in Quantum GIS or similar?


texts: ModelBuilder - use in line variable in expression

texts:
I have an in line variable of two fields (field1_field2) in modelbuilder.  Using the make feature layer tool, I would like to make an expression that queries field2.  Can an expression in the make feature layer be used to parse and query the in line variable.  I've tried the query below and some variations of it but not having any luck.

EVENTKEY = %Value%.split("_")[-1]

EVENTKEY is the field I'm querying and %Value% is the in line variable consisting of field1_field2


texts: Quick add vertices when creating feature (QGIS 3.2)

texts:
I am creating a new feature (polygon) in QGIS 3.2.2 and would like to easily add the vertices as I delineate a boundary.  Currently, I can do this one by one as the cursor snaps to the vertex of the features in the layer I am referencing.  But there are many vertices and (ideally) I would like to quickly go over the vertices and if they are within snapping distance they are automatically added to the new feature.

Is this possible?


texts: Why is Mosaic or Mosaic To New Image causing a 1-pixel gap between adjacent images?

texts:
I have a series of adjacent PNG images, and I'm trying to merge/mosaic them into a single seamless image. The images have clean "cuts" between them, with no overlaps.

However, when I use either the Mosaic or Mosaic To New Image tools, a 1-pixel gap has been inserted between some sets of images (but not between all sets of images):

Identifying within the gap shows Nodata. Here is the line in Python which creates the Mosaic:

arcpy.MosaicToNewRaster_management(list_of_images,outputImage,"#","8_BIT_UNSIGNED","#",1,"LAST","FIRST")

Does anyone have any suggestions as to what I need to check/change to avoid this gap?

Thanks


texts: Options for getting data from SQL Server into Leaflet app?

texts:
I have polygon data in SQL Server 2012 that I need to get onto a Leaflet application - we're not talking thousands of polygons, but simplified U.S. counties where only one state will be displayed at a time. The database gets updated every 5 minutes and the map needs to reflect that. I know there are several ways to get the data to geojson/json and put that onto the map (such as the jQuery $.getJSON method). We've thought about creating a web service to pull the data and do the conversion to geojson, then get it onto the map with Ajax. The concern there is that this sounds slow and CPU-intensive. Thoughts on that?

I'd lke to know if there are any methods that would allow me to more directly query SQL Server (from within the application perhaps) and get the polygon data into my Leaflet app, keeping in mind that it needs to be refreshed every 5 minutes.


texts: Export a shapefile attribute for points contained in that file for very large datasets with QGIS?

texts:
I've imported a multi-shape shapefile with ~300.000 shapes (445MB shapefile) and a CSV file with ~2.500.000 point coordinates into QGIS (100MB CSV file). I converted the CSV file to shapefile as well, resulting in a (~150MB shapefile)

My multi-shape file has a field named CD_GEOCODI with a unique code for each shape within the file. My point data file has 3 fields: lat, long and var1.

It is easy to visually spot points within the specific shapes. I would like to add the CD_GEOCODI variable from my multi-shape file to the export of the points point_data, so the resulting data will now have 4 fields: lat,long,var1,CD_GEOCODI in which the last field, CD_GEOCODI is the CD_GEOCODI variable from the shape which the point is within.

I've tried four methods:

Processing > Toolbox > Saga > Vector Point Tools > Add polygon attributes to points
Vector > Data Management Tools > Join by location
Install Point Sampling Tool plugin > Plugins > Analyses > Point Sampling Tool
Install MMQGIS plugin > MMQGIS > Combine > Spatial Join

All of them seem to work with small datasets but all of them freeze and make unpratical slow progress (based on the velocity of temp files growth in size, since QGIS is freezed) with my huge dataset.

My environment is: Win10 x64, qGIS 2.18.9, 6th generation i7-6700HQ, 8GB RAM, 1TB HD 5200

How can I achieve this goal with with those huge datasets in QGIS?


texts:
When deploying ArcGIS Portal, we have the ability to add a draw widget to our map app using Web AppBuilder.  However, when users draw on the map, and they try to share this markup with other users, the drawings do not show up.

Is there a way for our users to markup the map using the draw widget, or similar functionality, and be able to share that with other users?


texts: Which projection is suitable for the measurements of this great circle?

texts:
I have drawn a circle on google earth pro in such a way that its radius starts from earth center 0,0 (lat/long) and ends at 0,90 (lat/long). This circle divide the earth into two equall halfs as shown in image. It is in line shape of length 39,999 km.

Then I have converted this circle in ArcGIS format and visualize there, it seem to be in a rectangle shape, obviously due to projection . When I select web mercator or other projected coordinate system of world they all give too much length as compared to its actual length 39,999 km.

My questions are:

What method or suitable projection should I adopt in order to calculate its length on ArcMap?
Length of this circle on google earth is 39,999 km which is different from the earth circumference which is about 40,075 km, why?


texts:
I have two line dataset that represent the same roads:

First one has very precise geometry
Second one has an approximative geometry, but very interesting attribute.

What I want to achieve is to have the geometry of the first one and the attribute of the second one. I wish to find an automated way to do it. There are many roads and it would takes decades to do it by hand.

I prepared some examples:

At first I have this data:

And I want to transfer the good attributes to the good geometry:

To obtain something like that:

I had a look at the transfer attribute in ArcMap, however, it seems that it can't split the destination geometry in order to represent all the originals attributes values.
I probably need something that first conflate the bad geometry on the good one and then split it according to the attributes.


texts: Measure distances between specific features in QGIS for Desktop

texts:
I'm completely new using PyQGIS and developing plugins.

Currently I have a layer with the X and Y coordinates of electric distribution poles, where each pole is a feature. What I need to do is measure the distance between all the features in the map, but it must be done pole by pole. By example, take the first pole and get all the distances between the first pole and all the others, then, take the second pole and repeat until the end of the dataset.

Any ideas or references for doing this?

As a reference, I was doing the same thing in Matlab with a dataset of 20 features with something like this:

for i=1:20
  for j=1:20
      if i~=j
        distance = euclidean(P1(i), P1(j));

Where P1(i) and P1(j) are the X-Y coordinates of two different points (Each one Written as a complex number) and "euclidean" is a function that calculates the distance between these two points.


texts: Choosing projection / CRS for short-distance based analysis?

texts:
I'm having a bit of meltdown about projections and coordinate reference systems, and was hoping for a few pointers.

I have 3 GIS projects, at 3 sites in Europe (the largest of which measures about 10x30km) all using similar data.

My data layers include point features (from GPS, in lat lon), derived buffers around these points, manually created polygon features (representing landuse) which I want to convert to raster for analysis and some externally derived rasters reprenting imagery, elevation, etc. I will also have some GPS points representing animal locations.

I'm currently using ArcMAP, but also trying to migrate to R and/or QGIS (in preparation for losing my Arc license). I guess my problem isn't software-specific anyway.

I want to measure distance accurately (e.g. 1000m buffer around point features; nearest-neighbour distance) as well as area (point density; percent cover of different landuse within a buffer). Currently I'm working in WGS84 as this is the system my GPS points come in. I'm wondering whether this is the best projection (though as I understand it, it is not a true projection) to be using for these small-scale analyses. I've read that it might be better to use a site-specific transverse Mecator centred on each site. Will it make much difference? (I also want to make some nice figures for presentation, so it would be good if e.g. 1km buffers were actually round rather than eliptical - maybe this is a different point though)

I realise this sort of question has been asked a million times before. I'm not being lazy, but I currently have massively restricted web access so can't afford to figure out the answer myself.


texts: Using the Sketch Widget for ArcGIS JavaScript API, can you manually toggle between the reshape and transform?

texts:
I am using the sketch widget for ArcGIS JavaScript: https://developers.arcgis.com/javascript/latest/api-reference/esri-widgets-Sketch.html

I am building a custom toolbar with additional tools and would like to emulate the behavior of the "reshape" and "transform" buttons in the default toolbar. I don't see a way to do this through the Sketch Widget API.

Is this possible?


texts: How to calculate distances in a point sequence?

texts:
I'm looking for help calculating the distance between a sequence of points that are in a single Shapefile in QGIS.  Below is what my data look like and a blank distance column I've added to show how I'd like the distance to look.  I want to know the distance between point 1 and 2, 2 and 3, etc. I'd like distance to be in meters or Km, but currently my Shapefile is in a projection whose units are decimal degrees.

ID  LAT         LON         TIME        DISTANCE
1   10.08527    124.59833   21:24:37    0
2   10.08523    124.59830   21:25:07    ?
3   10.08526    124.59832   21:25:37    ?
4   10.08526    124.59831   21:26:07    ?

A number of people have asked similar questions, but none quite get at what I want to do.  This post is close, but this is in PostGIS, not QGIS  Calculate distances between series of points in postgis

This post got me part of the way there, but since I'm new to QGIS the answer doesn't provide enough detail for me.  For example, once I've installed the GRASS plugin, I think I need to save/import the Shapefile with my series of GPS points into GRASS so I can use the v.distance module, but I don't know how to do that.
 QGIS calculate distance of point along a line

Is the GRASS v.distance module the only way to go?  Or is there a more straight forward way?  If v.distance is the only way would someone be able to point me to or explain more step-by-step how to do this?


texts: QGIS character spacing in label inserted in print layout

texts: I normally insert labels directly through Layer Properties, which allows one to vary character spacing. However, I need to insert labels in Print Layout for a particular project and cannot find a way to manipulate character spacing in Print Layout. Can it be done?


texts: ArcGIS Web AppBuilder - Links to related sites within map

texts:
I am creating a map application using the Web AppBuilder in ArcGIS Online. The underlying map was created in ArcMap.
I know it is possible create a link to a webpage in the pop-up boxes, but is it possible to create a url that would take the user to a related point in the same map?
Example: I have created a map of historic locations within a town. However, the current school used to be located on a different site and I want users to be able to easily navigate from one location to the other.  I have point data for both the current and historic location.  If I were to create a section called 'Related locations' within the pop-up box, could I generate a link within the current school site's pop-up so that when someone clicks on it, it would take them to the point for school's previous location?


texts: Why does the Flow Direction tool output a raster with only one value?

texts:
I have a DEM with values that range from 0 to 1040. I'd like to create a layer with the flow directions by means of the Flow Direction tool in the Hydrology section of the ArcToolbox. However, when I try to do so, it creates a raster with only one value: zero. I'd like to have a raster with the values 1, 2, 4, 8, 16, 32, 64, and 128.
Does anyone know what I'm doing wrong, and how my mistakes can be repaired? Here are some screenshots of what I tried:


texts: Rating Randomness of Electromagnetic Conductivity (EC) Data Collected?

texts: So we have electromagnetic conductivity (EC) data that we have been collecting over some forest ground which logs an EC value every couple of seconds.  We like to see the data have nice distinct patterns to it so when you add symbology to the data with Class Breaks you get are nice clear and distinct areas of varying EC.  Some times we get data that is much more random in natural and doesn't have the nice distinct zones and areas of EC but rather has more of just random values across the area we collected data (i.e. big swings in EC values across small areas in spatially).  I would like to find a way to statistically rate how random the data is.  If there are nice distinct zones, it should give a rating of cleaner data.  If the data is much more random without good distinct areas of EC and lack of more consistent changes across the field, then it would give a low rating of clean data and high randomness of data.  I know my terminology is bad here but is there something in ArcMap I could use to measure the randomness of the collected data?


texts: can't initialize projection information in GRASS from R

texts:
I want to perform a small calculation in GRASS using the r.watershed-function. I convert my raster DEM to a SpatialGridDataFrame to import it into GRASS. However when initializing GRASS the projection information of the DEM-SpatialGridDataFrame is not kept/imported (i.e. projection(streams_g_buf) # [1] NA). Can somebody help me to spot the problem?

require(raster)
require(rgrass7)

## get Malawi DEM
mwi_dem <- getData("alt", country="MWI", mask=TRUE)

# reproject raster to UTM 36 S
mwidem_trans <- projectRaster(mwi_dem, crs=CRS("+init=epsg:32736"))
# as SpatialGridDataFrame for GRASS
mwidem_trans2 <- as(mwidem_trans, "SpatialGridDataFrame")

projection(mwidem_trans2) # [1] "+init=epsg:32736 +proj=utm +zone=36 +south +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"

## Initialize GRASS
initGRASS("/usr/lib/grass72",
          home=tempdir(), override=TRUE, SG=mwidem_trans2)

## Write SpatialGridDataFrame
writeRAST(mwidem_trans2, "DEM", zcol = 1, flags="overwrite")

# GRASS analyses
execGRASS("r.watershed", flags="overwrite",
          parameters=list(elevation="DEM", stream="streams_g", threshold=30))
execGRASS("r.buffer", flags="overwrite",
          parameters=list(input="streams_g", output="streams_g_buf",  distances=200))

## Read back into R
streams_g_buf <- raster(readRAST("streams_g_buf"))

projection(streams_g_buf) # [1] NA


texts: Processing multiple files simultanously using python

texts:
I have 11000 point shapefiles to run Near tool on against another point shapefile.
I am looking for a way to process 10 files at the same time.

Is there a way to process so many files using Python and arcpy in a resonable time considering that each file takes 5 min on average?


texts: Iterating over a raster using rasterio, and then saving to a shapefile

texts:
I'm trying to create a script to polygonalize a raster using rasterio. I want to only polygonalize pixel values that are not null.

From this SO answer, I derived this function:

def create_geometry(raster_input):
"""
Creates a shapely polygon from a raster image
:param raster_input: raster path
:return: shapely geometry
"""
mask = None
results = None
pixels = []

# with rasterio.drivers():
with rasterio.open(raster_input) as src:
    image = src.read(1)  # first band
    results = (
        {'properties': {'raster_val': v}, 'geometry': s}
        for i, (s, v)
        in enumerate(
        shapes(image, mask=mask, transform=src.affine)))
return results

That function is called by this function, which then saves to shapefile using the polygonalize_geometry() function:

def shape_builder(raster_input, output_path):
"Builds a polygon out of the pixel polygons"

# polygons = []
polygons = create_geometry(raster_input)

for polygon in polygons:
    polygon = shape(polygon[0]['geometry'])
    print(polygon)
    polygonalize_geometry(polygon, output_path)

I then convert it into shapely geometry, and then write a shapefile. However, the shapefile only contains one feature that encompasses one pixel (shown below). I believe that it is only operating on the first pixel it finds. The SO answer I linked above says that the function is a generator, but I haven't ever been able to figure out how to use them. I'm pretty sure that's the issue.

So the script is somewhat working.

The entire script is here.

Where did I go wrong?


texts: Multiple records for a single point feature in a shapefile

texts:
Is it possible to have a single point feature associated with multiple records in an attribute table?

For example, if I had a shapefile of type POINT containing 100 points indicating fire hydrants but I wanted the attribute table to contain records of all use of the fire hydrants, such that when I clicked on say fire hydrant #80, it selects all the records in the attribute table with use of fire hydrant #80.

I've not tried this yet, but as I understand it, the relationship is 1:1 (1 feature to 1 record). In the given example, if I had 100 features but the attribute table had 500 records, I believe that would be an invalid shapefile correct? The only way to resolve this is to "overlap" features, such that when it appears I am selecting a single feature, there are really, for example, 6 features representing fire hydrant #80 on top of one another. If possible, I would like to not do that kind of workaround.


texts: Getting nearest polygons in both sides of polyline with ArcPy?

texts:
I have polyline shapes representing roads and another shape layer with polygons that are representing buildings. I have to find some way with ArcPy to measure the distance from a point within the polyline road,usually in the middle, to the nearest buildings at both sides of the road. Graphically it would be something like this

In short, I'd need some way to use a kind of Near tool but with the functionality about measuring in the direction I want to.


texts: I need to create a shapefile with crossed streets. I can't split the lines to number each segment with the house number ranges.


texts: Is it possible to configure WFS layer to use GET requests?

texts:
Accordingly specification manipulation with WFS services can be perform using POST or GET (for example, http://demo.opengeo.org/geoserver/wfs?service=WFS&version=1.1.0&request=GetFeature&typeName=topp:states) requests. Is it possible to configure OpenLayers for using only GET requests (by default POST is used). Now my code looks like:

wfs = new OpenLayers.Layer.Vector('Example of WFS layer', {
    strategies: [new OpenLayers.Strategy.BBOX()],
    visibility: true,
    protocol: new OpenLayers.Protocol.WFS(wfs_property)
});


texts:
I'm wondering if anyone has any insight into the potential to auto-calculate the value of a field, either area or distance (ideally) or a formula (based on the values of say 2 other fields, or an existing field and another value).

Perhaps this is a plugin that presents a dialogue asking to what layer you want the auto-calculation to apply, and what formula or combination of fields you would like included in the calculation.

These tools are of great use to planners who are using a GIS as a drawing/design tool.

Thoughts?


texts: Using SLD to style according table column value

texts: How to style various points with different image according to value(text) in the column of table(layer)?


texts:
I have an ArcScene *.spk file, but no longer have ArcScene and only a Basic License for ArcGIS.

I have access to QGIS.

Is there any way to open / import the *.spk file to extract some of the features out of the file?


texts: Storing the GML application schema in postgres

texts:
I have to store the application schema for a GML file which is in Xsd format in PostgreSQL table.

I want to use the ISO standards for that.

From what I have searched I think ISO 19109 has to be used.

Can anyone suggest some thing


texts: Selecting point to obtain water basin leading to it

texts:
For the past few days, I have been working on building a watershed model including drainage basins and find that I am not looking for a global watershed and basin distribution, I am looking to generate those that lead to one specific point.

The idea is to select a point (either by click or creating a point feature or even a polygon feature) and defining the drainage basins as well as watersheds what will inevitable pass by that point.

Using QGIS, I'm wondering if there is a way to create such a tool or if there is a process/algorithm that will give me such a result.

My raw data is a DEM of the area of interest, treated to fill the voids and all.


texts: Extracting the geometry of the parts not intersecting from 2 shapes overlapping each other in postgresql?

texts: This seems simple enough, yet i do not find a clear awnser on that peculiar matter. What i am trying to do is to extract the parts of a shape that doesn't intersect with another when the 2 of them are overlapping. In short, i am looking for something like ST_intersection but in reverse.


texts: Filtering with setSubsetString using PyQGIS

texts:
I would like to create a layer with a filter on one of my columns (GAMME_DE_S). I saw that I can use setSubsetString but no way to use it.
uri = r"C:\impression_plan\FicheCana.shp"
layer = QgsVectorLayer(uri, "", "ogr")
resultat_layer=layer

layer.setSubsetString("GAMME_DE_S=VSR")


texts: Creating Voronoi for points that consider separated land or gaps in land

texts:
As can be seen in the picture, I've created a Voronoi for points in the UK. However, the Voronoi does not take into consideration the split in land in the south-west region of the country. Therefore, it does accurately display the areas of land nearest to the points. How do I get past this?
I have

UK Boundary vector
point data layer (text converted file, not shp)
Voronoi layer for points


texts: Customized Web-based for data cube visualization

texts:
I am stepping into web mapping for the first time. I wanted to develop software for animating multi-spectral, multi-temporal data (in the form of an orthorectified data cube).

User should be able to select parameters such as the band to visualize, the time-range or the time-step and animation speed parameters e.t.c.

The software should be web-based. I have data in the form of GeoTIFF and JP2 images Raster data from EO1 Hyperion satellite.

The result should look something like this https://glovis.usgs.gov/app.

I am not sure if we are using the right tech Stack:

Rasterio(python lib for accessing and manipulating spacial data, I have some experience with this lib)
OpenStreetMap as the base layer
React.js to provide frontend

I need help with the tech stack. I heard we can use qgis2web plugin to export my maps to Leaflet and GeoJSON, should I use the plugin instead?

This is the first time I am working with spacial data, so I am not quite sure on how to render TIFF images on the base map.


texts: OpenLayers 3 VectorTile Layer Reload style on zoom

texts:
I have a vectorlayer in my map and would like to change the style of features when the zoom level changes. Is this possible?

var vlayer = new ol.layer.VectorTile({
    renderMode: 'vector',
    preload: Infinity,
    source: new ol.source.VectorTile({
        attributions: '&copy; ' + 'OpenStreetMap contributors',
        format: new ol.format.MVT(),
        tileGrid: new ol.tilegrid.TileGrid({
            extent: ol.proj.get('EPSG:3857').getExtent(),
            resolutions: resolutions
        }),
        tilePixelRatio: 16,
        tileUrlFunction: tileUrlFunction,
        projection: 'EPSG:3857'
    })
    , style: stylefunction()
});

var map = new ol.Map({
    target: 'map',
    layers: [
        vlayer
    ],
    view: new ol.View({
        center: [0,0],
        zoom: 15
    })
});


texts: QGIS Attribute table - How to merge the values of $area field of several elements by the same attribute (point ID)?

texts:
I have a land use vector layer with polygons representing forest, farming, urban aggregates, etc. Each polygon has a field with its area. I was trying to intersect it with buffers of different diameters (50m, 100m, 500m) around points and then obtain the area of each land use inside that buffer in QGIS 3.2.2.

The problem is the intersect tool in Vector -> Geoprocessing Tools does not change the field Area on the resulting intersected layer, rather it persists with the value of the entire polygon, even though only a fraction of it is inside the buffer.

I tried several solutions but none worked to solve this, so I decided on a workaround. On the intersected layer, I added a field $area which shows the area in map units of each cut polygon and then another column with the percentage, which is the map units of that cut polygon dividing by the overall area of the circular buffer. I won't get the actual area of each land use but at least a percentage of occupation which is good enough.

Now the problem is that for each point there are several fragments of different polygons pertaining to the same land use. I wish to combine for each point the values of the same land use.

On this screenshot, you can see highlighted there are 3 pieces of land use "forest" on the 100m buffer around the same point 24. I wish to merge and sum these percentages and end up with a total percentage for each land use for each point.

Better yet if it were possible to create another table where the elements (lines) are merged by their land use, thus avoiding repeated values and providing a single element and respective total percentage value per land use type per point. On the example above, all three highlighted elements of land use Áreas florestais of point 24 would be collapsed into a single element with a single value for Point = 24, land use = Áreas florestais and Perc = sum(0.068+0.037+0.147)


texts: Zooming In and Out using ArcPy with ArcGIS Pro?

texts:
I have code that zooms a layout to a bookmark and then exports that layout to PDF. This step helps me produce maps on a large scale because I have about 40 bookmarks. However, if someone asks me to create the maps again but zoomed in or out more than that creates a problem because then I have to change all the bookmarks. I'm hoping there is some kind of code that would zoom to the bookmark and then I could tell it to zoom in or out a certain amount from there.

Here is the code I have currently:

import arcpy

p = arcpy.mp.ArcGISProject(r'C:\arcGIS_Shared\Python\CenterHeatMaps.aprx')
lyt = p.listLayouts("Layout_King")[0]
mf = lyt.listElements("MAPFRAME_ELEMENT", "Map*")[0]
bkmks = mf.map.listBookmarks()
m = p.listMaps()[0]
OnLyrList = m.listLayers('Topographic')
OnNxtList = m.listLayers('Open_GoHealth_Centers')
CenterList = m.listLayers("CASum*")

#Loop through center layers that begin with 'SumWithin'
for center in CenterList:

   #Loop through all bookmarks until finding a match for the active center
    for bkmk in bkmks:
        if bkmk.name == center.name[19:]:

            #Turn off all layers
               lyrList = m.listLayers()
               for lyr in lyrList:
                   lyr.visible = False

            #Turn on desired layers
               for lyr in OnLyrList:
                   lyr.visible = True
               for lyr in OnNxtList:
                   lyr.visible = True
               center.visible = True

            #Zoom to active bookmark
               mf.zoomToBookmark(bkmk)

            #Change Layout Title
               for elm in lyt.listElements("TEXT_ELEMENT"):
                   elm.text = "Annual_CA Summary by Block 2017 - " + center.name[19:]

            #Export to PDF
               lyt.exportToPDF(r"C:\arcGIS_Shared\Python\Export\CenterHeatMaps_CA" + "\\Annual_CA" + center.name[15:] + ".pdf")
               print(center.name[15:])

Is what I'm asking possible?
I haven't been able to find anything online.


texts: Generating Dynamic Legend with OpenLayers and GeoServer?

texts:
I am using OpenLayers 3 and GeoServer to create a map with multiple layers. I am able to generate legend for each layer from GeoServer using  GetLegendGraphic and display them using the following:

I tried to control the legend visibility based on the visibility of its corresponding layer. I searched various sources and found How to add a legend and toggle the legend with a layer in Openlayers? which did not work for my case.

Now I have the following problem:

When a Map layer is made invisible, the legend is still displayed. How, can I control the legend visibility dynamic i.e. make the legend visible only when its corresponding layer is selected and make the legend invisible when its corresponding layer is invisible.

Is there a working example?


texts: ArcPy Select Layer by Attribute for Date and Time

texts:
I have a list of trips with a field, "Trip Start Timestamp" thats type is date where it appears on the attribute table as M/D/YYYY HH:MM:SS AM/PM
I need group all trips with the same trip start timestamp, for example 06/01/2019 12:00:00 AM. I am running the arcpy select layer by attributes as shown below:

arcpy.SelectLayerByAttribute_management('TestSample_10_Pickups', 'NEW_SELECTION', 'Trip_Start_Timestamp = "2019-06-01 12:00:00 AM"')
Traceback (most recent call last):
  File "<string>", line 1, in <module>
  File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\management.py", line 8123, in SelectLayerByAttribute
    raise e
  File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\management.py", line 8120, in SelectLayerByAttribute
    retval = convertArcObjectToPythonObject(gp.SelectLayerByAttribute_management(*gp_fixargs((in_layer_or_view, selection_type, where_clause, invert_where_clause), True)))
  File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\geoprocessing\_base.py", line 506, in <lambda>
    return lambda *args: val(*gp_fixargs(args, True))
arcgisscripting.ExecuteError: ERROR 000358: Invalid expression
Failed to execute (SelectLayerByAttribute).

I am completely lost as to what this error is and if it has to do with the formatting of the date and time.
Does anyone know how to select layer attributes by date AND time?


texts: When do I need to use a proxy with OpenLayers2?

texts:
When exactly do I need to use a proxy with openlayers? I've done a lot of reading on this but can't seem to find a definitive answer. Wondering if this is the route of my problems.

For instance, if I have one server running Geoserver (WFS) on port 8080 and Apache (serving OpenLayers) on port 80 would I need to use a proxy?

If the servers are physically separate would I need to use proxy?


texts: Layered PDF from QGIS - Layers can't be turned on/off

texts:
I was following this guide on producing a layered PDF from QGIS, but somewhere along the way I seem to miss something. Exporting layers via the composer into a svg works fine, also loading the svg into scribus, splitting the groups and shifting the polygons to a new layer works.

The exported PDF even contains layers, but these don't disappear when being turned off.
For reasons of reproducibility, here is also the scribus file.

Perhaps my question should be attached to Can QGIS preserve layers when exporting to PDF? from 2015.

Trying this under Win7, using QGis 2.18.11, Scribus 1.5.2 and Acrobat Reader XI.


texts: Can't import OSM data to PostGIS, permission denied

texts:
I use Postgres 9.1 on Xubuntu 11.10. I've created a db 'aus' and user 'aus', granted ALL on the db to the user. Now I try to import the data, and it outputs this error:

$ osm2pgsql --create --database aus --username aus --prefix planet australia.osm -S /usr/share/osm2pgsql/default.style
osm2pgsql SVN version 0.70.5

Using projection SRS 900913 (Spherical Mercator)
Setting up table: planet_point
NOTICE:  table "planet_point" does not exist, skipping
NOTICE:  table "planet_point_tmp" does not exist, skipping
SELECT AddGeometryColumn('planet_point', 'way', 900913, 'POINT', 2 );
failed: ERROR:  permission denied for relation geometry_columns
CONTEXT:  SQL statement "DELETE FROM geometry_columns WHERE
    f_table_catalog = '' AND f_table_schema = 'public' AND f_table_name = 'planet_point' AND f_geometry_column = 'way'"
PL/pgSQL function "addgeometrycolumn" line 133 at EXECUTE statement
SQL statement "SELECT AddGeometryColumn('','',$1,$2,$3,$4,$5)"
PL/pgSQL function "addgeometrycolumn" line 5 at SQL statement

What's wrong with my database?

edit: the solution:

$ sudo su -l postgres
$ psql -d aus
# grant all on database aus to aus;
# alter database aus owner to aus;
^D
^D
$ osm2pgsql ...
here we go...


texts: Display an moveable image layer in openlayers 3

texts:
so im trying since a month now to solve the above mentioned problem.

The idea is to display an image of a city map and drag it with the mouse or move it with a slider to the desired location.

That's my code so far:

var img_extent = null;
var overlay_image_layer;
var img_view_extent;
var img_source;
var map_view_extent;
var multiplier = 100;
var img_heigth;
var img_width;

var upload = document.getElementById('imageUpload');
upload.addEventListener('change', dateiupload);

function dateiupload(evt) {
    var _URL = window.URL || window.webkitURL;
    var dateien = evt.target.files;
    var uploadDatei = dateien[0];
    var reader = new FileReader();

    var senddata = new Object();
    senddata.name = uploadDatei.name;
    senddata.date = uploadDatei.lastModified;
    senddata.size = uploadDatei.size;
    senddata.type = uploadDatei.type;

    reader.onload = function (theFileData) {
        senddata.fileData = theFileData.target.result;
    };

    objectURL = _URL.createObjectURL(uploadDatei);

    var image_html = new Image();
    image_html.src = objectURL;

    img_width = 0;
    img_heigth = 0;

    image_html.onload = function () {
        img_heigth = image_html.height;
        img_width = image_html.width;

        // get base image position
        // [minx, miny, maxx, maxy]
        map_view_extent = map.getView().calculateExtent(map.getSize());
        img_view_extent = [map_view_extent[0], map_view_extent[1], map_view_extent[0] + img_width * multiplier, map_view_extent[1] + img_heigth * multiplier];
        console.log(img_view_extent);

        img_source = new ol.source.ImageStatic({
            //projection: ol.proj.get('EPSG:3857'),
            url: objectURL,
            imageExtent: img_view_extent
        });

        overlay_image_layer =
            new ol.layer.Image({
                zIndex: 1,
                source: img_source
            });
        overlay_image_layer.addEventListener('change', function () {
            console.log('changed');
            map.render();
        });

        map.addLayer(overlay_image_layer);
    };
}

function change_on_x(event) {
    var temp = [img_view_extent[0] + (parseInt(event.value)), img_view_extent[1], img_view_extent[2] + (parseInt(event.value)), img_view_extent[3]]
    console.log(temp);
    img_source.set('imageExtent', temp);
    overlay_image_layer.setMap(map);
    overlay_image_layer.changed();

}

function change_on_x is called by an <input type="range">.

Can someone help me out here? :)

My map look like this:

var map =
new ol.Map({
    target: 'map',
    view: new ol.View({
        center: [1492062.75385, 6895234.55805],
        zoom: 11
    }),
    layers: [
        //base_map_group_layer,
        //overlay_image_layer,
        OSM
    ],
    interactions: ol.interaction.defaults().extend([
        new ol.interaction.DragRotateAndZoom()
    ])
});


texts: Get Road name from coordinates to be used in Overpass api

texts:
I have following Overpass api quer that takes the name of the road and give back entire road coordinates

http://overpass-api.de/api/interpreter?data=[out:json][timeout:25];(node["name"="%@"];way["name"="%@"];relation["name"="%@"];);(node["ref"="%@"];way["ref"="%@"];relation["ref"="%@"];);out body;>;out skel qt;

I am using following query to get the name of road from coordinates

http://nominatim.openstreetmap.org/reverse?format=json&lat=%f&lon=%f&zoom=18&addressdetails=1

 NSDictionary *dict = [response JSONValue];
 //NSString *placeId = [dict objectForKey:@"place_id"];
 //NSString *osmId   = [dict objectForKey:@"osm_id"];

 NSDictionary *roadDict = [dict objectForKey:@"address"];
 NSString *roadName = [roadDict objectForKey:@"road"];

Currently i have an example of Roadname @"John F Foran Freeway" using coordinates

[[CLLocation alloc] initWithLatitude:37.731267 longitude:-122.423580];

When i use this roadname in my overpass api. I get no results. Why its not detecting the road? Sometimes on few roadnames it work but mostly it don't recognize the road name. What is wrong?

http://overpass-api.de/api/interpreter?data=[out:json][timeout:25];(node["name"="John F Foran Freeway"];way["name"="John F Foran Freeway"];relation["name"="John F Foran Freeway"];);(node["ref"="John F Foran Freeway"];way["ref"="John F Foran Freeway"];relation["ref"="John F Foran Freeway"];);out body;>;out skel qt;


texts: Measure Tool Error on GIS Flex Application

texts:
On my organization's Flex Application, we use a Measure Tool that is provided by the Geometry Service from this website: https://tasks.arcgisonline.com/arcgis/rest/services/Geometry/GeometryServer.

After drawing a polygon for square footage or acreage measurement, and/or a line distance measurement, this error will pop up:

[FaultEvent fault=[RPC Fault faultString="HTTP request error" faultCode="Server.Error.Request" faultDetail="Error: [IOErrorEvent type="ioError" bubbles=false cancelable=false eventPhase=2 text="Error #2032"]. URL: http://tasks.arcgisonline.com/ArcGIS/rest/services/Geometry/GeometryServer/lengths?polylines=%5B%7B%22paths%22%3A%5B%5B%5B1586837%2E9747440305%2C310141%2E3985628667%5D%2C%5B1617457%2E4437312013%2C287527%2E401429212%5D%5D%5D%7D%5D&f=json&sr=102685"] messageId="9BFCC0A8-2A2C-0213-D8AD-02E8928637F5" type="fault" bubbles=false cancelable=true eventPhase=2]

Any thoughts on why this error occurs?  Could the Geometry Service be getting too many requests?


texts: Why am I getting this error when using count points in polygon?

texts:
I'm using QGIS 3.16.1 and trying to work with the 'count points in polygon' feature.  When I try to run it I get this error:
*QGIS version: 3.16.1-Hannover
QGIS code revision: b381a90dca
Qt version: 5.11.2
GDAL version: 3.1.4
GEOS version: 3.8.1-CAPI-1.13.3
PROJ version: Rel. 6.3.2, May 1st, 2020
Processing algorithm…
Algorithm 'Count points in polygon' starting…
Input parameters:
{ 'CLASSFIELD' : '', 'FIELD' : 'NUMPOINTS', 'OUTPUT' : 'TEMPORARY_OUTPUT', 'POINTS' : 'C:\\Users\\EVG\\Documents\\GIS FILES\\School\\Final Project\\Coyote_Final.shp|layername=Coyote_Final', 'POLYGONS' : 'C:/Users/EVG/Documents/GIS FILES/School/Final Project/shp_bdry_mn_city_township_unorg (4)/city_township_unorg.shp', 'WEIGHT' : '' }
Feature (15) from “city_township_unorg” has invalid geometry. Please fix the geometry or change the Processing setting to the “Ignore invalid input features” option.
Execution failed after 0.03 seconds
Loading resulting layers
Algorithm 'Count points in polygon' finished*

I'm not sure what the invalid geometry refers to.  I zoomed to the attribute on the map and it looked fine to me.  Then I tried deleting feature 15 from the attribute table since I didn't need the data anyways, ran it again, and got the same error.  How do I change the processing settings to "ignore invalid input features?"  I'm not seeing a setting to change it in the 'count points in polygon' menu.


texts: Output PostGIS geometry data as centroid or polygon depending on zoom level (GeoServer)

texts: This is what I want: when I zoom out, cities of the country should appear as points and at a certain level when I zoom in, I want them to appear as polygons (showing the area of the city and boundaries). How can I do that? I am thinking of creating two tables "cities_points" and "cities_polygons" in my PostGIS database. Any idea? (actually not only the cities, but it's the same for other components of my map).


texts: I've just tried following a simple tutorial to get familiar with QGIS. I'm trying to get a label to show up on a map without success. I followed the tutorial twice with the same result. In the Layer Styling window I've tried changing settings but do not see the label on the map. What am I not seeing?


texts: Need Advice - Displaying very large amounts of data on Leaflet js from database

texts:
I'm researching the best way to display a very large weather dataset (850mb +/-) for folks to use on their mobile and desktop devices. This huge dataset (1 point every 12km across the USA) that I'd like to show on a leftlet map will be using OpenStreet map layers. This database will be updated daily if not hourly with fresh data.

I'm looking for thoughts on what type of database to use, how to best display (speed) this data and any other advice you professionals have to share. There will be around 12 layers users could switch on an off.

I have been looking into GeoServer for example, and think that might be a good route.  I am not interested in using Google fusion tables and would like to run this on my own server as the data is already produced locally.

Also, the data is currently available is in the NetCDF format. Looks like GeoServer can use this.


texts: Calculating Distances Between Geometry Columns of Differing SRIDs

texts:
I am attempting a relatively simple task but have been running into trouble. I have two tables, one which has lat and lon points which were geocoded using google (SRID 3857) and another table which contains block boundaries obtained from the Census (SRID 4269). Below are the steps I have taken.

After converting the google map points to geometry. I attempted to make things easier for myself by calculating the centroids of the block boundaries with the following command

  CREATE TABLE mytable AS SELECT var1, var2, ST_Centroid(geom) as geom FROM oldtable;

Next I attempted to convert the centroids in the new table to SRID 3857 using the following

  CREATE TABLE block3857 AS SELECT var1, var2, ST_Transform(geom, 3857) as geom FROM mytable;

After creating a geospatial index on both the new census table and the table from google maps I attempted to merge the google maps geometry with the census centroids as follows:

          CREATE TABLE merged AS
            SELECT a.var1,b.var2, b.var3
            FROM googlemapstable as a,block3857 as b
            WHERE ST_DWithin(a.geom,b.geom,4828.032)

To my surprise the merge yielded no results. I plotted some of the points on top of each other in R and found that many of them were overlapping and should have merged. My next step to troubleshoot was to look at the distances between the points. I took only a small subset of the points in one city and found that ST_Distance was yielding results of almost 1,000,000 units.

I took a look at the geometry columns' underlying lat and lon points and found that the google map points looked normal for the US, but the census lat and lon points were different than i have ever seen, with the whole number portion of the point being 8 digits long. I went back and found that this change in the lat/lon occured AFTER converting to SRID 3857. Does anyone know if the ST_Distance function makes use of the lat and lon points to do the calculation and if this is what is causing my problem? Does anyone have any solutions to what is going on here?


texts: Dotted or dashed lines are drawn as solid in the output file

texts:
With QGIS 2.12, I am displaying polygon outlines as dotted lines.  The Print Composer screen correctly displays those outlines as dots. Here's a screenshot from the Print Composer window, showing dotted lines:

However, when I Print Composer > Export to Image (or pdf) the output file displays solid lines, not dots.  Here's a screenshot from the output pdf.  Note the solid lines:

The same problem happens if I select dashed lines.  I'd like to have dots or dashes in my output media instead of solid lines.

Is this a bug, or am I doing something wrong?


texts: Visualising a raster PostgreSQL table in QGIS

texts:
I have imported a raster tif image into PostgreSQL with this command:

raster2pgsql -s 32643 -I -M filepath.tif -F -t 100x100 public.databassename > filepath.sql

And imported the output SQL file inside the PostgreSQL database running this:

psql -U postgres -d databasename -f filepath.sql

After connecting to this database in QGIS, I am not able to visualise the tif image due to missing geometry content. Please help me out to visualise the tif raster imagery in QGIS.


texts: How to define image map with pins in Javascript?

texts: I have image map in .png format. All I need to do is to add a sidebar where users click on for example "restaurants" and pins with restaurants will show on an image map. When user mouseover a pin, it will show some more informations. I would like to do it in Javascript, but not with Google Maps. I have really small amount of pins for four different categories. What is the best way to do it in Javascript?


texts:
I am now working with VIIRS/NPP Active Fires by using python gdal. But I cannot read the data inside the files.
Filename = "NPP_AVAF_L2.A2012019.0600.P1_03110.2014057125956.hdf"
Here is some information that I can get by using gdal.

gdal.RasterCount = 0

And gdalinfo:

gdalinfo NPP_AVAF_L2.A2012019.0600.P1_03110.2014057125956.hdf
Driver: HDF4/Hierarchical Data Format Release 4
Files: NPP_AVAF_L2.A2012019.0600.P1_03110.2014057125956.hdf
Size is 512, 512
Coordinate System is `'
Metadata:
  AlgorithmType=OPS
  Beginning_Time_IET=[1.7056441e+15]
  BeginningTime=060027.600000Z
  DayNightFlag=Day
  EastBoundingCoord=123.866
  Ending_Time_IET=[1.7056444e+15]
  EndingTime=060609.000000Z
  EndTime=2012-01-19 06:06:09.000
  HDFEOSVersion=HDFEOS_V2.17
InputPointer=NPP_GRCMAE_L1.A2012019.0555.P1_03110.2014057115525.hdf,NPP_GRCMAE_L1.A2012019.0600.P1_03110.2014057115623.hdf,NPP_GRCMAE_L1.A2012019.0605.P1_03110.2014057115525.hdf
  InstrumentShortname=VIIRS
  LocalGranuleID=NPP_AVAF_L2.A2012019.0600.P1_03110.2014057125956.hdf
  LongName=VIIRS/NPP Active Fires 5-Min L2 Swath ARP 750m
  LPEATE_AlgorithmVersion=NPP_PRVAF 1.5.07.01
  LUTs_used=VIIRS-AF-EDR-AC-Int_v1.5.06.02_LP
  NorthBoundingCoord=27.8258
  Number_Fire_Pixels=256
  NumSCEA_RDR_TimeSegments=[18]
  NumSci_RDR_TimeSegments=[4]
  PGE_EndTime=2012-01-19 06:05:00.000
  PGE_Name=PGE330
  PGE_StartTime=2012-01-19 06:00:00.000
  PGEVersion=P2.3.0
  Platform_Short_Name=NPP
  ProcessingEnvironment=Linux minion5609 2.6.18-371.1.2.el5 #1 SMP Tue Oct 22 12:51:53 EDT 2013 x86_64 x86_64 x86_64 GNU/Linux
  ProcessVersion=P1_03110
  ProductionTime=2014-02-26 12:59:56.000
  ProxyDataType=Operational Data
  Resolution=Imagery
  SatelliteInstrument=NPP_OPS
  ShortName=NPP_AVAF_L2
  SouthBoundingCoord=3.93342
  StartTime=2012-01-19 06:00:27.600
  Unagg_DayNightFlag=TS 0: Day; TS 1: Day; TS 2: Day; TS 3: Day
  WestBoundingCoord=90.5542
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0,  512.0)
Upper Right (  512.0,    0.0)
Lower Right (  512.0,  512.0)
Center      (  256.0,  256.0)

This is what I can see from HDFViewer on same file:

I know that the file carries information of active fire points, but I cannot read the data inside into an array.
From gdalinfo above, I can see Number_Fire_Pixels=256, but how can I get lat,lon of these points?

How can I read the data by using python?

UPDATE: Here's the link to the file that I'm working with.


texts: How do I autopopulate a field in Collector?

texts:
I am wondering how I can use location to automatically generate a field in my attribute table with the name of a geographic region.
I want to:

create a point on the Collector app
automatically populate the field for geographic region (polygon) that the point lies in.

Basically, how do I avoid having to do post processing on my attribute table using a spatial join?
My GIS experience is fairly limited...


texts: Converting ECEF coordinates into geodetic LLA?

texts:
I am trying to convert ECEF coordinates into geodetic LLA.

This is the coordinates:

X = 1199816.63380
Y = -4812201.6785
Z = 4016140.62425

I checked the coordinates using this and it returns correct values:

Lat: 39.111677542
Long: -76
Alt: 12068.6198

But using my code, I get these values:

Lat: 179.79495170382071
Long: -75.99999999999999
Alt: NaN

And this is my code:

const RADIUS_EQUATOR = 6378137.0;
const RADIUS_POLES = 6356752.3;
const e = Math.sqrt((RADIUS_EQUATOR*RADIUS_EQUATOR / RADIUS_POLES*RADIUS_POLES) / RADIUS_EQUATOR*RADIUS_EQUATOR);
const eprime = Math.sqrt((RADIUS_EQUATOR*RADIUS_EQUATOR - RADIUS_POLES*RADIUS_POLES) / RADIUS_POLES*RADIUS_POLES);

function radianToDegree(r) {
  return r * (180 / Math.PI);
}

function ECRtoLLA(ECR) {
  let X = ECR.x;
  let Y = ECR.y;
  let Z = ECR.z;

  let a = RADIUS_EQUATOR;
  let b = RADIUS_POLES;

  //Auxiliary values first
  let p = Math.sqrt(X*X + Y*Y);
  let theta = Math.atan((Z*a) / (p*b));

  let sintheta = Math.sin(theta);
  let costheta = Math.cos(theta);

  let num = Z + eprime * eprime * b * sintheta * sintheta * sintheta;
  let denom = p - e * e * a * costheta * costheta * costheta;

  //Now calculate LLA
  let latitude  = Math.atan2(num, denom);
  let longitude = Math.atan2(Y, X);
  let N = getN(latitude, a);
  let altitude  = (p / Math.cos(latitude)) - N;

  if (X < 0 && Y < 0) {
    longitude = longitude - Math.PI;
  }

  if (X < 0 && Y > 0) {
    longitude = longitude + Math.PI;
  }

  latitude = radianToDegree(latitude);
  longitude = radianToDegree(longitude);

  return {
    latitude: latitude,
    longitude: longitude,
    altitude: altitude
  };
}

console.info(ECRtoLLA({
  x: 1199816.63380,
  y: -4812201.6785,
  z: 4016140.62425
}));

Any pointer or help?


texts:
I'm looking for a way to quickly transform simple lines into wavy, zigzag lines and so on in QGIS.
I already tried to use the method reported by @christoph in a question posted in this forum (Drawing wavy, wiggly lines in QGIS?) but I can't make it work.
In ArcGIS this is a default option, I would like to make it so easy to apply it to my projects as well.

the results I'd like to achieve are these:


texts: Keeping symbol in Map but not showing it in Legend using ArcMap

texts:
I am trying to find a solution to not display the name and the symbol on the Legend on my ArcMap Layout but I do want to keep the symbol in the Layout .

Here is the screenshot

I have ArcGIS Desktop 10.5.1 here


texts: Multiple mxds with the same Data Driven Pages properties

texts: I have 20 different mxds that display data for 20 different sites.  Is there a way to bring Data Driven Pages properties from one mxd to another so the new mxd will link to different tables and bring in the new attributes so I don't have to insert all of the dynamic text again?


texts: Location based service with normal image

texts:
We have an image for a park which has been rendered like this:

Now we want to add some geographic service based on the data: building an Android Application with the following features:

Display the map with drag/move/zoom operation.
Add overlays like markers,lines on the map.
Show the location of the user(by gps or network service).
Calculate best routes between two points.

Is this possible?


texts: Swapped coordinates when adding SHP created in QGIS to ArcMap

texts: I've created a shapefile layer in QGIS with contours in national coordinate system - 5514. Everything seems okay but when I add the layer to ArcMap, the whole layer is rotated and placed in totally different area. The problem is that the same point in QGIS has coordinates -540933,-1231037 and in ArcMap numbers are swapped and without minuses to 1231037, 540933. I don't even know how this can happen and what to do with that. Have anybody dealt with the same problem?


texts: Using Construct Points on all lines in shapefile using ArcPy?

texts:
I need to use the editing tool 'Construct Points' on several shapefiles that each have ~1300 polylines. I need to do this to create points every 1km along each of the polylines, which represent the orbits of a satellite. I've tried several other methods such as densify but none of them solve my problem as well as 'Construct Points'. The only problem is that you can only edit one line at a time, which would obviously take too long to do manually.

So I'm wondering, is there a way to implement 'Construct Points' in python so that it loops through all polylines and applies the tool with the same parameters?

I feel like it should be possible, but I've never used editing tools in arcpy so I'm not sure.

The output of 'construct points' is a point shapefile for the points that were constructed along the line used as an input. I'm wondering if I could do it in python such that it constructs the points on the line, then appends the output points to a larger point shapefile.

The following is pseudo code of what I have in mind, although the syntax is of course not correct. If anyone could shed any light, that would be great. I really would like to know if this is even possible with arcpy before I spend time trying to pursue it.

# using just one shapefile as an example (pseudo code, incorrect syntax):

orbitshapefile = (shapefile with ~1300 polylines)
pointshapefile = (empty shapefile to store points)

StartEdit(orbitshapefile) # open editing on the shapefile
for orbitline in orbitshapefile:
    editor.ConstructPoints(orbitline, 1km, ..(other parameters).., outputpoints) # call editing function
    pointshapefile.append(outputpoints) # append output points to 'big' shapefile

I just did construct points manually on just 2 lines in editor and found that using the same point shapefile on the second line that was used with the first adds those points to that shapefile, doesn't overwrite them. However, I need an indication of which line these points came from. There is an 'Id' attribute in the resulting table,  but everything is 0. Would there be a way to carry over the line's FID to all the points that came from it (alt. I could do a separate script after that uses intersect b/w line and points to come up with a lineID attribute). An updated version of the above code:

for orbitline in orbitshapefile:
    lineID = orbitline.FID
    editor.ConstructPoints(orbitline, 1km, ..(other parameters).., outputpoints)
    outputpoints[pts that were just created in this iteration].SetAttribute(lineIDField = lineID)

Some additional thoughts: I'm not sure how to subset the points that were just created-- maybe I could subset by saying something like: points in outputpoints s.t. lineID == 0 (because all points that were created in previous iterations will have their lineID attribute set to something other than 0)

Any thoughts on this rationale?

Can anyone AT LEAST tell me if using an editor tool iteratively in python is even possible?


texts: Reclassify says: 'Unique values cannot be obtained from raster, using values from dialog'

texts:
I am using ArcGIS 10.2. Whenever I reclassify a raster layer I get the following message:

  Unique values cannot be obtained from raster, using values from dialog

This does not prevent from proceeding with the reclassification. I found that if I convert the data to integer prior to reclassification, I do not receive this message.

However, I do not want to convert all my data to integer, especially since I am using data which is represented by values of less than one, hence needing the additional step of multiplying the rasters by 1000 before integrating them.

What does this error message mean? Can I just ignore it? Does it affect the accuracy of the reclassifications?


texts:
I need to show label of each feature based on attribute field value. I know that I can use ${attribute} notation for this purpose. But this approach works only if I make StyleMap and pass it to vector layer constructor. But in this case I can not modify style of each features separately. How I can make layer with labels and with not null value of style property of each feature?

My current config:

var custom_style = {
    fill: true,
    fillColor: "#feb24c",
    fillOpacity: 0.2,
    strokeColor: '#3182bd',
    strokeWidth: 1
};

var l = new OpenLayers.Layer.Vector("Coverage", {
    strategies: [new OpenLayers.Strategy.BBOX()],
    protocol: new OpenLayers.Protocol.HTTP({
        url: "http://***:8088/collect",
        format: new ol.Format.GeoJSON()
    }),
    style: custom_style
});

If I add label: ${attribute} to custom_style object, then text ${attribute} will be shown on map, not value of attribute.


texts: 'File is too big' error when uploading QGIS Plugin

texts:
I have written a simple QGIS plugin which zooms the map canvas to a user specified UK postcode. I would like to share the plugin but the csv file, which contains stripped down OS Code-Point open data is just over 40mb. This is too large to upload to the qgis repository and impractical to download and share.

I have also tried storing the postcodes in a JSON file and pickle file, both of which are also to large a file size. Is there any alternative way of storing the postcodes which may result in a smaller file size?

Thanks


texts: How to serve OpenStreetMap Tiles locally using ArcGIS Server 10?

texts:
We want to serve OpenStreetMap Tiles locally using ArcGis Server 10 to use as base maps in our ESRI Silverlight and Flex GIS clients.

We've downloaded our regions OSM file and created PNGs using Maperitive for Zoom 1 thru 18 and tested with Leaflet.  Also used a script to generate PGW files for each PNG.

Is this goal possible, should we use raster datasets or something different ?


texts: How I can make a child table appear in the parent table of attributes in QGIS?

texts:
I used "relate" to connect a many to one child table to a parent table (addresses to land parcels, one parcel can have several parcels) and I see it in the "form view but it is not displayed in the regular view in the attribute table, the same way a joined table would be (I don't see additional columns).
Is it normal?
How can I change this?
From what I understand "join" only work with one-to-one relationships.


texts: QGIS find steepest path on a DEM [Path of a drop of water]

texts: Is there an algorithm or plugin in QGIS to find the steepest path on a DEM given a point (the path followed by a drop of water)?


texts: Changing default IDE for editing Validation code in ArcToolbox?

texts: Do you know a way that I can change the application that opens the validation scripts in ArcGIS 10.2? At the moment hitting Edit... opens up Notepad, which ain't that friendly an  IDE.


texts: Container 'Control visibility by expression' for relations

texts:
I have two relations in a project on QGIS 3.10. They have the same field ('fid') of the parent but different fields ('topo_part';'topo_part_loc') of the same child. I want to see these relations in the attribute form of the layer parent. Actually, I want to see only the relation that relate some geometry childs because one relation keeps out the other one.

So I've made two container in the Drag and Drop Designer. Each of them contains a different relation of the layer. Now, i'm looking for an expression to control visibility.

I've tried many expressions but none of them work.

Where is the mistake? Any other solution?


texts: Varying distance cross section along the river

texts: I have used the script to draw the profile with this gis.stackexchange.com/a/50841/8104. This script is perfect. Can someone provide me idea, if in my case spacing is not same. I have points with different spacing. I want to use those points instead of fixed spacing. I would be thankful if you provide some guidelines


texts:
I'm having trouble to calculate the cellsize of Esri ASCII data.
I would like to use meteorological data from JMA following data with QGIS3 software which is used
in my workplace.
http://database.rish.kyoto-u.ac.jp/arch/jmadata/data/jma-radar/synthetic/original/2014/08/09/Z__C_RJTD_20140809000000_RDR_JMAGPV__grib2.tar
※Z__C_RJTD_20140809000000_RDR_JMAGPV_Ggis1km_Prr10lv_ANAL_grib2.bin
This data's mesh size is 1km square.
ncols 2560
nrows 3360
xllcorner 118⇒degree
yllcorner 48⇒degree
catalog(written in Japanese) https://www.data.jma.go.jp/add/suishin/catalogue/format/ObdObs001_format.pdf
But I don't have any idea how to decide (calculate) correct cellsize parameter.
Could you give me a hand?


texts:
I am trying to open an empty project using PyQGIS, i.e. doing this:

All the answers I can find relate to already existing projects. i.e. they load a file like this:

project = QgsProject.instance()
project.read('file.qgs')

I've tried passing no parameters/empty string to the read method, but with no effect. Looking at the methods that QgsProject.instance() provides, it doesn't seem like there is another method that would do this either.


texts: Can I use Representations (ArcGIS 10) to enclose a group of polygons with an outline without altering data?

texts:
ESRI describes representations as follows: A feature class representation is a property of a feature class that allows you to
  specify and store a series of rules that dictate the way the features
  in the feature class should be drawn.

I've used a few simple representations in the past to correct the appearance of some line features as well as complex outlines of polygons features with varying levels of success. I am still trying to fully understand the usefulness of representations. So far it seems like too much effort with little reward.

Right now I have a situation where I have a group of polygons around which i would like to show an outline. I know I can merge them but the group has several types and is always changing and I would like to make this outline dynamic; if a polygon is added or removed, the outline will change accordingly. I basically want to blow all the inner lines away without having to alter the data in any way.  In terms of symbology, it must remain a polygon outline with no fill.

Can this be done with representations? I was thinking of perhaps using the "Enclosing Polygon Effect" with "Close Path" method but that doesn't do anything. If I can't use representations is there another dynamic way of doing this without having to alter data or using masks or overlaying layer symbologies?

Below are images of what the layer looks like and what I want it look like.

Here is what the layer looks like (yellow outline):

Here is what I want it to look like:


texts: Updating initial feature class within for loop using ArcPy?

texts:
I have a forested area, that suffered different disturbances over multiple years. Each year, disturbance data was collected on different forest extend and some disturbances can overlap (beetle infestation or clear-cut).

My goal is to clean and summarize the dataset, i.e. unify the extent and allow the forest area to be disturbed only once over my time period.
Simply, I wish to clip the disturbance data to the available forest area in current year, and update the extend of forest cover after every disturbance.

I am trying to apply arcpy to automate the process.

Pseudocode:

read forest area & list yearly disturbances
create for loop to loop through list of disturbances

clip disturbance by available forest
update available forest by erasing the disturbance from current forest cover
do the same for next disturbance

I thought this could be fairy simple to do. However, my script either return empty dataset, erase the initial forest or says "that output and input data equals".

My desired output is the set of cleaned (not overlapping disturbances limited on forest) disturbances. I am trying to automate the process, as I have 3 different disturbances types over multiple years, together 18 disturbances featureclasses.

How to update the forest polygon and spatially restrict every disturbance?

import arcpy
import os

# Allow files to overwrite
arcpy.env.overwriteOutput = True

# Define directories
path = "U:\\project\\analyzed"
arcpy.env.workspace = os.path.join(path, "analyzed.gdb")
outWD = os.path.join(path, "cleaned.gdb")
aoi = "aoi_fin_03"

# Read all disturbance data
disturbLs = arcpy.ListFeatureClasses()

# Intersect available forest with disturbances in each year,
# update forest cover
for dist in disturbLs:
    # Clip the disturbance to forest extend
    outClip = os.path.join(outWD, dist)
    arcpy.Clip_analysis(dist, aoi, outClip)

    # Update available forest cover by erasing current disturbance
    updAoi = "in_memory\layer" + aoi
    arcpy.CopyFeatures_management(aoi, updAoi)
    arcpy.Erase_analysis(updAoi, outClip, aoi)

#Delete in memory files to avoid the lock
arcpy.Delete_management("in_memory")


texts: Identifying road segments completely cut off by polygon using QGIS?

texts:
I have a road network (black line) and I want to identify all segments that are completely cut off by my flood polygon (blue polygon).

I have clipped the road network to get all roads within the road network and I have been able to identify some road segments that are cut off by the flood zone polygon (green lines and pink lines), but not all of them. I did this through a combination of the 'intersection' tool and the 'Spatial Query' plugin.

I'm using QGIS for this.

I used the 'Topology Checker' tool to identify all 'dangles' but it highlights road segments that are, for example, just dead ends but are connected to the wider road network. I guess I'm looking for 'double dangles'?

I'm basically looking for road segments that are completely disconnected from the larger road network but they can be connected to one another. 'Dangle islands'?

If an agent is placed on a 'dangle island' and is not able to access the full network, then that road network is a 'dangle island'!


texts: How to stop select interaction without unselecting all features

texts:
I use interaction:

var select = new ol.interaction.Select({
    layers: [my_layer],
    style: new ol.style.Style({
            fill: new ol.style.Fill({
             color: '#A00000'
         }),
             stroke: new ol.style.Stroke({
             color: '#000000'
         })
            })});

for selecting features on my_layer:

 var map = new ol.Map({
        target: 'map',
        interactions: ol.interaction.defaults().extend([select]),
        layers: [
          new ol.layer.Tile({
            source: new ol.source.OSM()
          }),
          my_layer,
          vector]});

Now, I want to measure distance between points using sniped code from:
http://openlayers.org/en/v3.1.0/examples/measure.html
To avoid selecting features during measure process I have to remove select interaction from the map

map.removeInteraction(select);

But this removes selection from the selected features too.
How to avoid this problem about unselection after removeInteraction(select)?
Or maybe can I add any condition in this interaction to stop selecting features during measuring?


texts: How to Make automatic classes with algorithm in qgis modeler?

texts:
I have 1 column of numeric data in one attribute table of a layer made of points.

I want to regroup this data in classes (style) with the Qgis modeler and each class must have a specific lenght that i previously calculated with the modeler.

Which algorithm do i have to use in the modeler in order to create these classes automatically ?


texts: Is d.mon xN deprecated in GRASS in favor of something in the wxGUI?

texts:
I'm trying to learn command-line GRASS (6.4.3) using the KyngChaos binaries on a Mac (10.8.4) from Neteler and Mitasova's book, Open Source GIS: A GRASS GIS Approach (3rd ed.). I've learned from this site that some of the commands in the book have been deprecated with the introduction of the wxGUI. Is d.mon one of those commands that has been deprecated? I've found that I can still invoke a X11 graphics monitor using d.mon, but is it preferable to use the graphics monitor in the wxGUI? If so, then can I access the wxGUI's graphics monitor from a terminal session? I'm finding that the command console in the wxGUI behaves somewhat unpredictably compared to a terminal (e.g. I can't cd to change working directory).

Also, am I likely to continue running into inconsistencies between the book and the current version of the software as I work my way through the book?


texts: Mean function from rts package, returns empty dataset, R

texts:
I am using the rts package in R and specifically the mean function, over a a gridded climatic dataset with many NA values.
When I use apply.yearly(my_rts_object,mean) all I get is an empty raster time series.

Does anyone knows if a large number of NAs can introduce such a problem to the rts mean function?
I did not find any infos regarding NAs in rts documentation and I would really like to know if anyone has faced a similar issue (other functions like min and max work fine and return values).


texts: Why missing address IDs in postgresql/postgis query results?

texts:
I have two tables "ad" and "st" in my postgres DB. Table ad has 87 rows while st has 16,060 rows. Following code works perfectly to select st.geom (streets' geometry) within 50 meter radius with respect to every address locations.

SELECT
      ad.gid As address_id,
      ad.geom As address,
      st.geom As streets
FROM st
JOIN ad ON ST_DWithin(ad.geom, st.geom, 50.0)
Order by address_id

However, what I noticed that in the query results some address IDs are missing (for example; 5, 51, 56 etc.) which is not desired as I needed to search nearest streets for all 87 rows of ad table. Here is the screen shot of query results.

Using this command, I cross checked that address_id_5 does exist in the respective column.

Select ad.gid from ad;

Can some one please help me to understand this ambiguity? There are seven address locations missing out of 87. I have tried the same on other data as well and the problem persists. Both tables contain spatial indexes already.


texts: Retrieve pixel value using SNAP-engine and Python for specific geographic coordinates

texts:
I'm trying to retrieve pixel values in Sentinel-3 data for specific geographic coordinates using the SNAP-Engine.
Until now, I transform the geographic coordinates to the corresponding X and Y with the following code using these links:

https://stackoverflow.com/questions/51194477/converting-geographical-coordinates-long-lat-to-pixel-position-x-y-for-a-se
https://senbox.atlassian.net/wiki/spaces/SNAP/pages/19300362/How+to+use+the+SNAP+API+from+Python)

import sys
sys.path.append('/path/to/.snap/snap-python')
from snappy import ProductIO
import numpy as np
from snappy import GeoPos

p = ProductIO.readProduct('/path/to/file.xml')
iwv = p.getBand('IWV')

g = p.getSceneGeoCoding()

def XY_from_LatLon(ProductSceneGeoCoding, latitude, longitude):
    #From Latitude, Longitude satellite image, get the x, y position in image
    pixelPos = ProductSceneGeoCoding.getPixelPos(GeoPos(latitude, longitude), None)
    x = pixelPos.getX()
    y = pixelPos.getY()
    if str(x)=='nan':
        raise ValueError('Latitude or Longitude out of this product')
     else:
         return x, y

How can I retrieve the data from a pixel using image coordinates (X, Y) in SNAP-Engine ?

After some more searching in SNAP-Engine docs, I found the method: getPixelFloat() that returns the pixel value given x, y.
The code above was updated to the following :

import sys
sys.path.append('/path/to/.snap/snap-python')
from snappy import ProductIO
import numpy as np
from snappy import GeoPos

p = ProductIO.readProduct('/path/to/file.xml')
iwv = p.getBand('IWV')

g = p.getSceneGeoCoding()

def XY_from_LatLon(ProductSceneGeoCoding, latitude, longitude):
    #From Latitude, Longitude satellite image, get the x, y position in image
    pixelPos = ProductSceneGeoCoding.getPixelPos(GeoPos(latitude, longitude), None)
    x = pixelPos.getX()
    y = pixelPos.getY()
    if str(x)=='nan':
        raise ValueError('Latitude or Longitude out of this product')
     else:
         return x, y

x, y = XY_from_LatLon(sg, <lat>, <lon>)
print(iwv.getPixelFloat(x, y))

but it raises the following error:

Traceback (most recent call last):
File "snap.py", line 28, in <module>
print(iwv.getPixelInt(x, y))
RuntimeError: no matching Java method overloads found


texts: Script to either rotate existing text element or insert graphic element

texts:
I have about 300 maps that contain a text element that needs to be rotated.  I am pretty new to python but I have heard how great it is.  This seemed like the perfect task for python but I have had no luck finding any examples on how to do this.  The text is "Figure 1-1" and it just needs to be rotated 90 degrees and moved to a set location in the corner.  If rotating is not an option I thought I could just manually rotate the text and save it as a graphic and then write a script that would insert the graphic in the desired location.  This also doesnt seem possible.

Any help would be much appreciated.


texts: How to call ESRISignAddIn Silently (no UI)?

texts:
I am creating a build workflow and I ran into some problems.

When I call ESRISignAddIn.exe from the command prompt it is throwing up a dialog saying "Add-In successfully signed". I was expecting it to report that back within the cmd dialog and not open up a new messagebox.

Have anyone else seen this or have worked out how to sign an add-in "silently"?

The /s switch does not work.

This issue occurs with the 10.2 version of the exe (latest at time) and earlier versions.

"c:\Program Files (x86)\Common Files\ArcGIS\bin\ESRISignAddIn.exe" "C:\temp\my.esriAddIn" /c:"C:\temp\newCert.pfx" /p:myPass


texts: Connecting Geodesign and attributes of two layers using QGIS?

texts:
I'm kind of new in QGIS and I wonder if there is a function or a tool to connect Geodesign and attribute of two layers, for example:

I have a line vector layer named as (MT) with the attribute below

And a point layer named (HH) with the attribute below

So when I draw a point on the (MT) line layer the field "Route" in the point layer should take the info. of the field "Route Type" in the (MT) line layer as shown below


texts: Getting HTML-Values in a QGIS-form by pyQT5

texts:
I have created a form for point-objects, that is delivering some values (getting by a database-connection) in the "html-area" in the form.

So far the form works well and I´m getting the values that I need:

But since hours I try to get these values by PyQT5.

No problem to get the "NON"-HTMl-values like:

bottom_elevation = myDialog.findChild(QDoubleSpinBox, "BOTTOM_ELEVATION")

But how can I get by Python the entries in the html-section?


texts: Modifying straight line segments to follow geodesic / great circle paths using QGIS?

texts:
I have a kml file showing a dodecahedron stretched across the surfce of a globe, which I downloaded here: http://montalk.net/cgi-bin/coordinates.py

I would like to re-project this dodecahedron so I can create tiles for a world globe. However, the file just consists of straight line segments between the vertices, so the lines no longer follow their geodesic / great circle paths across the surface of the map when re-projected onto a flat map.

I run QGIS (and Google Earth).

How should I go about modifying this file so that the edges of the pentagons follow the correct great circle arc path across the surface of the map whatever projection they are in?


texts: How can I use Distance Matrix tool in QGIS with metric distance for global data?

texts:
I have 2 point layers spread out all across the globe. I would like to use the Distance Matrix tool, and would like to get the distances in metres or kilometres. I'm aware that I have to change the projection to a UTM projection.

But since my data is across all longitudes, which projection should I use? I tried a few (EPSG 3857) to no avail; I still get the distances in decimal degrees.


texts:
I copied rasters across into my geodatabase for a client but when they open the raster dataset in ArcMap they appear in tiles as shown in the image below.
When I drag and drop individual rasters over they appear as they should but when adding the dataset I had a projection error message. So I added the correct projection to the dataset but I still have the same problem when adding the dataset as a whole to ArcMap.

What can I do to fix this?


texts: ArcGis: Count categorical variable per nuts3 by overlapping point-file and polygon-file

texts:
Download from this website the NUTS_RG_60M_2006 files (shapefile) and the csv-file T_DischargePoints_West.csv (source) and save them on your computer. Open them through ArcGis and plot xy-points. You get this result:

The cvs-file contains a number of variables. There is one categorical variable (dcpWaterBodyType). I now would like to know how many points from each category of dcpWaterBodyType there are in each nuts3-polygon. If possible, I also would like to determine automatically what the percentage of points of WaterBodyTypeA, WaterBodyTypeB, WaterBodyTypeC... is in one polygon and this for all the polygons.

Can somebody explain me how to do this in ArcGis or R?


texts: Problem Visualizing 3D Terrain Model - USGS DEM

texts:
I recently downloaded a GeoTIFF 1 arc second from the USGS National Map Viewer.  I added the raster to ArcScene and attempted to visualize the data in 3D.  I set the Base Height to Floating on a Custon Surface, changed the Resolution to that of the Original Surface (3612 rows and columns; 0.000278 cell size).  However, when it's done processing, my DEM just looks a mess; there are giant spikes all through the data.  I've tried downloading DEMS from other areas and repeating the process but I get the same result.

What am I doing wrong?  I suspect perhaps there are Null values but I'm not sure.  Any help would be greatly appreciated.


texts: QGIS processing algorithm - add saved geopackage layer to map

texts:
I am still new to Qgis and Python and most likely this is really simple but however I do not manage  to do the following:
I am trying to write a python processing algorithm for the processing toolbox. It uses various native processing algorithms (like buffer) and these produce several output layers that I want to save to a single geopackage to make the results more organized. Because of this I can not use standard ProcessingParameterVectorDestination, because then you can only save one layer per gpkg.
The Layers are successfully saved to my geopackage as shown below, but I do not manage to display them in the map directly from my processing toolbox script.
Steps:

let my users choose a geopackage to save all the produced layers

       self.addParameter(
            QgsProcessingParameterFileDestination(
                self.OUTPUTFILE,
                self.tr('Output Geopackage'),
                'Geopackage (*.gpkg)'
            )
        )

and
destination = self.parameterAsFileOutput(parameters, self.OUTPUTFILE, context)

run several processing scripts and use processing.runAndLoadResults() to add the layers like in this Question
Example:

        zaun = processing.runAndLoadResults('native:buffer',
            {
                'INPUT': parameters[self.INPUT],
                'DISTANCE': zaunabstand,
                'SEGMENTS': 5,
                'DISSOLVE': False,
                'END_CAP_STYLE': 1,
                'JOIN_STYLE': 1,
                'MITER_LIMIT': 2,
                'OUTPUT': 'ogr:dbname=\''+ destination +'\' table=\"zaun\" (geom)'
            },
            context=context,
            feedback=feedback)['OUTPUT']

When I use the python console to check if the layers are there, it seems like it worked but I cant see them on the map:

layers = QgsProject.instance().mapLayers()
print(layers)
{'OSM_Standard_efcb35f0_1b35_4f6c_ac58_14682871cc31': <QgsRasterLayer: 'OSM Standard' (wms)>, 'Umring_46c1f9a4_f251_44f1_a1d6_cad5aa3ad847': <QgsVectorLayer: 'Umring' (ogr)>, 'test_a7d3776d_c9f7_44ef_ae2b_26422643c173': <QgsVectorLayer: 'zaun' (ogr)>, 'test_e7c8812e_eb24_4497_a942_1cd185c82306': <QgsVectorLayer: 'tische' (ogr)>}

Most likely this is really trivial but how can I show the layers in the TOC and on the map programmatically from my toolbox script?
I also tried to add them with QgsProject.instance().addMapLayer() but this seems not to work out of a processing algorithm.
I wonder if it is possible to use QgsProcessingOutputVectorLayer to output the result of my buffer but I have not idea how.


texts:
We frequently have to update excel sheets that contain reports as number, images, and pie charts and import it in our data view in ArcMap (sometimes 13 times a day). We do it by copying the row/column from excel and paste in in ArcMap Layout view . We want to automate this task by linking excel to Layout view, using Insert--> Object and grab that excel file, you can save so much time and and step this way. The problem is when we do this, the excel file gets messed up. The alignment of the texts changes and some texts don't show up and so many other problems.

Is there any other way to do this and not lose the format from excel?

We have ArcMap 10.1.


texts: Symbolising One-to-many relationships using ArcGIS Desktop?

texts:
I have a polygon feature class which has various extent areas with a name which is unique

I also have a table (in a file geodatabase) that has the extent names with a score value.  But there can be many score values per extent.

I want to visually represent the score values either by colour symbology or pie/bar chart symbology in ArcGIS Desktop (I have an ArcInfo license too).

It's one of these questions that get bounced round the Esri forums quite a bit, so I was wondering if anyone here had any suggestions or any work arounds that they have done in the past?


texts:
Just getting acquainted with GeoDjango. Have Basic Tutorial app created with sample world data. In trying out the admin, I appreciate seeing the polygon on the map and a couple basic tools for adding a new polygon and a tool for moving existing points.

Question:

What options are there for actually DELETING existing points?

Are there other tools for easily creating/editing map polygons. I am targeting Google Maps API. My GeoDjango App will be the back end for generating Google Maps. Looking for options for creating/modifying boundaries/polygons to be maintained by GeoDjango.

Thank you :)


texts: Changing color interval of raster layer in QGIS?

texts:
I have a map with a raster layer that I took from Natural Earth website.

My client wanted me to transfer the sea color so it will look black, and the landmass into gray.


texts: Fetching query from query layer using ArcObjects and C#?

texts:
As shown in the image below, I need to fetch query mentioned within query section of a layer using ArcObjects in C#.

            IDocument doc = ArcMap.Application.Document;
            IMxDocument mxDoc = doc as IMxDocument;
            IMap map = mxDoc.FocusMap;
            ILayer mapLayer = map.get_Layer(0);
            IFeatureClass c = (mapLayer as IFeatureLayer).FeatureClass;


texts:
Curious if the following code has a blatant error that has been overlooked with respect to rotating a single polygon feature.

The *.geometry().rotate(...) returns a '0', indicating it was successful - yet the display and data-source do not depict any rotation; even after reload / refresh of layer, legend, and canvas.

Either the problem is in syntax, or maybe (because this is a Dev release) a not yet fully operational function?

inputShapeRef.startEditing()
featureSearch = QgsFeatureRequest().setFilterFid(int(inputShapeID))
featureRef = inputShapeRef.getFeatures(featureSearch)
for feature in featureRef:
    feature.geometry().rotate(self.MaxRotation,self.baseCentroid)
inputShapeRef.commitChanges()

UPDATE:
Single line of code was missing, and the following is the updated script.  Sucessfully rotates the shape to the desired angle; quite rapidly!

    inputShapeRef.startEditing()
    featureSearch = QgsFeatureRequest().setFilterFid(int(inputShapeID))
    featureRef = inputShapeRef.getFeatures(featureSearch)
    for feature in featureRef:
        feature.geometry().rotate(self.MaxRotation,self.baseCentroid)
        inputShapeRef.changeGeometry(feature.id(),feature.geometry())  # <=== + = success
    inputShapeRef.commitChanges()

FINAL THOUGHTS / QUESTIONS:

Does anyone know if this functionality carries over to an entire layer, the layers extent, or selection of multiple features?
The SAGA translate processing function fulfills this need - but is data-thirsty. The
low data-overhead and rapid return of this new single feature rotation would be
highly beneficial for an entire layer or selection of multiple features therein.

For example, say a point array needs to be rotated in sync with an associated bounding polygon - rather then rotate each point individually, is there a means of rotating the entire array simultaneously without creating a new dataset?


texts: Uninstalling/reinstalling QGIS plugin where uninstall is disabled?

texts:
Plugin Raster->TerrainAnalysis->slope which had initially worked, started showing what another questioner referred to as a "crazy histogram":
Slope analysis from NED data gives crazy histogram?

Before finding his question and the answer, I endeavored to uninstall and reinstall the plugin.

I found/still find (QGIS 2.0.1, now upgraded to QGIS 2.4, on Ubuntu 14.04), the plugin "uninstall" button grayed out/disabled.

It's working, but I am worried about having a plugin I installed that does not offer an option to be uninstalled.

Failing anything else I would be tempted to use the path information to rename the .so library and see if it disappeared.


texts: OpenStreetMap - assembling a hierarchy for Russian Federation

texts:
I'm currently working with OSM/Nominatim and getting ready to create a hierarchical tree of Russia's administrative districts starting from district, and ending with city. As the PostgreSQL from Geografik being parsed into the database, and has taken quite a while, I'm wondering if there are any other ways to gather such hierarchy.

One that I found was Overpass API, but its query language is quite cryptic.

The other was Nominatim which I'm also not yet 100% how to operate with.

Are there are instruments out there that would allow me to extract specific administrative level(s) (3 to 6) data for Russia quickly?

Update

Having read about Nominatim configuration, I'm importing a much lesser dataset for administrative boundaries and places through settings/import-admin.style.


texts: How to select multiple features using advanced search SQL in QGIS

texts:
I have a layer of planning units covering a large area (there are over 16000 PU's).  Three of the columns in the attribute table hold information about which planning units were included in the solutions from a conservation planning programme (Marxan).  So the columns have an ID and then values ranging from 0 to 100 (the number of solutions that PU was included in).

What I am trying to do is to select the PU's from each column that were in all 100 solutions (so have a value of 100) - BUT only the PU's that match i.e. only PU's that have a value of 100 in ALL three columns (the columns represent different run scenarios).

I have tried SQL searches such as "ssoln_1 = 100 WHERE ssoln_2 = 100 AND ssoln_3 = 100" but my SQL is very basic and I am having no luck.  Any suggestions would be hugely appreciated - I hope it all makes sense.


texts: SpatialPointsDataFrame too slow FOR LOOP

texts:
I've run into a small problem where my for loop has been running for 3.5 hrs+. I have a Geodataset in SQL Server with a field called Shape that stores co-ordinates in geometry data type. Firstly, I connect R to my DB via ODBC  and retrieve the information I want (also converting the Shape column to something readable)

library(sp)
library(rgeos)
library(RODBC)
ch<-odbcConnect("SpatialAnalysis", rows_at_time=1)
df<-sqlQuery(ch, "select OBJECTID, LOT_NO, Shape.STAsText() as WKT FROM SRC_PLI_QLD")
cnt<-sqlQuery(ch, "select count(OBJECTID) from SRC_PLI_QLD")

This has 2.5 million points. I now read them into a SpatialPointsDataFrame, reading the first element first.

point.sp <- SpatialPointsDataFrame(readWKT(df$WKT[1]),
                                   data=data.frame(OBJECTID=df$OBJECTID[1], LOT_NO=df$LOT_NO[1]))

Now I read the rest of the elements. AND this is where the problem is. It has been 3.5 hrs and still running.

for (n in 2:as.integer(cnt)) {
  point.sp <- rbind(point.sp,
                    SpatialPointsDataFrame(readWKT(df$WKT[n]),
                                           data.frame(OBJECTID=df$OBJECTID[n], LOT_NO=df$LOT_NO[n])))
}

What is the problem in the above mentioned for loop? Is there another way I can do this?


texts: How can I use custom style from MapBox with leaflet L.vectorGrid.protobuf

texts:
I have a custom style, that I made in MapBox. I don't have it in a digital version.

How can I use it with leaflet?

From here I understand how to add it. But I still need to apply styles.
From this example I also understand, how to use my custom MapBox styles with MapBox GL.

But how can I use the custom styles with leaflet VectorGrid protobuf module?

Btw, I was very impressed with how much faster VectorGrid handles multiple markers or polygons and how fast it can change individual styles on 1000nds of objects


texts: Geo- referencing Projected coordinated raster with GCS - ArcGIS

texts:
I have raster displaying Lat/ Long grid, but map information shows that, it s been mapped in Gauss Krugger projetion/ Beijing 1954. So, technically we need to geo reference the map in Gauss Krugger/ Beigjing 1954 Projectec Coordinated System in ArcGIS.

So, basically what I need to do is, georeference raster by input of Lat/Long, but need to update the georeferencing by projected coordinated system mentioned above.

In Global mapper, we can load the raster and enter Lat/Long available in the map, straight away it will be recalculated to Projected coordinated system.

Can we get anything same in ArcGIS?