texts:
I can't seem to get a basic basemap background showing from any examples for contextily v 1.0rc2. I'm on python 3.8.0, 32 bit, windows 10. My packages and version are listed at the bottom.

It was supremely difficult to get it installed with all its dependencies but that may just be my setup. I haven't included the code as I'm using 3 examples, urls below.  The error partially initialized module 'contextily' indicates the package may not have installed properly.

I get the following errors:

basemap, extent = ctx.bounds2img(xmin, ymin, xmax, ymax, zoom=zoom, url=url) AttributeError: partially initialized module 'contextily' has no
attribute 'bounds2img' (most likely due to a circular import)```

From example:
https://share.cocalc.com/share/56cef78b3fda5e0338c2da71b222e6b275a9193e/contextily.ipynb?viewer=share

attribute 'add_basemap' (most likely due to a circular import)```

From this example:
https://geopandas.readthedocs.io/en/latest/gallery/plotting_basemap_background.html?highlight=basemap#adding-a-background-map-to-plots

attribute 'Place' (most likely due to a circular import)```

From this example  (from the actual developers of contextily?):
https://github.com/darribas/contextily/blob/master/examples/plot_map.py

Can anyone tell what the problem is?

Or happy to try an alternative package if it's more supported.

I managed to get a GDAL loaded raster file into a matplotlib figure using imshow and color it - I just don't know how I would get a tile layer with a URL like OSM into the figure.

Installed packages:

['-yinstaller==3.5', 'affine==2.3.0', 'alabaster==0.7.12', 'altgraph==0.16.1', 'arrow==0.15.6', 'attrs==19.3.0', 'babel==2.8.0', 'cartopy==0.18.0', 'certifi==2019.11.28', 'cffi==1.14.0', 'chardet==3.0.4', 'click-plugins==1.1.1', 'click==7.1.2', 'cligj==0.5.0', 'colorama==0.4.3', 'contextily==1.0rc2', 'cryptography==2.8', 'cycler==0.10.0', 'decorator==4.4.2', 'deprecated==1.2.10', 'descartes==1.1.0', 'docutils==0.16', 'fiona==1.8.13', 'flask==1.1.2', 'future==0.18.2', 'gdal==3.0.4', 'geographiclib==1.50', 'geojson==2.5.0', 'geopandas==0.7.0', 'geopy==1.22.0', 'geos==0.2.2', 'gitdb==4.0.5', 'gitpython==3.1.3', 'idna==2.9', 'imagesize==1.2.0', 'itsdangerous==1.1.0', 'jinja2==2.11.2', 'joblib==0.15.1', 'kiwisolver==1.1.0', 'lxml==4.5.1', 'markupsafe==1.1.1', 'matplotlib==3.2.0', 'mercantile==1.1.4', 'munch==2.5.0', 'networkx==2.4', 'numpy==1.18.1', 'ogr==0.12.1', 'osgeo-easy==0.0.1.dev16', 'osgeo4w-installer==1.4.0', 'osmnx==0.11.4', 'overpass==0.7', 'overpy==0.4', 'packaging==20.4', 'pandas==1.0.1', 'pefile==2019.4.18', 'pillow==7.0.0', 'pip==19.3.1', 'proj==0.2.0', 'psycopg2==2.8.4', 'pycparser==2.20', 'pygithub==1.51', 'pygments==2.6.1', 'pyinstaller==4.0.dev0+9dd34bdfba', 'pyjwt==1.7.1', 'pyopenssl==19.1.0', 'pyparsing==2.4.6', 'pyproj==2.5.0', 'pyqt5-sip==12.7.0', 'pyqt5==5.13.2', 'pyshp==2.1.0', 'python-dateutil==2.8.1', 'python-gitlab==2.2.0', 'pytz==2019.3', 'pywin32-ctypes==0.2.0', 'pyyaml==5.3.1', 'rasterio==1.1.4', 'reportlab==3.5.34', 'requests==2.23.0', 'rtree==0.9.4', 'setuptools==44.0.0', 'shapely==1.7.0', 'six==1.14.0', 'smmap==3.0.4', 'snowballstemmer==2.0.0', 'snuggs==1.4.7', 'sphinx==3.0.4', 'sphinxcontrib-applehelp==1.0.2', 'sphinxcontrib-devhelp==1.0.2', 'sphinxcontrib-htmlhelp==1.0.3', 'sphinxcontrib-jsmath==1.0.1', 'sphinxcontrib-qthelp==1.0.3', 'sphinxcontrib-serializinghtml==1.1.4', 'sqlalchemy==1.3.15', 'urllib3==1.25.8', 'werkzeug==1.0.1', 'wrapt==1.12.1', 'xlrd==1.2.0']


texts: I have downloaded a 100Km map from OS, but the problem is that it is far larger than I need and hence takes a long time to load every refresh in QGIS. Is there a way to cut down the size of the accompanying shape files to a more manageable size of, say, 10Km square?


texts: Exporting composer with 2 maps in QGIS: wrong extent in one map after export

texts:
After creating this map in the QGIS 2.18.2 map composer

(note that the red dot in the upper left map is shown in India),

I then export it to a PNG image file. However, the overview map in the upper right corner seems to have changed its extent and zoom during the export process, as the red dot is now located in Iran.

Any ideas why and how to solve this?


texts: Openlayers ajax loading works on FF but doesn't work on Chrome. Any idea why?

texts:
I'm using this ajax loading panel
http://thisiscontext.com/tools/jQuery-showLoading

And I have the following code :

jQuery('#map').showLoading();
var request = OpenLayers.Request.POST({
        url: "service.asmx/DeleteStopPoint",
        data: "{'TripId':'" + currentTrip + "','PointId':'" + feature.attributes.PointId + "'}",
        async: false,
        headers: { "Content-Type": "application/json; charset=utf-8" },
        callback: refreshMap
    });
jQuery('#map').hideLoading();

what happens is that in FF I do see the loading panel appearing before the request is made and hidden once its over... but in chrome it doesn't happen. it looks like the panel appears and immiedtly disappear (because if I comment out the hideLoading function , it appears AFTER the POST)

any idea why ?


texts: Applying pixel reliability mask to time series MOD13Q1

texts:
I am working with a 23-band GEOTiff in which each band is a 16-day MVC MODIS NDVI image.  The MODIS .hdfs include a "pixel reliability" image with the following values:

     [-1]: Fill/No Data-Not Processed.
     [0]:  Good data      - Use with confidence
     [1]:  Marginal data  - Useful, but look at other QA information
     [2]:  Snow/Ice       - Target covered with snow/ice
     [3]:  Cloudy         - Target not visible, covered with cloud

I want to apply a quality mask using ENVI to elimate pixels with values of -1, 2 and 3.  I know how to do this for a single image using the "Build Mask" function in ENVI.  If I use the 23-band NDVI.tiff file as my input and use the 23-band pixel reliability file as my mask, will this mask each NDVI image with its respective pixel reliability (i.e. January 1st NDVI will be masked with the January 1st pixel reliability mask), or will the "Build Mask" perform some other operation?

If I am thinking about this in the wrong way, is there a better way to apply the quality mask to the time-series stack?

Many thanks!

Emily


texts: Restricting zoom levels based on Bing Map availability?

texts:
My site is built on ArcGIS Server 10 with the JavaScript API 2.3, and uses Bing Maps as the background mapping. Although the site is aimed at NSW, Australia we allow the user to zoom anywhere in the world, at any of the supported scales (ie, the tiling scheme as specified in the ArcGIS Help File.)

Bing Maps road/satellite datasets are available at varying scale ranges, with populated areas supporting larger scales than empty rural areas. When zooming somewhere with no map/satellite data, a "broken camera" icon is shown. We'd like to avoid showing this ugly icon.

How can we restrict the maximum zoom scale where the map/satellite is not available, but allow full zoom-in where it is available?

Google Maps does this quite well - as you move around the world you can see the Zoom Slider expand/contract depending on what datasets are available at that location. They presumably have metadata showing the presence/absence of data, so we seek emulate this when using  Bing Maps, with no access to the metadata.


texts: Empty shapefile after project tool in for loop?

texts:
I am running a for loop to process a years worth of data. The entire script runs fine for the first run through, but when it runs again (when date is 1/2/2013), I get the following error: 'NoneType object has no attribute 'XMin'.

Before the project tool is run the 1/2/2013 shapefile looks fine. After the project tool is run, the new shp is empty. If I start the script at date = dt(2013, 1, 2) the script runs fine. So there is nothing wrong with the shapefile.

I am not sure if something is stuck in the memory? I added these 4 lines of code at the end of the loop but this doesn't change anything.

    sr = None
    srProj = None
    shpFile = None
    shpFileProj = None

#

date = dt(2013, 1, 1)

for i in range(365):

    print (date)
    fileFolder = date.strftime('%Y') + os.sep + date.strftime('%m') + os.sep + date.strftime('%d')
    fileName = "1day_observed_shape_" + date.strftime('%Y') + date.strftime('%m') + date.strftime('%d') + ".tar.gz"
    filePath = "N:\\zip\\Precip\\" + fileFolder + os.sep + fileName
    print filePath
    date += timedelta(days=1)

    # Copy the tar file to the extracted folder and to the archive folder on the desktop
    shutil.copy(filePath, cfg.exportPath)
    shutil.copy(filePath, "C:\\Users\\alou\\Desktop\\precip_tar\\2013")

    os.chdir(cfg.exportPath)
    #tarFile = tarfile.open(fileString[30:], 'r:gz')
    tarFile = tarfile.open(fileName, 'r:gz')
    tarFile.extractall('.')

    # List the shapefile
    for shpFile in glob.glob("*.shp"):
        print shpFile

    ################## creating grid START ###################
    # Delete prj file
    os.remove(shpFile[:-3] + "prj")

    logger.info("Define the projection")
    # Define projection
    sr = "C:/load_assist_files/precip/HRAP_GCS.prj"
    arcpy.DefineProjection_management(shpFile, sr)

    logger.info("Reproject the data")
    # Reproject the layer
    srProj = "C:/load_assist_files/precip/HRAP_Projected.prj"
    shpFile = "C:/extracted/" + shpFile
    shpFileProj = shpFile[:-4] + "_proj.shp"
    arcpy.Project_management(shpFile, shpFileProj, srProj)

    # Set local variables point to raster
    valField = "Globvalue"
    outRaster = "C:/extracted/outRaster"
    assignmentType = "MAXIMUM"
    priorityField = ""
    cellSize = 1

    # Set environment settings point to raster
    arcpy.env.outputCoordinateSystem = shpFileProj
    arcpy.env.extent = shpFileProj
    print shpFileProj
    print arcpy.env.extent
    XMin = int(arcpy.env.extent.XMin)
    YMin = int(arcpy.env.extent.YMin)
    XMax = int(arcpy.env.extent.XMax)
    YMax = int(arcpy.env.extent.YMax)


texts: The difference when converting WGS84 to EPSG:3857 in MapInfo and GDAL

texts:
Please help me solve the problem with the coordinate shift.
I have two TAB files, with coordinate systems - one in meters, the other in degrees.

GDAL says they are:

file A: LOCAL_CS["Nonearth",UNIT["Meter",1.0]]
file B: GEOGCS["unnamed",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563],TOWGS84[0,0,0,0,0,0,0]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]]

An online conversion service says that:

file A: WGS 84 / Pseudo-Mercator (EPSG:3857)
file B: WGS 84 (EPSG:4326)

In file A, I have a polygon, in file B there is a point with coordinates in degrees:

13.685983, 6.276257

I want to convert it to the coordinate system of file A (in meters).
When converting to all online services using both GDAL and openlayers, it translates to the point:

1523516.658565399004146, 700071.198322909418494

And only in MapInfo does it correctly translate to a point:

1521811.2, 697887.6

Question: Why do such differences arise? Could this be due to the fact that
"degree" in GDAL 0.0174532925199433, not 0.017453292519943295?


texts: How can I get location data (country,city, etc.) by giving latitude and longtitude using PostGIS?

texts:
Currently I'm using the following query :

SELECT ST_SetSRID(ST_MakePoint(x, y),4326);

and it returns me location in numbers something like :

010100000000000000008012300000000000001230

but I wanted to get location data as text instead.


texts: Export summary statistics to other format than Excel

texts:
Is there a tool from the modelbuilder that is able to export a summary statistics table to .csv instead of .xls?
As I need to standardize the solution in modelbuilder, so preferably via modelbuilder.

ArcMap 10.4.1 - Basic Version Concurrent use


texts: Covering shapefile to raster using rasterize in R,return same result

texts:
When i using the blow code to change 'thre' property of shp to raster:

r1 <- rasterize(city.1, r,'thre',fun='min')
r2 <- rasterize(city.1, r,'thre',fun='max')

here, city.1 is shp, thre is the property that to be transformed, r is the raster layer

I got same result no matter what the fun is - any hints on solving the problem?


texts: Overlaying two shapefiles and projecting one of them from Geographic to Projected Coordinate System?

texts:
I am trying to project a shapefile which is in Geographic Coordinate System i.e. Geocentric Datum of Australia 1994 to a Projected coordinate system i.e. Map Grid of Australia Zone 55 using GDA 94 datum, also called as GDA 1994 MGA Zone 55 I think.

I have several shapefiles and this is the only one in a different coordinate system.

I tried to change the projection in Arc Catalog - right clicking the layer going into properties and selecting the projected coordinate system I want but when I actually add the data into ArcMap with other layers already there, I get an error saying:

  Warning, inconsistent extent! One or more of the added layers has an
  extent that is not consistent with the associated spatial reference
  information. Re-projecting the data in such a layer may lead to
  unexpected behaviour

The layer gets added but it doesn't overlap other layers.


texts: Projection (Datum?) slightly off for plotting Geopandas

texts:
I'm trying to plot some polygons I obtained from a 2002 Census in Chile.
I'm loading them to geopandas, converting them to epsg='4326' and then sending them to Folium (which is a wrapper for Leaflet).

polygons = gpd.read_file('lim_polygon.shp')
gjson = polygons.to_crs(epsg='4326').to_json()
my_map = folium.Map(location=[-33.45, -70.61], zoom_start=12)

The real limits should be that red dotted line you see to the left of the polygons. This tells me the to_crs conversion is almost right, but off by about two or three blocks.

polygons.crs

gives

{'ellps': 'intl',
 'no_defs': True,
 'proj': 'utm',
 'south': True,
 'units': 'm',
 'zone': 19}

If I export the converted json and try it in geojsonio I get the same result, so I don't think it's an issue with Leaflet or Folium.

If it is that there's a different Datum.

How do I fix this?

Here are the shapefiles in case you want to look at them

Does 'ellps': 'intl' in the crs dictionary define a Datum?

proj -ld in the terminal gives me

nzgd49 intl         towgs84=59.47,-5.04,187.44,0.47,-0.1,1.024,-4.5993
                          New Zealand Geodetic Datum 1949

gdalsrsinfo lim_polygon.shp, run in the terminal, gives me

PROJ.4 : '+proj=utm +zone=19 +south +ellps=intl +units=m +no_defs '

OGC WKT :
PROJCS["User_Defined_UTM_Zone_19S",
    GEOGCS["GCS_User_Defined",
        DATUM["User_Defined",
            SPHEROID["User_Defined_Spheroid",6378388.0,297.0000000284015]],
        PRIMEM["Greenwich",0.0],
        UNIT["Degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["False_Easting",500000.0],
    PARAMETER["False_Northing",10000000.0],
    PARAMETER["Central_Meridian",-69.0],
    PARAMETER["Scale_Factor",0.9996],
    PARAMETER["Latitude_Of_Origin",0.0],
    UNIT["Meter",1.0]]

Is Geopandas using that info to convert?

I read some documents that suggest that the data from Chile might be using either SIRGAS or PSAD56, but at this point I don't know how to convert datums (or declare a new crs in geopandas in case it was misspecified) to test if this would solve the issue.


texts: SRID, Shape2SQL and SQL Server 2008r2 query

texts:
I'm using Shape2SQL to load the Ordinance Survey's Meridian 2 map into SQL Server 2008r2 (which is working very smoothly and painlessly).

According to this page (https://alastaira.wordpress.com/2011/02/16/loading-ordnance-survey-open-data-into-sql-server-2008/)

  Ordnance Survey data is defined using the National Grid of Great
  Britain, so always make sure that you load your data into a column of
  the geometry datatype, using SRID 27700.

I'd not heard of SRID before so haven't been touching that value when loading the SHP files. The roadnode.shp file is currently loading and I've noticed that SRID has been set at 4326.

Once the map files are all in SQL Server is it a simple matter of changing the value in the SRID (or whatever it is called) column to 27700, or will that alteration affect other data in the table?


texts: How to get started on a GIS Programming Career?

texts:
I am already working as GIS Consultant in UK.  Being an engineering consultancy (not IT or GIS based) company there wasn't any need of GIS development in the organisation.

But due to some recent project requirements my company is trying to promote us towards GIS development.

So my question is that,  which route should I progress my career, like:
Python - C# - .Net
OR
Java?

Quicker advice will be much appreciated as I am in India for my 3 weeks and have an opportunity to join some basic programming courses.  Recent visit to the local computer institute, C# and .Net was the main stress.


texts: Populating Access (*.mdb) table with text from Text Elements to inventory MXDs using ArcPy/Python?

texts:
I am trying to write a piece of code to populate an Access mdb table with information from text elements in an mxd (text elements all have "names"). The purpose of this is to create a catalog of my mxd's and what they hold without having to manually open up every mxd. Note that the "names" of the text elements are all the same in every mxd. So what I will end up with is a script which grabs the text from a text element and inserts it into the Access table in a specific place. The Access table will be 1 record per mxd.

The version of ArcGIS I am using is 10.2.2.

Note that I am relativity new to python so I do apologize if I have not explained this well and will happily clarify anything.


texts: Identify groups of stacked geometries in PostGIS

texts: I need to create a column in my PostGIS table with values identifying different groups of stacked/duplicate geometries, such that one group of stacked geometries might get a value of 1, another group a value of 2, and so on. I've come across this post but it indicates how to delete duplicates while saving one, but I do not want to delete the duplicates. I also do not just want to identify duplicates with a boolean value, as this would not distinguish between groups.


texts: Using ArcGIS model to batch export set of layout view as JPEGs?

texts:
I have a large number of layout views I would like to export as JPEG maps using ArcMap 10.0.  The layout view is based upon a user defined grid.  Normally what I would do is manually select a grid, zoom to the selected grid and export the layout as a JPEG.  But in this case there are over 1000 to do.

Is it possible to write an ArcMap model that will allow me to Select a Grid (table record) and then Zoom to the selected grid and then export the JPEG, which I can run as a batch?

Using ModelBuilder I imagine I would need to add the "Select" tool and then associate that with a "zoom to selected" script.  This is the part I fall down at.  While I have a script that does this, I cannot associate it with my selected record.

My grid is a polygon feature dataset (either shapefile or gdb) - basically a user defined grid for a layout at a particular scale. I really just need to zoom to each grid, which will have a numeric label 1 - 1000, in a field called ID, then export the layout of that map grid.


texts: Sync-ing features from ArcGIS Online to file geodatabase?

texts:
Using ArcGIS Desktop Advanced 10.1 SP1, ArcGIS Online for Organizations

I'm new to ArcGIS Online and we just purchased an ArcGIS Online Organizational account... I'm interested in taking a feature layer (eg. fire hydrants) and publishing it on ArcGIS Online (which I'm able to do...) and allowing users/public edit the layer by adding and updating hydrants.  What I'd like to do then on a regular basis is to sync or download the hydrant layer back to my file geodatabase that I have on my laptop.  The only way I've been able to figure out to get data back to my geodatabase is to download the layer from ArcGIS via shapefile and copy/replace the features in my geodatabase... wondering if there is a less cumbersome way to go about this...

I've referenced Downloading feature class from ArcGIS Online web map? but when I download the .pkinfo as @artwork21 suggested, it only seems to reference my feature layer on ArcGIS Online, not download anything to a geodatabase as he indicated...


texts: How to calculate frontages of city parcel data in ArcMap 10.0?

texts: I need to calculate the length of polygon frontages (along streets) from our city's parcel (i.e. lot) data in order to do an analysis on lot frontages in various zoning districts. How can I do this?


texts: Is there any site which provide classified image download for soil,vegetation or water resources?

texts:
(its particularly for indian states,i need.....,but kindly quote if any other info you have)
the raster data may also work for me in above aspects.


texts: Proj4 definition for WGS84 coordinates multiplied by 1000000

texts: For some reason, a vendor found it convenient to move from a Gauss–Krüger coordinate system (EPSG:31466) to something like WGS84 but then multiplied by one million. Hence lon = 6090805, lat = 50777584, should become 6.090805, 50.777584. Is there already a proj4 definition to transform this "coordinate system"? If not, could anyone help with a projection string?


texts: Adding features to existing layer through CSV file in PyQGIS

texts:
I have an existing QGIS project that was prepared with the QGIS GUI application. It contains several layers, one which I need to update with pyQGIS. The project data is stored in an sqlite database. Since I inherited the project from my predecessor without proper handover, I have to do a lot of guess work on how things are set up. My plan is now to do the following:

Prepare CSV file with the relevant feature information
Remove existing features from the relevant layer (currently 526 features)
Add features from my CSV file to the layer (640 features to be added)

However steps 2 and 3 are going wrong. This is my code:

import os
from qgis.core import *
from PyQt5 import *
from PyQt5.QtSvg import *
from PyQt5.Qt import *
import pandas as pd

print("loading map with projects...")
app = QgsApplication([], True)
app.setPrefixPath(r"C:/Program Files/QGIS 3.10/apps/qgis", True)
app.initQgis()
Project = QgsProject.instance()
Project.read("filename.qgz")

layers = Project.mapLayers()

print("updating projects in map...")

lay = layers["layername"] # choosing the relevant layer that needs updating

lay.startEditing()

fields = lay.fields()

print("number of features: {}".format(lay.featureCount()))

listOfIds = [feat.id() for feat in lay.getFeatures()]
lay.deleteFeatures( listOfIds )
print("Number of features to be removed: {}".format(len(listOfIds)))
print("number of features: {}".format(lay.featureCount()))

lay.triggerRepaint()
projects = pd.read_csv("projects.csv")

featureList = []
feature_names = ["pk","x","y", "ProjectName", "Category", "Capacity"]
for index, row in projects.iterrows():
    x,y = [row["x"],row["y"]]
    feature = QgsFeature(fields)
    geom = QgsGeometry.fromPointXY(QgsPointXY(x, y))
    feature.setGeometry(geom)

    for f in feature_names:
        feature[f] = row[f]

lay.addFeatures(featureList)

lay.endEditCommand()
lay.triggerRepaint()

print("number of features: {}".format(lay.featureCount()))

I tried various alternative commands for updating the features. However when I render the map it is empty. Also when I open the updated project in the QGIS application it appears empty.

The code returns the following:

loading map with projects...
updating projects in map...
number of features: 526
Number of features to be removed: 526
number of features: 0
number of features: 0

So the layers are apparently not properly added. Also when running the code again, it again removes 526 features, despite deleting them in the previous run and adding 640 features. When opening the underlying SQLite database it seems that everything in there was updated properly.

How can I properly update the features in my layer, i.e. remove the old ones and add new ones?


texts:
Is there any extension to get all upstream sections of a polyline shapefile which represent an hydrographic network in QGIS?

It could work like this:

you choose a specific branch and a direction
then the algorithm will select all the sections recursively until you have all the branches.


texts: Unable to stop ArcGIS Server automatic cache generation

texts:
I'm running ArcGIS Server 10.1 on a Windows box on Amazon, via remote desktop.

I've set up some services on which I have enabled caching using the Bing/Google/Esri scheme. I explicitly set these to Update Cache Manually.

For some reason, ArcGIS Server is changing these services to Update Cache Automatically. If I change it back to Manual, then immediately re-open the dialog, it has reverted to Automatic. If I choose View Cache Status, it's shown as "Tile Generation is in progress". Cancelling the tile generation has no effect - it immediately restarts.

I am creating these services by first creating *.sd files, in which I specify the cache levels and that I wish for Manual cache generation, then creating ArcGIS Server services from the *.sd files.

Why are these services changing irrevocably to Automatic cache updating?

Thanks


texts: How to republish Hosted Feature Layers on ArcGIS Online without losing the data?

texts:
I am desgning an App on ArcGIS Online for deploying tactical plans over the terrain. My workflow is next:

Design dataset and feature layers using ArcMap 10.2, using our Enterprise Geodatabase (Oracle)
Create a Map whith layers and simbology in ArcMap 10.2
Publish this map as a hosted feature layer in ArcGIS Online, using our enterprise account. The service is configured as Feautre Layer service, to allow editing to AGOL user
Go to ArcGIS Online and author a Map with this hosted feautre layer and an App with WebApp Builder
User can now start to work with the app, in a testing phase. User are now generating testing plans, and identifing new requirments a functionality needed

This new requirments, are some like, add new layers, or new fields to the existing layers, or new subtypes....

So with this needs, we continue developing this App, using this workflow:

Edit the Layers and Maps with ArcMap to incorporate new requirments, updating our local entreprise Geodatabase
Publish the service again, modifing the existing one. (Modify option of the publishing process)

The problem we have is that with this workflow we are losing all the data generated in ArcGIS Online Hosted Feature Layers. So the question is: What is the correct workflow we must follow to assure that AGOL Hosted Feature layers data generated by users is not lost?


texts: Buffer neighbouring polygons without overlap using QGIS

texts:
Is it possible to create a buffer around neighbouring polygons, so that the new polygons do not overlap?

Preferably in QGIS but any other tool will do.

So instead of the first result I would like to obtain something similar to the second:


texts: Creating polygon from DXF file using QGIS?

texts:
I want to create a polygon from a DXF file. I used v.in.dxf to convert the DXF file to GRASS vector map format. when I run the polygonize algorithm, it gives me the following error message(here I used "a.dxf" as the input)

Algorithm Polygonize starting...
Processing lines...
Noding lines...
Polygonizing...

No polygons were created!.
I have created a script that converts DXF to shapefile. In the script, first, DXF is converted to GRASS vector map format, then run the v.clean algorithm. After that polygonize algorithm is run, finally v.dissolve algorithm is run. My tool works correctly for most of the DXF files. But when I inserted the above-mentioned "a.dxf" as the input to my script, it gives "'NoneType' object has no attribute 'getitem' See log for more details" error.

When I use the same file in ArcGIS, it creates a polygon.

How should I overcome from this problem?


texts:
I have been trying to get the Disclaimer widget found here to function with a relatively normal CMV application. I have reviewed the several other posts on Stack Exchange troubleshooting this widget but have had no luck.

I've copied all the files to my project directory in what I think is the correct location and added the requisite configuration to viewer.js.

When I load my page, the disclaimer does not appear and I see the following errors in Firebug:

The code I am using is copied directly from the developer samples, so I must either have a file in wrong location or have something configured incorrectly. When I comment out the disclaimer widget config in viewer.js, the errors disappear in Firebug.


texts: Mapinfo 15.2.2 not displaying all polygons from a MSSQL Spatial DB table

texts:
I am trying to open a fairly large table stored on our MSSQL Spatial DB - the table has around 22,000 polygons.

The table opens perfectly on Mapinfo 12.5 (both 64bit and non-64bit), but for some strange reason only a few hundred polygons display on Mapinfo 15, however when you browse the table, all 22,000 are present.

Seems like the tabular data is present, but no spatial data i.e. cannot run queries where the missing polygons are.

I have done some testing since: When a table is pushed onto the server using EasyLoader via Mapinfo 12.5, the error appears on Mapinfo 15, and vice versa.

Any idea as to why the polygons do not display?


texts: Showing polyline corrdinates in text input

texts:
This is my code:
map.on("click", function(e) {
    var latlng = e.latlng;
    document.getElementById("lat1").value = latlng.lat;
    document.getElementById("lng1").value = latlng.lng;
    console.log(latlng);
    });

How can I show all polyline coordinates, displaying in console.log, into text input?


texts:
I have a point database marking different locations (presence of an animal) on a grid of 250 m x 250 m over a certain period of time. Each one of these points have a date.
By intersecting the grid file with the point shapefile I can calculate how many points are in each square of my grid. What I'm trying to figure out is how to calculate (and show) for how many days the animal was present in each square of my grid. Is it possible? (I don't have the pro-license)


texts: DeleteField_management gives ERROR 000308 Invalid field type?

texts:
I am attempting to delete a list of fields in a feature layer but am getting an error:

Execute error: Failed to execute. Parameters are not valid. ERROR 000308 Invalid Field type

Here is a code snippet:

arcpy.MakeFeatureLayer_management(roads, tmprds)

arcpy.DeleteField_management(tmprds, ["FIELD_NAME"])

I have listed the fields in my feature layer using:

fclass = tmprds
fieldnames = [f.name for f in arcpy.ListFields(fclass)]
print (fieldnames)

and I get a list of layers formatted like:

[u'FIELD1', u'FIELD2, u'FIELD3']

I have tried using both IN_MEMORY and real locations which don't seem to make a difference.

Where am I going wrong?

Thanks


texts:
I have a file Trees.style that I'm trying to load into ArcGIS 10.0  I can't get it to be recognized.

Can somebody point me to detailed (read bonehead) loading directions?


texts:
I want to connect my GeoNetwork catalog to QGIS by way of the CSW protocol (metasearch extension). The first step has been done in fact I have connected QGIS to GeoNetwork by copying this link http://localhost:8080/geonetwork/srv/en/csw?  and then the information of my catalog appear in the metasearch  windows.

But the problem  is when I try the search by map extent an error windows appears

Search error: 'Raised exception while searching metadata : java.lang.IllegalArgumentException: Error when parsing spatial filter (version: 1.1.0):<ogc:Filter xmlns:ogc="http://www.opengis.net/ogc">
  <ogc:BBOX>
    <ogc:PropertyName>ows:BoundingBox</ogc:PropertyName>
    <gml311:Envelope xmlns:gml311="http://www.opengis.net/gml" srsName="urn:ogc:def:crs:OGC:1.3:CRS84">
      <gml311:lowerCorner>-5 -4.5</gml311:lowerCorner>
      <gml311:upperCorner>29.5 30</gml311:upperCorner>
    </gml311:Envelope>
  </ogc:BBOX>
</ogc:Filter>.
Error is: java.lang.RuntimeException: Parsing failed for Envelope: java.lang.RuntimeException: Could not create crs: urn:ogc:def:crs:OGC:1.3:CRS84'


texts: Creating new shapefile in QGIS from "Save selection as..."?

texts:
When selecting a feature from census boundary file and using "save selection as" to create a new shapefile, I get the following error:

"Export to vector file failed.
Error: creation of layer failed (OGR error:Failed to open Shapefile `brent.shp'.
)" boundary file attached

I have re-installed QGIS 1.8 several times, but that didn't help.
Also, I cannot seem to create new vector files from joins using "Save selection as..."

Has anyone found a solution to these problems?


texts: Correlating Worldclim data with vif function in R

texts:
I want to corelate 19 variables from Worldclim V2.0 using the vif test. There's no ready code like for V1 so I wrote a code below but I'm facing one error. I'm new to R.
library (raster)
library(usdm)

setwd("C:/Users/Areej/Documents")
bioFinal<- c("wc2.0_bio_30s_01.tif", "wc2.0_bio_30s_02.tif",
             "wc2.0_bio_30s_03.tif", "wc2.0_bio_30s_04.tif",
             "wc2.0_bio_30s_05.tif","wc2.0_bio_30s_06.tif",
             "wc2.0_bio_30s_07.tif", "wc2.0_bio_30s_11.tif",
             "wc2.0_bio_30s_12.tif", "wc2.0_bio_30s_16.tif",
             "wc2.0_bio_30s_17.tif", "wc2.0_bio_30s_08.tif",
             "wc2.0_bio_30s_09.tif", "wc2.0_bio_30s_10.tif",
             "wc2.0_bio_30s_13.tif", "wc2.0_bio_30s_14.tif",
             "wc2.0_bio_30s_15.tif", "wc2.0_bio_30s_18.tif",
             "wc2.0_bio_30s_19.tif")
bioFinal

#check collinearity structure - try two methods vifstep and vifcor
# 1 no cor, 2-3 low, 4-5 moderate, >5 high
library(usdm)
v1 = vifstep(bioFinal, th=10 )
v2 = vifcor(bioFinal, th=0.9)
v1
v2

when I run it I have this error for both vifstep and vifcor :
  Error in (function (classes, fdef, mtable)  :
  unable to find an inherited method for function ‘vifstep’ for signature ‘"character"’

How can I fix the error, and can use a zip file in the code instead of listing all these layers?


texts: Forcing browser to cache map tiles and re-use them in OpenLayers

texts:
OpenLayers introduced client side caching in Version 2.11 and 2.12, along side their performance improvements for other map functions. Needless to say, these functions are intertwined and rely on each other so I don't think I can just pull out the caching code and use it with version 2.10.

Is there any way to force the client browser to cache the map tiles and re-use them on subsequent requests?

I am hoping this can be done for my WMS TileCache tiles retrieved from my ms4w Apache server.

Here is a live example on my host: http://www.spatialanalysis.ca/cachemytiles/ Check out the network panel for my WMS layer and then check the OSM base layer. The OSM base layer returns "from cache" tiles once they've been loaded once. I have a feeling this is an Apache setup option.


texts: System requirement specification for Big Data Postgesql, PostGIS, and GeoServer based service

texts:
I want to start a GIS portal for around 3000-5000 visitors per day. MY database will be of an approximate size of 10 Terrabytes. The maps and layers are moderately complex, satellite-map overlay is supported. I am basically trying to project maps on the satellite images so people using my system can figure out the route they want to take. I will be using Geoserver, POST GIS, POSTGRE SQL for this purpose. Kindly, let me know what might be the hardware and/or software requirements to start such a project. Also kindly let me know the number and types of computing hardware that I may require for a smooth operation.
Thanks in advance.


texts: ERROR 000358: Invalid expression "Parcelles_class.OBJECTID" = 1 OR "output_table1.src_OBJECTID" = 1 Failed to execute (Select)

texts:
I'm executing this code. I don't know if i'm missing something,

I'm using Arcpy, Python, Arcgis 10.1

whereclause1 = '"Parcelles_class.OBJECTID" = 1'
whereclause2 = '"output_table1.src_OBJECTID" = 1'
arcpy.Select_analysis("Join_Output_3","target_neighbor_parcel",whereclause1 + " OR " + whereclause2)

I got this error while i'm executing it on PYTHON window:

  ERROR 001143: Background server threw an exception.

And this while i execute it in IDLE :

  ERROR 000358: Expression non valide Parcelles_class.OBJECTID = 1 OR
  output_table1.src_OBJECTID = 1 Ã‰chec de lâ€™exÃ©cution de (Select).


texts: Capturing points within 1km of polygon boundary

texts: I have a polygon and a set of datapoints that fall both inside and outside the polygon.  The points each have XYZ data that I want in the output.  I want to select only those datapoints that fall within 1km of the polygon edge/boundary.  Any ideas on what Arc toolbox feature to use?  I am operating with ArcMap 10.3.1


texts: Getting currently clicked layer from MapLayerTable

texts:
I've added a MapLayerTable to my application, and I'm loading multiple layers using shape files successfully. Now I would like to select features on the map of the selected layer. For this, I would need to find the currently clicked layer.

I've already used the following code to check whether the layers are active and selected.

for (Layer layer : mapPane.getMapContent().layers()) {
    if (layer.isSelected()) {
        System.out.println("Layer selected : " + layer.getTitle());
    }
}

But what I need is a method to select only one layer. E.g.: clicked layer, right clicked layer. I hope there would be a method to select one layer.


texts: Using PROJ.4 to convert from 3D Cartesian ECEF coordinates to Lat/Long coordinates with inverse Icosahedral Snyder Equal Area (ISEA) projection

texts:
I have a set of X, Y, Z coordinates in the Earth Centered Earth Fixed (ECEF) geocentric system. I need to convert them from ECEF coordinates into Latitude/Longitude coordinates using an inverse Icosahedral Snyder Equal Area (ISEA) projection. I see that this projection is available as part of the PROJ.4 library (http://proj4.org/projections/isea.html). However, I am not fully sure how to use this to transform coordinates between the two systems.

Any guidance on how to use the ISEA projection (both forward and inverse) with PROJ.4 will greatly help, along with any sample code.


texts: How to Iterate through files in a folder QGIS, using python

texts:
I am using QGIS 3.16, and I am trying to run the “service area from layer” function on approximately 400 layers (200,000 starting points in total), I am also due to be repeating this at a range of distances so looking to use python. Also, I cannot merge together as there are 1,000s km of road and just not feasible to set it running for next month.
I have been trying to run the “service are from layer” function via python and I can get it to accept my road and node layer, which have to have the same name (but are in different folders), and then output this to a different folder. The issue is that the table contains a row for each service area, but there is no geometry written. Below is the code that I have been working on but reading about I can see that there is something relating to FeatureSinks but not sure how I go about including this.
I have tried, but I think my Python skills are just a little too rusty.
My code is below.
import os
import processing
import glob
import string

#Folder where my road networks are stored
Select_folderRoads='//Mac/Home/FINAL/March2021/r1/'
#Folder where the starting points/nodes are stored
Select_folderSDH='//Mac/Home/FINAL/March2021/d1/'
#Folder where the output is to be stored
Save_results='//Mac/Home/FINAL/March2021/path250m/'
os.chdir(Select_folderRoads)

for lyr in glob.glob("*.sqlite"):

        #Individual road network layer
        inlayer = Select_folderRoads +lyr
        #Ouptut filename
        outname = Save_results+lyr
        #Layer containing start points/nodes
        startname = Select_folderSDH+lyr
        inlay = lyr[:-7]
        inlayc = inlay.upper()

        testinlayer = str(inlayer+'|'+"layername="+inlayc)

        processing.run("native:serviceareafromlayer",
            {'LAYERS':  inlayer,
            'DEFAULT_DIRECTION': 2,
            'DEFAULT_SPEED': 50.0,
            'DIRECTION_FIELD': '',
            'INCLUDE_BOUNDS': False,
            'INPUT': testinlayer,
            'OUTPUT_LINES':outname,
            'START_POINTS':startname,
            'STRATEGY': 0,
            'TOLERANCE': 0.0,
            'TRAVEL_COST': 250.0,
            'VALUE_BACKWARD': '',
            'VALUE_BOTH': '',
            'VALUE_FORWARD': ''})


texts: How do I change 'Coordinate Display' for a project, programatically?

texts:
From the user interface in QGIS, I can change the coordinate display (unit) for a project by going to:

Project > Properties... > General (tab) and under the Coordinate Display section can change the setting for Display coordinates using to Decimal degrees.

How do I do this with pyqgis?

Things tried so far:

Searching the results of QSettings().allKeys() but there are 748 - so I may have missed it.
QSettings().setValue('QgsCollapsibleGroupBox/QgsProjectPropertiesBase/mCoordinateDisplayGroup/collapsed', True) allows me to collapse and expand the group box that contains the setting - useless, but so close!
looking at QSettings().value('qgis/measure/displayunits'), but making changes in the GUI doesn't seem to change the value of this, nor does guessing with QSettings().setValue('qgis/measure/displayunits', 'degrees')

Also tried setting the destination CRS to one that uses decimal degrees, as per AWGIS's suggestion, but again it's made no difference:

QgsMapCanvas().setDestinationCrs(QgsCoordinateReferenceSystem('EPSG:4326'))
QgsMapSettings().setDestinationCrs(QgsCoordinateReferenceSystem('EPSG:4326'))

Bonus points if you can point me to a reference where I could have looked this up myself!


texts: Getting Error on GeoServer :Error obtaining new data store and Exception

texts:
Here is i added GeoServer Log

WARN [data.store] - Error obtaining new data store

Caused by: java.lang.NumberFormatException: For input string: "9090/mapd"
        at java.lang.NumberFormatException.forInputString(NumberFormatException.java:65)
        at java.lang.Integer.parseInt(Integer.java:580)
        at java.lang.Integer.valueOf(Integer.java:766)
        at com.mapd.jdbc.MapDConnection.(MapDConnection.java:80)
        at com.mapd.jdbc.MapDDriver.connect(MapDDriver.java:52)
        at org.apache.commons.dbcp.DriverConnectionFactory.createConnection(DriverConnectionFactory.java:38)
        at org.apache.commons.dbcp.PoolableConnectionFactory.makeObject(PoolableConnectionFactory.java:582)
        at org.apache.commons.dbcp.BasicDataSource.validateConnectionFactory(BasicDataSource.java:1556)
        at org.apache.commons.dbcp.BasicDataSource.createPoolableConnectionFactory(BasicDataSource.java:1545)
        at org.apache.commons.dbcp.BasicDataSource.createDataSource(BasicDataSource.java:1388)
        at org.apache.commons.dbcp.BasicDataSource.getConnection(BasicDataSource.java:1044)
        at org.geotools.data.jdbc.datasource.AbstractManageableDataSource.getConnection(AbstractManageableDataSource.java:48)
        at org.geotools.jdbc.JDBCDataStore.createConnection(JDBCDataStore.java:2081)
        at org.geotools.jdbc.JDBCDataStore.createTypeNames(JDBCDataStore.java:992)
        at org.geotools.data.store.ContentDataStore.getTypeNames(ContentDataStore.java:308)
        at org.vfny.geoserver.util.DataStoreUtils.getDataAccess(DataStoreUtils.java:97)
        at org.geoserver.catalog.ResourcePool.getDataStore(ResourcePool.java:638)
        ... 119 more

so how to solve this error and Exception?


texts: Getting geometry network data from File Geodatabase

texts:
I have a GDB file with the power distribution utility data set. I was able to get data (File Geodatabase Feature Class) from the GDB file by using QGIS's OpenFileGDB. The remaining information which I need is represented by File Geodatabase Geometry Network. The Geometry network should contain connectivity information of the features.
I have tried QGIS OpenFileGDB and also ESRI FileGDB but none of them gave me the required information.
Does anybody know how to extract this data from the GDB file?


texts:
I'm using the GDAL .Net library to open a PBF file I've downloaded from the GeoFabrik website. The file opens fine in QGIS.
However when I call the Gdal.Open method on the file, it throws an exception stating:

System.ApplicationException: '`south-africa-and-lesotho-latest.osm.pbf' not recognized as a supported file format.'

My code:
GdalConfiguration.ConfigureGdal();

var data = Gdal.Open("C:/south-africa-and-lesotho-latest.osm.pbf", OSGeo.GDAL.Access.GA_ReadOnly);

I find this weird since, when the ConfigureGdal method is called, it prints out a list of drivers to the output window, one of them being OSM PBF:

GDAL 180: PGDUMP-PostgreSQL SQL dump
GDAL 181: OSM-OpenStreetMap XML and PBF
GDAL 182: GPSBabel-GPSBabel

I've made some changes to the GdalConfiguration class that gets created when installing the GDAL Nuget package. These changes include adjusting the directories to set the correct gdal data and driver directory, as well as setting the osminfo.ini config file.
The problem is, I don't actually know what those directories should contain to know for certain I've chosen the right ones.


texts: Configure mapconfig.xml for null tag_value

texts:
I am using osm2pgrouting version 2.3.6. When using osm2pgrouting for .osm file, I created tables ways. However, when I compared it with table planet_osm_line, it lost a lot of line.

I determined that all lines having "highway" = null were not generated in ways table. I think the root cause is all lines that configured in mapconfig.xml file will be generated. But in mapconfig.xml file, it only allows tag_value to be defines specific text, not null. Ex: tag_value name="motorway" id="101"

So, my question is:

Is my determination correct or not?
If correct, how can I define value null in mapconfig.xml so that the final ways table will not lose data.


texts: Pan viewer window changes scale first time

texts:
I have written a simple tool to pan a viewer window (if present) to the centre of the map view, which reduces the number of clicks required and mouse travel to Pan Viewer to Map Location, it works mostly as anticipated but if the scale of the map viewer is set to a round number (in this case 24000) the scale is expanded on the first use but not on 3rd or subsequent uses:
 Scale set before running.
 Scale after running first time.
 Scale for further instances.

The code:

IPoint CentreOfScreen = new PointClass();
IPoint CentreOfInset = new PointClass();

// derive the centre of the map from the active view extent
IEnvelope pMapExtent = ArcMap.Document.ActiveView.Extent;
CentreOfScreen.PutCoords((pMapExtent.XMax + pMapExtent.XMin) / 2,(pMapExtent.YMax + pMapExtent.YMin) / 2);

// get the application windows, could be of a few
// types but we'll filter the one we want later
IApplicationWindows pAppWin = (IApplicationWindows)ArcMap.Application;
ISet pDatWindows = pAppWin.DataWindows;

// reset the enumeration and get the first object so that it's not null
// Leaving pDatWin unset will mean the loop will not progress
pDatWindows.Reset();
Object pDatWin = pDatWindows.Next();

while (pDatWin != null)
{
    if (pDatWin is IMapInsetWindow)
    {
        // this is the type of data window we want, not an overview window
        IMapInsetWindow pInsWin = (IMapInsetWindow)pDatWin;
        IMapInset pIns = pInsWin.MapInset;

        double pPreScale = 0;
        if (pIns.UsingZoomScale)
        {
            pPreScale = pIns.ZoomScale; // is set correctly
        }

        // derive the centre of the map viewer from its visible bounds
        IEnvelope pExt = pIns.VisibleBounds;
        CentreOfInset.PutCoords((pExt.XMax + pExt.XMin) / 2, (pExt.YMax + pExt.YMin) / 2);

        // find how far to move the edges of the extent by finding the displacement
        // of the centre of screen to the centre of the map viewer
        double dX = CentreOfScreen.X - CentreOfInset.X;
        double dY = CentreOfScreen.Y - CentreOfInset.Y;

        // apply the offset to a new envelope
        IEnvelope pNewExt = new EnvelopeClass();
        pNewExt.PutCoords(pExt.XMin + dX,   pExt.YMin + dY,
                        pExt.XMax + dX, pExt.YMax + dY);

        // assert the new extent on the inset
        pIns.VisibleBounds = pNewExt;

        if (pIns.UsingZoomScale)
        {
            pIns.ZoomScale = pPreScale; // again set correctly
        }
        pInsWin.MapInset = pIns; // reset the inset window map inset

        // refresh and yield until the inset is redrawn
        pInsWin.Refresh();
        System.Windows.Forms.Application.DoEvents();
    }
    pDatWin = pDatWindows.Next(); // go to the next data window
}

Reasserting the scale seems to have no effect; The math is fairly simple so it's unlikely I'm making a mistake there, if there were a mistake in the math I would expect the extent to keep growing on 2nd and subsequent iterations.

Has anyone who has worked with a viewer window in ArcObjects that can establish why the scale changes and knows how to keep the original scale?


texts: Cutting all polygons with each other - mega slicing

texts:
Draw a big circle.  Then a smaller one on top and then a smaller one and a smaller one.

You see a stack of circles like a bulls eye.  If you had drawn these in QGIS the other way you would only see the largest one as it hides the others.

Is there a way in a single layer to slice every feature with every other?

Basically i have contours and i only want the band between each polygon rather than the stack of progressively bigger polys like the bulls eye example.

I have perfect data where all polygons are closed and sit perfectly inside each other.


texts: DSAS cannot find database table for baseline layer

texts:
I am trying to calculate change statistics using a baseline and shorelines for 5 years in DSAS VERSION 4.4. The default parameters, transect casting stages go well and there is no trouble in detecting baselines and shorelines, but while calculating change statistics, the error message:

  unable to find database table for baseline layer 'baseline'

pops up as shown in the screenshot image below.

Could anybody guide me on this?


texts: How to capture real time weather data of points using OpenWeatherMap

texts:
I want to set points on a map layer with autofilled data in some fields in the attribute table, which already works fine with date/time and username/observer.
Is there also a way to capture real time weather data every time a point is set on the layer? (E.g. temperature, humidity and sight)
I have access to an API Key of OpenWeatherMap.


texts: How to convert Lat & Long from Onemap.sg to Google Maps

texts:
I have a set of address for Singapore, I would like to take the Lat & Long from onemap.sg (SYV21: X, Y / Easting, Northing)  and then store it in my databse , convert into (WGS84: Longitude, Latitude) used by google maps and then display on Google Maps.

I have looked at this -- http://sgonemap.codeplex.com/ but could not get much out of it.

It will be really nice if someone can help me on this.


texts:
My question is if it would be possible somehow to upload the data taken in the field using QFIELD in the cloud. I would like to know if there is any similar workflow as with, for example in ArcGIS, Collector/Survey123 upgrading the data directly into an ArcGIS Online.
I have seen they are going to release something called QFIELDCLOUD, so will this approach be similar to ArcGIS Online?


texts: Extracting time series of one pixel from the SMAP soil moisture (soil moisture in %) product

texts:
Currently I am writing on my script to extract a time series of SMAP pixel values of the band 'smp', however, my studyarea is rather small and hence the script I found online in this community from AADR
Extracting Soil moisture value using Google Earth Engine?
doesn't fit entirely for my question, since it develops on the mean of several pixels.
The script I am working with:
   var table = ee.FeatureCollection("users/marinnvum/DownloadNDVI"),
    geometry = /* color: #d63000 */ee.Geometry.Polygon(
        [[[16.62138331919588, 48.21312568441956],
          [16.619409213360917, 48.20875013496642],
          [16.62061084299959, 48.2084355385735],
          [16.625159869488847, 48.2077777398707],
          [16.62558902293123, 48.209007530577786],
          [16.624473223981035, 48.21020869299206],
          [16.625288615521562, 48.211409827232984]]]);

// Dataset

var imageCollection = ee.ImageCollection('NASA_USDA/HSL/SMAP_soil_moisture')
                .filterBounds(geometry)
                .filterDate('2018-01-01','2020-09-23');

// Select band -> smp = surface soil moisture (%)

var soilmoisture = imageCollection.select('smp');

// gget scale

var scale = imageCollection.first().projection().nominalScale();

// Take pixel // I know that here is the problem, presumably with ee.Reducer.count() --> I only need one pixel from the time series which has 324 elements

var imageCollection = imageCollection.map(function(image){
  return ee.Image(image.setMulti(image.reduceRegion(ee.Reducer.count(), geometry, scale)));
});

// Get the images from the dictionary and print it. --> I also have the feeling that aggregate_array doesn't give the relevant data output

var means = imageCollection.aggregate_array('smp');
print('Bodenfeuchte in %', means);

// put in a chart

print(ui.Chart.image.seriesByRegion(imageCollection, geometry, ee.Reducer.count(), 'smp', scale));

Map.centerObject(geometry, 10);
Map.addLayer(soilmoisture, {min: 0, max: 1}, 'soilmoister imageCollection');

I replaced ee.Reducer.count() already with ee.Reducer.toList in order to have the relevant values as a list and then, export them. But eventually my values I am looking for are NULL.
Interestingly if I //Create a list of image objects.
   var imageList = imageCollection.toList(324);
print('imageList', imageList);

I get provided with all relevant images, however, smp or ssm or whatever band you look for, has no values at all if you check them:
   print('imageList'+imageList.length)


texts:
I am Very New to GML format.

first of all, I have came to the assumption that a schema file will display gml in a particular style as defined by the schema, similar to the way i load a .qml on a shapefile in Qgis.

If this is somewhat correct my problem is;

I can load GML into qgis (as vector layer) and see point and line layers and also data layers with no geometry that refer to said points and lines.

I can view GML text and see the schema location/details.

I also have schema as an .xsd

How do I view a GML as more than just points and lines or an xml doc?


texts: Calculate distance for point to point based on value in field of first point

texts:
I am trying to figure out the best way to calculate the distance for each student from the school they attend. We have 2 shapefiles - the student locations with a field identifying the school they attend, and the school point file.

The best we could think of was to use the Proximity toolbox – but then we would need to do each set of students for each school individually I think.


texts:
I need a California roads network suitable for building a network dataset and routing to 120 locations. It needs to at least cover the Northern half of California. 2010 Census roads need alot of cleanup or I would try using those. Any ideas for places to start?

Thanks,


texts: Adjusting position of individual graduated symbols in ArcMap?

texts:
The house symbols in the map are overlapping each other.

Is there any tool in ArcMap that I can manually move symbols individually?


texts:
Running algorithm via window operation for output, it indulge pixel shift.

This pixel shift is 2 rows down and 1 column right.

CODE:

def GetToUTM(Data, pixel_sz, Zone):

gdal.Open('\\myData\\asterDem.tif')
Past_ref = osr.SpatialReference()
Past_ref.ImportFromWkt(Data.GetProjectionRef())
newref = osr.SpatialReference()
newref.SetUTM(abs(Zone),Zone > 0)
transform = osr.CoordinateTransformation(Past_ref, newref)
geotransform = Data.GetGeoTransform()
nx, ny = Data.RasterXSize, Data.RasterYSize
print nx,ny
(ulx, uly, ulz) = transform.TransformPoint(geotransform[0], geotransform[3])
(lrx, lry, lrz) = transform.TransformPoint(geotransform[0]+ geotransform[1] * nx, geotransform[3] +
                                           geotransform[5] * ny)
getfile = gdal.GetDriverByName('GTiff')
dataout = getfile.Create('PyAster.tif', int((lrx-ulx)/pixel_sz), int((uly-lry)/pixel_sz), 1, gdal.GDT_Float64)
newGeoTrnsm = (ulx, pixel_sz, geotransform[2], uly, geotransform[4], -pixel_sz)
dataout.SetGeoTransform(newGeoTrnsm)
dataout.SetProjection(newref.ExportToWkt())
gdal.ReprojectImage(Data, dataout, Past_ref.ExportToWkt(), newref.ExportToWkt(), gdal.GRA_NearestNeighbour)

dataout = gdal.Open('\\untitled\\name_sake.tif')'

The Image is getting re-projected and pixel issue is resolved under my threshold.

Now that original image has value of range from 100-1211(float64), re-projected one report it from 0- 1190.43(float64) which is incorrect.

How can I make the necessary correction?


texts: Removed features on openlayers 3 is re-rendered on zoom in/out

texts:
In map I have Vector with Features inside, upon selecting some other options those features should be changed, therefore, I grab the map, select vector that I need, grab source from that map, and then I can remove feature. There are multiple ways that can be done:

this.map.getLayers().a[i].getSource().removeFeature(this.map.getLayers().a[i].getSource().B[2]);

or

this.map.getLayers().a[i].getSource().clear()

Where i is iteration through all layers.

In both cases, selected feature/features are removed from map. BUT if I try to zoom in/zoom out, that features that I have removed are re-rendered again. Is that some sort of cache? How can I remove them permanently?


texts: .IMG file conversion to a QGIS readable format

texts:
I recently found QGIS and enjoy being able to play with it, most importantly it can convert almost anything to .mbtiles, a format that is readable by foreflight on iOS.

Now I have an OLD .IMG file that we use in our GPS to follow a certain powerline. My knowledge of gps data is limited and that file is the only copy I have of that

I would like to read that file in QGIS, but I can't find anything online, only that... well I can't

I have been able to open the file in GPSMapEdit, I can see the line, I can see the poles and everything.

I have been able to convert parts of it to a .kml and then upload it to QGIS but a lot of information from the original file is not there

I've thought about saving it to a .mp file from MapEdit but still.. no luck..

TL:DR I want to convert a Garmin IMG to a readable format for QGIS


texts: Add help tab in plugins from the toolbox

texts:
Is it possible to add in QGIS 3 (Python Script) a help tab on the plugins created for the toolbox as it was in QGIS 2?

Looking in the API I found the parameters helpString that will be deprecated in 3.4, shortHelpString which presents the short help on the right side of the main plugin window, and helpUrl that creates a button with an external link.

Is there any other parameter that enables the use of a help tab?

I would like to create this tab to make it easier to read the plugin's manual even without internet connection and without opening external files. Being able to display a local HTML would be a great option.


texts: Moving a shape around the earth while maintaining the proper projection

texts:
What is the transformation needed to be taken on every data point of a shape to
rotate that shape around the globe.

I'm trying to figure out the math behind this kind of transformation so I can create something similar to this using leaflet and d3 https://gmaps-samples.googlecode.com/svn/trunk/poly/puzzledrag.html

Edit: To be clear, I am looking to achieve this using leaflet, D3, and javascript not desktop GIS software. My end product will a function that takes a point (lat,lng) shifts it by an amount (shift_lat,shift_lng) and returns (new_lat,new_lng)


texts: How can I get timestamps from Garmin's .gdb files?

texts:
I have got the some Gps data in the form of Garmin's .gdb files. I have converted it to CSV files using GPSBabel.

The CSV output contains Latitude, Longitude, and some kind of free form text remarks. I do not see any timestamps.

Do  Garmin's .gdb files contain timestamps?

This is the first time, I'm working with this file format, so my understanding of the format/data is spotty at best.

Where can I find more information about this format? How can I get the TimeStamps of the Gps points?


texts: Normally, when you use a hillshade as a layover in ArcGIS, you set the transparency at 80% or so (in the Layer Properties, display settings). As a result, the underlying layers, for example land use, will display with a grey haze, even the plain areas. How do I get a Hillshade where plain areas (middle grey tones) will display as transparent?


texts: Is there a source for historical shapefiles and tiles? Like the ancient Greece/Rome, medieval Europe and Japan? Ideally with multiple snapshots during different periods in their development, since they've gone through quite a bit of restructuring over their lifespans.


texts: Calculate row differences in seconds of a datetime column with ArcGIS

texts:
I have a point feature class in ArcGIS desktop. There is a timestamp column named TS in datetime format in the attribute table.

I need to calculate and write time differences of rows from previous rows to a new column named RowDiff. Results will be in "seconds". RowDiff datatype is short integer.

I've found these links but I have no idea about implementing the functions:
https://desktop.arcgis.com/en/arcmap/latest/analyze/arcpy-time/esritimedelta-class.htm
https://www.programiz.com/python-programming/datetime

I've found the following code but results are in string format like 0:00:06.00005 and useless for me. What I need is row differences in seconds and an integer data type, not decimals.

import arcpy

from arcpy import da

gpsdata = arcpy.GetParameterAsText(0)

#set firstRow variable to true
firstRow = True

with arcpy.da.UpdateCursor(gpsdata, ["OBJECTID", "TS", "D_T"]) as cursor:

    for row in cursor:
        gap2 = row[1]

        #check if its the first row
        if firstRow == True:
            gap1 = gap2
            #set firstrow to false
            firstRow = False
            continue
        row[2] = str(gap2 - gap1)
        gap1 = gap2
        cursor.updateRow(row)

Sample results should be as follows:


texts: Looking for a more effcient process for retrieving information

texts:
Using ArcGIS Desktop Standard 10.2.

Need to retrieve specific information from multiple (10+) feature classes that are adjacent to/overlap/ or in close proximity to specific areas of interest (polygons).

Currently using 'Identify Features' from the basic toolbar, drawing a polygon around the area of interest and then filtering through the attributes of the features that are returned. Then the specific information that is needed is written down on a template form and later summarized in a formal report.

This is very tedious and time consuming considering the large backlog of these areas of interest. Need a tool or script that would autorun these areas of interest against a select group of feature classes and output a report with the information I am looking for. Tried the 'Intersect' tool but the license level only allows 2 input feature classes. Do not have access to Spatial Analyst tools with a Standard license.

Does such a tool or script exist? Is there a more efficient alternative?


texts: Extract raster pixel values from shapefile in QGIS

texts:
I have been encountering difficulties trying to extract the raster pixel values from a polygon shapefile that has been assigned the same CRS.

The problem under inspection is how to get a list or a table with coordinates and pixel values from a polygon shapefile over a singleband raster for further investigation of values.

What is the easiest and most efficient way to perform this otherwise seemingly very simple task?


texts:
Are the new(ish) NADCON 5.0 transformations available in ArcGIS?

If so, at what versions?  If not, can they be added?

Specifically, I am looking to go from NAD83(2011) to NAD27.

I am wondering if I need to try and get the NADCON 5.0, uses Esri's built-in two step transformation in 10.4.1 (via WGS84 ITRF08 or is it just a dummy?), or just setting the data to NAD83 and using the base NADCON.

The two step that is available to me in ArcMap now is WGS_1984_(ITRF08)_To_NAD_1983_2011 + NAD_19827_To_WGS_1984_79_CONUS

How do these two transformations relate to NADCON or NADCON 5?


texts: OpenLayers 3 displays WMTS tiles slightly downscaled

texts:
I'm having trouble with a WMTS source in my OpenLayers 3-driven web application.

The map is shown and the layer is rendered, but the 256x256 px tiles are rendered slightly downscaled. For instance, the section covered by this tile

looks like this when screen grabbed from the map

(216x216 px)

The source in question is set up using WMTS' getCapabilities, so there's minimal manual labour going on:

function getSource(getCapabilitiesUrl) {
    var parser = new ol.format.WMTSCapabilities();
    return fetch(getCapabilitiesUrl)
        .then(function (resp) {
            return resp.text();
        })
        .then(function (xml) {
            return parser.read(xml);
        })
        .then(function (wmtsCap) {
            return new ol.source.WMTS(ol.source.WMTS.optionsFromCapabilities(wmtsCap, {
                layer: 'my-layer',
                matrixSet: 'EPSG:32633'
            }));
        });
}

The projection is defined like this:

proj4.defs('EPSG:32633', '+proj=utm +zone=33 +ellps=WGS84 +datum=WGS84 +units=m +no_defs');

Can anyone point me in the right direction as to how I can get OpenLayers to display the tiles at the right scale (i.e. unscaled)?

The problem exists in a clean environment without any CSS styling, so this should not be the source of the problem. I'd make a live demo for you to try, but the geoserver is unfortunately not available to the public.

EDIT

This is the (anonymised) request for the 256x256 tile shown at the top of the post:
https://example.com/geowebcache/service/wmts?layer=my-layer&style=&tilematrixset=32633&Service=WMTS&Request=GetTile&Version=1.0.0&Format=image%2Fpng&TileMatrix=32633%3A5&TileCol=9&TileRow=12


texts: GDAL and Python - Don't print GDAL Error messages

texts:
According to the GDAL Python Documentation, GDAL functions do not raise exceptions or errors but do print errors to the console via stdout and stderr.

I have some code in my script that calls GetStatistics on a GDAL raster band.  If the band only contains nodata cells (which is possible because it is a tile of a larger dataset) the GetStatistics call fails with an error, but the script continues.  I would prefer not to see the error messages that are still printed to the console when this happens as it interferes with other print statements that I have written into my own script.

Is it possible to turn off these messages?


texts:
I want to perform a line density analysis for my line data, based on an attribute.

The heatmap option (kernel density analysis) seems to be only possible for point data.

I have tried to convert my line data into points, and create a heatmap, but the attributes for the weight were lost on the conversion.

This is a line density analysis performed with ArcGIS. I was looking to do the same with QGIS.

Is it possible to perform a similar line density analysis in QGIS?


texts: How to show Turn Restriction data from OS MasterMap Highways on the map?

texts:
I want to obtain the location and details of all turn restrictions from the OS MasterMap Highway Network datasets I have.

I tried loading the Turn Restriction dataset into QGIS but it does not read the location of the features, as the dataset itself does not contain any geometry.

I know that all Turn Restriction features are related by GML-ID to the features in the Road Link dataset, but I was told that this could not be done within QGIS as the software does not understand the links in the data.

I tried translating both datasets from GZ format to SHP using Astun Technology Loader, hoping the Turn Restriction output would contain geometry, but it failed, exporting only a DBF file (see below). So I guess translating to SHP won't help?

So How can I link those two datasets and get the turn restriction shown on the map in QGIS? Preferably methods without much coding?


texts:
I am looking to update row values within specific fields of a polygon geopackage using an if statement--similar to how Esri users would use an Update Cursor . For example, say there is a text field called campsite with values ranging from 1W - 10W. I would like to change the value 10W to 20.

There is an example provided by @gene which seems to write a single value for each row in the attributes.

for feature in fiona.open(a shapefile): # = for row in cursor:
    feature['properties']['Elevation']= 15 # = row[1] = 15

However, when I implement the solution using a geopackage the attributes do not change. Here is my attempt:

import fiona

gpkg = r'/path/to/campsites.gpkg'

for feature in fiona.open(gpkg): # = for row in cursor:
    if feature['properties']['Campsite'] == '10W':
        feature['properties']['Campsite'] = '20'

How can I update specific row values within a polygon geopackage using if/else control flow?


texts: How to match MODIS reflectance files (e.g.,MOD09GQ) with correct cloud mask (e.g., MOD35_L2)

texts:
I hope to be working with MODIS data, mostly the V6 250m reflectance files (e.g., MOD09GQ). The naming convention of these files makes it easy to identify which .hdf covers my area of interest (i.e., MOD09GQ.A2016341.h18v04.006.2016343073521, where h18v04 is the corresponding tile location).

I will need to supplement the reflectance data with other datasets, such as the cloud cover product MOD35_L2. The naming convention for these .hdf files do not include the tile location, but instead a 24hour data start time (i.e., MOD035_l2.A2016341.1125.006.2016341174916, where 1125 is the data capture time).

How do I know which supplimentary .hdf file to download, that corresponds to the location of my gridded data?

My guess is that the gridded reflectance products metadata indicate a capture time for each corner. Does this metadata time correlate with the associated supplementary MODIS products?

Will I therefore need to repeat the process for the other three corners as well, to cover the whole area?

Is there an easier way to identify which suplolimentary .hdr files I need?


texts: Auto delete surrounding objects from different layers

texts:
I'm quite new to GIS technology so pardon me if I make any mistakes or use the wrong terms here.

So, I got a map, with point layer A, B, C.

A is the main point layer. I would like to make something like this but got no idea how I'm gonna do this.

When I start an edit session, click on a point in layer A, delete it, all surrounding object within a radius of r meters in layer B and C would be deleted as well.


texts: QGIS: How to invert the values in a layout legend?

texts:
At 3.16.3 I have added this legend to a map layout.  Notice how the lowest value is at the top, with values increasing towards the bottom:

I would like to instead have the lowest value at the bottom, with values increasing towards the top, maintaining their associated color patches (blue = lowest value, pink = highest value).  I used MS Paint to flip the image above to show the desired effect (of course this flipping also flips the text so that it is unreadable):

The QGIS documentation does not discuss such re-ordering the legend.
Nor have I found the ability to re-order in the Layer Styling panel, where the lowest value always is placed at the top:

Note that I'm not talking about inverting the color ramp.  Rather, I want to invert the value order.


texts: Why does ArcGIS Server 10.1 change path variable to g_ESRI_variable_1?

texts:
I've had a problem with my geoprocessing service finding the data it needs for execution. I found out what's causing this is ArcGIS's assignment of a new variable "g_ESRI_variable_1" in the new copy of the code that lives on the server. However when I try to edit the code back to its original state it says "Permission Denied". I have already given rw permissions to my publisher server on the file server hosting the data and registered the folders in Server Manager.

Here are the top few lines of my original script:

import os, regUtil, regQueries, urllib, sys
from arcpy import GetParameterAsText

outpath = GetParameterAsText(1)
root_path = r'\\SERVERNAME\\Data\tiger\TigerData'

And here are the top few lines of the copy in the arcgisserver folder automatically made upon creation of the geoprocessing package:

import os, arcpy
g_ESRI_variable_1 = os.path.join(arcpy.env.packageWorkspace,u'tigerdata')
# Esri end of added variables

# Esri start of added imports
import sys, os, arcpy
# Esri end of added imports

import os, regUtil, regQueries, urllib, sys
from arcpy import GetParameterAsText

outpath = GetParameterAsText(1)
root_path = g_ESRI_variable_1

It seems very weird that ArcGIS would go in and change a variable automatically. Why does this happen and how would I either prevent this from happening or edit the script back to the way I wrote it?


texts: SLD styling rule for first and last geometry feature in PostGIS table

texts:
Writing a SLD file to style a Points layer which I publish as WMS with GeoServer, is it possible to specify a rule that is only applied to the first and last features in the underlying PostGIS table?

The table is sorted by means of that the point features represent a road track which has a start- and endpoint. My goal is to use a specific style for the start- and endpoint, different from how all the other points are styled. However, the table's length is dynamic so that the first feature will always have featureid==1 but the last one is unknown in advance and thus cannot be hard-coded.

Is there any possibility in SLD to filter for the first and last feature in a column? Naturally, the features of interest have the lowest and highest featureIDs out of all features so something like min(featureid) and max(featureid) would also work, if available.

Access to the database table could also be provided if it is possible in any way to filter for the result of a SQL query like SELECT min(featureID).

This part of the documentation refers to functions called startPoint and endPoint but apparently these can only style the start- and endpoint of a linestring, not the first and last point out of a column of points.


texts:
I have coordinates in QND95 (Qatar National Grid) and I need to convert it to WGS84, specifically EPSG:4326. There are already many tools available for this, but the problem is I want to use the conversion in my mobile application.
So basically I'd like to know how to convert QND95 to WGS84, any formula to transform them, which I can use in my application.


texts: MapServer - Show Polygon Layer with non-transparent border and transparent background

texts:
In MapServer, I want to show a polygon WMS layer (watersheds). I want to show the boundaries (outlines) of the polygons as a non-transparent thick line. And I want to show the inside (fill) of the polygons with 50% transparency. Is this possible?

I tried to set OPACITY 50, but this changes the transparency of both the outline and the fill.
The second option I tried was to set TYPE line, but this shows the fill as complete transparent and the GetFeatureInfo no longer works by clicking inside the polygon.

What I would like to get, is:
The outlines of my polygons have 0% transparency.
The fill of my polygons has 50% transparency.


texts: How to define output format of raster calculator in QGIS?

texts:
Is it possible to define the format of my target raster in QGIS raster calculator?

I have an Int16 raster which becomes float after processing the following command:

myraster@1>0

As I'd like to use GDAL sieve (which only works with Integer rasters), it would be great not to change the raster format by hand again.


texts:
How can I create something like this as a layer in QGIS?
I currently have one large Grid set with ID'#'s, and I'm wanting to create nested layers of 5 x 5 grids 5 levels down. (so for example grid A is 62500m x 62500m, and within grid A I want 25 tiles of 12500m x 1500m, and in each of those I want 25 tiles of 500m x 500m, and in each of those I want 25 tiles of 100m x 100m) with an ID specific to the tile it's embedded within (ex: 2_3_4_5_A)


texts: Memory issue loading raster into numpy array

texts:
I am trying to open a historical imagery tiff dataset into a numpy array so I can do some image processing on it afterwards. However, it gives me 'unspecified error' when trying to run the (very) short code. Is this a memory issue? Or am I missing a crucial step? I have put a screenshot of what I'm working with below.
import numpy as np
import arcpy
from arcpy import env
import arcpy.sa as sa
arcpy.CheckOutExtension("Spatial")
env.workspace = r"Z:\adit\Desktop\2020-2021\Remote Sensing\project\data\1942_Aerials_Central_Sonoma_County"
env.overwriteOutput = 1

historical_imagery = "1942_Central_Sonoma_County.tif"
historical_imagery_rast = sa.Raster(historical_imagery)

#print historical_imagery_rast.mean
#print historical_imagery_rast.standardDeviation

hist_image_raster = arcpy.RasterToNumPyArray(historical_imagery_rast)


texts: How to know why a pgRouting function returns "FAIL"?

texts: I am getting a "FAIL" message back when executing pgr_createTopology on a particular table, but the reason is not obvious. Is there any way to know precisely why the function is failing? Perhaps some sort of pgRouting log?


texts: Alternatives to ArcGIS Image Server Extension to serve Mosaic Dataset?

texts:
I am wondering if anyone has any resources for alternatives to ESRI's ArcGIS Image Extension for Server?

We are looking to serve JP2000 images that have been consolidated using the Mosaic Dataset Method in ArcGIS.

any thoughts?


texts: Count the number of raster cells of each land use type in neighborhood

texts:
I have a raster about different land use types(5 types),each raster cell has a categorical value(1-5). ArcGIS Focal Statistics can calculate mean, range, and even variety of each processing cell, but I want to count the number of  each land use type for the processing cell in the neighborhood area.

For example,the processing cell in pic has 3 land use type1, 2 type2, 2 type3, 1 type4 and 1 type5. (Neighborhood: 3*3 rectangle),how can I do this?


texts: Defining CRS with tifffile - export .tif

texts:
I am trying to figure out how to define/keep my images CRS throughout a splitting script I have written. I tried adding metadata={'init': 'epsg:3857'} when calling tifffile.imsave(dir_out, image), but when I view the split images in QGIS, the CRS is still undefined.
Is there something I am missing here?
import os
import tifffile

def tile(filename, dir_in, dir_out):
    name, ext = os.path.splitext(filename)
    img = tifffile.imread(os.path.join(dir_in, filename))

    windowsize_r = 10000
    windowsize_c = 10000

    i = 0
    for r in range(0, img.shape[0] - windowsize_r, windowsize_r):
        for c in range(0, img.shape[1] - windowsize_c, windowsize_c):
            window = img[r:r + windowsize_r, c:c + windowsize_c]
            out = os.path.join(dir_out, f'{name}_{i}_{ext}')
            tifffile.imsave(out, window, metadata={'init': 'epsg:3857'})
            i = i + 1

tile('Orthomosaic_export_MonFeb01193821460106.tif',
     r'D:\ortho',
     r'D:\model_images')

It seems that when I am slicing my image window = img[r:r + windowsize_r, c:c + windowsize_c] that extent is what is being saved in place of my image. When in fact I want the extent from the image.


texts:
Could someone please point me in the right direction for creating a geotiff DEM with Geotools that from gridded X,Y,Z coordinates.

I am currently doing this with:

GridCoverage2D coverage = new GridCoverageFactory().create("DEM", dem, WORLD);

Values that are integer,integer,double for x,y,z however, after re-projecting data the coordinates can be double,double,double for x,y,z. GridCoverage2D will only accept integers for x,y. I can round these to integers however I am not sure that is the correct approach.

How can I use the in built functions to create a raster dem without losing precision?


texts: Problem with convert georeferenced jpg to geotiff

texts:
I'm beginner with GIS tools. I have jpg file with GPS data and I need reference this file to Microstation as georeferrenced raster.

I tried use function Raster -> Conversion -> Translate with defined target SRS but without success. When recieved geotiff is referenced to Microstation he is not visibile. It is visible when SRS is set in Microstation model but problem is that raster lands on 0,0 coordinates. I enabled GDALTools plugin.

Thanks in advance for help. Regards, Krzysztof


texts:
I am compiling QGIS in Ubuntu 20 with Ninja. Compilation build succesful. Now I have made some changes in UI file and .cpp file and want to recompile again like we do in Visual Studio.
How to do this?


texts: Selecting features with NULL Geometry from Shapefiles using QGIS?

texts:
I have been given a Shapefile that contains records with NULL geometry, along with records with non-NULL geometry. This situation appears to cause problems in QGIS (1.8 and 2.0) when editing the geometries.

I have been able to remove the NULL records by using the Marquee Select tool in the Map pane to select all the features in the dataset and saving the selection as a new shapefile, what I would like to know is how  to identify and delete the NULL geometry records in situ in QGIS. I have been unable to do this using the 'Select By Expression' functionality in QGIS 2.0 with expressions such as the following:

 $area = 0

 $geometry IS NULL

 geomToWKT(  $geometry ) =  ''

I have uploaded a sample TAB file containing five records, two of which (ID 3 and 4) have null geometries. The format is a MapInfo TAB file as I am having problems creating a Shapefile that preserves the NULL geometries! (I am not sure of the provenance of the Shapefile that sparked the question to begin with)


texts: How to Merge, Buffer, and then Intersect multiple features based off of the Buffer feature in ModelBuilder

texts:
I have created what I thought to be a simple ModelBuilder model and the final output (.shp format) is always an empty. As you will see in the ModelBuilder image below:

Merged 3 track logs now called TrackMerge
Buffered the newly created merged file, TrackMerge, with a 500 meter buffer. It is now called BufferOutput.
I performed an intersect of the BufferOutput.shp with 3 separate polygon shapefiles that are in a folder or Workspace called WeedTracksMergeScript. I did this by using the Iterate Feature tool. I did this by selecting the Workspace called, WeedTracksMergeScript. I gave it a Wildcard name of *Infest because all of the file names that I wanted to be iterated began with Infest.
The output below is simply called Outputs.shp

I've tried numerous renditions of this and have not been able to make the final Intersection of the 3 additional shapefiles with the shapefile BufferOutput.shp. I've tried skipping the Iterate feature tool and feeding the 3 shapefiles directly into the Intersect tool to intersect with BufferOutput.shp and the output has also been empty.

I know the BufferOutput.shp intersects with the polygon shapefiles in the particular file or Workspace called WeedTracksMergeScript. I have visually inspected this to make sure these shapefiles intersect the BufferOutput.shp.

Can anybody point out what I might be doing wrong?


texts: Is there a way to export online ESRI basemap layers to raster?

texts: I have been working on polygonizing the vegetation areas in some counties from imagery. When I have the mosaic raster locally I can easily get the polygons using Image Classification Toolbar. But when it comes to using online ESRI imagery I don't know what to do. I was wondering if someone could suggest me a way to help in saving these basemap as raster locally. When I right lick on the layer it doesn't give me the Export option.


texts: Moving cursor backward i.e. opposite of row = rows.next()?

texts:
Can anyone tell me if there is a way to go up one row in python cursor?

i.e. opposite of row = rows.next()


texts: QGIS Server: Project name shown as layer. No other layers are shown

texts:
I have installed QGIS Server on Ubuntu 14.04

I do get a valid WMS response trying to reach

http://localhost/cgi-bin/qgis_mapserv.fcgi?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetCapabilities

Next step I've done is trying to add a *.qgs map file and loading it back in QGIS 2.18.

There should be 2 layers served as WMS in that particular project, however I only get to see one called identified by 0 which is also my project title...

How come my layers are not showing up ?

 <Layer queryable="1">
<Name>gisoost.be</Name>
<Title>Gisoost.be</Title>
<CRS>CRS:84</CRS>
<CRS>EPSG:31370</CRS>
<EX_GeographicBoundingBox>
<westBoundLongitude>3.32657</westBoundLongitude>
<eastBoundLongitude>3.53157</eastBoundLongitude>
<southBoundLatitude>51.007</southBoundLatitude>
<northBoundLatitude>51.1356</northBoundLatitude>
</EX_GeographicBoundingBox>
<BoundingBox CRS="EPSG:31370" maxx="91244.8" minx="77053.8" maxy="203167" miny="189038"/>
</Layer>


texts: Where can I find free Bathymetric data for Indian Ocean?

texts:
I am looking for bathymetric data, to place in a map of India for visualization only.

I don't really need very high resolution data, since the the map will be at a very small scale: Maybe 1:1Million or 1:5,00,000.

Which are the sources that I can freely use?


texts: How to convert NAV4 coordinate system to WGS84?

texts:
I would like to know whether it is possible to convert the coordinates given in a NAV4 coordinates system into WGS84. I could not find out what NAV4 is and if it is actually a standardized CRS.

Here is an excerpt from the XML document which contains the coordinates:

<stop stopID="2309" name="Rasthaus" locality="Unterwittelsbach">
  <place selected="0">
    <placeID>36</placeID>
    <placeName>Unterwittelsbach</placeName>
    <OMC>9771113</OMC>
  </place>
  <coords>
    <coordSystem>NAV4</coordSystem>
    <longitude>4437685</longitude>
    <latitude>788825</latitude>
  </coords>
</stop>

Is there a EPSG reference which uniquely identifies NAV4? I could not see any helpful information in the document which help to identify the CRS used.


texts: QGIS 3.2 using a USGS Topo for a basemap

texts:
On QGIS 3.2 I add, through the XYZ Tiles tab in the browser, this URL:

http://services.arcgisonline.com/ArcGIS/rest/services/USA_Topo_Maps/MapServer/tile/${z}/${y}/${x}.jpg

All that is displayed when it is added as a layer is a page of postage stamp size (5x5) pictures of the entire world. Moving or zooming simply refreshes the same view of a 5x5 squares of the world. I've tried this a number of ways and all FAIL. Any suggestions would be helpful.

The goal is to have USGS Topo Maps as a base map.


texts: How to define disturbed/clipped contours elevation for Volume Polygon?

texts:
I want to calculate surface volume for different layer on TIN surface.
I am working on an area in Norway which is 1 km * 1 km. The aim is to calculate surface volume of different imaginary elevation.

Because it is a clipped area of a big Island,some contours have been clipped as well ! Now the problem is those areas that have been clipped considering as zero or plane surface, although they are almost high hills and have around 140 m height !

Is there any way to prevent this error or fix the disturbed contours??

I attached the picture for become more clear!


texts: Label point layer with attribute of surrounding line features (QGIS 3.x)

texts:
I would like to label the points of my point layer with the name attributes of the surrounding line features. The result should be a list, with the names of all lines next to the point (f.e. 10 map units).
In the screenshot you can see, which result I do expect.


texts: Populate many-to-many intermediate table betwen two non spatial tables in ArcGIS

texts:
Is it possible to populate a many-to-many intermediate table between two non spatial tables in ArcGIS?

The intermediate table shows in ArcMap, but the origin and destination tables do not show in Attributes.

[EDIT]
Following the recommendations of @Ben S Nadler and the instructions in the link provided:

  Open the Attribute dialog box by clicking the Attributes button on
  the Editor toolbar. The selected features will appear in the left
  panel of the Attribute dialog box.

the Attribute dialog box does not show any tables and the intermediate table can not be populated.

The non spatial table to table many to many relation was finally simulated through creating the intermediate table and its fields manually and two one to many relations with it.

However, the relationship is therefore not working properly:

the intermediate table can only be filled manually or relation by
relation and,
using the relation icon in the attribute tables of the
origins only brings you to the intermediate table and therefore
require one more step.


texts: Changing the legend element anchor point using ArcPy

texts:
Is it possible to changed the legend element anchor point using ArcPy?
This is the default anchor point:

and i want this anchor point location:

I know how to changed the legend element position x,y. But if the anchor point isn't in the right place, then i will get unwanted result.


texts: Getting all information about ways from python overpass library?

texts:
I'm having problems getting started with the overpass QL python library. I want to get a collection of ways that match the search criteria, but I can't seem to get anything but the points. I need to know which points belong with which ways so that I can draw them properly.

$ ipython
In [1]: import overpass
In [2]: api = overpass.API()
In [3]: ways_featureset = api.get('way(42.819,-73.881,42.820,-73.880);>;')
In [4]: ways_featureset['features'][0]
Out[4]: {"geometry": {"coordinates": [-73.8925346, 42.8155621], "type": "Point"}, "id": 42001122, "properties": {"highway": "traffic_signals"}, "type": "Feature"}
In [5]: non_points = [way for way in ways_featureset['features'] if way['geometry']['type'] != 'Point']
In [6]: non_points
Out[6]: []

I think that I need to do something similar to the "out skel" call, but you aren't able to do that in the library according to the docs (because the library handles appending out for you).


texts: Converting SDE.ST_GEOMETRY polyline to SDO_GEOMETRY with SQL: Why are there so many nulls?


texts: Calculate a field using text from another field separated with comma

texts:
I have an address field containing the following text: "64,ramkoti,hyderabad,500070"

I need to copy hyderabad to another field.

What would be the code used in the field calculator tool to achieve this?


texts: How to find stream order in QGIS 2.18.12?

texts: In QGIS 2.18.12. I read that using these two commands r.stream.order or v.stream.order one can find stream order. But I am unable to find them in my processing toolbox?


texts: Fulltext Search with Spatial Search with pagination

texts:
We tried to use postgresql+postgis to save our spatial data with points/lines/polygons.

Though the postgresql have a built-in full-text search engine, we have to use lucene for some reasons, and though the Lucene have a spatial engine, it does not support the line/polygon interception. :(

However we need the full-text search combined with the spatial search and do the pagination.

How do you handle this kind of problem?


texts:
I'm having trouble rendering my WFS features in projection EPSG:4326, the below code works when I switch projection to EPSG:3857 (and change the center), it works fine. I've checked the response from the server and its coming back as EPSG:4326 format.  The projection in my Geoserver for this layer is EPSG:4326.  Am I missing a setting somewhere in the openlayer code?

import WFS from 'ol/format/WFS';
import GML from 'ol/format/GML';
import Tile from 'ol/layer/Tile';
import Vector from 'ol/layer/Vector';
import {bbox as bboxStrategy} from 'ol/loadingstrategy.js';
import Select from 'ol/interaction/Select';
import DragPan from 'ol/interaction/DragPan';
import Modify from 'ol/interaction/Modify';
import {Draw} from 'ol/interaction';
import VectorSource from 'ol/source/Vector';

import MouseWheelZoom from 'ol/interaction/MouseWheelZoom';
import {click, pointerMove, altKeyOnly} from 'ol/events/condition.js';
import Style from 'ol/style/Style';
import Stroke from 'ol/style/Stroke';
import Snap from 'ol/interaction/Snap';
import {Map, View, Feature} from 'ol';
import OSM from 'ol/source/OSM';
import {fromLonLat} from 'ol/proj';

var formatWFS = new WFS( {
}
);

var formatGML = new GML({
    featureNS: 'http://wamp-dev/us',
    featureType: 'crop_other',
    srsName: 'EPSG:4326'
});

var xs = new XMLSerializer();

var sourceWFS = new VectorSource({
    loader: function (extent) {
        $.ajax('http://192.168.5.106:91/ows', {
            type: 'GET',
            data: {
                service: 'WFS',
                version: '1.1.0',
                request: 'GetFeature',
                typename: 'US:crop_other',
                srsname: 'EPSG:4326',
                bbox: extent.join(',') + ',EPSG:4326'
            }
        }).done(function (response) {
            console.log(response);
            sourceWFS.addFeatures(formatWFS.readFeatures(response,
                {dataProjection: 'EPSG:4326', featureProjection:EPSG:4326'})
        );
        });
    },
    strategy: bboxStrategy
});

var layerWFS = new Vector({
    source: sourceWFS
});

var interaction;

var interactionSelectPointerMove = new Select({
    condition: pointerMove
});

var interactionSelect = new Select({
    style: new Style({
        stroke: new Stroke({
            color: '#FF2828'
        })
    })
});

var interactionSnap = new Snap({
    source: layerWFS.getSource()
});

var map = new Map({
    target: 'map',
    controls: [],
    interactions: [
        interactionSelectPointerMove,
        new MouseWheelZoom(),
        new DragPan()
    ],
    layers: [
        new Tile({
            source: new OSM({
                url: 'https://cartodb-basemaps-{a-d}.global.ssl.fastly.net/dark_nolabels/{z}/{x}/{y}.png',
                opaque: false,
                attributions: []
            })
        }),
        layerWFS
    ],
    view: new View({
        projection: 'EPSG:4326',
        center: [-119.934411,36.627444],
        //center: [-13350873.59,4401244.07],
        zoom: 15
    })
});


texts:
QGIS 1.8.0.
Scipy installed.
Statsmodels installed.
Using Sextante to try to calculate MCP.

When I run the MCP calculation from the Toolbox, a layer is added to the file, but the layer doesn't contain any information.

The input layer contains a set of points and the field I select from that layer contains OS grid references.

Does anyone have any ideas what I'm doing wrong? Any help much appreciated.


texts: How to get pixel size of Data Frame in arcMap layout?

texts:
The data frame where the map diagram will be presented in, appeared on the screen not in the actual size which is part of A4 page template, I want to get the actual size by getting the pixel size and convert it to centimeters.

How to get the pixel size of Data Frame?


texts:
I want to calculate azimuth between two points so I did it this way like described here :
How to calculate the angle at which two lines intersect in PostGIS?

DROP TABLE IF EXISTS angletmp;
CREATE TABLE angletmp AS
WITH
intersections AS
(
    SELECT
        ST_Transform
        (
            (ST_Dump
            (
                ST_Intersection
                (
                    a.geom,
                    b.geom
                )
            )).geom,
            3785
        )
        AS intersection,
        a.geom as geom1,
        b.geom as geom2,
        a.id as id1,
        b.id as id2
    FROM
        pipe a,
        pipe b
    WHERE
        ST_Intersects
        (
            a.geom,
            b.geom
        )
),
buffers AS
(
    SELECT
        intersections.intersection,
        ST_ExteriorRing
        (
            ST_Buffer
            (
                intersections.intersection
                , 10
            )
        )
        AS extring,
        ST_Transform(intersections.geom1,3785) AS geom1,
        ST_Transform(intersections.geom2,3785) AS geom2,
        intersections.id1,
        intersections.id2
    FROM
        intersections
),
points AS
(
    SELECT
        ST_GeometryN
        (
            ST_Intersection
            (
                buffers.extring,
                buffers.geom1
            )
            , 1
        )
        AS point1,
        ST_GeometryN
        (
            ST_Intersection
            (
                buffers.extring,
                buffers.geom2
            )
            , 1
        )
        AS point2,
        buffers.intersection,
        buffers.extring,
        buffers.geom1,
        buffers.geom2,
        buffers.id1,
        buffers.id2
    FROM
        buffers
)
SELECT
    points.point1,
    points.point2,
    points.extring,
    points.geom1,
    points.geom2,
    abs
    (
        round
        (
            degrees
            (
                ST_Azimuth
                (
                    points.point2,
                    points.intersection
                )
                -
                ST_Azimuth
                (
                    points.point1,
                    points.intersection
                )
            )::decimal % 180.0
            ,2
        )
    )
    AS angle,
    points.id1 AS id1,
    points.id2 AS id2
FROM
    points
LIMIT 1000
;
ALTER TABLE angletmp ADD COLUMN id serial;

but I get this error :
ERREUR:  Argument must be POINT geometries

I can't understand the origin of the error and I'm wondering if you could help me solve it


texts: I have Shapely geometries where latitude and longitude where mapped to X and Y the wrong way around. Latitude is mapped as X, longitude as Y. How can I fix the coordinates so that the current X and Y switch places?


texts: How to merge vector buildings height layer with raster DEM to get solar maps in r.sun in Grass GIS?

texts:
Welcome, this is my first post.
I have vector layer with buildings heights and vector/raster DEM file.
I need to create solar maps with sun radiance values on the roof (and on surface, ofc). I am trying merge buildings with DEM, but I dont know how excatly to do it. I converted vector buildings to raster, made pixels outside buildings to null, and tried merge this file in 'r.mapcalc'. I think I should use 'if' function, but I dont know to how. Could anybody help me with it? It is possible to make flat surface in DEM under buidings in order to place buildings there?