texts: Python Syntax Error 000989 when running Python Snippet of successful geoprocess?

texts:
I am having an issue with a simple python snippet.

Here's the background: I wanted to include an additional step in a nightly update script that replaces characters in a field that simply do not play nice with a web application that allows users to query the data.  These characters included angle brackets ("<" and ">") which I was trying to replace with "Less Than" and "Greater Than".

So when I was testing the Field Calculation, I finally got it to work the way I wanted to, and as I have done with so many other snippets, I simply copied the Python Snippet of the geoprocess, and pasted into my script.

Below is the snippet:

arcpy.CalculateField_management("path/to/my/layer","field","RemoveBad( !field!)","PYTHON_9.3","def RemoveBad(x):/n   x = x.replace(\"<=\", \" Less Than or Equal To \" )/n   x = x.replace(\">=\", \" Greater Than or Equal To \")/n   x = x.replace(\"<\", \" Less Than \")/n   x = x.replace(\">\", \" Greater Than \")/n   x = x.replace(\"=\", \" Equal To \")/n   x = x.strip()/n   return x")

Now in my script, I do have a log that reports geoprocessing errors, so when I run the script, I am getting the following:

Executing: CalculateField path\to\my\layer field RemoveBad( !field!) PYTHON_9.3 "def RemoveBad(x):/n   x = x.replace('<=', ' Less Than or Equal To ' )/n   x = x.replace('>=', ' Greater Than or Equal To ')/n   x = x.replace('<', ' Less Than ')/n   x = x.replace('>', ' Greater Than ')/n   x = x.replace('=', ' Equal To ')/n   x = x.strip()/n   return x"
Start Time: Mon Aug 25 09:35:32 2014
Failed to execute. Parameters are not valid.
ERROR 000989: Python syntax error: Parsing error SyntaxError: invalid syntax (line 1)
Failed to execute (CalculateField).
Failed at Mon Aug 25 09:35:32 2014 (Elapsed Time: 0.00 seconds)

So I started combing through the script but I cannot find the issue, which leads me to believe it has something to do with the angle brackets.

To reiterate, this python snippet is from a successful Field Calculator geoprocess.

Per Alex's feedback I created the following variable:

RemoveBad = """def RemoveBad(x):
   x = x.replace(\"<=\", \" Less Than or Equal To \" )
   x = x.replace(\">=\", \" Greater Than or Equal To \")
   x = x.replace(\"<\", \" Less Than \")
   x = x.replace(\">\", \" Greater Than \")
   x = x.replace(\"=\", \" Equal To \")
   x = x.strip()
   return x"""

Then I updated the snippet:

 arcpy.CalculateField_management("path/to/my/layer","field","RemoveBad( !field!)","PYTHON_9.3",RemoveBad)

The error is now gone and script is working as expected.


texts: Merging polygon and point Shapefile in ArcGIS for Desktop?

texts:
I have two shapefiles that I want to merge.
One is of polygons, the other is a point file.

The polygons are missing information in the attribute table that the point file has. In general, the two line up pretty well, but others are not so exact.
What would be the most efficient way of merging these two shapefiles to have all of the information I need?

More information: I'm working on a historic mapping of my University's campus. The current edition has a point file for the buildings around campus, with years built, razed, building name, etc.

I want to instead use a polygon file and have building footprints. I contacted the university, and was given a shapefile that has the footprints.

However, none of the building have any information associated with them, including the name or building number.
My goal is to be able to name the building footprints in the most efficient way possible.

The problem is, there are no matching fields, and the only matching information I have is the actual location.


texts: How does GIS detect the error causing field in "ExecuteError: ERROR 001156"

texts:
I always get all types of errors, and would like to know if i can use the information in these error messages generated by arcpy to fix the issue and repeat the step where it occurred.

for example in batch merging:

   layerList[shp1,shp2,shp3,shp4]
   arcpy.Merge_management(layerList,"C:\\Data\\title.shp")

and I get the following error:

   ExecuteError: ERROR 001156: Failed on input OID 0, could not write value 'Kanata' to output field ErrorCuasingField
Failed to execute (Merge).

How can I get the ErrorCausingField in my script, so i can maybe delete it before repeating the merging process.

Also I have a another question addressing this issue specifically: (Batch merging field mapping).


texts: Creating selection record using ArcPy and Python script tool?

texts:
I want to create a tool that list records of a field in the selection part and move some of records to the right pane.

How can i create the UI using ArcGIS and a Python script tool in a standard toolbox (TBX)?

I can create the selection part but not the sort and move buttons.


texts: Postgis: WGS84 Lat Lon Point to Gauss-Kruger 3

texts:
How do I convert a WGS84 lat lon point in PostGIS to Guass Kruger 3 (31467)?

Somehow, I do not seem to get it right. I tried the usual way for converting projections:

st_transform(st_setsrid(st_makepoint(50, 8), 2346), 31467)

However, this gives an out of bounds error. I assume PostGIS expects the values to be Cartesian coordinates.


texts: MS SQL Server's STDistance output differs from Google Earth/Maps straight line distance

texts:
I am testing SQL Server STDistance function in the following way:

DECLARE @g geography;
DECLARE @h geography;

SET @g = geography::STGeomFromText('POINT(53.553813 9.99158)', 4326);
SET @h = geography::STGeomFromText('POINT(52.5234051 13.4113)', 4326);
SELECT @g.STDistance(@h);

These points represent, respectively, Hamburg and Berlin.

The STDistance call above returns 394612.166639773 (aprox. 394km). However, on Google Earth/Maps, when I calculate the straight line between the two points above, I get aprox. 254km.

If I understood correctly, both SQL Server SRID 4326 and Google Earth represent/use the WGS84 reference system, so I don't know what else could be causing this difference.

Any ideas?


texts: Create shapefile of bookmark extents and/or automatically print bookmark extents to PDF?

texts:
Is it possible to export the extents of all bookmarks in an active map document to a shapefile containing the polygons of the shape and/or to automate the print to PDF of all the bookmark extents?

The code for doing this with extents at How do I insert a DATAFRAME_ELEMENT value into a dynamically generated shapefile containing its extent?


texts: Why Can't ArcPy Search the Newly Created Feature

texts:
I'm running the following tool in ArcGIS 10. I have an ArcPy Script like this:

import arcpy
from arcpy import env
env.workspace=r"C:/sde/lcimrl.sde"
fc="SDE.ATTACHMENTS\SDE.attachment_points"
fcRows = arcpy.SearchCursor(fc, "OBJECTID = 11365")
for row in fcRows:
   print "I found one"

The feature with OBJECTID 11365 is just created.

It can find the feature if I run it as standalone script. But when I run it in the ArcMap's interactive python window or import it into a ArcTool toolbox and then run it, it cannot find the feature.

I have to force the refresh of the SDE connection in the ArcCatalog, and then I can get the point. Or I have to wait for a while, like several hours, then it can find the point.

I tried arcpy.ClearWorkspaceCache_management() and arcpy.refreshTOC(). They do not help.

Any help is appreciated.


texts:
I am building a Choropleth map in D3.JS.

In order to do that, I have downloaded a map from Eurostat.

From the GeoJSON file, I want to extract the country for the visualization, in my case Portugal.

As I am planning to use OGR2OGR, I have installed the GDAL package on Anaconda and in order to extract the map of Portugal, I am using the following code:

ogr2ogr -where FID='PT' pt.geojson CNTR_RG_01M_2016_3857.json

But it is triggering an error.

  org2org -where FID="PT" pt.geojson CNTR_RG_01M_2016_3857.json
   SyntaxError: invalid syntax

How do I use org2org on Anaconda?


texts:
After installation, I got the error messages below when I click the Spatialite (V.7) Plugins:

An error has occurred while executing Python code:

Traceback (most recent call last):
  File "C:\Users\uqipurna/.qgis2/python/plugins\QspatiaLite\qspatialite.py", line 68, in run
    dlg = QspatiaLiteApp(self.iface)
  File "C:\Users\uqipurna/.qgis2/python/plugins\QspatiaLite\QspatiaLiteApp.py", line 80, in __init__
    self.changeDb(max(0,self.connections['idx']))
  File "C:\Users\uqipurna/.qgis2/python/plugins\QspatiaLite\QspatiaLiteApp.py", line 146, in changeDb
    self.drawTreeView()
  File "C:\Users\uqipurna/.qgis2/python/plugins\QspatiaLite\QspatiaLiteApp.py", line 151, in drawTreeView
    model=self.db.createTreeModel()
  File "C:\Users\uqipurna/.qgis2/python/plugins\QspatiaLite\Classes\Database.py", line 107, in createTreeModel
    self.getTables()
  File "C:\Users\uqipurna/.qgis2/python/plugins\QspatiaLite\Classes\Database.py", line 81, in getTables
    header,data=self.executeQuery(query)
  File "C:\Users\uqipurna/.qgis2/python/plugins\QspatiaLite\Classes\Database.py", line 46, in executeQuery
    cursor.execute(query,params)
DatabaseError: file is encrypted or is not a database

What is going on with my installation?


texts: Modifying PostGIS database using geoserver

texts:
I want to build a mobile app using GeoServer and PostGIS database and I am little confused about how user will interact with database.

I want end users to be able to upload new geospatial data to the database (not from any map view but from a formula). In order to be able to upload data, user shall first sign in (or log in).
Each data will be associated to the user id.

Do I need to create a web service to handle accounts creation and database modification (wrapping SQL query) or can I directly use GeoServer?


texts: How did QGIS fix a shapefile's character encoding, and can I do this on the command line with ogr2ogr?

texts:
Solved: I'm not sure why, but it works fine if I specify --config SHAPE_ENCODING "ISO-8859-1" in my ogr2ogr command. For whatever reason, this allowed ogr2ogr to go from the original TigerLine shapefile, as-is, directly into a MySQL database. My final ogr2ogr instruction, in full:

ogr2ogr --config SHAPE_ENCODING "ISO-8859-1" -f "MySql" MySql:"usa_basemap,host=127.0.0.1,user=myUser,password=myPass,port=3306" -lco engine=MYISAM "C:/path/to/data/tl_2012_us_county.shp" -nln county -nlt "geometry" -s_srs EPSG:4269 -t_srs EPSG:3857

Background: I was finally able to use ogr2ogr to import the TigerLine 2012 USA Counties shapefile into MySQL, but only after using QGIS to save the shapefile with a new encoding value (UTF-8).

I'm wondering: What can I do differently, right inside the commandline with ogr2ogr, to avoid having to involved QGIS in the workflow? (Nothing against QGIS, it's just that I want to be able to automate this..)

Here's what I did before involving QGIS:

1) Changed the MySQL max_allowed_packet size as is warned might be necessary on the OGR MySQL driver page.

In MySQL Workbench, I set the new parameter value like this:

SET GLOBAL max_allowed_packet=1073741824;

2) Next, I tried to import it into MySQL with ogr2ogr like this:

ogr2ogr -f "MySql" MySql:"basemap,host=127.0.0.1,user=myUser,password=myPass,port=3306" -lco engine=MYISAM "C:/path/to/data/tl_2012_us_county.shp" -nln county -nlt "geometry" -s_srs EPSG:4269 -t_srs EPSG:3857

But that bombs with an error that spills outside of the buffer, it ends as shown below, note the appearance of the special chars, which are lost.

[....hugely....truncated....console....output....]
322122,-7366573.8710948359 2062754.2910927597))',1) , '72', '045', '01804502', '
72045', 'Comer�o', 'Comer�o Municipio', '13', 'H1', 'G4020', '490', '41980', 'A'
,       73557129,         319735, '+18.2250402', '-066.2202984')

It's worth noting that after hitting this brick wall and Googling around, I tried the tips suggested by Mateusz Łoskot in this discussion; to restate the salient portion:

  *Windows Command Prompt can work with UTF-8 characters if you change codepage to UTF-8:

  0) Open new prompt (cmd.exe)

  1) Change font to Lucida Concole

  3) chcp 65001*

But this didn't overcome the obstacle.

Finally, wondering if QGIS could help, I used QGIS to export the shapefile, specifying UTF-8 as the character encoding, then the ogr2ogr command I provided above worked. Bingo. Now I have "Comerío Municipio" in my database..

..but what gives?

So my question: Is QGIS doing something I can add to my ogr2ogr instruction that automatically resolves the special characters issues I encountered so I can automate this task in the future and not have to manually involve QGIS?


texts: Problem with "Select by Expression" in QGIS

texts:
I am using QGIS 3.10 for Mac and facing some problems with the "Select by Expression" tool. The function I'm using is not selecting any features, even though it should. I have tested it with a few other functions, and the tool works just fine––this particular function seems to be hairy, and I'm not sure why.

I need to select features in the layer that have a certain population density. As such, I am using the following function in "Search by expression":

"TotPop_" / "sq_mi" >= 1000

I'm not sure what's wrong with the function, because it's not any selecting features at all, even though I manually went through and found a couple that apply! (The table is, unfortunately, too long for me to sift through them all manually.)

Does anyone have ideas about what I might be doing wrong?


texts: Installing pip packages from bat file and OSGEO Shell?

texts:
After installing QGIS (that from the 2.14 version has pip installed by default) I want to automatically install some pip packages within a bat file.

One method to do that is writing a file.txt with all the needed packages, open the OSGEO Shell as admin, cd in the directory of the file.txt and launch python -m pip -r file.txt

How do I make a bat file that automatically runs the Shell ad admin without asking for the password (that is the same behavior of right click and run as administrator)?


texts: Why doesn't the Kernel Density Estimate function appear in AniMove plugin?

texts: can anyone tell me why the Kernel Density Estimate function does not appear under my AniMove plugin in QGIS? I can only perform the MCP. I run Windows OS.


texts: PyQGIS raster interpolation from vector layer

texts:
I want to generate raster from line vector layer using PyQGIS. I was following this post:
How to call the interpolation plugin from the python console?

but I have a problem with:

layer_data.InterpolationAttribute =0

It works, but when I try to set different index it still takes values from first column. It looks like 0 index is default and can't be modified. I am working on Ubuntu 16.04.

Maybe this is a reason?


texts: Intersection problem: too many elements?

texts:
I have this code where an IGeometry shape (4-5 circles not filled) should be intersected with another layer's features

                // Convert from IGeometry to ADF Geometry Polygon
                ESRI.ArcGIS.Geometry.IPointCollection com_pointcollection =
                        (ESRI.ArcGIS.Geometry.IPointCollection)multipointPoints2;
                ESRI.ArcGIS.ADF.Web.Geometry.Point[] new_adf_points =
                                        ESRI.ArcGIS.ADF.Web.DataSources.ArcGISServer.Local.Converter.FromIPointCollection
                    (com_pointcollection);
                ESRI.ArcGIS.ADF.Web.Geometry.PointCollection new_adf_pointcollection = new
                    ESRI.ArcGIS.ADF.Web.Geometry.PointCollection();
                for (int i = 0; i < new_adf_points.Length - 1; i++)
                {
                    new_adf_pointcollection.Add(new_adf_points[i]);
                }

                ESRI.ArcGIS.ADF.Web.Geometry.Ring new_adf_path = new
                    ESRI.ArcGIS.ADF.Web.Geometry.Ring();
                new_adf_path.Points = new_adf_pointcollection;
                ESRI.ArcGIS.ADF.Web.Geometry.RingCollection new_adf_pathcollection = new
                    ESRI.ArcGIS.ADF.Web.Geometry.RingCollection();
                new_adf_pathcollection.Add(new_adf_path);
                ESRI.ArcGIS.ADF.Web.Geometry.Polygon new_adf_polygon = new
                    ESRI.ArcGIS.ADF.Web.Geometry.Polygon();
                new_adf_polygon.Rings = new_adf_pathcollection;

                // Now get the intersect query

                ESRI.ArcGIS.ADF.Web.SpatialFilter spatialFilter =
                    new ESRI.ArcGIS.ADF.Web.SpatialFilter();
                spatialFilter.Geometry = new_adf_polygon;
                spatialFilter.ReturnADFGeometries = false;
                spatialFilter.MaxRecords = 10000;

                System.Data.DataSet queryResultsDataSet = new System.Data.DataSet();
                for (int index = 0; index < layerIDs.Length; index++)
                {
                    System.Data.DataTable queryResultsDataTable =
                  queryFunctionality.Query(queryFunctionalityq.Name, ConvertFromLayerNameToLayerId(m_Map, m_ResourceMapItem, m_firstLayerName), spatialFilter);
                    queryResultsDataTable.TableName = gisResource.Name + "_" + layerNames[index];
                    queryResultsDataSet.Tables.Add(queryResultsDataTable);
                }

I'm converting the IGeometry (which I've drawn on the screen and I know it's correct) to a web ADF Polygon and then using this polygon for an intersection query

However some features which are clearly not intersecting the circles' borders are being included in the result

To better explain the situation I did this screenshot:

Why is that? The red squared features are all different so it's not a problem of "same feature"

Update: with less features it's the same, if they intersect the circle line they get colored, but the others get colored too.. simply adding them seems like telling the system "I'm intersected too"


texts: Adding date to label in QGIS Atlas composer?

texts:
I'm trying to add the current date into a label for use in a composer for colleagues to run on a daily basis. The label is currently;

Bins to be delivered: [% "Count" %]

All I need to add is the current date, but when I make the label this;

Bins to be delivered: [% "Count" %] as at $now()

nothing happens.

Do I need to enter something in the 2 end brackets to bring the date in? Or is it something else I need to do?


texts:
I have to make a symbology based on a unique color value so I add a plus column to the attribute table with the number of color values. I went to the layer properties - style - data defined properties and trying to change the outline color, but it is so complicated. I have 5 different value numbers so if the value number is 2, it can be red, and so on... What codes should I use?

The map should look like this with different sizes and colors based on the unique value


texts: ArcPy Raster to Polyline outputting empty shapefiles and a 999999 error?

texts:
I am running a lengthy script to post-process raster channel paths into polylines so I can measure their length using Add Geometry. Script runs well until I go to convert the rasters to lines using the Raster to Polyline tool in arcpy, shapefiles are output when the error is overwritten but attribute tables are empty and they can't be drawn in arcmap. Input Raster is integer type, data has coordinate system, channel paths are 1's and all other cells are nodata. Shapefiles are missing the optional .sbn and .sbx supporting files.

I have tried running the Raster To Polyline tool from within ArcMap and it works fine.

Workarounds accepted. Script is long so a snippet where failure occurs is below.

arcpy.env.workspace = CPC_rename  #sets workspace
CPC_lines = CPC_rename + '\\' + folder7  # sets output workspace
rasters = arcpy.ListRasters('*', 'TIF')  # Lists Rasters
for x in range(0, len(rasters)):  # set up a for loop
    linesOut = os.path.join(CPC_lines) + '\\' + rasters[x] + '.shp' # setup the output file
    arcpy.SetSeverityLevel(0)  # override the 999999 error so it runs all the way through
    arcpy.RasterToPolyline_conversion(rasters[x], linesOut , "ZERO", 0, "NO_SIMPLIFY", "Value")  # apply the arcpy tool
    print 'something' # print something to confirm its running


texts: What is the word for "Bar" progress and write it in python script?

texts: I am a beginner in python and have written a few scripts and I am wondering if there is a term for python to call and show the progress of the bar telling you from start to finish like you see in geoprocessing at the lower right corner in ArcMap ?


texts: Exporting image from app-side using Google Earth Engine

texts:
I've built an app that displays an image. I want the user of the app to be able to download the data from the viewport in a GeoTIFF format for further use in GIS.

I've got my export function running on the developer side and can download the data via the Tasks pane as per instructions. But is there a way to open this download up for the external user? I'm starting to fear not....

Here is the relevant code that creates/adds a button to execute the export function.

// Export data function
var exportData = function(range) {
  print("Exporting data...")
  Export.image.toDrive({
    image: ndvi,
    description: 'imageToDriveExample',
    folder: "GEE_Exports",
    fileNamePrefix: "NDVI",
    dimensions: "10000x10000",
    scale: 30,
    shardSize: 100,
    fileDimensions: 1000,
    fileFormat: "GeoTIFF"
    //maxPixels: 1e12
  });
}
//...add update map button to panel
var exportDataButton = new ui.Button({
  label: 'Export Data',
  //style: {stretch: 'horizontal'}
});
// Register the function to the button click event.
exportDataButton.onClick(exportData);
// Add to button panel
buttPanel.add(exportDataButton);

If this is possible, will I need a cloud account?


texts: Display Feedback "Stutters" When Moving Features

texts:
I have a custom tool for moving features that implements IDisplayFeedback and geometry-specific CoClasses to show the user where they are moving the feature.

    IMovePolygonFeedback polygonFeedback;
    bool feedbackStarted;

    public override void OnMouseDown(int Button, int Shift, int X, int Y)
    {
        feedbackStarted = false;

        //A BUNCH OF STUFF HERE TO SELECT A FEATURE ON MOUSEDOWN

        geometry = feature.Shape;

        polygonFeedback = new MovePolygonFeedbackClass();
        polygonFeedback.Display = doc.ActiveView.ScreenDisplay;
        polygonFeedback.Start(geometry as IPolygon, doc.CurrentLocation);
        feedbackStarted = true;
    }

    public override void OnMouseMove(int Button, int Shift, int X, int Y)
    {
        IPoint point = doc.ActiveView.ScreenDisplay.DisplayTransformation.ToMapPoint(X, Y);
        if (polygonFeedback != null && feedbackStarted)
            polygonFeedback.MoveTo(point);
    }

What sometimes happens when the mousedown and mousemove happen quickly is that the feedback display stutters for a second, leaving a (what I assume to be) graphic, then continues on with a new graphic moving around with the mouse until mouseup.  Looks kinda like this:

http://i.imgur.com/99fj8.png

At mouseup I change the feature shape to the new position and store the feature.  I union the original feature extent with the new feature extent and pass that envelope to a partial refresh.  Problem is sometimes parts of the "stutter" aren't in that envelope so I get dangly bits:

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

I'm trying to avoid straight refreshes but these dangly bits are not pretty and don't go away until a refresh happens (forced by user or move/pan, etc.)

Questions:

Does anybody know why the "stutter" is happening and how can I prevent it?
Is there a way to partially refresh the dangly bits?
What is the feedback?  Is it a graphic?  Can I trap it during the onmove to add it to the partial refresh envelope?
Is there an infinitely better way to do this that I'm mising?


texts:
I did a dGPS survey last summer. Unfortunately, we didn't have a base station. So I did my survey in post-processed Kinematic PPK mode. Because we borrowed the equipment, we don't have access to post-processing software. So here I am with Trimble files (filename.t01) from our survey that need to be post-processed to get accurate positions. Fortunately there are online services that can do the post-processing for you. After converting the files to RINEX format (filename.o) I was able to successfully have them processed by an online service at Natural resources Canada.

However, the output is providing me with a track and not with the individual points from the survey. These were marked during the survey and are marked in the RINEX file with

          *** Start of Occupation ***

Is there any way to efficiently get out the GPS positions for these individual markers?


texts:
Occasional requirement to visualise arbitrary vector data in QGIS, adding features from Python.
[understanding that they cannot be overlayed with any data on known projections]
How can I specify Unknown / non-earth CRS for the layer?
crs = QgsCoordinateReferenceSystem(??)


texts: Match centroids from one dataset with polygons of another

texts:
I want to make a data frame showing which voting precincts are in each school zone. To do this, I want to see which precincts' center points are in which school zones.
library(sf)

precincts <- st_read("voting_precincts.shp")
precincts$centroid <- st_centroid(
  x = precincts
)

elemZones <- read_sf("school_zones.shp")

So although I have the centroid coordinates for precincts, I do not know how to "join" them with the school zones data.


texts: Conversion of workspace to ordinary .tab

texts:
Does anyone know a way for converting a wokspace (*wor., *.mvs) to an ordinary MapInfo table (.tab)?

I have worked in workspace with a high number of GPS-field registrations. I have made changes (edit shape etc) to the registered entities in the workspace-mode. However, I discovered that these changes are stored in the workspace-mode only - not in the original tabs. Layer control shows a list of the original tabs but when I open the original tabs in a new window, the changes/edit does not appear.

Therefore I would be very happy of there is a way to convert the whole workspace to a new .tab that may be opened in the MapInfo .tab-mode.


texts: Using Shapefile with Annotations in QGIS?

texts:
In my project I downloaded some map polygons from the Open Government portal of Canada.

One of the files contains annotations for provinces, territories, water bodies, and more.

When I import it into my QGIS project I only see polygons filled with a solid color (see map figure). But I believe it would be possible to fill the polygon with text, that is also provided in the attribute table.

I got as far as displaying the text at specific X and Y coordinates, but it did not fit into the given polygons.

Maybe I could make use of the other files that were provided for my "base_anno.*" file?

Overall I think it should be possible to bring the text into the boxes provided by the shapefile. I just don't know where and how I could do that


texts: How to get image pop ups in QGIS using oracle blob field with image?

texts:
I have a point feature in oracle database which has as attribute a blob field that have image inside for every point feature taken in field.

Is there anyway that i can decode and display it on a popup display in Qgis?

I can do hyperlink or popup display if the image files are outside database as file and their path as attribute but i want to know if they can be read directly from oracle blob field.


texts: ArcGIS Server Featurelayer: Error on adding features via JavaScript API after switch from Oracle db to PostGIS

texts:
We have a webapp which consumes services published on ArcGIS Server. All database relevant operations are done with ESRI ArcGIS abilities.
We now switched from Oracle 11g to PostgreSQL 9.4 with PostGIS 2.5. Although everything is done with Esri ArcGIS abilities it is necessary to make changes on the code (this may be bad design of the app or maybe it's unavoidable... doesn't matter at the moment). I made these changes as far as I can see and almost everything is working fine, except adding/editing features via the exposed feature layers.

When trying to add a feature ArcGIS Server will report an error like this in the log:

  Error: Underlying DBMS-Error[Error: Geometry has Z dimension but column does not ({featureclass})::SQLSTATE=22023].

The table is defined with Point and another table with Polygon. So it is correct that the column does not support Z dimension. But I'm using the ArcGIS JavaScript API 3 which doesn't support Z coordinates at all according to the point documentation. The featurelayers also don't support Z and M values.

When looking with Chrome developer tools what data are sent to the server I get something like this:

[{
    "geometry": {
        "x": 6359.682509019719,
        "y": 7720831.9973746035,
        "spatialReference": {
            "wkid": 3857
        }
    },
    "attributes": {
       <mylist of attributes>
    }
}]

The geometry looks fine and has nothing with Z. I checked the attributes multiple times that they have the correct spelling, having valid values and are all lower case (what is default for identifiers in PostgreSQL).

When using this example (of course with the attributes) on the FeatureServer REST website and adding the feature via "addFeatures" I get following JSON feedback:

{
 "error": {
 "code": 400,
 "message": "Unable to complete operation.",
 "details": [
   "Rowbuffer creation failed."
  ]
 }
}

This is some more information but I still don't know what does cause this problem. In the server logs I will again have the error message above.

What could cause such an error message when using featurelayers in combination with PostGIS?

Edit:

The problem only occurs when adding new features to the layer. If I'm deleting or editing a feature everything works fine and the changes are passed to the database without errors.


texts: Add Layers with Leaflet from Different Spatial Data Frames in R

texts:
Very new to R's Leaflet package and have two spatial data frames.

1) An outline of a city
2) Its major highways

I wish to show both of these on the same map. I get how to do this with ggmap or tmap. I'm little confused when it comes to Leaflet.

Here is some code that takes public data:

toronto - an outline of the city
all_roads - Toronto's road network

I then view each separately with Leaflet. How do I view multiple layers from different spatial data frames on the same map?

library(sf)
library(leaflet)

# Read in Toronto boundaries as a spatial data frame
# WGS84 latitude and longitude
# Data: https://www.toronto.ca/city-government/data-research-maps/open-data/open-data-catalogue/#4a30ca3e-70f0-bf79-1e76-459ebda8458c
toronto <- st_read("/Users.../torontoBoundary_wgs84/citygcs_regional_mun_wgs84.shp")

# Read in street data
# (WGS84 projection)
# https://portal0.cf.opendata.inter.sandbox-toronto.ca/dataset/toronto-centreline-tcl-2/
all_roads <- st_read("/Users.../centreline_wgs84/CENTRELINE_WGS84.shp")

# Just highways
xway <- filter(all_roads, FCODE_DESC == "Expressway")

map <- leaflet(data = hood_shp, options = leafletOptions(minZoom = 11, maxZoom = 13))
map %>% addPolygons(color = "#aaaaaa", weight = 1.0, smoothFactor = 0.5,
                    opacity = 1.0, fillOpacity = 0.7)

highway <- leaflet(data = xway)
highway %>% addPolygons(color = "brown", weight = 1.0, smoothFactor = 0.5,
                        opacity = 1.0, fillOpacity = 0.7)


texts:
I am very new to GIS stuff, I am looking at a dataset displayed in WMS, I have attached the image of what it looks like. I want to change the feature colour of Mine and ERW to Blue and Pink respectively. I the link for the WMS is https://data.opendevelopmentcambodia.net/dataset/mineerw-casualties-2005--2013/resource/fda1f89d-6ba6-465f-b832-94c2c99c18f7
How would I go about changing the feature colours? There are a lot of different file types to download so I am very confused.


texts: In gdal how to save 3 dimensional array as stacked raster image?

texts:
I have constructed a 3-dimensional numpy array by stacking multiple 2-dimentional rasters.

>>> print arr_list_stack.shape
(30L, 370L, 400L)

As the arr_list_stack.shape shows that this 3-dimensional numpy array has 30 bands, by using gdal I want to save this array as raster image which has all the 30 bands.

I have a method which can save 2-dimentional arrays

def CreateGeoTiff(Name, Array, driver, NDV,
                  xsize, ysize, GeoT, Projection, DataType):
    Array[np.isnan(Array)] = NDV
    DataSet = driver.Create(NewFileName, xsize, ysize, 1, DataType)
    DataSet.SetGeoTransform(GeoT)
    DataSet.SetProjection( Projection.ExportToWkt() )
    DataSet.GetRasterBand(1).WriteArray( Array )
    DataSet.GetRasterBand(1).SetNoDataValue(NDV)
    return NewFileName

but CreateGeoTiff() fails to save 3-dimentional arrays and returns the error

  File "C:\Users\Documents\Quadium\SPI\Code\testing_stack_arr.py", line 39, in CreateGeoTiff
    DataSet.GetRasterBand(1).WriteArray( Array )
  File "C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\osgeo\gdal.py", line 1235, in WriteArray
    callback_data = callback_data )
  File "C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\osgeo\gdal_array.py", line 359, in BandWriteArray
    raise ValueError("expected array of dim 2")

What changes should I make in CreateGeoTiff() method that will enable me to save 3-dimentional array.


texts:
I am making a map using layers stored in a PostgreSQL (PostGIS) database. I am using GeoServer to successfully load layers with Leaflet in WFS format, but I wonder how to make SQL request on it through my JavaScript code.

How can I do that ?
Or more generally, how can I can simply connect to a PG database within a JS script ?


texts: Relocating georeferenced image, shapefiles and grid using Spatial Adjustment tool

texts:
I somehow georeferenced an image about 120 miles east of where I should have georeferenced it. This should be a simple fix, except that I then proceeded to digitize that extremely detailed image (read: time consuming) into two shapefiles and also gridded that data all based on the incorrect location.
A friend recommended to just re-georeference the original image and use the Spatial Adjustment tool to simply move the shapefiles and grid over to the right location.
The problem is: I've never used the Spatial Adjustment tool and I'm having trouble figuring out how to do this properly. There are three options: edge snap, rubbersheet and transform.
Can anyone help or suggest an easier way to achieve this without having to re-digitize everything?


texts: Set edit widget to read-only via PyQGIS for QGIS 2.8

texts:
I am trying to set an editor widget to read-only using PyQGIS. This would be the equivalent of setting the widget to Text edit and unchecking editable via the GUI.

I am using QGIS 2.8.3 and am have set the widget in this way:

myLayer.setEditorWidgetV2(fieldNameIDX, 'TextEdit')

I notice in the API for QGIS 2.15 that there is a new QgsEditFormConfig class with a readOnly method, but I believe I can't use it because of my QGIS version and I cannot immediately update.


texts:
Problem with adding icon to layer control.

CSS:

.fa-b:before {
        font-family: Courier;
        font-weight: bold;
        content: 'B';
    }

JS:

var geojsonLayer = new L.GeoJSON.AJAX("https", {
        pointToLayer: function(feature,yx){
          var marker = new L.marker(yx, {
            icon: L.AwesomeMarkers.icon({
              icon: 'fa-b',
              markerColor: 'cadetblue',
              prefix: 'fa'})});
          return marker;}
      });

var clusterGroup = new L.markerClusterGroup();
      geojsonLayer.on('data:loaded', function() {
      clusterGroup.addLayer(geojsonLayer);
      map.addLayer(clusterGroup);
      })

var groupedOverlays = {
        "Group": {
          '<icon class = fa-b > B': clusterGroup
        }}


texts: Forcing Python Script Tool Parameter Input Field to Begin at Specific Folder?

texts:
I'm setting up a Script Tool (Arc 10.2.2) where the first parameter will be a CSV file chosen by the user, always from the same folder.

When the user clicks the first parameter's windows navigator drop down (to select the desired .csv) I want the ensuing catalog window to begin at the folder I specify.

I've already tried setting the folder path in the "Default" field of the Parameter Properties Window but that didn't work.  This made the path appear in the Parameter field, but when the drop down is clicked that path is ignored and the user is shown an entirely different location.

I'm quite sure the proper validation code will achieve my desired functionality, I just don't know what that is.  I have read many ESRI web pages on the topic but have yet to find the proper coding.

How do you force a Tool Parameter Input Field to Begin at a Specific Folder?

I seek the proper function and which module it should be in (i.e. which of the standard 3 validation modules: _init_(self) , initializeParameters(self) , or updateParameters(self).).


texts:
I'm a recent Computer Science graduate who prides himself on his knowledge of software development and problem solving, but has no idea how GIS systems work. (Read: I'm a programmer.)

I've been charged with developing (ie. coding) a new webmap from scratch.

I can program in any language that has bindings to COM (Microsoft Component Object Model), so I've been writing in PHP.


texts: Using GDAL tool rasterize with QGIS so that points represent center of grid

texts:
I want to convert point data into a grid and used the gdal tool rasterize. The tool works but unfortunately the result looks like this.

The converted point does not represent the center of the grid but the topleft corner. I need it to be in the center. I couldn't find any settings in the tool nor Explanation in the GDAL-Documentation to accomplish that.

Does anyone have suggestions?

I rasterrized the points into raster and NOT raster to point.

If I upload the data directly as raster I get them correct but that doesn't always work. In the picture the underlying raster image is the correct one.


texts: QGIS dissolve is changing the shape of polygons rather than just combining them

texts:
QGIS novice here.  I have a shapefile for US census PUMA (public use microdata areas).  One of the attributes in my shapefile is the metropolitan area of the PUMA.  I want to combine PUMA's which are located in the same metropolitan area.  To do this, I am using the dissolve command, and choosing the metropolitan area attribute as the dissolve field.  For most of the United States, things work as expected, but there are a few areas in which there is weird behaviour.  For example, instead of combining the Chicago pumas into a big polygon, QGIS redraws a Chicago polygon of slightly different size.  See the following pictures to see what I mean:

Is this a QGIS bug, or is there something I am doing wrong?  Is there anything I might try to fix this behaviour?


texts: Export WebMap using PrintTask Failing Silverlight API

texts:
I am using a custom tool to export a webmap to PDF. The tool works successfully through ArcGIS 10.1. The problem crops up when consuming the tool through my Silverlight application using PrintTask. The tool works in an intermittent fashion i.e. sometimes it works, and sometimes it gives an error.

I turned to ArcGIS Manager log and extracted the following error:

Traceback (most recent call last): File     "C:\arcgisserver\directories\arcgissystem\arcgisinput\Toolbox\ExportPDSitePlanPart.GPServer\extracted\v101\printingtoolbox\ExportPDSitePlan.py", line 48, in result = arcpy.mapping.ConvertWebMapToMapDocument(Web_Map_as_JSON, templateMxd) File "c:\program files\arcgis\server\arcpy\arcpy\mapping.py", line 234, in ConvertWebMapToMapDocument return webmap_tuple(*convertArcObjectToPythonObject(arcgisscripting._convertWebMapToMapDocument(webmap_json, template_mxd, notes_gdb, extra_conversion_options))) RuntimeError: Layer "PD Working Layer": Failed to create layer from service at http://10.58.161.70/arcgis/rest/serv...eatureServer/0. Layer "PD Working Layer": Failed to create layer from service at http://10.58.161.70/arcgis/rest/serv...eatureServer/0. Failed to execute (ExportSP). Failed to execute (ExportSP).

I can confirm that I can access the service from the same server where the tool is running. (i.e. 10.58.161.70).

Can someone point me in the right direction?

Thanks


texts: Join County Data by state and county name without FIPS - join by two attributes?

texts: I have data by county that I want to join with census data but when I do the join by county name it messes up because some states have the same county names and I don't know how to get around this without manually adding FIPS to all the counties. Is there a way to do a join by saying to make sure the State and County match, not just the County name?


texts: Creating geopackage layers as output with QGIS batch processing autofill option

texts:
Is it possible to use the Calculate by Expression tool to output geopackage layers. For example, with the batch option for the polygon to line function.  With the following expressions I can create individual shapefiles or individual geopackage files but not layers in a geopackage file;
This creates a separate geopackage file for each layer;@project_folder ||'/' || @INPUT  || '_line'

This creates a separate shapefile for each layer;
@project_folder ||'/' || @INPUT  || '_line.shp'
However, I can't find a way to output layers to the one Geopackage file. The Polygon to Line function seems to default to individual geopackage files if no extension is given.  The following creates output - not individual layers but separate geopackage files;
@project_folder ||'/test.gpkg layername=' || @INPUT ||'_line
The file names end up like this;


texts:
I need to display polygons as different colors based on a value in a separate data table that I've imported from Excel. The polygon attribute table and data table both have a field with the polygon ID#, which I'm trying to use to join the 2 tables.

There are multiple entries for each individual polygon in the data table (each month over 14 years has a value). I've tried doing a one-to-many join based on the polygon ID (called 'toparea' below), but only the first set of values is appended to the shapefile (ie, only Jan 2000 will join to the polygon attribute table). I need to display the 'new_cpue' values for Jan through Dec for 2000-2014. Here is what I'm getting when I try to join:

I created the shapefile with the polygons manually from a set of coordinates, and it has all the underlying spatial information. There is no spatial information in the second table I imported except for the ID# of the polygon.

Relating the 2 tables works fine. When I relate I can select the polygon and get the information contained in the related data table, but I can't display the information for the polygon from the related data table, which is ultimately what I need to do.


texts: Abbreviate Road Labels in QGIS only when required

texts: I see that it is possible to create an abbreviations list, by using 'Apply label text substitutes' in the Labels Properties, but this abbreviates ALL labels. Is there a way to use the abbreviations only when required, ie. to fit labels on shorter roads and not abbreviate when there is enough room to show the full label?


texts: GDAL's -a_ullr switch - cell centers or outer edges?

texts: I am using gdal_translate to extract a NetCDF out to tiff. In the gdal_translate documentation, there is an option -a_ullr to define the upper-left and lower-right corners (bounds) of the grid. Are these the cell centers of those grid cells, or the outer edges of the grid itself?


texts:
I setup an sde connection for Oracle database, where multiple databases are sitting. All databases have same owner, sharing same password but identical user (e.g., BUI, EICU, FRNSW and LPI). While configuring the sde connection it is showing all databases (featuredatasets) in one big list, rather than individual connections.
My previous look was

and now it is showing as

Hope I am missing some properties. How I can parse individual databases ?


texts: Is it possible to add a stationary pop-up window to a QGIS point?

texts: I'm a QGIS newbie, and, am just learning where to find things in QGIS.  I'm using version 3.8.  When a user clicks on a "point" in QGIS I want a window to pop-up that can contain text, a photo, and, an HTML pointer to a website.  Essentially, I want to create the equivalent of a "PlaceMark" in Google Earth and in Google Maps. I've learned how to use the "map tip" feature in QGIS already (i.e. when I hover over the icon, text and data appear. 'Map tip' is nice but I want user's to be able to click on the QGIS point and have a window open that contains static info and ability to have user re-direct to a web site.  Please let me know where I should begin my education as to how to do this, or, let me know if it is even possible.


texts: Getting points within a range from a postgres database?

texts:
ok so this is query i'm making which returns no result:

SELECT
    building.gid
FROM
    building
WHERE
    ST_Dwithin(building.geom, ST_GeomFromText('POINT(long lat)',4326), 400)

If I run this in QGIS, I get a syntax error, but if I run this on pgAdmin, i get no error, but simply 0 rows....

my long and lat are in degrees.

what could possibly go wrong?


texts: Selecting the largest from sets of duplicates in QGIS

texts: I have a poly line layer where each feature has a Geometric Length and numerical ID. Some of these IDs are duplicates. I would like to be able to select the longest line from each set of duplicates to create a new layer with just the unique IDs.


texts: OL3 How to reproject vector layer on fly (ol 3.14)

texts:
I have a map in OL 3.14 with a WMTS in a specific epsg.
I add a WMS layer in an other projection and it's ok, OL reproject right.
But with a vector Layer (geoJson, in the same projection as WMS), Openlayers don't display my feauture.

For my local projection, I use proj4.

My code :

var projectionL = ol.proj.get('EPSG:2154');
var projectionCC = ol.proj.get('EPSG:3946');

var ortho2013Extend= [1372000.0,5207000.0,1403020.0000000002 , 5238000.0];
//11 niveau de tuilage
var matrixIds = new Array(11);
//donc 11 niveaux resolution à gérer
var resolutions = new Array(11);

var size = ol.extent.getWidth(ortho2013Extend) / 302;
for (var z = 0; z < 11; ++z) {
  // generate resolutions and matrixIds arrays for this WMTS
  resolutions[z] = size / Math.pow(2, z);
  matrixIds[z] = z;
}
var ortho2013=new ol.layer.Tile({
  opacity: 0.7,
  extent:ortho2013Extend,
  source: new ol.source.WMTS({
    url: 'https://sigar.agglo-larochelle.fr/CeramikServer/rest/wmts/cdalr/ortho2013?',
    layer: 'ortho2013',
    matrixSet: projectionCC,
    format: 'image/jpeg',
    transparente:true,
    projection: projectionCC,
    tileGrid: new ol.tilegrid.WMTS({
      origin: ol.extent.getTopLeft(ortho2013Extend),
      resolutions: resolutions,
      matrixIds: matrixIds,
      tileSize:302
    }),
    style: 'default'
  })
})
var feuxSonoresLayer = new ol.layer.Vector({
   source: new ol.source.Vector({
        url: 'feux_sonores_cc46.geojson',
        format: new ol.format.GeoJSON(),
        projection:projectionL
    })
});
map = new ol.Map({
    layers:[ortho2013, pluLayer, wmsHotSpotWifi, feuxSonoresLayer],
    target: 'map',
    //aucun controle
    //controls:[],
    view: new ol.View({
        projection: projectionCC,
        center: [1377000,5227000],
        extent: ortho2013Extend,
        zoom:14
    })
});

I have found some examples but in old version (3.8) and it used ol.source.GeoJSON, but this source doesn't exist now.

How to do ?


texts: Setting ModelBuilder to ignore error and continue?

texts:
I am running a couple of models to firstly, reproject a directory of raster images and secondly to clip them to a polygon boundary.  When running the clip the model exits if a raster falls just outside the clip polygon.

Is there a way that I can change a setting to have ModelBuilder note the error and move on without quitting?

I do not understand python so please, no suggestions about how it could solve the problem.  If there is no setting I can change I will try other approaches.


texts: It seems that there isn't a plugin to perform the right-hand-rule on a layer or selected features. Can this be done by using some Python code? Is there a function that works in the same way as ST_ForceRHR?


texts:
I'm currently designing a map for methamphetamine labs in my state.  Ideally what I need is the center lat/lon coordinates for each county in my state or the lat/lon coordinates for each county seat.

Does anybody know where I could find such a dataset?

I'm looking for coordinates for counties in the state of Indiana.


texts: Calculating values from feature service without adding field to service

texts:
I am accessing an ArcGIS Online feature service where each day the values in the table are accumulating rather than displaying that day's data. This is not my feature service so I do not have editing access.
For example, in the example below the value for Day 2 is 5 and the value for Day 3 is 15.
Day 1 10
Day 2 15
Day 3 30
Is there a way to calculate these values in a ArcGIS Online web map or application without writing into the service?


texts:
I suppose you could interpolate from the Census data, but is there a data feed that has accurate DAILY total population of the U.S.?

The finer the granularity, the better (so if there was daily populations by state or county that would be great).


texts:
I have GIS data on individual events (with timestamps), and I am running a panel data regression, with the outcome variable being the number of events by district by month. My independent variables also include individual coordinates in each district.

I wonder if the right way to go forward is to merge the different layers in QGIS (district shapefiles, event data etc) and export the resulting attribute table as CSV so that I can run the regression in R. If yes, I am struggling a bit with what exactly to do. Is there a way in which I can export the data directly as a panel, which includes the number of events for each district each month? Or am I on completely the wrong track?


texts: Identifying "long and narrow" polygons in with PostGIS

texts:
I have a set of polygons representing large areas, say city neighborhoods. I want to identify the large overlapping areas between them.
But there's a problem: sometimes these polygons will overlap along their perimeters (because they were drawn with little precision). This will generate long and narrow overlaps that I do not care about.
But other times there will be big overlaps of robust polygons, meaning large areas where a neighborhood's polygon overlaps another. I want to select only these.
See the picture below of just the overlaps. Imagine I wanted to select only the blue polygon in the lower right corner.

I could look at areas, but sometimes the narrow ones are so long they end up having areas as large as the blue polygon. I've tried to do a ratio of area / perimeter, but that has also yielded mixed results.
I've even tried using ST_MinimumClearance, but sometimes the large areas will have a narrow part attached to it, or two very close vertices.
Any ideas of other approaches?

In the end what worked best for me was using a negative buffer, as suggested by @Cyril and @FGreg below.
I used something like:
ST_Area(ST_Buffer(geom, -10)) as neg_buffer_area

In my case, units were meters, so 10 m negative buffer.
For narrow polygons, this area returned zero (also, the geometry would be empty). Then I used this column to filter out the narrow polygons.


texts: How to draw a polygon between a number of polygons automatically

texts:
The weight “empty” area which is bounded by the orange boundary should represent a land use “Roads” so I should fill it with polygons. Now my question is how to fill these empty areas with polygons automatically?

I could use Autocomplete Polygons to fill the gaps, but it’s time consuming since I have more than 150 Urban master plans! Is there any method to do this more quickly and more precisely?


texts: Why is arcpy.StackProfile_3d -> arcpy.SaveGraph_management creating empty graph?

texts:
In a python script I'm trying to save a graph to a .png file with graph output from arcpy.StackProfile_3d but the .png file shows an empty graph. The table created from arcpy.StackProfile_3d shows values. Does anybody know why the graph is still empty?

import arcpy
# Set Local Variables
profileLine = "C:/temp/data.gdb/profile"
profileTargets = "C:/temp/data.gdb/raster"
profileTable = "C:/temp/value.dbf"
graphName = "Sample Graph"
outGraph = "C:/temp/graph.png"

# Execute StackProfile
arcpy.StackProfile_3d(profileLine, profileTargets, profileTable, graphName)

# # Execute SaveGraph
arcpy.SaveGraph_management(graphName, outGraph, "MAINTAIN_ASPECT_RATIO")


texts: Several IF statements with Indentation error

texts:
I wrote a very long "if statement" script but something doesn't work. It returns several Indentation errors but I cannot understand where...could be something else?

from PyQt4.QtCore import *
from PyQt4.QtGui import *
import processing
#define layer
layer_sezcens=QgsVectorLayer(r'C:\Users\federica\Desktop\POLIMI\DOTTORATO\GIS MAPPING\Maggiolina_map\sezcens2.shx',"sezcens2","ogr")
#add new field
layer_sezcens.startEditing()
field_mainAge=layer_sezcens.dataProvider().addAttributes([QgsField("MainAge",QVariant.Int)])
layer_sezcens.commitChanges()
#modify new field
layer_sezcens.startEditing()
field_mainAge=layer_sezcens.dataProvider().fieldNameIndex('MainAge')
field_minDiffer=layer_sezcens.dataProvider().fieldNameIndex('MinDiffer')
field_FlatsClass= layer_sezcens.dataProvider().fieldNameIndex('maxFlatsC')
print field_mainAge,field_minDiffer,field_FlatsClass
for feat in layer_sezcens.getFeatures():
        if feat['E3']==NULL:
        expr_mainAge=0
    else:
        if feat('MinDiffer')>0:
            if feat['MinDiffer']==feat['MaxFlatsC']-feat['E8']:
                expr_mainAge=1
            elif feat['MinDiffer']==feat['MaxFlatsC']-feat['E9']:
                expr_mainAge=2
            elif feat['MinDiffer']==feat['MaxFlatsC']-feat['E10']:
                expr_mainAge=3
            elif feat['MinDiffer']==feat['MaxFlatsC']-feat['E11']:
                expr_mainAge=4
            elif feat['MinScarto']==feat['MaxFlatsC']-feat['E12']:
                expr_mainAge=5
            elif feat['MinDiffer']==feat['MaxFlatsC']-feat['E13']:
                expr_mainAge=6
            elif feat['MinDiffer']==feat['MaxFlatsC']-feat['E14']:
                expr_mainAge=7
            elif feat['MinDiffer']==feat['MaxFlatsC']-feat['E15']:
                expr_mainAge=8
            elif feat['MinDiffer']==feat['MaxFlatsC']-feat['E16']:
                expr_mainAge=9
        elif feat('MinDiffer')=0:
            list_mainAge=[feat['E8'],feat['E9'],feat['E10'],feat['E11'],feat['E12'],feat['E13'],feat['E14'],feat['E15'],feat['E16']]
            max_AgeList=max(list_mainAge)
            if feat['E8']=max_AgeList:
                expr_mainAge= 1
            elif feat['E9']=max_AgeList:
                expr_mainAge=2
            elif feat['E10']=max_AgeList:
                expr_mainAge=3
            elif feat['E11']=max_AgeList:
                expr_mainAge=4
            elif feat['E12']=max_AgeList:
                expr_mainAge=5
            elif feat['E13']=max_AgeList:
                expr_mainAge=6
            elif feat['E14']=max_AgeList:
                expr_mainAge=7
            elif feat['E15']=max_AgeList:
                expr_mainAge=8
            elif feat['E16']=max_AgeList:
                expr_mainAge=9
    print expr_mainAge
    layer_sezcens.changeAttributeValue(feat.id(), field_mainAge, expr_mainAge)
layer_sezcens.commitChanges()


texts:
I have two raster data layers from cropscape for 2014 and 2015 that have been clipped to Fresno County. I want to look at the change in almonds between these two years (if and where almonds were added or taken out). I used the difference function in the Image Analysis window to obtain the temporary difference/change layer. What I'm unsure how to do, is make this layer useful and classify it so I can see areas of positive/negative change and no change.

I watched a couple video tutorials explaining how to perform change analysis using landsat data, but I'm not sure how to apply this to my cropscape data. For instance, to use the remap function to classify the difference data, I need to know what definition parameters to enter, but I'm not sure of what these would be.

Does anyone have any suggestions or other functions that would be better to use?


texts: Adding data from table to grid in ArcGIS for Desktop?

texts:
I have a table with 720 rows and 156 columns. I have a grid with 720 rows and 156 columns. I would like to assign the values from the table to the corresponding values in the grid.

I am a beginner to ArcGIS.


texts: What is extra file created in satellite image data with format IGWX and XML file

texts:
Generally XML in raster file contain image metadata with its detail information regarding coverage,image source, coordinate extent etc.But here, when I open this image in my software it placed in unusual place with random placement.
While I tried open same image after deletion of IGWX and XML file it correctly placed.
What is reason of this back drop file.
Attached snapshot shows before files generated.


texts: GeoDjango raster field display on Leaflet

texts:
In the new Django version support raster field and I follow this manual for GeoDjango raster field and this question for add raster I have built a simple demo.
Here the code :
settings.py
postgis database backend
'django.contrib.gis.db.backends.postgis'

models.py
from django.db import models
from django.contrib.gis.db import models

class Elevation(models.Model):
    name = models.CharField(max_length=100)
    rast = models.RasterField()

All works fine without error. But I have one BIG question where I don't know how to complete it.
How do I display this raster in some web map like Leaflet?


texts: Calculating sequential numbers into sorted table not dependent on OBJECTID?

texts:
I am having a problem with adding sequential numbers to the new column based on other sorted value. I would like to use 'Sort by' column to be able to generate numbers in the 'Required number'column.

I have a script listed below (Field calculator) that allows me to generate numbers in 'Script used' column. Unfortunately, the numbers start from 0 instead of 11.

mxd = arcpy.mapping.MapDocument("CURRENT")
lr=arcpy.mapping.ListLayers(mxd)[0]

tbl=arcpy.da.TableToNumPyArray(lr,("objectid","MapOrderNum"))#MapOrderNum is a field that you want to use to sorb by
bs=sorted(tbl,key=lambda x: x[1])
aDict={}
for i,row in enumerate(bs):
 aDict[row[0]]=i
def sortSeq(objectid):
 return aDict[objectid]

sortSeq ( !OBJECTID!)


texts: How to set up QGIS + georeferencer for easiest georeferencing

texts:
The goal is to make converting scanned paper maps to map tiles as easy as possible.

I've done a few, but the slow part is manually entering coordinates of points. What I'd like to do is for each point on the scanned paper map, click on a satellite image (ie a Google Map) to automatically set its coordinates. I know this should be possible with QGIS and the georeferencer plugin (I'm using 1.7.2) but...

If I set QGIS to use WGS84 projection, I haven't found a way to display Google Maps in QGIS so that the georeferencer plugin tool will take coordinates from it.

If I set QGIS to use EPSG:900913, I can use the Openlayers plugin to add a Google Maps layer, and the georeferencer plugin tool happily takes coordinates from it. But gdal2tiles won't process the resulting geotiff file that georeferencer produces.

So currently I set QGIS and georeferencer to use WGS84, type in lat/lon pairs manually, and then gdal2tiles processes the resulting geotiff successfully.

Has anyone found a workflow that avoids the problems I described above?

Thanks & regards, David


texts: QGIS unable to export OSM topology to SpatiaLite

texts:
QGIS 2.18.2
I am downlaoding OSM data for a QGIS project.
The first 2 steps work fine: Vector/OpenStreetMap/Download data and Import topology from XML.
However when I try to execute step 3 - Export topology to SpatiaLite - I get an error "Failed to export OSM data: Unable to create table.  Duplicate column name FIXME."  The error occurs for all 3 export types, point, line and polygon.


texts: Extract areas with certain pixel values from raster

texts:
I have land cover data from GEE (https://developers.google.com/earth-engine/datasets/catalog/ESA_GLOBCOVER_L4_200901_200912_V2_3) for specific countries. I would like to extract the values from the raster for certain regions. Specifically, I'd like to have a SpatialPolygonsDataFrame with one variable for each land cover class containing the area (absolute or relative) of that region covered in a specific class. My first idea was to transform the raster into a SpatialPolygons object and then calculate the intersection, but from what I read transforming rasters is quite inefficient and functions like extract should be used. How can I achieve what I want without transforming the raster?

edit:
Sorry, I feel like title and text don't really match. My original idea was to convert the raster, then create a subset only containing areas with a certain values, then calculating intersections and distances. Is it possible to do this with the raster too?


texts: Disappearing polygons when zooming in Feature Cache extent

texts:
I'm currently using ArcMap (v 10.6.1.) and an enterprise geodatabase (v 10.6.0) on SQL Server 2017.
I encountered a problem, when I build Feature Cache on full map extent and zoom in to my polygon feature. The polygon geometry is being cut and reduced, despite the fact that I do not use any of editing tools at that moment. Polygons are cut out of the shape of zoomed in extent (1st and 2nd picture). This happens while Feature Cache is working. When I zoom out until the point when Feature Cache isn't working anymore polygon appears to be back (3rd picture).

I tried to install some ArcGIS 10.6.1. patches:

ArcGIS (Desktop, Engine) 10.6.1 General Update Patch
ArcGIS (Desktop, Engine) Text Performance Patch
ArcGIS Desktop 10.6.1 Exporting Data From a Layer Patch
ArcGIS Desktop Dialog Initialization Performance Patch
ArcGIS Desktop TLS Patch ArcGIS 10.6.1 (Desktop, Engine, Server) Buffering Degenerated Polygon Patch

But none of them helped.
Has anyone encountered a similar problem? Is there any solution currently available?


texts: Point buffer not round - jagged circle like thing

texts:
Any ideas why running a point buffer (50 and 150 meters) would result in this shape?

ArcGIS 10.2 - All projections match


texts: How to save a QGIS workspace with database link for offline use?

texts: Is it possible to save current workspace of QGIS for offline usage like in mapinfo? My workspace has link to database. and when I reopen workspace unlink in mapinfo it can't open it.


texts: Bing and Esri satellitary basemaps are misaligned by about 2 meters

texts:
Did you ever notice a shift (2 meters) between satellitary basemaps from Bing and Esri?

The problem basically arose when GIS technicians have to deal with drawings made in CAD using Bing basemaps as reference that do not overlap properely with the Esri basemaps.

And this is particularly important when one wants to share a map through ArcGIS Online, which obviously uses Esri basemaps (I know it is possible to use Bing maps with ArcGIS Online, but this is not the point).

As far as I know, both basemaps have the same Web Mercator reference system.

Is there any way to deal with this other than shifting the whole data by a delta x, y in your opinion?


texts:
I have a custom map which I get from a mbtiles file with the bounds -13,25,49,60 and zoom from 4 to 10.The display div is dynamic, I mean that it depends on window size and the extent restriction works but still tries to load images around the map if the display div is bigger.
How can I limit the zooming out below 4?
Here is my code:

var map;
var fromProjection = new OpenLayers.Projection("EPSG:4326");   // Transform from WGS 1984
var toProjection   = new OpenLayers.Projection("EPSG:900913"); // to Spherical Mercator Projection
var extent= new OpenLayers.Bounds(-13,25,49,60).transform( fromProjection, toProjection);
var mbTiles = new OpenLayers.Layer.TMS("base", "base.php", {
    getURL: mbtilesURL,
    type:'png',
    transitionEffect: "resize",
    isBaseLayer: true,

});
function mbtilesURL (bounds) {
    var db = "img/map/base.mbtiles";
    var res = this.map.getResolution();
    var x = Math.round ((bounds.left - this.maxExtent.left) / (res * this.tileSize.w));
    var y = Math.round ((this.maxExtent.top - bounds.top) / (res * this.tileSize.h));
    var z = this.map.getZoom();
    return this.url+"?db="+db+"&z="+z+"&x="+x+"&y="+((1 << z) - y - 1);
}
map = new OpenLayers.Map('map', {
    projection: new OpenLayers.Projection("EPSG:900913"),
    displayProjection: new OpenLayers.Projection("EPSG:4326"),
    layers: [mbTiles],
    restrictedExtent : extent,
    numZoomLevels: 11,
    controls: [
        new OpenLayers.Control.Navigation({
            dragPanOptions: {
                enableKinetic: true
            }
        }),
    new OpenLayers.Control.PanZoom(),
    new OpenLayers.Control.Attribution()
    ]
});
var position = new OpenLayers.LonLat(9.4679,37.6831).transform( fromProjection, toProjection);
var zoom = 4;
map.setCenter(position, zoom );
map.zoomToMaxExtent();
var switcherControl = new OpenLayers.Control.LayerSwitcher();
map.addControl(switcherControl);
switcherControl.maximizeControl();

I made it work! Had to add :

   resolutions: [4891.969809,2445.984905,1222.992452,611.4962262,305.7481131,152.8740565]

and to add to z+4 to get to right zoom level. I used this site http://simonmercier.net/blog/?p=88 to get the right resolutions to use.


texts: Return random string from a list of strings? Field Calculator, Python

texts:
I am having an issue with producing the correct syntax to place random strings from a list of min into fields with the field calculator.

I have tried the following code:

import random

ran = ['kid', 'many', 'love', 'play']

random.shuffle(ran)

for i in ran:

print i

I am only getting one value? I would like to have a random values throughout the entire field.


texts: QGIS 3.16.5 v.distance returns UnicodeDecodeError

texts:
I'm trying to run v.distance on QGIS but get the following message/s:
2021-04-14T17:45:55     CRITICAL    Traceback (most recent call last):
  File "C:/PROGRA~1/QGIS3~1.16/apps/qgis-ltr/./python/plugins\processing\algs\grass7\Grass7Algorithm.py", line 434, in processAlgorithm
  Grass7Utils.executeGrass(self.commands, feedback, self.outputCommands)
  File "C:/PROGRA~1/QGIS3~1.16/apps/qgis-ltr/./python/plugins\processing\algs\grass7\Grass7Utils.py", line 398, in executeGrass
  for line in iter(proc.stdout.readline, ''):
  File "C:\PROGRA~1\QGIS3~1.16\apps\Python37\lib\encodings\cp1255.py", line 23, in decode
  return codecs.charmap_decode(input,self.errors,decoding_table)[0]
 UnicodeDecodeError: 'charmap' codec can't decode byte 0x8c in position 113: character maps to

How do I fix it?


texts: Python adding grouplayer to muliple maps

texts:
I have the following code which is supposed to add a one group layer to each of the maps. In this case the group layer is always the same.

Desired output:

Map1 = 1 group layer added
Map2 = 1 group layer added
Map3 = 1 group layer added

However this is the result I am getting:

Map1 = 1 group layer added
Map2 = 2 group layers added
Map3 = 3 group layers added

Can someone tell me why/how the results are being compounded?

Note: I used to have this code within a "for" loop. Split it out so could try to figure out why not working correctly.

    SDESurveyed = (r"M:\projects\SDE Surveyed Layers")

    #This grouplayer is blank. Used because Python currently cannot create a grouplayer on it's own
    ClientLayer = arcpy.mapping.Layer(r"F:\TOOLS\PythonTools\ClientSurveyedData.lyr")

    #Look for group layer the matches Client name
    for filename in os.listdir(SDESurveyed):
        if fnmatch.fnmatch(filename[:-4], Client):
            surveyedFile = SDESurveyed + "\\" + filename
            #Create a layer file to the path
            SDElayer = arcpy.mapping.Layer(surveyedFile)

    #SDEkeys = ["AC_","CN_","AM_"]
   # for item in mxd:
        #if item in SDEkeys:
            #imxd = mxd[item]
    imxd = mxd["CN_"]
    print imxd
    df = arcpy.mapping.ListDataFrames(imxd,"Project Area")[0]
    for lyr in arcpy.mapping.ListLayers(imxd, "*", df):
        if lyr.name == "ClientSurveyedData":
            arcpy.mapping.RemoveLayer(df, lyr)
            print("PREVIOUS CLIENT DATA LAYER REMOVED...")
    print("ADDING NEW CLIENT DATA LAYER...")

    targetClientLayer = arcpy.mapping.ListLayers(imxd, "COMPANY SURVEYED DATA", df)[0]

    # Add the layer to the layer group
    arcpy.mapping.AddLayerToGroup(df, targetClientLayer, ClientLayer, "BOTTOM")

    #Target Layer to where layer will open
    targetGroupLayer = arcpy.mapping.ListLayers(imxd, "ClientSurveyedData", df)[0]

    # Add the layer to the layer group
    arcpy.mapping.AddLayerToGroup(df, targetGroupLayer, SDElayer, "BOTTOM")
    print "LOOPED"
    imxd.save()
    #del targetClientLayer
    #del targetGroupLayer

    imxd = mxd["AM_"]
    print imxd
    df = arcpy.mapping.ListDataFrames(imxd,"Project Area")[0]
    for lyr in arcpy.mapping.ListLayers(imxd, "*", df):
        if lyr.name == "ClientSurveyedData":
            arcpy.mapping.RemoveLayer(df, lyr)
            print("PREVIOUS CLIENT DATA LAYER REMOVED...")
    print("ADDING NEW CLIENT DATA LAYER...")

    targetClientLayer = arcpy.mapping.ListLayers(imxd, "COMPANY SURVEYED DATA", df)[0]

    # Add the layer to the layer group
    arcpy.mapping.AddLayerToGroup(df, targetClientLayer, ClientLayer, "BOTTOM")

    #Target Layer to where layer will open
    targetGroupLayer = arcpy.mapping.ListLayers(imxd, "ClientSurveyedData", df)[0]

    # Add the layer to the layer group
    arcpy.mapping.AddLayerToGroup(df, targetGroupLayer, SDElayer, "BOTTOM")
    print "LOOPED"
    imxd.save()

    imxd = mxd["AC_"]
    print imxd
    df = arcpy.mapping.ListDataFrames(imxd,"Project Area")[0]
    for lyr in arcpy.mapping.ListLayers(imxd, "*", df):
        if lyr.name == "ClientSurveyedData":
            arcpy.mapping.RemoveLayer(df, lyr)
            print("PREVIOUS CLIENT DATA LAYER REMOVED...")
    print("ADDING NEW CLIENT DATA LAYER...")

    targetClientLayer = arcpy.mapping.ListLayers(imxd, "COMPANY SURVEYED DATA", df)[0]

    # Add the layer to the layer group
    arcpy.mapping.AddLayerToGroup(df, targetClientLayer, ClientLayer, "BOTTOM")

    #Target Layer to where layer will open
    targetGroupLayer = arcpy.mapping.ListLayers(imxd, "ClientSurveyedData", df)[0]

    # Add the layer to the layer group
    arcpy.mapping.AddLayerToGroup(df, targetGroupLayer, SDElayer, "BOTTOM")
    print "LOOPED"
    imxd.save()


texts:
I recently migrated some map services from Geoserver to QGIS Server and up till now everything works quite nice. I only miss the possibility to cache tiles, mostly used by Openlayers.

Is there any alternative to GeoWebCache for QGIS Server? I see I could use its standalone version but want to avoid Java and Tomcat on the server. So what would be the best alternatives?

I would prefer any open source alternative written in PHP or Python. So far I had a look on TileStache, but before giving it a shot I wanted to hear your recommendations.


texts: How can I fill empty raster cells with values from another raster in QGIS?

texts:
I have 2 rasters:

1 that has many gaps in between pixels:

And one that is more filled:

I want to fill the empty cells from raster 1 with the cells from raster 2.  Cells that are filled in raster 1 should be kept.

How can I do this using QGIS?


texts: Update, Insert Cursors don't work on Features that have annotation?

texts:
I've been working on a script to update my feature classes by grabbing new data from one feature class, and applying it to another feature class.  It does this with several cursors, basically get value from the new data and set value to the old data.  It works just fine on feature classes without feature-linked annotation.

However, if the feature class has feature-linked annotation it will run through the script and not save the updates or appends.  Is this a known limitation, am I doing something wrong, or is there a work around to handle this?

My other post shows the layout of my code for reference.


texts: Converting vertical datum from WGS84 to EGM96 in ArcGIS Desktop?

texts:
I am comparing two different Digital Elevation Models

SRTM (Datum: EGM96)
Cartosat (Datum: WGS84)

How to convert WGS84 to EGM96 in ArcGIS Desktop?


texts: Select features using an expression, QGIS

texts:
I've a problem when using Select features using an expression. I've an attribute table with many features: It's a shapefile from a previous CSV table. The data are localities for water analysis, where several parameters appear. But I only want few of those locations.

What I do is:in Attribute table,  I open Select by Expression, Then in Fields and Values, I select my column of interest "Parameters", and from Load Values- All Unique button, I choose the names I want to pick : alfa- hexaclorociclohexano, beta-hexaclorociclohexano, gamma-hexaclorociclohexano.

The expression I use is:

"Parametro"  =  'alfa-Hexaclorociclohexano'  +   'beta-Hexaclorociclohexano'  +  'delta-Hexaclorociclohexano'  +  'gamma-Hexaclorociclohexano'

But I can't select my values...  If I just pick every value independently, i.e.: 'alfa-Hexaclorociclohexano', they are selected. Then, what my problem can be ?, should I not use the + symbol ?  I also tried with the And symbol, but it's not working either.  Should I try another option like the Query Builder ?

Select by Expression seemed the easier option to me but I'm not able to continue.


texts: Assigning values to field based on values in other field using Python Parser of ArcGIS Field Calculator?

texts:
I have a table with two fields A and B, and i want fill field B with "some values" from field A:

FieldA:
123456_CC;
123_D;
1234567_AB;
123456768_F

I just need the "text" that is associated with the value of the field

FieldB:
CC;
D;
AB;
F

What expression to create this?


texts:
I have a 30m raster dataset for the continental US, and would like to divide it into manageable chunks to convert it to vector, etc.

What's the best way to accomplish this? I don't want any overlap within the tiles, and don't care about the shape of the resulting rasters. The Split Raster tool in ArcGIS would seem the natural choice, but it requires that I specify a resampling method. I don't want to resample my data at all-- I want this the same data as the original raster, just in smaller pieces.

Happy to pursue python-based solutions if need be, but hoping there's some out of the box solution that doesn't involve transforming the data.


texts: Turning off high resolution imagery in ArcMap World Imagery basemap?

texts:
I was helping a patron find imagery in an area using ArcMap's World Imagery as a basemap, and found a troubling issue. We just needed low resolution imagery for locating study sites, as we were reloading the map often. However, on opening the context menu and attempting to remove layers from rendering, we were greeted with an error, "This type of layer cannot be turned off". Screenshots follow:

I've tried changing visual scale range, but that is impossible with the layer.

Is there any way to force the layers off or remove the high resolution imagery?


texts: Can't load requested DLL The specified module could not be found

texts:
I want to import cartopy and geopandas,it keeps showing me:
Can't load requested DLL: D:\Applications\python-3.8.2\Lib\site-packages\osgeo\gdalplugins\ogr_FileGDB.dll 126: The specified module could not be found.
But question is the DLL is right there. I have tried to change the PATH to D:\\Applications\\python-3.8.2\\Lib\\site-packages\\osgeo\\gdalplugins with double slash, and it still won't work.
My gdal was installed with pip install gdal and then added the path to environment variables
And I also tried to reinstall gdal but nothing changed.
So what should I do?


texts: Kriging metal concentration with irregular set of sample points

texts:
I need to perform some interpolation to map a metal concentration around a platform. The distribution of my stations is not regular (see figure below where stations are the green points labelled with the relative metal concentration).

What is the correct procedure/method for kriging? Are there any other specific settings, techniques, or values to look at when dealing with my hotspots and kriging?


texts: Analytical Hierarchy Process (AHP) tool: Input layers do not have same grid extent

texts:
I am using this tool to create a suitability map for landfill site selection. However, when I gave the input raster layers, the .csv file for comparison table and run the tool, an error message saying- "Unable to execute algorithm - Input layers do not have same grid extent"

I infer from this that the size of all input layers must be same, but do not know if that can be done or not. And if yes, then how?


texts: Adding large amount of data to web map using ArcGIS API for JavaScript?

texts:
I've got problem, I'm creating map based application, and now I have to add to map about 100 000 graphic points, from coordinates.

I was trying in many ways to do it, but always process of adding points takes about 7-10 seconds, or my web browser just hang on.

is it even possible to add this amount of data to map without having problems with performance?

I'm using ArcGIS API for JavaScript with ArcGIS Server 9.3.1.


texts: LULC Classification of Sentinel-2 Using XGBoost in R

texts:
I am trying to do Land cover classification using xgboost in R language.
These are the following steps:
> library(raster)
Loading required package: sp
Warning message:
R graphics engine version 14 is not supported by this version of RStudio. The Plots tab will be disabled until a newer version of RStudio is installed.
> shp <- shapefile("C:\\Users\\Admin\\Desktop\\Sat\\Training Samples_1.shp")
> ras <- stack("C:\\Users\\Admin\\Desktop\\Sat\\Final TIFF\\subset.tiff")
> ras
class      : RasterStack
dimensions : 6735, 11679, 78658065, 7  (nrow, ncol, ncell, nlayers)
resolution : 10, 10  (x, y)
extent     : 512504.8, 629294.8, 1276355, 1343705  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=43 +datum=WGS84 +units=m +no_defs
names      : subset.1, subset.2, subset.3, subset.4, subset.5, subset.6, subset.7

> shp
class       : SpatialPolygonsDataFrame
features    : 4
extent      : 515066.7, 704552.5, 1265814, 1399478  (xmin, xmax, ymin, ymax)
crs         : +proj=utm +zone=43 +datum=WGS84 +units=m +no_defs
variables   : 6
names       :   Classname, Classvalue, RED, GREEN, BLUE,  Count
min values  : Barren Land,          1,   0,   112,  106, 275174
max values  :       Water,          4, 162,    91,   75, 893831
> library(xgboost)
> vals <- extract(ras,shp)
> classes <- as.numeric(as.factor(shp@data$class)) - 1
>
> train <- data.matrix(vals)
>
> classes <- as.numeric(as.factor(shp@data$class)) - 1
>
> xgb <- xgboost(data = train,
+                label = classes,
+                eta = 0.1,
+                max_depth = 6,
+                nround=100,
+                objective = "multi:softmax",
+                num_class = length(unique(classes)),
+                nthread = 3)
Error in xgb.DMatrix(data, label = label, missing = missing) :
  REAL() can only be applied to a 'numeric', not a 'list'

I get the Error:
Error in xgb.DMatrix(data, label = label, missing = missing) :
REAL() can only be applied to a 'numeric', not a 'list'.
I found the tutorial for this here


texts: intermittent failures with SelectByLocation in an ArcPy addin

texts:
ArcGIS 10.2, Windows 10

this is the line that fails:

arcpy.SelectLayerByLocation_management(WRS_layer, "INTERSECT", gridLayer, "", "NEW_SELECTION")

The WRS_layer is in a file gdb, the gridLayer is in an SDE gdb. Both are polygon feature classes.  SDE layers are being edited by a user-initiated edit session, but in this addin button there is no arcpy.da.Editor.  The SelectByLocation will work fine for a while, and then starts giving an error message:

Traceback (most recent call last):
File "C:\Users\...\ESRI\Desktop10.2\AssemblyCache\{1E...D8E}\NORDIS_Landsat.py", line 432, in onClick
   arcpy.SelectLayerByLocation_management(WRS_layer, "INTERSECT", progressLayerName, "", "NEW_SELECTION")
File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\management.py", line 6812, in SelectLayerByLocation
    raise e
ExecuteError: ERROR 999999: Error executing function.
General function failure
General function failure [WRS2_Nordis]
Failed to execute (SelectLayerByLocation).

Here is the referenced line in management.py

6806      from arcpy.geoprocessing._base import gp, gp_fixargs
6807    from arcpy.arcobjects.arcobjectconversion import convertArcObjectToPythonObject
6808    try:
6809        retval = convertArcObjectToPythonObject(gp.SelectLayerByLocation_management(*gp_fixargs((in_layer, overlap_type, select_features, search_distance, selection_type), True)))
6810        return retval
6811    except Exception, e:
6812        raise e

Arcmap doesn't crash, but only a restart of ArcMap solves the problem. I don't think this is a code problem, but I hope that someone else has experienced something similar and can suggest a cause and/or solution.

EDIT: as requested, here is the top part of the button's code immed after def onClick(self)

mxd = arcpy.mapping.MapDocument("CURRENT")
# SET SOME VARIABLES
Error = False
WRS_layer          = "WRS2_Nordis"
WRS_fields         = ["PR"]
username           = getpass.getuser().lower()
rastersPresent     = []
rastersMissing     = []

# 0. VERIFY THAT THE REQUIRED LAYERS ARE PRESENT IN THE MXD.
print
layersPresent = []
for lyr in arcpy.mapping.ListLayers(mxd):
    layersPresent.append(lyr.name)
for lyrReqd in (WRS_layer, progressLayerName):
    if not lyrReqd in layersPresent:
        Error = True
        print "ERROR: cannot find layer  ", lyrReqd, "  in the TOC."

# 1. GET THE SELECTED GRID CELL
d = arcpy.Describe(progressLayerName)
if len(d.FIDSet) == 0:
    ProgressGrid_Select = 0
else:
    ProgressGrid_Select = len(d.FIDSet.split(";"))
if ProgressGrid_Select <> 1:
    print " You have", str(ProgressGrid_Select), " ProgressGrid cells selected  -- should be exactly 1."
    Error = True

# 2. GET THE INTERSECTING WRS SCENES
# make a WRS list of lists, with PathRow

if not Error:
    # SelectLayerByLocation HONOURS THE SELECTION SET IN in_layer
    print 'DEBUG: attempting SelectLayerByLocation_management(WRS_layer, "INTERSECT", progressLayerName, "", "NEW_SELECTION")'
    arcpy.SelectLayerByLocation_management(WRS_layer, "INTERSECT", progressLayerName, "", "NEW_SELECTION")
    #  make list of WRS scenes, with PR (data type Long)
    WRS_scenes = []
    # SearchCursor (in_table, field_names, {where_clause}, {spatial_reference}, {explode_to_points}, {sql_clause})
    with arcpy.da.SearchCursor(WRS_layer,WRS_fields) as cursor:
        for row in cursor:
            WRS_scenes.append(row[0])

    if len(WRS_scenes) == 0:
        print "No WRS_scenes intersect with selected Progress Grid cell."
        Error = True
    else:
        print "Found " + str(len(WRS_scenes)) + " WRS scenes for selected progress cell:"
        WRS_scenes.sort()
        prs = ''
        for pr in WRS_scenes:
            prs += ("  " + str(pr))
        print prs
    print


texts: I have a SRTM DEM and need to find whether it has voids or gaps in the data.


texts: Geoserver Vector tiles hollow issue with layer

texts:
I am using vector tiles from PostGIS, Geoserver for the world layer,
I have multiple layers for the global level,
I have seen multiple holes and missing tiles of some boundaries in the layer of vector tiles.
I have tried multiple methods to resolve the issue still not found the solution.

geometry st_isvalid,
geometry topology.

WMS layer with the same geometry.
no hollow.

same layer for vector tiles hole in Layer

Is this related to the multiPolygons/geometryCollection issue or something else?


texts: Field Collection App for Geoserver or Any OS GIS Server

texts:
I am looking for a mobile app/interface for field data collection and mobile editing of GIS data that relies/syncs back to some type of GIS server or database. I am aware of Fulcrum and AmigoCloud, but those are third party. I have yet to roll out a GIS server simply for this reason, I still have to use Collector/Fulcrum/etc to collect the data and have no way to edit the data via mobile. In my initial search for open source GIS server solutions I assumed that this functionality already existed, but I have not yet found a good solution.

Some workarounds I have looked at:

Google pre-filled Form using gps location from a leaflet map
Memento database for Android for data collection/synced with google drive
Various other android gps collection apps
hosting a uMap installation which allows editing of features
Appsheet for Google Sheets (comes close but has its limits)
Mapcentia hosted server (they are supposed to be working on an Android field data collection tool)

Even if the field collection tool could sync to cloud storage that would get me closer.


texts: Selecting more than one point using QGIS

texts:
I am trying to convert the monasteries layers from:
https://harvard-cga.maps.arcgis.com/apps/View/index.html?appid=b4dad589ee4b40e39a93bb0b026e4ee8
into a csv table. For that I am using the WMS link:
https://ags.cga.harvard.edu/arcgis/services/darmc/roman/MapServer/WMSServer
I do that via Layer> add Layer> add WMS/WMTS layer..., then I click new and paster the WMSServer link. I am then able to see all the monastic houses points, so up to here everything all right. Now, comes my problem: Although when I click on one dot I can see the associated info, I don't know how to click more than one spot. I have tried to use the identy features by radius, but then I cannot gather the points.


texts: How do you specify the width of a column that has been extruded from a point in arcscene?

texts: I've extruded a set of points based on a z variable but I want to increase their width so they're easier to see.


texts: Is there a way to automatically run the Workflow Manager Postinstallation (WMXPostInstallation.exe)?

texts: I'm looking for a way to script the WMX post-install and, for example, provide the sde connection info as parameters.  I don't see any command-line parameters for WMXPostInstallation.exe--are there any undocumented command line parameters.  Or is there different approach to automating post-installs?


texts: Query WMS for available layers using QGIS 3

texts: I believe we can do it using ArcGIS and I want to know if it is also possible, to query objects in a WMS layer using QGIS 3, my research so far have been unsuccessful. The identification tool does not work, nor the one for consulting the attribute table.


texts:
I'm working on my own geoportal. I have prepared a php page which users can describe their own SLD configuration. The SLD file will be stored in GeoServer style folder. But when I tried to view maps I couldn't see the changes on map. But if I go to the styles editor and click to the Submit button in my new SLD files, I can see the new SLD configuration on my map. I think Submit button is making some another configuration that ı don't know.

Are there any way refreshing the style edits?


texts: Extracting Lat/Lng from Shapefile using OGR2OGR/GDAL

texts:
I'm trying to extract lattitude/longitude values for the streets in this shapefile:
http://gispub02.sfgov.org/website/sfshare/catalog/stclines_streets.html

I've seen it reommended in many places that I use this command to do so:

(Below does work! Spelling error fixed.)

  ogr2ogr -f CSV output.csv STCLINES_STREETS.shp -lco GEOMETRY=AS_XYZ

However, when I do so, I get a CSV containing a lot of useful information about the streets -- names, neighborhoods, etc -- but no actual lat/long values.

What am I doing wrong? I'm a GDAL noob so am sure I'm missing something obvious.


texts: How to create a heatmap with numbers using Leaflet map?

texts:
I wish to create a heatmap like this one at http://www.trulia.com/local/san-francisco-ca/tiles:1|points:1_crime

With the number changes at different zoom level - the number can appear in a mouse-hover window instead of directly showing on the map.

Any idea how to make this with Leaflet map?

EDIT:

I have a csv file with three columns:

    lat  lng  count

Any tips how to easily load it for
Leaflet.markercluster?


texts: Is it possible to expose Bing Maps as a WMS endpoing via GeoWebCache?

texts:
Is it possible, using GeoWebCache or some other software, to expose Bing Maps as a WMS endpoint? I'm trying to incorporate Bing Maps into custom built software. The software is designed to communicate with WMS web services. Being able to access Bing Maps (Google Maps or some other mapping API would also work) via WMS would be great, since it would require virtually no changes to our software.

I've noticed that GeoWebCache can communicate with Bing Maps, but it appears that its web service exposes the map using the same protocol as Bing itself.


texts: Connecting to ArcSDE without ArcPy license

texts: Currently, I do not have a license for ArcMap and I am wondering if it is possible to get data from an ArcSDE without ArcPy. I have been able to leverage GDAL to read / write to a file geodatabase but am curious if anyone has accomplished this before.


texts: Cut, rotate and interpolate raster heightmap with QGIS

texts:
I have some raster heightmap (example, GeoTIFF). I need to cut a rectangle at a specified angle, rotate this rectangle to a horizontal matrix and interpolate it.

I'll try to explain in more detail (using a combination of QGIS & GIMP).

1) Open map.

2) Create rectangle.
I create a new vector layer with polygon type and insert 4 points with the mouse. Maybe I can insert a rectangle as it is easier and more accurate?

3) Clipper.
Menu: Raster -> Extraction -> Clipper.

4) I got a rectangular piece of map, but do not know what to do next.

Next steps I will do in GIMP. (because I don't know how to show it in QGIS)

5) Rotate the map horizontally.

6) Set grid with a given accuracy and consider the average value of the height, which fell in the grid cell.

7) Save new XY matrix to some file.

How can I do it only with QGIS?
How can I do it with only GDAL and Proj4?


texts:
I'm new to all this systems and world !
I have some points represented in the XY coordinate in Tunisia and I'm looking for a library to transform those points.
for example : x = -92336.468 / y = 99043.202 and the origin of the XY coordinate is : lat: 33.8869 ; long: 9.5375
when I use this website converter from XY to Lat-Long I get the adequate result when testing the results on google map .
What's the algo of the conversion or if there is a ready Library that would be awesome.


texts: Update field based on values in other fields using search and update cursors

texts:
I'm a Python beginner, I'm doing a risk analysis with a point feature, and I'm trying to write a script that allow me to update a field (proximity) based on different numeric values within different fields.

Each field (field_1, field_2, field_3, field_4) has just one possible value ( respectively 1, 2, 3, 4) and in each row there is recorded only one of these different values.
Example: If in a row there's 1 in field_1  it's impossible to find other values in the other fields.

I would like to put everything in one field (proximity) to indicate the risk level.

I tried to run the following script to test only with field_1, but it's not updating, it's just printing "Nothing for each row".

with arcpy.da.UpdateCursor("point", "proximity") as cursor:
for row in cursor:
    with arcpy.da.SearchCursor("point", "field_1") as scursor:
        for row in scursor:
            if row[0] == 1:
                cursor.updateRow(row)
        else:
            print ("Nothing")


texts: Alternatives to joining attributes in QGIS?

texts:
In QGIS 2.18, I can join attributes stored in a csv-sheet to a shapefile, which works fine. However, when saving the respective layer, the joined attributes are not stored.

Is it possible to add attributes from a separate table permanently to a shapefile?


texts: Generating a token for ArcGIS REST Services in VBA

texts:
I am trying to generate a token in Excel VBA to add data to a table on our server.

I am using the following code to generate my token:

Dim tokenHTTP As Object
Set tokenHTTP = CreateObject("MSXML2.ServerXMLHTTP")
Dim tokenJSON As String
tokenJSON = "[{""f"" : ""json"", ""username"" : ""username"", ""password"" : ""password"", ""client"" : ""requestip"", ""expiration"" : 1}]"
Dim tokenURL As String
tokenURL = "https://myserver/arcgis/tokens/generateToken"
tokenHTTP.Open "POST", tokenURL, False
tokenHTTP.setRequestHeader "Content-type", "application/json"
tokenHTTP.send (tokenJSON)

The response is a link to the token usage page.

Can anyone help me with this?


texts: Couldn't load plugin 'processing' Ubuntu 14.04.5 Import Error

texts:
I've upgraded QGIS to 2.18 in our Ubuntu 14.05.5 recently, in the process we have some dependency problems with PostGIS but solved updating sources repositories (using http://qgis.org/ubuntugis & http://qgis.org/ubuntugis-nightly-release trusty).
But we still have a problem loading 'processing' plugin in QGIS 2.18.0 start up, getting the following error message box, because undefined symbol:

GDALTranslateOptionsFree causes an Import Error:
Couldn't load plugin 'processing'

ImportError: /usr/local/lib/python2.7/dist-packages/GDAL-2.1.0-py2.7-linux-x86_64.egg/osgeo/_gdal.so: undefined symbol: GDALTranslateOptionsFree
Traceback (most recent call last):
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 308, in loadPlugin
    __import__(packageName)
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 607, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
  File "/usr/share/qgis/python/plugins/processing/__init__.py", line 28, in
    from processing.tools.dataobjects import *          # NOQA
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 607, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
  File "/usr/share/qgis/python/plugins/processing/tools/dataobjects.py", line 37, in
    from processing.algs.gdal.GdalUtils import GdalUtils
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 607, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
  File "/usr/share/qgis/python/plugins/processing/algs/gdal/GdalUtils.py", line 32, in
    from osgeo import gdal
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 607, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
  File "/usr/local/lib/python2.7/dist-packages/GDAL-2.1.0-py2.7-linux-x86_64.egg/osgeo/__init__.py", line 21, in
    _gdal = swig_import_helper()
  File "/usr/local/lib/python2.7/dist-packages/GDAL-2.1.0-py2.7-linux-x86_64.egg/osgeo/__init__.py", line 17, in swig_import_helper
    _mod = imp.load_module('_gdal', fp, pathname, description)
ImportError: /usr/local/lib/python2.7/dist-packages/GDAL-2.1.0-py2.7-linux-x86_64.egg/osgeo/_gdal.so: undefined symbol: GDALTranslateOptionsFree

Python version: 2.7.6 (default, Oct 26 2016, 20:33:43) [GCC 4.8.4]
QGIS version: 2.18.0 Las Palmas, d8f9d2e

Python Path:
/usr/share/qgis/python
/home/tximikel/.qgis2/python
/home/tximikel/.qgis2/python/plugins
/usr/share/qgis/python/plugins
/usr/local/lib/python2.7/dist-packages/GDAL-2.1.0-py2.7-linux-x86_64.egg
/usr/lib/python2.7
/usr/lib/python2.7/plat-x86_64-linux-gnu
/usr/lib/python2.7/lib-tk
/usr/lib/python2.7/lib-old
/usr/lib/python2.7/lib-dynload
/usr/local/lib/python2.7/dist-packages
/usr/lib/python2.7/dist-packages
/usr/lib/python2.7/dist-packages/PILcompat
/usr/lib/python2.7/dist-packages/gtk-2.0
/usr/lib/pymodules/python2.7
/usr/lib/python2.7/dist-packages/ubuntu-sso-client
/usr/lib/python2.7/dist-packages/wx-2.8-gtk2-unicode
/home/tximikel/.qgis2//python


texts: Problem in cell dimension assignment from "save as" in QGIS

texts:
I am trying to save one raster in QGIS by changing the cell dimensions. I select 30 x 30 as cell resolution, but then the output raster has a resolution of 29.9985 x 29.9985.

How can I assign exactly 30 as cell resolution?


texts: Segmenting a sequence of GPS readings into coherent 'GPS Tracks'

texts:
Consider the problem of wanting to extract semantically coherent 'tracks' from a sequence of GPS readings.

In this context, I mean 'semantically coherent' to be arbitrary and subjective. For example, consider the Moves App: it's a mobile app which sits on a smartphone and reads the GPS position at a somewhat constant rate. Afterwards, it generates a 'timeline' of sorts from the readings into understandable segments, such as:

What are the methods/algorithms available to create such segments?

Currently, I use a naive ad hoc algorithm: 'replay' the measurements, adding them to the current segment, start a new segment if next position is too far away in time and/or space from the end of the current segment. I have moderate success with this approach, but this makes it too sensitive to outlier points and is not very prone to parallelisation in this form.


texts:
I am trying to set snap mode within QGIS plugin.

So far, I have tried:

QgsSnappingUtils().setSnapToMapMode(QgsSnappingUtils.SnapAdvanced)

and

QgsSnappingUtils().setSnapToMapMode(2)

without success.
QGIS Snapping options (Settings>>Snapping Options...) continue with:
 Snapping mode = Current layer

Any hint? (Using v 2.14.0)


texts: Calculate mean of raster excluding values in ArcMap?

texts:
I want to calculate the mean value of a specific region a raster data layer (being a GeoTIFF file), using the Zonal Statistics as Table tool in ArcMap.

However, some of the values are faulty and need to be excluded in the calculation.

How can I do this?


texts: FME - Get parts of lines within polygons and transfer polygon ID to lines in result

texts:
I've got a set of polygons (my investigation areas - identified by IDs) and a set of lines that intersect these polygons.
My task is to only get the lines within the polygons and give them the containing polygon's ID. I tried to do this with FME and used a Clipper, which gets the work done on the geometry part: It gives me the parts of the lines that are inside the polygons, but I don't know how to transfer the said ID, which I need to sort and group the lines by in my result.