texts: Deriving a value from a delimited text layer and displaying that value in a separate layer

texts:
Using QGIS 2.18.14

I have a delimited text layer imported, this contains coordinates that essentially plots a linear grid across the country. Each coordinate in the file denotes a corner of the grid and each coordinate has an value attributed to it.

I can visualise the coordinates with their attributed value but what I want to do is use the lower left-hand coordiates attributed value and display it in the middle of the corresponding square.

As follows:

How do I achieve this in QGIS?


texts: Adding Multiple raster datasets to mosaic dataset (for temporal change) in ArcMap?

texts:
I have a number of raster datasets containing lake water temperature data. I am looking to make these time enabled as to be able to scroll through them. I have converted them to .TIFF format and have tried using the "Add Raster to Mosaic Dataset" Tool. However, I am told to find my .TIFF files while "workspace" is selected. I cannot find these files unless I change "workspace" to "raster dataset"

Am I missing a particular setting in the tool? OR Am I missing a particular step?


texts:
I am trying to create relationship classes in a Feature Dataset. The Relationship Class created successfully, I want when I update a record in the Origin table the Destination's would be updated, but it didn't.

Any solution OR procedure?


texts: CARTO Mobile SDK.UWP (C#) Get always get current map's south-west MapPos and north-east MapPos

texts: How can I always get the current map's south-west MapPos and north-east MapPos when I zoom in/zoom out the map or drag the map? I want to create a dynamic bounding box for Carto Map.


texts: Finding raster calculator input in QGIS?

texts: I am using the raster calculator in QGIS 2.18.11, and I want to find the equation/input that went in the "raster calculator expression" box for rasters that have already been created with the raster calculator.


texts: How to deal with direction in pgr_ksp (K shortest path)?

texts:
I am working on an Android app for navigation with PostGIS and pgrouting. I have implemented pgr_dijkstra and it's working fine. I am not considering direction in routing at this point and just want routes from my source to destination.

I am facing a problem when I am implementing k-shortest path (pgr_ksp) for finding alternative paths from source to destination. I guess pgr_ksp always consider direction when calculating the route. So when I give it a particular source and target value in the query, its gives me nothing and the reason behind this is at some point on the route between that source and target, some road segments share starting or ending point, that is why it can't generate route.

Kindly tell me if there is any possibility that direction not considered here just like Dijkstra

In this picture A is my source, B is my destination. Ignoring the direction of travel, ksp should give me shortest path i.e. node 1, 2, 5, 6

But it is giving me 1, 2, 3, 4, 5, 6 which is not shortest path. startpoint and endpoint is also mentioned here with their road segment color

Kindly suggest me a solution.


texts: Displaying different attributes within layer with individual colors in QGIS 3?

texts: I am trying to display different soil types within a geological soil layer with individual colors per each type of soil. As of now it is displaying all as the same color.


texts: How to automatically convert line features into polygons using QGIS/PostGIS?

texts:
I have a master table of line features and a number of tables that I use to join (in the SQL sense) to produce various layers.  Some of theses layers consist of boundaries of administrative zones (think counties ;)  Using SQL I can either show the outline of all the zones or subsets of the zones.

Now I want to have a layer that consist of a set of polygons, one for each zone -- so I can use it with the atlas generator.

I would prefer a solution that can be automated, preferably in an sql query.

Some points:

Currently the line segments end points don't match exactly (how would one get them to?)
I know about ST_BuildArea but I think it expects the lines to form perfect closed polygons.  When I try and use it I get back n empty polygons where n is the number of lines in the input.

What I am looking for is something like ST_BuildArea which will allow some tolerance in the line end points.

Another way of coping with this would be to have a method of standardising endpoints of lines given some tolerance  Does such a thing exist?  Before I go and start re-inventing the wheel.  I figure it should be a fairly easy plugin to write.


texts: Ensure model runs from a specific layer?

texts:
I have some models for various routines but I really want them to run on a specific layer. I don't want to be able to run the model with any input layer as if done so by accident it can mean crashing out of QGIS due to the issue of having my model run on too big a dataset.

I would just like the model to only run if a specific named layer is selected. For example if the layer I want to use as my input is "map" then I don't want the model to run if any other input layer is chosen.

I have no python knowledge but I'm sure something could be done.


texts: Is there a way to export labels in QGIS?

texts: I've used the v2 labeling engine to place some nice looking road labels. Now I'd like to export those labels as a point shapefile for use in another program. Is there any way to accomplish this?


texts: Adding particular WMTS layer makes ArcMap experience Serious Application Error?

texts:
I am trying to access a WMTS using this URL from ArcMap 10.3.
I was able to get the list of layers (see picture below), and the service is shown in the Catalog Window.

Anyway, when I try to import it (by clicking and dragging a layer into the Table Of Contents), ArcMap always crashes. This seems not to happen with ArcGIS Pro and other GIS softwares (e.g. QGIS), while ArcMap 10.4 crashes as well. As far as I know, the WMTS I am trying to access is a static .xml file, which is not necessarily following the OGC standards. So the main problem could be some non-standardized lines in the xml, but this is just a guess, as I am not an expert of such services.

So I am asking if from the link I provided someone is able to understand what could be the reason of the problem (e.g. some non-standardized lines that show evident errors of interoperability), or at least if the problem occurs also on other machines other mine.

I eventually add that the Tile Matrix Sets shows two "subgroups" (see the picture):

GoogleMapsCompatible
WGS84


texts: Is there something like postgis for mariaDB?

texts:
Someone told me that PostGIS exists. Currently I'm storing latitude and longitude information un a MariaDB, but I want to know if there is something similar for MariaDB.

Unfortunately, migrating to PostgreSQL is not an option.


texts: ArcMap 10.1 Crashing when using field calculator

texts: I keep getting a serious error message when I try to use the field calculator. I've joined a .csv to my SDE feature class and am trying to set one of the columns equal to one of the columns from the .csv. But it keeps crashing! Is this happening to anyone else?


texts: Converting the local coordinates to GPS (lat, long)

texts:
I am having the XYZ points in some local coordinate frame. Each point I have, can also be represented as the distance from the origin and the azimuth angle it makes with one of the horizontal axis.

I also have the GPS coordinates(lat,long) and the heading corresponding to the (0,0,0) point of the local coordinate frame.

I need to convert all the points in the lat,long.

The approach I can think of is to use the Destination point given distance and bearing from start point formula from this link, as I have the:

starting point(long,lat corresponding to the origin)
distance between the origin and the actual point being converted to long,lat

But I am not sure what will be the bearing required in the formula. Does it relate to the azimuth, I have for each point?


texts: Caching aerials in a different coordinate system

texts:
I have a strange problem to deal with caching of the aerials. All of our Aerials are stored in SDE in State Plane WKID-2237 (NAD_1983_StatePlane_Florida_West_FIPS_0902_Feet). But we are planning to publish all our map services in WGS_1984_Web_Mercator Auxiliary Sphere (WKID-102100) to match up with the ArcGIS online/Bing maps.

What is the best way of caching the aerials in a different coordinate system than that of its source in storage? Transforming it and then storing it is being restricted- partly due to space limitations.

Can I do this by just changing the coordinate system in the Map Document to WKID 102100 and caching it without any Transformations? Are there any pitfalls with this approach?


texts: Seeking web page or small install GIS client to test WMS URL

texts:
I am using JMP 9, which supports WMS maps (as a client), but it doesn't seem to be working now.
Is there a web based, or small install GIS client on which I can test a WMS URL?
To clarify, I was looking for something that would help me find out the correct querystring parameters to send, as the service I was hoping to use did not specify them.


texts: How to group elements in QGIS legend side by side?

texts:
I would like to edit the legend in my map so that the entries of two polygons will appear together, i.e. the two entries in the legend should be displayed appear as "Study area and possible study area" as two elements side by side.
Is there a possibility to group the symbols in such a way?

Study Area and possible study area

 _ _   _ _
|_ _| |_ _|


texts: Open a numpy array as a GDAL raster in memory

texts:
I have a numpy array which I would like to open as a raster in GDAL. Currently the only way I can find to do this is to save the raster as an image to the file system and then open it again. I was hoping that there would be a way that I could just create an in-memory raster from a numpy array and then set its transformation and coordinate system manually.

The closest thing I can find is here: https://gis.stackexchange.com/a/37431/48798

And came up with this solution:

        nrows, ncols = self.paths_array.shape

        # Explanation for this mess: https://rasterio.readthedocs.io/en/latest/topics/migrating-to-v1.html
        top_left = (self.ref_img.meta["transform"][2], self.ref_img.meta["transform"][5])
        geotransform = (self.ref_img.meta["transform"][2],
                        self.ref_img.meta["transform"][0],
                        self.ref_img.meta["transform"][1],
                        self.ref_img.meta["transform"][5],
                        self.ref_img.meta["transform"][3],
                        self.ref_img.meta["transform"][4])

        driver = gdal.GetDriverByName("GTiff")
        raster = driver.Create("myraster", ncols, nrows, 1, gdal.GDT_Int16)
        raster.SetGeoTransform(geotransform)
        srs = osr.SpatialReference()
        srs.ImportFromEPSG(self.srid)  # My SRID
        raster.SetProjection(srs.ExportToWkt())
        raster.GetRasterBand(1).WriteArray(self.paths_array) # My Numpy array with data

By the way, I realize that it is sort of odd that I am taking the result of a rasterio raster and turning it into a gdal raster. The reason for this is that there is a very good algorithm that we found that is written using the Python GDAL bindings, and so we have to switch from rasterio to pure GDAL for this part

But it seems to be written for using the GTiff driver, which will create the file on disk. I thought I could avoid this by using the Memory driver:

    nrows, ncols = self.paths_array.shape

    # Explanation for this mess: https://rasterio.readthedocs.io/en/latest/topics/migrating-to-v1.html
    top_left = (self.ref_img.meta["transform"][2], self.ref_img.meta["transform"][5])
    geotransform = (self.ref_img.meta["transform"][2],
                    self.ref_img.meta["transform"][0],
                    self.ref_img.meta["transform"][1],
                    self.ref_img.meta["transform"][5],
                    self.ref_img.meta["transform"][3],
                    self.ref_img.meta["transform"][4])

    driver = gdal.GetDriverByName("Memory")
    raster = driver.Create("myraster", ncols, nrows, 1, gdal.GDT_Int16)
    raster.SetGeoTransform(geotransform)
    srs = osr.SpatialReference()
    srs.ImportFromEPSG(self.srid)
    raster.SetProjection(srs.ExportToWkt())
    raster.GetRasterBand(1).WriteArray(self.paths_array)

However, no luck. When I use the Memory as input to gdal.GetDriverByName(), then the result of driver.Create("myraster", ncols, nrows, 1, gdal.GDT_Int16) is None.


texts: Add metadata and coordinate system to tif file created with GDAL CreateCopy

texts:
I want to convert ascii grid data files from Deutscher Wetterdienst to GeoTIFF using GDAL in Python, but I have problems adding the correct coordinate system (since DWD uses their own RADOLAN grid its a bit problematic).

I'm using CreateCopy() which as far as I have read is using the gdal_translate command.
In gdal_translate I would just add the correct values with appending gdal_translate -a_srs="+proj=stere +lat_0=90.0 +lon_0=10.0 +lat_ts=60.0 +a=6370040 +b=6370040 +units=m" but I can't find a way to do it in the Python API.

Is there a way to add these values while translating, or adding them after creating the file using the dataset return of .CreateCopy()?


texts: how to start qgis with --configpath option on windows

texts:
how can I start QGIS on Windows 7 with the --configpath option?

I'm used to work with Linux where this procedure is very very simple (just write in a terminal qgis --configpath) with all the options.

Someone knows how to do the same thing on windows?


texts:
Newbie here. I've just loaded OSX Catalina on my iMac. I've installed GDAL 2.4 complete from the Kyngchaos site. The install works great for all GDAL commands except for the python ones (gdal_polygonize.py for example). Then I get this error:

-bash: /Library/Frameworks/GDAL.framework/Programs/gdal_polygonize.py:/Library/Frameworks/Python.framework/Versions/3.6/bin/python3: bad interpreter: No such file or directory

The python version I'm running is python 3.7.5

In my system there is no /Library/Frameworks/Python.frameworks folder. Instead, it seems to be here:
/usr/local/bin/python3.7

However, my PATH variable includes this folder:

/Library/Frameworks/GDAL.framework/Programs:/Library/Frameworks/PROJ.framework/Programs:/Library/Frameworks/SQLite3.framework/Programs:/Library/Frameworks/UnixImageIO.framework/Programs:/usr/local/bin:/usr/bin:/bin:/usr/sbin:/sbin

What am I doing wrong? how do I fix this?


texts: QGIS - exporting shapefile to convert to .JSON for use in PowerBI

texts:
I'm trying to export a shapefile I've been using in QGIS for use in PowerBI because PowerBI only supports .JSON extensions. So I've been converting them with mapshaper.org. This is my shapefile result in PowerBI:

I've tried a couple of things:

Exported my "custom" .SHP file as a .GeoJSON (I say custom because I deleted the geographies I did not need in QGIS Attributetable - maybe this is causing the problem?). I then converted it in mapshaper.org to .JSON. The shape appears fine in the mapshaper.org preview, but it's wonky when I upload it to PowerBI as per my attached image above.
Converted the original shapefile I downloaded from StatsCan. This has all geographies (most of which I don't use). When I load this into PowerBI nothing happens. It basically hangs the program. The file is roughly 92MB. Again, when I'm converting use mapshaper.org the map appears correct.

I know this a GIS forum and not PowerBI so I don't expect help on that end.


texts: Disappearing Pixels on the edge of Hillshade raster

texts:
When I've created Hillshade in ArcMap and received the output raster, I've noticed that output Hillshade has one pixel less than I have in the source DEM. It doesn't matter where I create it (ArcToolBox or arcpy). Aspect, Slope work fine without that artefacts.

I've tested this tool on ArcGIS Destop 10.4.1 and 10.6.1. May be somebody know how to solve this problem.

I've tried to do same in QGIS and discovered the option "Compute Edges" in Hillshade dialog (pic below) which solved the problem with edges. Is it possible to do same in ArcGIS (may be in script)? I need to processing a lot of rasters.


texts: Flooded roads and deprived communities: QGIS route analysis function to avoid flooded roads?

texts:
In QGIS I am analysing the flooding of road networks. I have created a vector line layer of all the stretches of roads in my study area where high risk flood extent areas intersect them. I now want to see which deprived areas would be most impacted by calculating the shortest route to essential services (avoiding these flooded roads).
I have polygon layers of regions of deprived communities and point layers of hospital locations.
Is there a plugin or function in QGIS that allows me to calculate the shortest route (and the time it takes, if possible) to the nearest hospital from the edge of a polygon, avoiding these flooded roads (line layer)?


texts: Example of RasterData class in ArcGIS for Flex

texts:
How would I use the RasterData class in com.esri.tasks.supportClasses for ArcGIS API for Flex? I haven't found any documentation for it outside the API and that doesn't give much info either. Could I construct a RasterData object and use it as a workaround for the limitations imposed upon a client's input parameters for a geoprocessor request to ArcGIS Server? For example,

var rasterOne:RasterData = new RasterData();
rasterOne.url = "http://example.com/rasterOne.tif";
rasterOne.format = tif;

var params:Object = new Object();
params.Input_Raster = rasterOne;

var gp:Geoprocessor = new Geoprocessor();
gp.url = "http://example.com/arcgis/rest/services/GPServer/ExampleGP";
gp.submitJob(params);

Assume that the GPService has one input parameter, a raster, and does a simple reclassify. Would this work, and would I even be able to publish the GP service with a raster as the input?


texts:
I am using an off-campus computer and need to add the Spatial Analyst to ArcMap in ArcGIS for Desktop.

Could someone please give step-by-step instructions on the order I need to install this?


texts: Determining if coordinate is within polygon?

texts:
I have a list of points describing the boundaries of Spain. I want to be able to tell whether a pair of lat,lon is within these boundaries. I have tried the following:

import shapefile
import matplotlib.pyplot as plt
from shapely.geometry import MultiPoint, Point, Polygon
from shapely.geometry.polygon import Polygon
sf = shapefile.Reader(r"\ESP_adm0.shp")
shapes = sf.shapes()
lat = []; lon = []
for i in range(len(shapes[0].points)):
    lon.append(shapes[0].points[i][0]);lat.append(shapes[0].points[i][1])

I know I am retrieving the points, because I'm able to plot and get the desired results:

plt.plot(lon,lat,'.', ms=0.1)

(plot in the link below)
plot result

I do the following to get the poitns into a polygon:

    coords = list(zip(lat,lon))
    spain_pol = Polygon(coords)

And then I use the contains function.

    spain_pol.contains(Point(0,42))
    spain_pol.contains(Point(42,0))

These both return false. In fact I haven't been able to get a single point I've tried to return a True.

I have tried all sorts of things, and I think I must be missing something fundamental.

Perhaps the facts Spain has islands and there's more than one polygon is the problem?


texts: Creating a fixed distance buffer from a virtual layer

texts:
First off, i'm trying to learn some python to automate scripts for basic tasks i carry out in QGIS (2.18).

Below I create a virtual layer (polyline) using the following script:

from qgis.core import QgsVectorLayer, QgsMapLayerRegistry
from qgis.utils import iface
from qgis.core import *
from qgis.utils import *
from qgis.gui import *

vlayer = QgsVectorLayer("?query=SELECT cm_code, cm_ndcode1 as 'nd_1', cm_ndcode2 as 'nd_2', COUNT() as 'Num CM', geometry FROM t_cheminement GROUP BY   cm_ndcode1, cm_ndcode2 HAVING   COUNT() > 1" , "CM_DUBLONS", "virtual" )

props = { 'width' : '3', 'color' : '255,0,0' } sl=QgsSymbolLayerV2Registry.instance().symbolLayerMetadata("SimpleLine").createSymbolLayer(props) s = QgsLineSymbolV2([sl])

vlayer.setRendererV2(QgsSingleSymbolRendererV2( s ))
QgsMapLayerRegistry.instance().addMapLayer(vlayer)

The above code works fine - but next I want to create a fixed distance buffer around the polylines.

I started off with this:

qgis_fixeddistancebuffer(input = CM_DUBLONS, output = 'c:/SIG_DEV_MODEL/buffer_test.shp', distance = 1000)

But it returns the error message below:

  Traceback (most recent call last):  File "< input >", line 1, in
   NameError: name 'qgis_fixeddistancebuffer' is not defined

I've also tried

processing.runalg('qgis:fixeddistancebuffer', CM_DUBLONS, 500, 6, dissolve, "C:/SIG/SIG_DEV_MODEL/CM_buf2.shp")

but likewise I get a similar error message "'processing' is not defined"

I'm probably missing something simple...

Otherwise I've since looked up the following:

Pyqgis Fixed Buffer Creation in Meters with Haversine

and

How to create dissolved buffer layer with pyqgis?

But it's not quite what I'm looking for, and I'm not sure how to adapt the script.


texts: Paste a feature from an existing shapefile into a row with null geometry in QGIS

texts: If I have an existing row in my PostGIS database with NULL geometry, is it possible to select a feature from another layer in QGIS, and paste that into the empty geometry column for the PostGIS row?


texts: I have tried to create field in TOC called with numbers with python. But, it has a problem: it couldn't do it because name's of fields are numbers. It is strange for me. When I do new field in shape file with mouse and put name field with number I haven't had problem.


texts: Need to find strings that contain the same word twice

texts:
I have attribute data with owner names. I need to select data that contains the last name twice.

For example, I may have an owner name that reads "BENNETT MCCARL & ARNETTE BENNETT".

I would like to select any rows in the attribute table that have a reoccurring last name such as the example above. Does anyone know how I can go about selecting that data?


texts: how to use "data defined override" in qgis 2.0?

texts:
How can I use the new "data defined override"-funktion in qgis 2.0? I did not find any further information about it in the handbook an somewere else.
I want to label my layer with different front styles for example:

if the field "type" = 'city' and 'town' -> the label have to be in size 12pt and bold

Thanks!


texts: Get geotiff Transform from point correspondences

texts:
I am trying to find a way to get the Geo Transformation (NW corner position, pixel size and rotation) for a given tiff image. I have a tiff which I know is using some coordinate reference system (an UTM projection using the WGS84 reference ellipsoid).

I would like to make it a geotiff using a software by associating, for example, point correspondences between the tiff image and a reference map. I would love to use something like Maptiler, but without the 10k*10k limitation on the image resolution.

For now, the solution I can think of is to downscale the input tiff, use Maptiler to get the GeoTransform, then use GDal to make this tiff a Geotiff by upscaling the pixelSizes by the inverse coefficient.

Are there any free alternatives to generate a geotiff by using correspondences?


texts: How to change page orientation and number of pages of print composer layout using pyqgis in QGIS 3.4?

texts:
I am writing a qgis layout script to export my results in a pdf. I am using the tutorial provided here but this tutorial do not teach how to turn the page orientation from landscape to portrait.

my code is like this, it works but creates with a landscape orientation

import os

p = QgsProject()
layout = QgsLayout(p)
layout.initializeDefaults()

# create a map and stretch it over the whole paper size
map = QgsLayoutItemMap(layout)
layout.addItem(map)

# picture
pic = QgsLayoutItemPicture(layout)
pic.setPicturePath("C:\\Users\\c.png")
pic.setLinkedMap(map)
pic.attemptMove(QgsLayoutPoint(1.4, 1.8, QgsUnitTypes.LayoutCentimeters))
pic.attemptResize(QgsLayoutSize(10, 10, QgsUnitTypes.LayoutCentimeters))
layout.addItem(pic)

# label
label = QgsLayoutItemLabel(layout)
label.setText("hello")
label.adjustSizeToText()
label.attemptMove(QgsLayoutPoint(15, 1, QgsUnitTypes.LayoutCentimeters))
layout.addItem(label)

# export the pdf
pdf_path = os.path.join(QgsProject.instance().homePath(), "output.pdf")
exporter = QgsLayoutExporter(layout)
exporter.exportToPdf(pdf_path, QgsLayoutExporter.PdfExportSettings())

How can i change the pages to portrait orientation. And how can i add new pages?


texts: Export AutoCad 2017 layer (with lines, blocks, and text) to GIS shapefile

texts: How can I convert (export?) a DWG block into a GIS shapefile so that the block's insertion point will become a point in the shapefile and the block's text attribute will carry over as an attribute in the shapefile?


texts: Filling one field based on another (same attribute table) with Python

texts:
I have one field with a list of integers. I want to fill another field with the sum of those integers. When I print the results on the Python console, everything is fine. But I don't know how to fill every cell: when I try return, all the features have the same result as the first one. I understand why that happens but I cannot figure out what is the right way. I've spent a while now looking for examples but no one corresponds exactly to what I'm looking for.
In the code below "ESPESSURA" is the field with the list of integers (ex: 40,67), attrs[13] the field I want to fill.

from qgis.core import *
from qgis.gui import *
from qgis.utils import iface
@qgsfunction(args=0,group='Custom')
def sum_list_str(value1,feature,parent):
    lyr = iface.activeLayer()
    features = lyr.getFeatures()
    idx = lyr.fields().lookupField('ESPESSURA')
    for feat in features:
        attrs = feat.attributes()
        attrs[13]=sum([int(i) for i in attrs[idx].split(',')])
        #return(attrs[13]) -> takes the first result and fills everything with it
        #print(attrs[13]) -> prints the values I'm looking for
        #lyr.updateFeature() -> QgsVectorLayer.updateFeature(): not enough arguments
        #lyr.updateFeatures() -> 'QgsVectorLayer' object has no attribute 'updateFeatures'
    #lyr.commitChanges()

(I'm pretty sure there are a lot of inelegant things in the code. Any comment is welcome.)


texts:
Using gdal_grid for spatial interpolation I have run across a problem: the nodata argument is ignored. The only solution I have found for this is to then use gdal_translate but it is rather cumbersome (see below). My bash scripting skills are very, very minimal. Can anyone suggest an improvement, like passing the outputted file directly to the next tool?

gdal_grid -zfield "D2012" \
-a invdist:power=2.0:smoothing=0.0:radius1=0.0:radius2=0.0:angle=0.0:max_points=0:min_points=0:nodata=-9999 \
-txe 647822 651672 -tye 7535500 7537800 -outsize 77 46 \
-of GTiff -ot Float64 -l probe2004_2012 probe2004_2012.shp ./Output/temp.tiff
gdal_translate -stats -a_nodata -9999 -a_srs EPSG:3006 ./Output/temp.tiff ./Output/2012_IDW_2_50.tiff
rm ./Output/temp.tiff


texts:
I have a health data set and an "urban stress". I'm trying to do bivariate analysis within QGIS using Joshua Stevens lessons at http://www.joshuastevens.net/cartography/make-a-bivariate-choropleth-map/.
As per Joshua's instructions, in QGIS the best way to do this is to write a CASE function to classify values into a preset organizational value.
For example:

CASE WHEN  "CrimeRate" > 4.0 THEN 3
    WHEN "CrimeRate" <= 4.0 AND "CrimeRate" > 2.0 THEN 2
    ELSE 1
END

This field calculator function worked perfectly.
The next example filled my column with NULL values:

CASE WHEN  "POOR345" > 52.9 THEN 'C'
    WHEN "POOR345" <= 52.9 AND "POOR345" > 40.3 THEN 'B'
    ELSE 'A'
END

Does anyone have an idea of why this is happening?


texts: How to create an annual NDVI from Landsat 8 imagery

texts:
I am relatively new to GIS.

I am studying wildlife (I have collar data that lasts a year) and my supervisor and I would like to create an annual NDVI map for my projected area. I was able to download a years with of imagery (although some does have cloud coverage) ranging from the year that I am working with.

Now, from what I understand - with one Landsat image (lets say representing a day out of a specific month). I choose Bands 2, 3, 4, and 5. I then use image analysis (in Arcgis) to composite those bands and calculate NDVI from there.

However, I would I do this with multiple landsat images that span across a year?
do I composite the bands and calculate the NDVI separately one by one and then merge the separate NDVI rasters?
Or do I simply choose bands 2, 3, 4, and 5 for all of the separate landsat images (so for example, bands 2-5 for april, 2-5 for may, etc) and then make one composite (representing all of the months) from there?

I tried looking for tutorials on how to composite and annual NDVI but to no avail.

I calculated my individual monthly NDVIs (had some cloud coverage but not directly on the mountain that I am working with). Every individual NDVI had values between -1 and 1 (looking good so far).
Read into cell statistics and decided to take all of these separate NDVIs and plug them in with the overlay being the mean.

Results: http://i.imgur.com/RXVS23Z.jpg (high being red and low being green).
Also, results with a separate raster (of the forests found in the mountains) overlaid and put on transparency for comparison.
http://i.imgur.com/aYRqSoc.jpg
The only thing that I am concerned about is how my values (from the cell statistics output) range from 0.53317 to - 0.232147 and not between 1 and -1


texts: A WMS layer does not seem to load correctly in projection

texts:
No idea how to ask this I am quite new to GIS.
I am trying to load this WMS service https://geoservices.informatievlaanderen.be/raadpleegdiensten/GRB-basiskaart-grijs/wms with EPSG:3812 but every time it is cut off halfway a province

The reason I need EPSG:3812 is that I want to put a layer on top off it and that layer only is available in EPSG:3812. If I add the above layer with EPSG:3856 for instance it loads perfectly but then the EPSG:3812 layer has jagged edges and I want as sharp an image as possible so I guess they both have to be the same CRS?


texts: Modelling stream in 3D based on known depths across it using ArcGIS Desktop?

texts:
I have an excel spreadsheet that contains the elevations across the profile of a stream. I was wondering how to model this in 3-d.
 Basically there are 6 monitoring stations that are within a stream. I know the elevation of each bank, the slope of each bank, and the depth across the stream at each location. The depth was measured every 20cm across the profile of the stream.
eg/ 0.35m across/depth 0.4m, 0.55m across/depth 0.95m, etc. I know that because there are only 6 stations, the result will probably be very blocky, but there are ways to fix that. Basically, I am just looking for somewhere to start. Like I said, I am pretty sure I have all of the pertinent information in an excel spreadsheet, just need to know how to get going.


texts:
I am trying to show a 4 panel map, each with a different dataset. So far I have this, but the maps do not seem to be linked and the layers are not loading.

Does anyone know what I am doing wrong?

    //////////NO2///////////
var NO2 = ee.ImageCollection('COPERNICUS/S5P/NRTI/L3_NO2')
            .select('tropospheric_NO2_column_number_density')
            .filterDate('2018-01-01', '2020-05-30')

// Make predefined data layers that can be selected.
var palettes = require('users/gena/packages:palettes');
var palette = palettes.colorbrewer.BrBG[11].reverse();

var NO2_band_viz = {
  min: 30,
  max: 300,
  palette: palette,
  opacity: 0.8
};

///population
var pop = ee.ImageCollection("WorldPop/GP/100m/pop")

/////nighttime lights/////
var dataset = ee.ImageCollection('NOAA/VIIRS/DNB/MONTHLY_V1/VCMSLCFG')
                  .filter(ee.Filter.date('2017-05-01', '2020-05-31'));
var nighttime = dataset.select('avg_rad');
var nighttimeVis = {min: 0.0, max: 60.0};

// Create a map for each visualization option.
var map1 = ui.Map()
map1.add(ui.Label('NO2'))
map1.addLayer(NO2, NO2_band_viz, 'NO2')
map1.setControlVisibility(true);

var map2 = ui.Map()
map2.add(ui.Label('Population'))
map2.addLayer(pop, {}, 'Population')
map2.setControlVisibility(true);

var map3 = ui.Map()
map3.add(ui.Label('Nighttime Lights'))
map3.addLayer(nighttime, nighttimeVis, 'Nighttime Lights')
map3.setControlVisibility(false);

var map4 = ui.Map()
map4.add(ui.Label('NO2'))
map4.addLayer(NO2, NO2_band_viz, 'NO2')
map4.setControlVisibility(false);

var linker = ui.Map.Linker(map1, map2, map3, map4);

// Enable zooming on the top-left map.
map1.setControlVisibility({zoomControl: true});

// Show the scale (e.g. '500m') on the bottom-right map.
map3.setControlVisibility({scaleControl: true});

// Create a grid of maps.
var mapGrid = ui.Panel(
    [
      ui.Panel([map1, map2], null, {stretch: 'both'}),
      ui.Panel([map3, map4], null, {stretch: 'both'})
    ],
    ui.Panel.Layout.Flow('horizontal'), {stretch: 'both'});

// Center the map at an interesting spot in Greece. All
// other maps will align themselves to this parent map.
map1.setCenter(21.2712, 38.4151, 12);

/*
 * Add a title and initialize
 */

// Create a title.
var title = ui.Label('Title', {
  stretch: 'horizontal',
  textAlign: 'center',
  fontWeight: 'bold',
  fontSize: '24px'
});

// Add the maps and title to the ui.root.
ui.root.widgets().reset([title, mapGrid]);
ui.root.setLayout(ui.Panel.Layout.Flow('vertical'));

Link to code here: https://code.earthengine.google.com/59c837c11e763ae298097a4a15911c08


texts: How to get distance values from nearest neighbor query for 5 closest features to multiple records in PostGIS?

texts:
I found this query on this blog https://geeohspatial.blogspot.com/2013/05/k-nearest-neighbor-search-in-postgis.html?showComment=1480796271673#c3250651227513965563 which returns the 5 closest records in one table to each record in another table:

SELECT l.gid, (SELECT ARRAY(SELECT c.city_name
  FROM bndry.us_cities c
  ORDER BY l.geom <#> c.the_geom LIMIT 5))
FROM offshore_meta.load_centers l;

This works great for my needs. The only other thing I would like it to be able to do is to return an array of the distance values in a second column for each item in the first array containing the names. I'm assuming I can use ST_Distance but am not sure how to apply it in this case. Can anyone help?


texts: QGIS - CRS for a fictional planet map (raster) - What am I doing wrong?

texts:
I have a procedurally generated DEM heightmap (1600px * 800px) showing an entire (fictional) planet in equirectangular projection (= equidistant cylindrical) that I uploaded in a QGIS raster layer with layer CRS set to [ESPG:4087- WGS 84]. My project CRS is set to [EPSG:4326 - WGS 84]. For the moment I consider that this fictional planet has the same dimensions as Earth (as to avoid having to create a custom CRS).
I want to set a CRS for this map in a way that:

The coordinates in the lower right corner show degrees going from -180°...+180° east-west, and +90°...-90° north-south.
The "measure line" tool shows accurate distances that vary with the latitude.

I've tried obtaining this result many ways but haven't succeeded yet (I've tried methods like writing a word file myself, write a custom CRS in proj4, using an external map reprojector, etc... ). The closest I've gotten is by georeferencing the four corners of my map onto coordinates (-180,90); (-180,-90); (180,90) & (180,-90) of the [ESPG:4326] CRS (so basically "pasting" this map on an Earth map):
==> This gives me coordinates in degrees like I wanted, and when I measure the length of the equator with the measure line tool, it's about 40000km, so it seems to work, right? Nope. When I measure the vertical sides at the poles, it also gives me 40000km... Since all meridians converge at the poles, the distance measured there should be zero.
So how could I solve this? One problem is that I don't know how the "measure line" tool actually works on the fundamental level:

Does it it just measure pixels and convert them linearly to a distance in meters? In that case I should use a CRS that conserves distances (thinking about the azimuthal equidistant projection, but I'm open to alternatives);
Or can it actually take into account the fact that the map layer is in an equirectangular projection and adapt the measured distances accordingly? I have no idea however how I would/should do this.

Another problem is that I'm quite new at QGIS and all this CRS stuff (project CRS vs. layer CRS, projected vs. geographic CRS, custom CRS's, ...) is confusing.

Note that the map looks like this in appropriate software:


texts: Print values of a field in arcgis using python


texts:
I've install Geoserver into Tomcat server 8.5, and since I've published new layer from PostGis data source, Geoserver show this message:

Server is currently busy. Please wait...

How can I stope it or fix it? On my server running Windows.


texts:
I have a point shape (ID, value) which I have converted to a path (points to path). Is it possible to get the lines between each point? The picture below shows what I mean. The colored path is what I want to get.

Is there a tool or workaround?


texts: Zooming to first queried feature using ArcObjects?

texts:
I'm attempting to query a feature class for all features that meet a certain criterion (in this case, the Status field is NULL). However, I'm not getting very far into the attempted method before ArcMap crashes completely. I have no reasonable idea on how to trace errors that cause a program crash, for starters, so I'm more confused than necessary, I expect. I'm building this in VS 2008 Express, for what it's worth.

 public void PerformAttributeQuery(ESRI.ArcGIS.Geodatabase.IFeatureClass fc)
    {
        ESRI.ArcGIS.Geodatabase.IQueryFilter queryFilter = new ESRI.ArcGIS.Geodatabase.QueryFilterClass();
        queryFilter.WhereClause = "[Status] IS NULL"; // create the where clause statement

Everything is fine for that first couple of statements, but this query causes a crash:

        // query the table passed into the function and use a cursor to hold the results
        ESRI.ArcGIS.Geodatabase.IFeatureCursor featurecursor = fc.Search(queryFilter, false);

It crashes before reaching the next line (I put in debug message boxes for testing), but I'll put it here for propriety:

        ESRI.ArcGIS.Geodatabase.IFeature feature = (ESRI.ArcGIS.Geodatabase.IFeature)featurecursor.NextFeature();
        ESRI.ArcGIS.ArcMapUI.IContentsView currentContentsView = ArcMap.Document.CurrentContentsView as IContentsView;
        ESRI.ArcGIS.Carto.IFeatureLayer featureLayer = (ESRI.ArcGIS.Carto.IFeatureLayer)currentContentsView.SelectedItem;
        ESRI.ArcGIS.Carto.IFeatureSelection fSelection = featureLayer as ESRI.ArcGIS.Carto.IFeatureSelection;
        fSelection.Add(feature);
        FindCommandAndExecute(ArcMap.Application as ESRI.ArcGIS.Framework.IApplication, "{AB073B49-DE5E-11D1-AA80-00C04FA37860}");
    }

I cannibalized the PerformAttributeQuery snippet to work with FeatureClasses instead of Tables (or so I had hoped), if you'd like a proper reference. My questions are:

How can I get better reporting of crash errors?
Where did I go wrong? Answer: I must shamefacedly admit that I used the wrong characters to denote the location of the feature table I was querying. Changing the [Status] in my queryFilter.WhereClause to \"Status\" quite literally solved all obvious problems.


texts: Creating a point corresponding to the maximum elevation of a DEM

texts:
I am wondering how can I create a point corresponding to the maximum elevation value of a DEM. I mean, from the DEM I can easily get the maximum elevation. For the sake of a viewshed analysis, I would like to create a point feature corresponding to the spot with maximum elevation. I guess it is simple, but I am having difficulties in finding the right procedure.

I am using arcgis 10.1 (with advanced license).


texts: CSS Spritesheets and picture marker symbols in the AGS Javascript API

texts:
I'm working on improving the load performance of a large ArcGIS Javascript application, and one area where we can have trouble is that there is a visible delay between putting a graphic with a PictureMarkerSymbol on the map and seeing the symbol load, particularly when there are a large number of features with various symbols.

Is there a way to pull picture marker symbols or something similar from a CSS spritesheet to get all of the images loaded as a group?


texts: Removing OUTER boundary of group of polygons?

texts:
I'm trying to figure out how to remove only the outer boundary of a group of polygons. For instance, imagine I have a layer containing the map of the United States, where each state is its own polygon/attribute. I want to create a new layer that consists only of the boundaries BETWEEN states. In other words, all of the coastline and borders would be gone, but the state grid would still be present (only the coastal or border states would be affected). Is there any way of doing this?

What I'm trying to do is replace an inaccurate boundary edge with a more accurate boundary edge. The problem is that my accurate boundary edge doesn't include the interior polygon edges. So I'm trying to take the interior polygon edges from one layer and combine it with the accurate boundary edge of the other. Unfortunately the answer to the question the way I asked it results in multiple lines that no longer include the original attribute data for the original polygons. I now have a feeling I'm going about it in entirely the wrong way!

I ended using a combination of what Nir and Midavalo recommended. I used Polygon to Line (instead of Feature to Line) to turn my "state" and boundary layer into lines, then selected by attribute (Left_ID = -1) to remove the outer boundary. Then I copied and pasted my more accurate boundary line into this layer, and used Feature to Polygon to convert the new combines lines layer into Polygons. Then since I only had about 12 or so polygons I manually entered the appropriate attribute data into the new polygons. Problem solved.


texts: Seeking USA CLU shapefiles for agricultural land?

texts:
I know Is geocomm.com permanently unavailable? was asked a few days back and is closed. That is the site I was hoping to buy from initially.

Data providers equivalent to USDA/NRCS/FSA Common Land Unit (CLU)? from 3 years also points to the geocomm site (which is now defunct).

Is it possible to find the 2008 FSA version of the CLU data which was the last version that was publicly available?


texts: Using remote sensing satellite for chlorophyll mapping?

texts:
Location of my research are reservoirs which have an area of 47 and 62 km2.

I need a better image resolution than Landsat(Pan).


texts: Extracting all raster values from buffer using R?


texts: What language do QGIS expressions use, what language should QGIS queries be written in?

texts:
This part of the docs says it supports only a small part of the SQL syntax but does that mean it's SQL based or am I better writing expressions in something else, MySQL?

This is similar to What expressions for QGIS field-calculator? but my query was not limited to the field calculator, it also applied to the filter and label options. I understand from the answer to that question that QGIS has it's own language for the field calculator is that the same for label by expression and filter by expression. I generally try SQL based syntax on these as a first approach with support from here when those fail.


texts: ArcGIS Select Tool: Do Not Create Empty Feature Class

texts:
I am creating an ArcGIS Model that imports a number of Bentley .dgn files to a geodatabase.  My workflow involves:

Running the CAD to Geodatabase Tool, which brings all the drawing features into a feature dataset.
Taking the feature types and then selecting features by the "layer" attribute and writing them to a feature class with an appropriate name.

My issue occurs at the Select tool.  It will create a feature class even if it did not find any features to place in the feature class.  This creates a whole bunch of empty feature classes since there are a number of layers which may or may not show up in any of my .dgns.  So, I would like to know if there is a way to keep my model from creating blank feature classes?


texts: Dividing Polygon into parts which have equal area using R


texts: Multiuser editing coupled with Google Street view?

texts:
I am currently looking for an interface with two essential features:

Ability for multiple users to edit data simultaneously
Ability for all users to use Google Street View within interface (As an add-in, for example)

My clients have to manually edit data for thousands of buildings based on Google Street View visual information. However they have minimal notions of GIS and I need to provide an editable interface where several users can log in and edit this data simultaneously.

I know ArcGIS Explorer has a Google Street View Add in, however does it support the multiuser environment that I am looking for?


texts: Openlayers: Parsing GeoJSON-properties results in a "data"-property in feature

texts:
I'm using OpenLayers 2.12.

I try to load a GeoJSON-response from server and (nearly) everything works fine: but the properties contained in the GeoJSON-response ("nummer" and "befund" below) were put in a data-attribute of the feature-objects, not in their attributes-property, which I would expect. As observed in firebug, in the attributes-property I only find the names of my properties and empty strings as values.

Did I miss something, some parameters when creating the formatter or something like that? It's ugly to copy the values in every feature back to attributes...

My parsing function:

function parseGeoJson(featuresText) {
    var formatter = new OpenLayers.Format.GeoJSON();
    var features = formatter.read(featuresText);
    if (features.length > 0) {
        editor.editLayer.addFeatures(features);
    }
}

My GeoJson:

{"type":"FeatureCollection",
"features":[
    {"type":"Feature",
     "properties":
          {"nummer":"1",
           "befund":"1.1.1"},
     "geometry":{
          "type":"Polygon",
          "coordinates":[
              [[15.37109375,81.11328125],[65.29296875,86.38671875],
             [36.11328125,58.7890625],[15.37109375,81.11328125]]]
      }
    }]
}


texts: Cut shapefile to raster in Python using GDAL/OGR

texts:
I'm using a huge shapefile containing the coastlines of the world:

And want to cut out polygons from the shapefile to fit to my geotif like so (the blue polygons represents the cut out ones):

I want to do this in order to mask the land afterwards (and maybe extent the covered area if it doesn't properly cover land).

How does one do this in Python using GDAL/OGR libraries (performance do matter)?

The shapefile consists of many polygons and I'm not sure whether merging all of this to one feature will speed up performance?

One could get the bounding box of the raster, but how to go from there I'm not sure of:

# Raster image to clip
InputImage = 'Geotif.tif'

# Polygon shapefile used to clip
Shapefile = 'worldcoastlines.shp'

rasterFormat = 'GTiff'
vectorFormat = 'ESRI Shapefile'

# Open datasets
raster = gdal.Open(InputImage, gdal.GA_ReadOnly)
wktproj = raster.GetProjectionRef()
sref = osr.SpatialReference()
sref.ImportFromWkt(wktproj)

ulx, xres, xskew, uly, yskew, yres  = raster.GetGeoTransform()
lrx = ulx + (raster.RasterXSize * xres)
lry = uly + (raster.RasterYSize * yres)


texts:
Is there a way to harvest data from the Eudat B2Share portal in GeoNetwork according to the CSW ISO19139 format?

Currently the only working protocol is from OAI/PMH with prefix oai_dc


texts: Strange exposure of a PostGIS layer in QGIS after importing OSM data into a database

texts:
I've imported a osm layer into my postgres/postgis database.

osm2pqsql -E 3857 -c -S /usr/share/osm2pgsql/default.style -H localhost -W -d landslide -U postgres oberpfalz-latest.osm.bz2

During the import there were several warnings like this:

Processing: Node(3261k 69.4k/s) Way(399k 7.27k/s) Relation(80 26.67/s) Standard exception processing way_id 33907: TopologyException: side location conflict at 1396469.99 6361915.7599999998

Accept of this, I think the import runs without any problems.

After that I added planet_osm_roads to my qgis-project (EPSG:3857).

The result is shown in the picture below. I tried to change the SRIDs, but no change.

Any ideas what this strange exposure comes from?


texts: Running network analyst for closest facility anaysis in ArcMap?

texts:
I’m trying to carry out closest facility analysis in ArcMap 10.2.2 on the network analysis layer created from the OS integrated transport and urban path networks. But when I click "solve" in the network analyst bar, to generate the routes between the incidents (postcode centroids) and facilities (green space access points), its failing to calculate routes for all my incidents.

An error window comes up once the analysis in complete listing hundreds of the incidents which it says are on traversable positions of the network e.g. "Location BS138RG in incidents is on a non-traversable network element position".

It repeats this message for 1715 of my incidents, then the list goes on to list 2716 postcodes which it says no facilities can be found for, I assume because of the incidents being on non-traversable parts of the network. All my access points and postcodes were originally located fine when I loaded their locations, shown from the image below. And routes have been calculated fine for about 85% of the incidents, but I can't work out why it's not working for them all, many of the facilities and locations now showing up as errors on the map.

I didn't set any restrictions in the analysis settings of the network analysis layer, so assumed this would mean that all parts of the network should be traversable. However, because I set drive time as the impedance, perhaps the problem's that in the speed profile for this impedance, the speeds are set to 0 mph along networks which aren't accessible by cars e.g. paths, bridleways. In which case, perhaps there's no way around some of the incidents being on non traversable parts of the network when calculating drive times, since cars are never going to be able to drive along these types of networks! Or perhaps the problem's that most of my facilities and incidents aren't actually situated on the network. But I thought selecting to "snap to position along network" when loading the locations of each of them, would prevent this being a problem.


texts:
I am new to remote sensing, so I would want to clarify my understanding of the meaning of Multi Temporal Images.

As far as I understand, multi temporal images are multiple images of the same  scene acquired at different times.

Is there more to their defintion, or are multitemporal images just images of a scene X at two different times, t1 and t2?


texts:
I have been using the zonal statistics tool to quantify the variability of slopes around archaeological sites.

I have noticed in some experiments that there are differences between the mean and the median of the values within the polygons that I use. We want to use the mean slope to calculate biomechanical indices in burials of each of the polygons.

Since there is a difference between the mean and the median of each polygon, it can indicate that there is no normal distribution in the slope values sampled. In that sense, using the mean and the standard deviation would be somewhat misleading.

Is there any way to calculate the interquartile range of the values of the pixels inside a polygon? I am using QGIS 3.6 on Windows 10 with 32 GB of RAM.


texts:
I am creating point symbolizer with external graphic, and it works great if I load an image/png, but when I try to load svg file I get gray squares on map

StyleFactory sf = CommonFactoryFinder.getStyleFactory();
StyleBuilder sb = new StyleBuilder();
Graphic graphic = sf.createDefaultGraphic();
PointSymbolizer pointSymbolizer = sf.createPointSymbolizer(graphic, "the_geom");
File file=null;
if(tacka.equals("Trigonometrijska tacka")){
    file = new File(PointStyle.class.getResource("/trigonometrijska.svg").getFile());
}else  if(tacka == "Poligonska tacka"){
    file = new File(PointStyle.class.getResource("/poligonska.svg").getFile());
}else if(tacka == "Linijska tacka"){
    file = new File(PointStyle.class.getResource("/linijska.svg").getFile());
}
else if(tacka == "Detaljna tacka"){
    i=1;
}else{
    i=1;
}
graphic.setSize(ff.literal(10));
ExternalGraphic external =
        sb.createExternalGraphic( "file:///" + file, "image/svg+xml");
graphic.graphicalSymbols().clear();
graphic.graphicalSymbols().add(external);

Mark mark1 = sf.getCircleMark();
mark1.setFill(null);
Stroke stroke = sf.createStroke(ff.literal(JSimpleStyleDialog.lineColor), ff.literal(JSimpleStyleDialog.lineWidth));
mark1.setStroke(stroke);
if(i==1){
    graphic.graphicalSymbols().clear();
    graphic.graphicalSymbols().add(mark1);
}

I suppose that something is wrong about format but can't figure out what.

<sld:Graphic><sld:ExternalGraphic>
<sld:OnlineResource xmlns:xlink="http://www.w3.org/1999/xlink" xlink:type="simple" xlink:href="file:///F:\workspace\AGGF_Geotools\target\classes\trigonometrijska.svg"/>
<sld:Format>image/svg+xml</sld:Format>
</sld:ExternalGraphic>
<sld:Size>10</sld:Size>
</sld:Graphic>


texts: GEE, reduceRegions vs reduceRegion, User memory limit exceeded

texts:
I am new in GEE.

I am trying to reduce a mean statistics over an image with reduceRegions.
I increased tileScale to 16 (max value possible), decreased the number of polygons, exported to asset for performing long running operations,
and the error 'User memory limit exceeded' kept appearing.

Than, I found (https://developers.google.com/earth-engine/debugging?hl=en#reduceregion-vs.-reduceregions) map a reduceRegion function divides work into smaller tasks. I set bestEffort as true and maxPixel parameter. But still the error keeps appearing. The scale is the only thing I would like to keep, 10 m of res.

The code works if I run it for two polygons but when the number increases I get the error. Do you have any advise?

Here is the link of my script:
https://code.earthengine.google.com/e1e80f9ca804651b6c1b3a3614d7e483


texts:
I have 365 rasters with the same extent, and I am trying to sum them together by month. I have the script below, which was working until I added the line out1 = Con(IsNull(raster), 0, raster) which is supposed to get rid of any nodata cells in each raster I am processing. Now my script keeps raising the exception:

NameError: name 'Con' is not defined.

I am using the Con function exactly as it is used in all of the ArcGIS help documentations, so I understand why it doesn't recognize Con.

Anyone have any ideas as to what is causing this error? My script is below.

#script to create scaled habitat maps

#import  modules
import arcpy

#check out spatial analyst extension
if arcpy.CheckExtension("Spatial") == "Available":
    arcpy.CheckOutExtension("Spatial")
else:
    arcpy.AddError("Spatial Analyst extension is unavailable. Please check to see if you have the Spatial Analyst License.")

from arcpy import env
rasters = "Z:/sharks/Models/Summer2013/RF/ALL/SB_RF/NSMPpresRFPdh/AR_NSMPpresRFPdh/outputRasters"
env.workspace = rasters
listRasters = arcpy.ListRasters()
listRasters.sort()
listRasters
janRasters = listRasters[0:30]
febRasters = listRasters[31:58]
marRasters = listRasters[59:89]
aprRasters = listRasters[90:119]
mayRasters = listRasters[120:150]
junRasters = listRasters[151:180]
julRasters = listRasters[181:211]
augRasters = listRasters[212:242]
sepRasters = listRasters[243:272]
octRasters = listRasters[273:303]
novRasters = listRasters[304:333]
decRasters = listRasters[334:364]

rasterList = [janRasters, febRasters, marRasters, aprRasters, mayRasters, junRasters, julRasters, augRasters, sepRasters, octRasters, novRasters, decRasters]
nameList = ["janHabitat.img", "febHabitat.img", "marHabitat.img", "aprHabitat.img", "mayHabitat.img", "junHabitat.img", "julHabitat.img", "augHabitat.img", "sepHabitat.img", "octHabitat.img", "novHabitat.img", "decHabitat.img"]

counter = 0
for list in rasterList:
    i = 0
    print ("processing list: %s" % (list,))
    for raster in list:
        print ("processing raster: %s" %(raster,))
        #convert nodata to zero
        out1 = Con(IsNull(raster), 0, raster)

        #sum rasters together
        if i == 0:
            outSUM = arcpy.Raster(raster)
            i += 1
        else:
            outSUM = outSUM + raster
            i += 1

    #save final output to the current workspace
    outSUM.save("Z:/sharks/Models/Summer2013/RF/ALL/SB_RF/NSMPpresRFPdh/AR_NSMPpresRFPdh/monthlyHabitat/" + nameList[counter])
    counter +=1


texts: Vector seems to be only a pixel of a raster after reprojection

texts:
I have the South America hydro1k data for slopes and I am trying to extract the data using point sampling tool with a grid in Brazil's map.
My problem is that, even though vector and raster seem to overlay -- and I have tried many CRS conversions in both the vector and the raster with no success -- when I use the point sampling tool, all the grids get the same value of band from the raster. It seems that the grid map corresponds to just one pixel in the raster. This probably has something to do with the projections, but I haven't been able to solve the problem. Here is the projection data for hydro1k raster

PROJCS["US National Atlas Equal Area",GEOGCS["Unspecified datum based upon the Clarke 1866 Authalic Sphere",DATUM["D_Sphere_Clarke_1866_Authalic",SPHEROID["Clarke_1866_Authalic_Sphere",6370997,0]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Lambert_Azimuthal_Equal_Area"],PARAMETER["latitude_of_origin",-15],PARAMETER["central_meridian",-60],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1]]

and my Brazil's map is in SIRGAS 2000.

This is what happens when I have the raster in the original projection and Brazil vector in SIRGAS2000.I am posting Brazil's map without the grid just for illustration purposes.

And this is what happens when I reproject Brazil vector to the raster projection.

That is, Brazil vector is much smaller than south america raster.With the sampling tool, all the grid points have the same value for the band.
Does anyone has any idea of how to fix this?


texts:
In Windows GDAL terminal (in Windows, GDAL has a specific terminal/prompt) I usually do loops and it works like a charm.

What I'm trying to do the same in a Shell terminal in Ubuntu with the following syntax:

for %i in (*.tif) do gdal_translate -of KEA %i %~ni_new.kea

Considering I'm in the right folder, I don't need to specify any path. I know Ubuntu (trusty in my case) has a different syntax but even after several tries, I can't figure out how it works in Ubuntu terminal.

Any ideas?


texts:
I do

from qgis.core import *
layer = QgsVectorLayer(r"D:\file.shp", "boundingBoxes", "ogr")
layer.featureCount()

outcome  is 428

layer.select(layer.pendingAllAttributesList())
layer.selectedFeatureCount()

outcome is 0

Why the results are not equal?


texts: How R-Trees, R*-Trees can be used to index the Road Network?

texts:
I am trying to understand the Tree Indexing method of the Road network. In R-Tree the Minimum Bounding Rectangle is used to store the data points, but how it can hold the line segments of the Map? Or each MBR will just hold the single edge of the roadmap? From a basic Spatial Join what I got is that the connected line segment's rectangle overlaps.
Still it seems a lot of confusing to me. Could anybody please help me with the simple yet easy explanation. A graphical visualization would be more helpful.


texts: Bands disappear after Raster Calculation ArcGIS

texts: I have a tif file containing 36 bands. After I applied raster calculation on the file the bands aggregate into one single band which I don't want because I need them for further calculation. People say that I can apply raster calculation on each band. I can do that but it is too troublesome .Is there any better way to keep all the bands as well as performing calculations on them?


texts:
Is there a way to close a geoprocessing operation when the dialog box has been hidden/pushed to the back of the windows? I ran a python process as VB Script and when the error message appeared I misclicked. The error popup was hidden behind the attributes table.

I finally terminated ArcMap.

Is there an easier solution?


texts: ArcSDESQLExecute to query and update from multiple SDE Databases

texts:
I am trying to use ArcSDESQLExecute to return the results of a join and then write the return query into a Feature Class located in a different database. Using Microsoft SQL Server and Arc10.2.2

Is there any way to use ArcSDESQLExecute to UPDATE a table located in one SDE database from the query results of tables/feature classes from another SDE_database? Would using a linked server work? Or is ArcSDESQLExecute limited to manipulating tables within the same database?


texts: Connecting to WMTS of GeoServer from QGIS?

texts:
I am trying to configure a WMTS using the version 2.5.2 of Geoserver. It seems to be quite straightforward.

The problem comes when I try to connect to the WMTS using QGIS 2.6.0. I have access to the Getcapabilities but when I add a layer the image is empty.


texts: display large amount landmarks in different scales in ArcGIS Server

texts:
I have large amount of Landmarks, I want to display on map using ArcGIS Server. I want to display these landmarks based on the scale, some thing like

 Scale 1:1000000 display 3-10 landmarks only
 Scale 2:500000 display 10-20 landmarks only
 Scale 3:100000 display 20-100 landmarks only
 Scale 4:5000 display 100-300 landmarks only

The displaying would be random landmarks or based on some filtering. its something like using QuadTree, or some thing like google. as far you zoom more landmarks appears.
It is also some thing like Maplex labeling (buffer labeling), but this time feature displaying not labeling.

How I can do that without to have multiple copy of layer with different scales and different queries


texts: Joining QGIS Property Parcel Mapping shapefile to an Excel file and displaying symbol in Parcel Mapping Shapefile

texts:
I am new to QGIS. I am a land surveyor and I have all of the property parcel mapping information in a shape file for my area.  I am trying to link/join an Excel file which has all of the survey files indexed which contain survey information and plans for various properties.
The parcel mapping file has  a field entitled PID.... the excel file has a file called PID and also a field called Plan No.
I am not sure of the correct terminology but I would like to join the two, and have a symbol displayed in the parcel mapping which would represent that I have a Plan No. for that parcel.
At this point I think I have joined the Parcel mapping to the Excel file, but do not know what to do next.


texts:
I am using QGIS 1.8.0 on Windows 7.
In nviz (via SEXTANTE), after selecting an altitude layer and the GRASS region extent, I hit OK; the progress bar quickly reaches 100% and then absolutely nothing happens. I tried loading other layers, including additional vector ones and changing the extent.
I have also reinstalled QGIS using OSGeo4W, cleaned the registries and reinstalled the plugin. What am I missing?


texts: set cells with specific values as NoData in afloating raster

texts:
I have a raster file which includes floating values. I need to set cells with specific values as Nodata. For Raster integer values, the setnull tool can be applied, is it any tool for floating point?

The image that I have problem with it is the outcome of linear spectral analysis of ENVI software. I imported the image from ENVI to ArcMap by using ASCII to raster tool.

I have tried the Con in the raster calculator of ArcMap but It did not change my image and I got my original image again as output. I need to set the values less than zero to Nodata.


texts: After two buffers operations (+/-), why is the result file bigger than the original?

texts:
First file ( .shp ) had 1000KB, after 2buffers ( buffers on polygon, 1st buffer +100,
2nd buffer -100 ( small things near polygon should be removed ))
in quantum gis file grown to 1500KB, and after the same 2 buffers in postgis file has 20000KB.
Anyone can tell me why?


texts: I have upgraded from ArcGIS 10.0 to ArcGIS 10.3.1 and I noticed that I lost my Folder Connections in the process.  Is there a way to recover my folder connections from ArcGIS 10.0?  I have all my local files on the c: drive backed up.


texts:
I am using ogr2ogr to load data from multiple GML files into PostGIS. As shown in the ogrinfo output below, each GML file has multiple layers.

    Had to open data source read-only.
INFO: Open of `D:\Projects\osmaps\data\vectorMapDistrict_unzipped\HP40.gml'
      using driver `GML' successful.
1: Ornament (Polygon)
2: SurfaceWater_Line (Line String)
3: TidalWater (Polygon)
4: SpotHeight (Point)
5: TidalBoundary (Line String)
6: FunctionalSite (Point)
7: SurfaceWater_Area (Polygon)
8: Road (Line String)
9: Foreshore (Polygon)
10: AdministrativeBoundary (Line String)
11: NamedPlace (Point)
12: Building (Polygon)

I want to load all layers except the Ornament layer. How do I get ogr2ogr to skip the specific layer?

My current script in Windows Powershell:

# Input gml directory
$gmlDir = "D:\Projects\osmaps\data\vectorMapDistrict_unzipped"
# Loop through gml files
Foreach($gmlFile in Get-ChildItem $gmlDir -filter *.gml) {
    #Load data
    ogr2ogr -f PostgreSQL PG:"dbname=osmaps host=localhost port=5432 active_schema=vectormapdistrict user=postgres" -append -progress --config PG_USE_COPY YES $gmlFile.fullname}


texts: How to determine the optimal raster resolution when interpolating a surface from points?

texts:
I'm looking for guidance on how to determine the desired output resolution, when interpolating a surface from irregularly-spaced point samples.

I have a series of boreholes taken across a city, and the density of the samples varies considerably - sometimes boreholes are located within ~5m of each other, while in other locations they are ~1km apart.

What techniques should I use to determine the applicable cell-size, when interpolating a surface from these points? Does the optimal cell-size depend on the interpolation method?

I'd like the highest resolution which is supported by the dataset. eg I assume that a 1m grid is no more accurate than a ~10m grid - but how do I determine that number?


texts: className of leaflet popup option does not apply to the intended DOM element

texts:
I have a problem similar to this post.  Basically, the CSS is not applying to the element I am trying to target, and I cannot see why not:
/* css file: */
.sitePopUp .leaflet-popup-content-wrapper{

    background: white !important;
    color: #0a10b9;
    box-shadow: 0 3px 14px rgba(0,0,0,0.4) !important;
    background-color: #61ead6 !important;
    border-radius: 0.3em;

}

.sitePopUp{
    background: white !important;
    color: #0a10b9;
    box-shadow: 0 3px 14px rgba(0,0,0,0.4) !important;
    background-color: #61ead6 !important;
    border-radius: 0.3em;

}

// Having both above styles affects the element begind the leaflet-popup-content-wrapper
// If one of these two styles are omitted, the style is not applied to anything

/* JS file */
var popup = L.popup()
                    .setContent(`Some content... `)
let customPopUpOptions  = {className: 'sitePopUp', 'max-width': '500'}
let customMarkerOptions = {'icon': cutsomIcon}
let marker              = L.marker( [latitude, longitude],  customMarkerOptions )
                                        .bindPopup(popup, customPopUpOptions)

The page source file (HTML) shows the sitePopUp class is added, but the element behind the popup is styled instead of the popup

Some content

    ×

Inspecting element shows:
element.style {
    opacity: 1;
    transform: translate3d(979px, 241px, 0px);
    bottom: -7px;
    left: -159px;
}
.leaflet-fade-anim .leaflet-map-pane .leaflet-popup {
    opacity: 1;
}
.leaflet-fade-anim .leaflet-popup {
    opacity: 0;
    -webkit-transition: opacity 0.2s linear;
    -moz-transition: opacity 0.2s linear;
    transition: opacity 0.2s linear;
}
.sitePopUp {
    background: white;
    color: #0a10b9;
    box-shadow: 0 3px 14px rgba(0,0,0,0.4);
    background-color: #61ead6;
    border-radius: 0.3em;
}

I expected that the first style where the decedent element is styled do the trick, so what is going on?
EDIT:
So, this is getting a bit confusing.  The picture below shows what I mean by applying the style to the element below.  HOWEVER, I tried to reproduce a static page on codesandbox.io and in there the style is applied correctly to the popup. see the link here


texts: How to remove unwanted columns from a polygon in R AND keep the result as a spatial data frame?

texts:
To simplify my problem, I'm gonna start at the source instead of explaining my full-map problem.

I have a shapefile of Ethiopia pulled from GADM (can be downloaded here). A lot of the attribute columns are unnecessary and just create extra clunkiness as I perform clips/joins/etc. with the shapefile. So I want to remove all but the first 4 columns. Doing this seems easy enough, and I tried and succeeded with this bit of code:

library(rgdal)

#read in shapefile
eth <- readOGR(dsn = "D:/Mapping-R/Ethiopia", layer = "ETH_adm0")

#subset data/remove columns
eth <- eth@data[, -(5:67)]

This does indeed leave me only with the first 4 columns, but it also makes the object a regular data frame rather than a spatial data frame ready for mapping. So how do I remove these unwanted columns while maintaining status of SPDF?


texts: Using topo to raster tool on temperature contours and getting error 010089?

texts:
I am using the Topo to Raster tool in ArcGIS Desktop 10.5 on a global temperature contour shapefile and keep getting this error:

  ERROR 010089: C:\Users\workspace is not a supported feature type.
  Failed to execute (TopoToRaster).

I just used this tool on a smaller contour file and it worked just fine. I have not changed any default parameters other than not enforcing the drainage.

Additionally, if anyone knows of a global sea surface temperature raster file that would eliminate my need for this tool.


texts: How does LayerLimit in WMS capabilities work in QGIS 2.18.11?

texts:
I set <LayerLimit>1</LayerLimit> in service metadata in capabilities file to restrict the number of layers that the client is permitted to request in a single GetMap request. But when I add a new WMS layer in QGIS, still I can select multiple layers simultaneously and multiple layers are requested coma separated ("LAYERS=layer1,layer2,layer3"). Does layer limit work differently? What is the correct way to restrict multiple layers request?

WMS version

<WMS_Capabilities xmlns:xlink="http://www.w3.org/1999/xlink" version="1.3.0">


texts: Network Dataset with one-way streets, ArcGIS 10.3

texts:
I want to create a Network Dataset with one-way streets in ArcGIS 10.3. I have an attribute "Direction of Traffic" with four statuses:

open in both directions
open in pos. direction
open in neg. direction
closed in both directions

If the latitude value of the endpoint is lower than the one of the startpoint, the status is 2 (pos). An example for a traffic circle (driving right):

the line from 12 o'clock to 9 o'clock is pos. (2)
the line from 9 o'clock to 6 o'clock is pos. (2)
the line from 6 o'clock to 3 o'clock is neg. (3)
the line from 3 o'clock to 12 o'clock is neg. (3)

How can I use this information to create a Network Dataset with one-way streets?


texts: Points to line accuracy concerns in ArcMap

texts:
I have a csv of a GPS trace with lat-longs that I want to convert into a line so I can measure the distance. My workflow was as follows:

1) Import csv into ArcMap 10.4

2) Convert latlongs to points and create a feature file

3) Use points to line tool in the Data Management toolbox to draw a line from the points (the points are in sequence and have unique IDs)

4) Measure the length of the line feature by calculating geometry.

The csv had 277 points; the calculated length of the line feature after step 4 was ~21 km. However, when I do the same calculation using other software (e.g. using geodist command in stata) or by manually measuring the distance between the points (using the ruler tool in ArcMap), this distance comes out to be ~18 km. I cannot figure out why this discrepancy is occurring. At some places, the points are overlapping (as shown in image below), but I would expect this to results in higher estimation of distance when measuring as distance between points compared to the length of the line, and not a lower one as is currently. Can someone provide insight as to why I am getting this discrepancy? As of now it seems the points measure (18 km) is the correct one and not the 21 km one. Why is the line length feature not showing the right results?


texts: QGIS - graphical modeler aka. processing toolbox - update only selected fields

texts:
I have two layers:

point layer representing some house addresses.
polygon layer representing some kind of territory.

What I'm trying to do is, create a model where I select the points inside the polygon using:

"select by location"
"extract selected feature" (+dependency),

Since the "select by location" can not be used in the modeler as an output. Then I am willing
to change only the selected points field "Satellite" that intersects with the polygon,  by changing it's value
for example with a string 'test'.
CASE 1:
I only managed to change this field in the New Layer: "TEST_LAYER" added in the field calculator.
However this is not my intention , I want to update the field "Satellite" in "Points in the Polygon"-Layer.
Extract by location example:

----->Resulting in:

Select by location + extract selected features example
https://i.stack.imgur.com/5g8tQ.jpg
Just trying the green colored function in the field calculator in the attribute table is working but not working in the modeler.
Execute a field calculation on just selected records from Processing Toolbox in QGIS 3?
CASE 2:
I managed to change the field "Satellite" to 'test' in the desired "Points in Polygon"-Layer but
without making use of the previous selection. The whole column will update wheter or not is selected. In this case I made use of "Update Field" which is being provided by the Networks  PlugIn.
https://i.stack.imgur.com/5BPtM.jpg
The first example 1) can not update the field because the input from the extract by location is not valid
The second example 2) is working changing the field in the "Points in the Polygon"-Layer but is not knowing anything about the selection. Dependency between the "Extract by location" and "Update field" is not helping.
Any ideas on how I can resolve this problem?


texts: How can I config/use MapCache/MapServer to serve my MBtiles as WMS service?

texts:
Currently I am using GeoServer for my Map Service as WMS.

I use TileMill for creating a customised map with styling and SHP files and then export as mbtiles.

I have tried TileServer-PHP & TileStream as temporary solution for z/x/y.

I wonder how I am able to setup a MapServer to consume MBtiles and serve out as WMS standard.

Please assist me on this as I am still quite new to GIS and have to use OpenSource technology.

Update

I use MapProxy to serve out the basemap into WMS OGC Standards


texts: Trouble adding and populating field on the fly in ESRI feature layer

texts:
I'm trying to build a mashup consisting of a Census tiger web feature layer, and some JSON. I want to make a map of the feature layer, then add a new field to the layer, and populate from the JSON. Here is the url to the feature layer:
https://tigerweb.geo.census.gov/arcgis/rest/services/TIGERweb/tigerWMS_ACS2015/MapServer/82

My JSON data is a simple array of objects:

   var data =  [{"P0010001":"4779736","P0050010":"185602","NAME":"Alabama","state":"01"},{"P0010001":"710231","P0050010":"39249","NAME":"Alaska","state":"02"},...]

I then start my webmap, etc:

require([stuff],function(stuff){
    var map = new Map...
    var featurelayer = new FeatureLayer etc
    map.addLayer(featurelayer);
});

Then below this I try to add a field to the layer, and populate with JSON. Here's what i've tried:

for(var i = 0; i < layer.graphics.length; i++)
                            {
                                for(var j = 0; j < new_data.length; j++)
                                {
                                    if(layer.graphics[i].attributes.STATEFP == new_data[j].state)
                                    {
                                        layer.graphics[i].attributes.NEWFIELD = new_data[j].P0010001;
                                    }
                                }
                            }

But this didn't work. I tried logging Object.keys(layer) and got a long list of keys, of which one is graphics, however, graphics doesn't seem to have any information with it.

I also tried calling map.redraw() after the loop, but that didn't help. There is an example similar to what I'm trying to do here: http://cartographicperspectives.org/index.php/journal/article/viewFile/cp76-hunt/1311 and they use arr.forEach, which is apparaenlty part of dojo?? I'm trying to accomplish the join using vanilla js.

So I also tried the layer.attr() method in the ESRI api. Since its a layer method I assume I need an array of graphics to use it. But when I call console.log(layer.graphics.length) i get 0.

So is it possible to programmatically add a field, and do a join such as this to a feature layer? I.e. add a field and populate it from JSON?


texts: Calculating bounding box coordinates based on center and radius

texts:
I need to calculate a bounding box for osm query. I get coordinates of center of the area (WGS84) and radius (meters). How can I do it?

I though about using some kind of projection, calculating projected box coordinates and then convert it back. Is it good idea? If so, which projection will suit the best to this task? I've read that miller's projection is a poor idea and it's better to use state plane projection. But on the other hand in my case bounding box doesn't need to be accurate. The most important for me is to calculate those coordinates quickly.

I would be grateful if you could show me some code example that do conversions essential for calculating bounding box.

I wanted to make this works on whole planet, but if there are better approximations I can assume that it should work in eastern Europe
Latitude: 48.342 - 55.279
Longitude: 13.645 - 25.071


texts: Converting from one Transverse Mercator based system to another produces inconsistencies?

texts:
As a civil engineer, the one semester of survey training we had leaves my knowledge quite lacking in the field of coordinate systems.

We're working on a project based in Gabon and we get information from two parties in two different coordinate systems. The one coordinate system is UTM WGS84 Zone 32N, which is convenient. The other coordinate system is unknown but we got the following information from the other company:

False easting: 500000
False northing: 500000
Central meridian: 11.5
Scale factor: 0.9996
Latitude of origin: 0.0
Linear unit: Meter (1.0)

I used Autodesk Civil 3D to create a new coordinate system based on the above information but I don't get accurate results. When I convert the survey from one company to Lat Long and then convert that back to the other company's system everything is far out.

After digging around on the internet I found that the Gabon Transverse Mercator System of 2010 uses a central meridian of 12 degrees so I changed mine accordingly.

After this conversion the two surveys are now only a few meters off, but it is still much too inaccurate for my comfort.

I suspect that you should be able to linearly convert between the two systems without any distortions seeing that both systems are Transverse Mercator based on WGS84.

Is my understanding of this correct, or is it harder than I expect it to be? Maybe I'm missing something big and the central meridian really is 11.5 degrees.

AndreJ's answer clears up some of my assumptions, which are clearly wrong. I am however not trying to just move the drawings from one origin to another, although this would have been possible if my assumptions were true.

I am trying to convert points from one system to lat-long coordinates. My expectation is that if I take those lat-longs and project them in my drawing in the other coordinate system then the point should be in the correct place. If the one survey is done in one of the older systems then I would completely understand that it is out by a few meters, though I expect it to be based on one of the newer systems.


texts: How can I retrieve a list of the most used GIS software?

texts: I would like to know how many user use a GIS software and compile a list of the most used GIS software. I thought to use the download counter of the sites but this is quite inaccurate how can I determine the number of users of proprietary software. Is there a way to find the number of licenses sell by e.g. ArcGIS?


texts:
Encode:

1- A way of encoding a FeatureType is :

SimpleFeatureType TYPE = DataUtilities.createType("Location", "geom:Point,name:String");
GML encode = new GML(Version.GML2);
encode.setBaseURL(new URL("http://localhost/"));
encode.encode(System.out, TYPE);

2- Another way of encoding:

 FeatureTypeTransformer t = new FeatureTypeTransformter();
 t.transform(TYPE,System.out);

Decode:

I have a xsd of a featuretype:

<?xml version="1.0" encoding="UTF-8"?>
<xs:complexType xmlns:xs="http://www.w3.org/2001/XMLSchema" xmlns="http://www.w3.org/2001/XMLSchema" name="jahanserver:test1_Type">
<xs:complexContent>
    <xs:extension base="gml:AbstractFeatureType">
        <xs:sequence>
            <xs:element name="the_geom" minOccurs="1" nillable="false" type="gml:GeometryAssociationType"/>
            <xs:element name="name" minOccurs="1" nillable="false">
                <xs:simpleType>
                    <xs:restriction base="xs:string">
                        <xs:maxLength value="2147483647"/>
                    </xs:restriction>
                </xs:simpleType>
            </xs:element>
        </xs:sequence>
    </xs:extension>
</xs:complexContent>
</xs:complexType>

And I want to decode it as a SimpleFeatureType. How can I do that?

Possible solution 1:
The GML2ParsingUtils Class may be helpful. It has a method named featuretype for parsing schema, but I don't know how to prepare the method arguments.

Possible solution 2:

GML gml = new GML(GML.Version.GML2);
gml.decodeSimpleFeatureType(url,name);

But I don't have any url. I have a HttpServletRequest containing a bare featuretype with complexType xml root.


texts: What is port used by GPS sensor that comes in tablet?

texts:
Wanted to see: what is the port used by the GPS sensor that comes in a tablet?
or I can get the GPS?

In a tablet with Windows 8.1
What I want is to get by vb code, the coordinates of my location, using the GPS sensor of the tablet.


texts: GDAL 3.0.4 invalid coordinate transformation result

texts:
I'm trying to transform coordinates from a modified EPSG:31255 CRS (false northing +y_0 = 0) to EPSG:31255 (false northing +y_0 = -5000000) and always get a wrong result using GDAL 3.0.4.
The transformation should remove 5000000 in y, but reverses xy as well.
With GDAL 2.4.0 everything works as expected.
Custom EPSG: +proj=tmerc +lat_0=0 +lon_0=13.3333333333333 +k=1 +x_0=0 +y_0=0 +ellps=bessel +towgs84=577.326,90.129,463.919,5.137,1.474,5.297,2.4232 +units=m +no_defs

EPSG 31255: +proj=tmerc +lat_0=0 +lon_0=13.3333333333333 +k=1 +x_0=0 +y_0=-5000000 +ellps=bessel +towgs84=577.326,90.129,463.919,5.137,1.474,5.297,2.4232 +units=m +no_defs

from osgeo import osr,ogr
wkt = 'POLYGON((32611 5320163,32611 5337745,54270 5337745,54270 5320163,32611 5320163))'
srs = '+proj=tmerc +lat_0=0 +lon_0=13.3333333333333 +k=1 +x_0=0 +y_0=0 +ellps=bessel +towgs84=577.326,90.129,463.919,5.137,1.474,5.297,2.4232 +units=m +no_defs'
trs = 31255

src = osr.SpatialReference()
src.ImportFromProj4(srs)
tgt = osr.SpatialReference()
tgt.ImportFromEPSG(trs)

transform = osr.CoordinateTransformation(src,tgt)
polygon = ogr.CreateGeometryFromWkt(wkt)
polygon.Transform(transform)
print(polygon.ExportToWkt())

Do I miss something here? Is there any Config Setting that prevents this behavior?


texts: Finding shortest path without barriers using the sf package in R

texts:
Say I have three ships near PNG in the pacific. I want to know the shortest path between all three, while not going aground. This can be done with the ArcGIS function "cost distance" and "shortest path", but I would like to do it in R with the cool new "sf" package. Alternativly using gDistance. Anyone know how to?

Example data:

library(marmap)
library(raster)
library(sf)

#Make three ship positions near PNG
pts <- data.frame(nr=1:3, x=c(145,145,150), y=c(-10,-1,-7)) %>%
  sf::st_as_sf(coords = c("x","y")) %>%
  sf::st_set_crs(4326)

#Get bathymetric data from the area
papoue <- getNOAA.bathy(lon1 = 140, lon2 = 155,
                        lat1 = -13, lat2 = 0, resolution = 4)

#Turn into raster and keep only depth lower than 0 meter (water)
r <- marmap::as.raster(papoue)
plot(r)

r@data@values[r@data@values>=0] <- NA
r@data@values[r@data@values<0] <- 1

plot(r)
plot(pts, add=T)


texts: How to add WFS layer from Geoserver to Openlayers2

texts:
Ok, I know that this question has been asked here many times and I apologize for asking it again but I've looked all answers and I still can't fix my code to work.

I would really appreciate if someone could tell me where I'm going wrong.
I'm trying to add from GeoServer vector layer with this attributes:

layer name:marine,
workspace:marine,
workspace url: "http://www.geoserver.org/marine"
epsg of layer on geoserver: 3765, epsg of map used in openlayers: 900913

This is part of my code that could help you see where I'm wrong:

<script>
        var geographic = new OpenLayers.Projection("EPSG:4326");
        var mercator = new OpenLayers.Projection("EPSG:900913");
        var world = new OpenLayers.Bounds(-180, -89, 180, 89).transform(
            geographic, mercator
        );
        var center = new OpenLayers.LonLat(16.00, 44.00).transform(
            geographic, mercator
        );
        var map = new OpenLayers.Map("map-id", {
            projection: mercator,
            displayProjection: geographic,
            units: "m",
            numZoomLevels: 16,
            maxExtent: world});
        map.addControl(new OpenLayers.Control.Permalink());
        map.addControl(new OpenLayers.Control.MousePosition());

        var osm=new OpenLayers.Layer.OSM();
        map.addLayer(osm);
        var bing = new OpenLayers.Layer.Bing({
            key: "AqTGBsziZHIJYYxgivLBf0hVdrAk9mWO5cQcb8Yux8sW5M8c8opEC2lZqKR1ZZXf",
            type: "Aerial",
        });
        map.addLayer(bing);

        var marine = new OpenLayers.Layer.Vector("marine", {
            strategies: [new OpenLayers.Strategy.Fixed()],
            protocol: new OpenLayers.Protocol.WFS({
                version: "1.1.0",
                url:"http://localhost:8080/geoserver/wfs",
                featureNS: "http://www.geoserver.org/marine",
                featureType: "marine",
                geometryName: "the_geom",
                srsName: "EPSG:900913"

            })
        });
        mapp.addLayer(marine);

Somebody mentioned on one topic that if I put my HTML file to geoserver-->data_dir-->www, it could help but it didn't.

I spent more then half day trying to figure out where I wrong and I would really appreciate help. I just can't understand why is adding WFS layer such a big problem for everyone??


texts:
I have a layer of polygons that represent two specific zones withing a city. I want to get a count of the number of parcel polygons in another layer that are within the two specific zones I am interested in. Some of the boundaries of the two layers touch each other. Is there a way to get a count of parcels in each zone in QGIS?

I've tried doing a join attributes by location using the zones as the target layer and the parcels as the join layer. I used contains with a precision of 2 to try an account for the boundaries touching. I have the attribute summary as sum and Keep all records. The outcome is a shapefile but there is no sum of the number of parcels in each zone even though I have a count column in the parcel attribute table.

Is this the best way of going about this? IS there a way I could select using an expression from either attribute table?

EDIT:
Here are the images of the polygons and the Join dialogue box for the options I used:


texts: I have a series of convex hulls that I want to overlay but I want to be able to see all of them at the same time. How do I do this so? I'm guessing making them slightly transparent will do the trick


texts: How Update one table row with combare with the other polygon table using their geom in SQL Server

texts:
a-> points
b -> Polygon

i want to update city trow in table a with compare b of geom in SQL server

this is the PostgreSQL Query

update a set city=b.name from b where st_contains(b.the_geom, a.the_geom)

how convert this to sql server or how update city trow in table a, on SQL server?


texts:
In QGIS I am trying to make a map with the Maidenhead Locator System overlayed. I am making these overlays with QGIS's create grid tool but the grids are not coming out right.

If I leave the map unprojected the bigger and smaller grid lines are parallel but any projection I use they are not parallel. Does anybody know how to fix this?


texts:
As I am trying to add raster to my mosaic dataset and were able to see them and now all of my raster are not showing.

Why did it could not build a footprints ? I have tried everthing to get them restore back...


texts: The graduated symbology doesn 't take into account the attribute I need in QGIS

texts:
I m working on a project on QGIS
I have a point shapefile vector that represents my store
In the attribute table, i have bunch of data such as size, unit per transaction and most importantly store turnover

When i want to graduate the symbology, the turnover isn't present in the list i can choose from

Does someone know why? and how to fix it?
Thanks


texts: Creating clusters with same number of points in each in CartoDB

texts: I have a map of about 4200 retailers and I would like to create clusters, with each cluster containing 10 retailers in CartoDB. How can I do this?


texts: Populate fields based on selected objects using the updatecursor

texts:
How do I get my script to only populate fields on selected objects while editing on the map?

# Importing arcpy
  import arcpy

# Set the workspace environment
  arcpy.env.workspace = r'C:\Users\aa2zz6\Desktop\TrackTraceDownload\UPDM2016.gdb'

  rows = arcpy.UpdateCursor('P_ExcessFlowValve')
  for row in rows:
    if row.ItemNumber == '100':
      row.Status = 'In Service'
      row.Owner = ''
      row.Diameter = '1'
      row.Capacity = ''
      rows.updateRow(row)
  del row
  del rows

  print "Updated Selected EFV Features on the map"


texts: I am running on MAC OSX El Capitan. My version of postgres is 9.4.5.0. My PostGIS version is 2.1.7. I cannot upgrade to PostGIS 2.2, because when I am trying to do, it pops "PostGIS requires SFCGAL >= 1.1.0 (found 1.0.5)". Than when I try to fix it by upgrading SFCGAL, it stops process in the steps of compiling "union.cpp" library. I have found on github, that the problem might be in the old version of Boost, so I have upgraded Boost, but the problem has preserved.


texts: I keep getting the same error when trying to clip a raster

texts:
Error: Cutline feature without a geometry.

Consistently get the above error when trying to clip a raster. I have looked at all the other questions relating to this and tried all the 'fixes' however I am still getting the same error. Anything I may have overlooked.

There are no errors when validating the polygon and both polygon and raster are in the same CRS - WGS84 / UTM 37S


texts:
When I use “ui.Chart.image.seriesByRegion” command in Google Earth Engine it required to define the spatial scale of the raster. Please, look at the below command, if I remove "250" it shows error. I did a search, but could not find clear solution for it.

I have tried to apply the solution given in Exporting images in native CRS and resolution, but it did not work as well "NDVI.projection is not a function".

Is there any possibility to keep native spatial scale of the image while producing chart?

var pts = ee.Geometry.Point([10.61134244284358, 50.00619843404531]);

var NDVI = ee.ImageCollection("MODIS/006/MOD13A1")
    .select(['NDVI'])
    .filterDate('2018-01-01', '2019-06-25');

print(NDVI);

var chart = ui.Chart.image.seriesByRegion(NDVI, pts, ee.Reducer.first(), 'NDVI', 250)
.setChartType('LineChart');

print(chart);


texts: Calculating relative share of area covered by water/percentage of polygon area that overlaps in ArcGIS for Desktop?

texts:
I need to calculate the share of the area within a 50km radius of major German city centers that is covered by water with either Arcmap 10.3.1 or Qgis 2.14.
My plan is to use the Global Lakes and Wetlands Database as a source , but if there is a dataset publicly available with higher resolution which would also include all the rivers I am very open to a recommendation.

I have already created a shapefile that represents the area of interest around each city in an EPSG:3035 projection. How can I now calculate the share of the area within each polygon covered by water? I have already tried to follow these instructions How to calculate a percentage of an area covered by polygons ArcGIS 10.1 but the intersect tool does not show any results while the tabulate intersect tool shows errors. What would be the suitable approach for this task?


texts: Split multiple polylines by distance in ArcGIS

texts:
Is there a method to split thousands of polylines by employing a tool like the "Editor" - "Split" method and using distance, in ArcGIS?

When I utilize the "split" - "distance" tool (fig1), it only splits one section leaving the remainder unsplit (fig2).

I have attempted the "COGO" tool - using "Proportion", and can utilize the "duplicate" button to assign the remainder of segments which seems like a viable option.  But I would like to find some way to automate this process.  Has anyone found a method to do this?  Or is Python (coding) the only way?

I am trying to split over 2000 lines into equal distance segments (fig 3) of ~600m. Each segment will have its angle measured, and will ultimately result in a rose diagram representing angular statistics.  The distance of each segment is important in my splitting regime (a long line will have a greater representation of a particular angle on my rose diagram than the neighboring shorter line segments).

I currently am using ArcGIS 10.2 with an advanced user license.


texts: Is setting raster CRS in QGIS georeferencer useless? Still have to choose the same for target SRS?

texts:
I'm trying to move out of ArcGIS into QGIS permanently for basic operations, so be patient please.
Here is the problem. I'm trying to georeference a topo map in EPSG:22523 (Corrego Alegre UTM). I did the regular procedure, which has worked before for WGS84 stuff, and is also shown in QGIS tutorials, including this one: http://www.qgistutorials.com/en/docs/georeferencing_basics.html

However the georeferencing just doesn't work. When I open the layer properties in the project, the CRS for the layer is actually WGS84 (4326). In the link above a user asks the same question, and his conclusion, after not getting an explanation, is that setting the layer CRS in georeferencer correctly is just useless. You still have to put that CRS in target SRS in transformation settings. However there isn't a single georeferencing tutorial where this is mentioned. In fact they all put target SRS different from layer CRS. Please help.


texts: There are no local edits to upload in ArcGIS Collector?

texts:
I edited a feature layer offline. I made sure that the Sync and Edit options were selected before I exported the service (and verified on ArcGIS online). I had no problems collecting the data, and could see there were 19 features that were not synced. When I tried to export it said "Push-only sync is enabled and there are no local edits to upload". I can see the local edits, and yes, push-only is enabled, but I am the only individual editing the layer so it doesn't matter. So why aren't the local edits being synced?

I am using an iPad Pro. I have a Bad Elf receiver linked to the iPad. I've used this system before, and am unsure what I did wrong!


texts:
I am working on a GIS project for Virginia and the data spans the entire state. There are two CRS available VA North (FIPS 4501) and VA South (FIPS 4502).

Is it possible to use just one of those CRS for the entire state, or do I have to split the data up into North and South?


texts:
I am trying to write a program which will render tiles from EPSG2180-based WMTS server. Given coordinates of a point (WSG84 or EPSG2180, I can convert from one to another) I want to render a tile containing it.

Each TileMatrix has info about TileWidth/Height, MatrixWidth/Height, TopLeftCorner and I understand what they mean. However, they don't seem sufficient to translate lat/lon to row/column. Apart from those each TileMatrix has mysterious ScaleDenominator value.

What is ScaleDenominator and what is the relationship between it and other values describing TileMatrix? How can I use it to translate lat/lon to tile's row/column?