texts: How do I make a map's background transparent in a layout?

texts:
I am trying to add an overview map to my layout. The main map is projected while the overview map is a 'globe' (image below). I want to make the white part of the overview map transparent.

I have made the background transparent in page properties and I have turned off background for both maps. There is no option to make the background of the main QGIS window transparent.

Is there anything else I can do other than export separately and add it as an image with transparency values?

EDIT: Forgot to mention, I have looked at How to export a map with a transparent background in QGIS? but it didn't work.


texts: How to run Processing qgis:extractbylocation?

texts:
I try to use the processing 'qgis:extractbylocation' in QGIS 2.8.2 from the QGIS algorithm provider.

import processing
# import the layers by a given uri
inputlayer = QgsVectorLayer(uriTr.uri(), "inputlayer", "postgres")
intersectlayer = QgsVectorLayer(uriOSM.uri(), "intersectlayer", "postgres")

# initialize a predicate instance
pred_intersects = processing.core.parameters.ParameterGeometryPredicate(enabledPredicates = 'intersects')
# run the processing
processing.runalg('qgis:extractbylocation', inputlayer, intersectlayer, pred_intersects, None)

However I get the following error:

Traceback (most recent call last):
    File "<input>", line 1, in <module>
    File "~/.qgis2/python/plugins/processing/tools/general.py", line 71, in runalg
        alg = Processing.runAlgorithm(algOrName, None, *args)
    File "~/.qgis2/python/plugins/processing/core/Processing.py", line 317, in runAlgorithm
        if not param.setValue(args[i]):
    File "~/.qgis2/python/plugins/processing/core/parameters.py", line 829, in setValue
        elif len(value) == 0:
AttributeError: ParameterGeometryPredicate instance has no attribute '__len__'

The algorithm needs the following inputs

processing.alghelp("qgis:extractbylocation")
     ALGORITHM: Extract by location
     INPUT <ParameterVector>
     INTERSECT <ParameterVector>
     PREDICATE <ParameterGeometryPredicate>
     OUTPUT <OutputVector>

ParameterGeometryPredicate is described here: https://fossies.org/dox/qgis-2.8.2/classprocessing_1_1core_1_1parameters_1_1ParameterGeometryPredicate.html

How do I need to initialize a correct instance of the ParameterGeometryPredicate class?


texts: Can I download the entire NRCS SSURGO database?

texts:
I am interested in downloading and using the entire USGS SSURGO database.

The California Soil Resource Lab provides instructions on working with these data in R, but it starts with the assumption that I have downloaded tables as .csv files - but doesn't say how to get these files, other than that the metadata is available as an Access template

Can I download the entire SSURGO database? Is it available in formats other than csv + Access?


texts: Very Inaccurate Geocodes in Canada (Quebec)

texts:
Edit: Found this page detailing the standard geocoding abbreviation types which may be of use when trying to match foreign street types. It appears that we'll have to create a much more robust version.

The problem is as listed in the title: Our current geocoding attempts (using a proprietary list of addressed and a proprietary streets network) in Quebec are coming up with abysmal match rates (40-50% overall). We have tried looking at problems with characters sets and postal code inconsistencies, but nothing seems to be nailing it down. Presuming these two factors are not the issue, what suggestions might be offered as possible pitfalls when geocoding in this area?

Important notes:

 1. Geocoding with ArcGIS 10

 2. Using FSALDUs as post code in "ABC XYZ" format.

Here are some example addresses.

    ARC_Street          ARC_ZIP     CITY                           STATE
1134 BRIDGE ST W        A7R 3L1     NAPANEE                         ON
75 RUE PRINCIPALE       G0A 1A0     SAINTE-CHRISTINE-D'AUVERGNE     QC
43 RUE PRINCIPALE       G0A 1A0     SAINTE-CHRISTINE-D'AUVERGNE     QC
325 RANG SAINT-JACQUES  G0A 1A0     SAINTE-CHRISTINE-D'AUVERGNE     QC
325 SAINT-JACQUES RANG  G0A 1A0     SAINTE-CHRISTINE-D'AUVERGNE     QC
125 RUE PRINCIPALE      G0A 1A0     SAINTE-CHRISTINE-D'AUVERGNE     QC
168 RANG SAINTE-ANNE N  G0A 1A0     SAINTE-CHRISTINE-D'AUVERGNE     QC
168 RANG SAINTE-ANNE    G0A 1A0     SAINTE-CHRISTINE-D'AUVERGNE     QC


texts: Rule Based styling with different forms and colors

texts:
I have a map about different locations of amphibians. I want that the colors illustrates the specific species and the form illustrates the stadium (young, adult). I know how to give different colors according to the species, but I didn't managed to also add different forms.
The problem is similar to the following one Setting symbology type and color based on data using QGIS but the answer there didn't helped me. I think it works by using rule based symbology, but I didn't managed to find out how I have to set my rules.


texts: Earth Engine - UI image selection linked maps, Split Panel

texts:
I'm trying to create a split panel where a linked map loads different images and wipe allows one to see clearly the difference between years.

But I'm stuck with a weird issue when placing the panels containing the image Selectors.

When I just add the left map, it works fine.  But when I add the right map, it takes both panels into the right map.

// Set a center and zoom level.
var center = {lon: -73.97, lat: 40.77, zoom: 13};

// Create two maps.
var leftMap = ui.Map(center);
var rightMap = ui.Map(center);

// Remove UI controls from both maps, but leave zoom control on the left map.
leftMap.setControlVisibility(false);
rightMap.setControlVisibility(false);
leftMap.setControlVisibility({zoomControl: true});

// Create a split panel with the two maps.
var splitPanel = ui.SplitPanel({
  firstPanel: leftMap,
  secondPanel: rightMap,
  orientation: 'horizontal',
  wipe: true
});

// Remove the default map from the root panel.
ui.root.clear();

// Add our split panel to the root panel.
ui.root.add(splitPanel);

////////////////////////
// Make the UI
////////////////////////

//////////
// Set the flow
/////////

var horizontalFlow = ui.Panel.Layout.flow('horizontal')
var verticalFlow = ui.Panel.Layout.flow('vertical')

// Year Select lists

var imageList = ["Image 1","Image 2","Image 3","Image 4","Image 5","Image 6"]

var imageSelect_1 = ui.Select(imageList,"Select Year 1 to Load")
var imageSelect_2 = ui.Select(imageList,"Select Year 2 to Load")

// query labels
var image1Label = ui.Label("Image 1: ")

var query1Panel = ui.Panel([imageSelect_1])

var image1Panel = ui.Panel({
                            widgets:[image1Label,query1Panel],
                            layout:verticalFlow,
                            style: {position: 'top-left'}
                            })

leftMap.add(image1Panel)

/// Image 2 panel
// query labels
var image2Label = ui.Label("Image 2: ")

var query2Panel = ui.Panel([imageSelect_2])

var image2Panel = ui.Panel({widgets:[image2Label,query2Panel],
                            layout:verticalFlow,
                            style: {position: 'top-right'}
                            })

// Error begins here.  if you comment the next line, left panel works
rightMap.add(image2Panel)

var linker = ui.Map.Linker([leftMap, rightMap]);

What am I doing wrong
I built this inspired by the Examples/User Interface/Split Panel


texts: Meaning of dash in house address number in USA?

texts:
I have a set of addresses with a dash between the numbers:

800 - 1100 HORD ST 75202
2400 - 2500 KITTRELL ST 75201

I thought this refers to houses from 800 to 1100 and from 2400 to 2500 on those streets. But there are other examples like these:

136-39 37 AVE
1221-102 CANYON ROCK CT

where that theory doesn't seem to fit.


texts: Error running ArcHydro tools - can't set locations or arguments correctly?

texts:
I am working on trying to do some batch watershed processing.

I have watersheds created using ArcHydro, and was trying to run the Longest Flow Path tool in a batch method, on different feature classes (watershed polygons) located in a feature dataset in a geodatabase, using the code below:

# import stuff
import arcpy, os, traceback, sys, time
import ArcHydroTools

RasterLoc = "C:\Entire_State\Souhegan_River"
VectorLoc = "C:\Entire_State\Souhegan_River\Souhegan_River.gdb"
ArcHydroTools.SetTargetLocations("HydroConfig","Layers",RasterLoc,VectorLoc)
arcpy.CheckOutExtension("Spatial")

arcpy.env.workspace="C://Entire_State//Souhegan_River//Souhegan_River.gdb//Watersheds_Separate"
Flow_Dir_Grid = "C://Entire_State//Souhegan_River//Layers//fdr"

feat_classes = arcpy.ListFeatureClasses()

for feat in feat_classes:

    try:
        in_name = arcpy.env.workspace + "//" +  str(feat)
        out_name = arcpy.env.workspace + "//" +  str(feat) + "_LFP"
        ArcHydroTools.LongestFlowPath(in_name,Flow_Dir_Grid,out_name)

    except:
        name = feat + "_LFP"    # Just here to inspect errors

arcpy.AddMessage("Done generating LFPs.")
arcpy.CheckInExtension("Spatial")

And I get this error:   "The default vector location () cannot be a directory. Modify the current geoprocessing workspace or set the Arc Hydro target locations. Failed to execute (LongestFlowPath)."

I don't really understand this error and how to fix it - would someone advise or provide an example of a working script that access this tool correctly?

Edit:  I have edited the script, to pursue a different path and got farther, with a different set of errors....This way it at least takes a while to fail, so it has me thinking it might be getting further along....?

Script:

# import stuff
import arcpy, os, traceback, sys, time
import ArcHydroTools

arcpy.env.overwriteOutput="TRUE"
arcpy.CheckOutExtension("Spatial")

Flow_Dir_Grid = r"C:/Entire_State/Souhegan_River/Layers/fdr"

RasterLoc = "C:\Entire_State\Souhegan_River\Layers"
VectorLoc = "C:\Entire_State\Souhegan_River\Souhegan_River.gdb"
ArcHydroTools.SetTargetLocations("HydroConfig","Layers",RasterLoc,VectorLoc)

feature_class = r'C:\Entire_State\Souhegan_River\Souhegan_River.gdb\Watershed_Worked_ProjSP83'
arcpy.MakeFeatureLayer_management(feature_class, "watersheds")

with arcpy.da.SearchCursor("watersheds", '*') as cursor:
    for row in cursor:
        try:
            id = row[3].encode('utf-8')
            arcpy.SelectLayerByAttribute_management("watersheds", "CLEAR_SELECTION")
            where_clause = """"{}" = {}""".format("HydroID", id)
            arcpy.SelectLayerByAttribute_management("watersheds", "NEW_SELECTION", where_clause)

            out_name = 'C:\Entire_State\Souhegan_River\Souhegan_River.gdb\Watersheds_Separate\\'+ str(row[3]) + "_LFP"

            ArcHydroTools.LongestFlowPath("watersheds",Flow_Dir_Grid,out_name)

        except:
            print ( "error for: " + str(row[3]) + "_LFP")

arcpy.AddMessage("Done generating LFPs.")
arcpy.CheckInExtension("Spatial")

Error:

"System.Runtime.InteropServices.COMException (0x80004005): Error HRESULT E_FAIL has been returned from a call to a COM component.
   at ESRI.ArcGIS.Geoprocessing.GeoProcessorClass.Execute(String Name, IVariantArray ipValues, ITrackCancel pTrackCancel)
   at ESRI.APWR.ArcHydroGPTools.ApConfigGP.CreateFeatureClassGP(String path, String layerName, esriGeometryType geomType, String& exMessage, String templateName, String hasZ, String hasM, ISpatialReference spatialRef, Boolean addtoMap) in C:\Builds\HydroSolutions\10.5_ArcHydro\Sources\ArcHydroTools\src\ArcHydroGPTools\ApConfigGP.vb:line 4253
Failed to execute (LongestFlowPath).

"


texts: When I add a grid to my map in print composer, the dotted line is rendered "on top of" my map labels making them hard to read. Is there a way to make the labels appear on top of the grid lines?


texts:
I want to load Shapefile in Flexviewer 2.2

Any sample would help a lot.


texts: Reading Esri File Geodatabase (.gdb) from .NET?

texts:
I have been reading various threads on reading (not writing) an ESRI file geodatabase (.gdb). What I need to do is process a feature class in a geodb and extract the vertices of each feature. I do not need a visual component as this will be just pure processing of data. Looking for an open-source framework I can readily use in my .NET project. So far I have looked at

DotSpatial: Libraries work great for shapefiles but currently no ESRI gdb support
MapWindow: Looks like it supports MSSQL, SQLLite, and Postgre geodbs but not ESRI geodbs
QGIS: Does have the capability to read ESRI geodbs but the API is currently only for Python and C++
ArcGIS: Using ArcObjects this is quite easy but would like this available for users that do not require checking out a license

None of these have worked for me.

How do I read a feature class from an ESRI file geodatabase purely for processing the records?


texts: From track+waypoints to elevation profile with waypoints

texts:
I have a track and waypoints (1) imported in QGIS. I want to calculate the closest point from every waypoint to the track (2). Then, I want a table with data that allows me to draw an elevation profile overplotting the waypoints in R. I should generate the following data:

Elevation profile from the track (already achieved using the complement "Profile Tool").
Data for waypoints which should include three columns from data calculated from the intersection at the closest point with the track: waypoint name, distance to the beginning of the the track, elevation.


texts: For a SDE FeatureClasses list all Roles that have been Granted any Privileges on it, and which Privileges each Role have on it

texts:
I have searched through the ArcObjects 10 SDK API and googled, but havnt been able to find what I am looking for, so thats the reason for me posting this question.

Q: How can I for a SDE FeatureClasses list all Roles that have been Granted any Privileges on it, and which Privileges each Role have on it?

The equivalent in normal Oracle would be something like this:

SELECT dtp.PRIVILEGE,
       dtp.grantable ADMIN,
       dtp.grantee,
       dtp.grantor,
       dbo.owner,
       dbo.object_type,
       dbo.object_name,
       dbo.status,
       '' column_name
  FROM ALL_TAB_PRIVS dtp, ALL_OBJECTS dbo
 WHERE dtp.TABLE_SCHEMA = dbo.owner AND dtp.table_name = dbo.object_name
       AND dbo.object_type IN ( 'TABLE', 'VIEW' )
UNION
SELECT dcp.PRIVILEGE,
       dcp.grantable ADMIN,
       dcp.grantee,
       dcp.grantor,
       dbo.owner,
       dbo.object_type,
       dbo.object_name,
       dbo.status,
       dcp.column_name
  FROM ALL_COL_PRIVS dcp, ALL_OBJECTS dbo
 WHERE dcp.TABLE_SCHEMA = dbo.owner AND dcp.table_name = dbo.object_name
       AND dbo.object_type IN ( 'TABLE', 'VIEW' )

The query would list out all object of the types in the list provided (in the example I have listed TABLE and VIEW) which privileges have been granted on them like SELECT, UPDATE etc... and to which Roles.

Related Q: SDE FeatureClass Granting and Revoking Privilege with ArcObjects


texts: Download OpenStreetMap data -> connection refused

texts:
When I try to download an OSM by means of the QGIS OpenStreetMap plug-in I always get:

Download failed. Connection refused.

When I try the same with the OSM downloader plug-in I get at least an error code:10061.

It doesn't matter if I select canvas, layer or manual, the result is always the same.

running QGIS 2.16.1, Windows 10 64bit


texts: Download hourly relative humidity using GLDAS using Google Earth Engine

texts: Does anyone know how I can download hourly GLDAS relative humidity (%) using Google Earth Engine? The only output that I am able to get is specific humidity at a specific location, and am having trouble converting this. Is there a way of just obtaining RH straight from GLDAS?


texts: Replace values in raster 1 with values from raster 2 in ArcMap 10.7.1

texts:
Raster 1 is a DEM of a study site.
Raster 2 is a smaller DEM of the corrected water depth of a river that runs through the site.
How do I replace the pixels within raster 1 with the corrected water depth pixels (raster 2)
They are both of the same resolution (0.2m)


texts: ArcGIS Online - floating pop-ups in web map

texts:
Is there a way to make pop-ups in ArcGIS Online floating or draggable?

In case of long text in pop-up, moving with a whole map is for map users impractical.


texts: Error with ROI Semi Automatic Classification Plugin QGIS 3.4

texts:
I have troubles for making ROI with SCP plugin in QGIS 3.4.
I dowloaded some landsat images with the plugin and trying to create ROI. To do it, I made a Bandset with DOS1 atmospheric correction, created a training file for making the ROI and display a 4-3-2 images with Landsat 7 corrected images. When I draw my polygon for ROI there's nothing anormal, but when I want to save those ROI, I get an error with those two messages :

Error [31]: Error calculating signature. Possibly ROI is too small

Information [9]: Data projections do not match. Reproject data to the same projection

My ROI aren't too small and all the bands and the training files are in the same projection so I don't know what's happening.


texts: Trying to perform A WFS CQL_FILTER for multiple TYPENAMES with different GEOM name identifiers

texts:
Using a WFS service https://geo.onroerenderfgoed.be/geoserver/wfs I want to retrieve the features with a CQL INTERSECT query.

I need to get results from 2 TYPENAME:
vioe_geoportaal:vast_htp,vioe_geoportaal:vast_dibe_relict

The problem is that the CQL query interacts with the geomtry object. In these 2 typenames they are defined different geom and geometry.

Request 1

Where the object is called geom is succesful according :

https://geo.onroerenderfgoed.be/geoserver/wfs?SERVICE=WFS&REQUEST=GetFeature&VERSION=1.1.0&CRS=EPSG:31370&OUTPUTFORMAT=application/json&TYPENAME=vioe_geoportaal:vast_htp&CQL_FILTER=INTERSECTS(geom, POLYGON((123664.03100000322 210489.5150000006,123705.24599999934 210489.5150000006,123705.24599999934 210587.86800000072,123664.03100000322 210587.86800000072,123664.03100000322 210489.5150000006)))

Request 2

Where the object is called geometry is succesful according :

https://geo.onroerenderfgoed.be/geoserver/wfs?SERVICE=WFS&REQUEST=GetFeature&VERSION=1.1.0&CRS=EPSG:31370&OUTPUTFORMAT=application/json&TYPENAME=vioe_geoportaal:vast_dibe_relict&CQL_FILTER=INTERSECTS(geometry, POLYGON((123664.03100000322 210489.5150000006,123705.24599999934 210489.5150000006,123705.24599999934 210587.86800000072,123664.03100000322 210587.86800000072,123664.03100000322 210489.5150000006)))

Now I want 1 request to encapsulate the 2 requests above.

I've tried concatenating them - but as a logical result of the difference in the name geom vs geometry

https://geo.onroerenderfgoed.be/geoserver/wfs?SERVICE=WFS&REQUEST=GetFeature&VERSION=1.1.0&CRS=EPSG:31370&OUTPUTFORMAT=application/json&TYPENAME=vioe_geoportaal:vast_htp,vioe_geoportaal:vast_dibe_relict&CQL_FILTER=INTERSECTS(geometry, POLYGON((123664.03100000322 210489.5150000006,123705.24599999934 210489.5150000006,123705.24599999934 210587.86800000072,123664.03100000322 210587.86800000072,123664.03100000322 210489.5150000006)))

But raises a logical error :

<ows:ExceptionReport xmlns:xs="http://www.w3.org/2001/XMLSchema" xmlns:ows="http://www.opengis.net/ows" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" version="1.0.0" xsi:schemaLocation="http://www.opengis.net/ows https://geo.onroerenderfgoed.be:443/geoserver/schemas/ows/1.0.0/owsExceptionReport.xsd">
    <ows:Exception exceptionCode="InvalidParameterValue">
        <ows:ExceptionText>Illegal property name: geometry</ows:ExceptionText>
    </ows:Exception>
</ows:ExceptionReport>

What can I do instead of launching 2 request?


texts: Attribute Indexes ArcSDE 9.3.1 for Attribute Joins on GUID field types

texts:
I create my joins on primary/foreign keys that are guids...I don't seem to be able to create an index on this field in my target table. After joining to this table it takes 4 minutes to open the source table.

Can indexes be built on guids fields..if not what work arounds do you guys see? Thanks.


texts: Performing QGIS modeler multiple raster calculations with data from external CSV?

texts:
I need to perform multiple calculations on a set of rasters based on values in a CSV metadata file.

I'm using the QGIS graphical modeler and can call a specific field in a CSV with the Table and Table Field inputs.

Since I have multiple rasters, and each raster has a corresponding value in the CSV, how do I call the specific row? (similar to a lookup for instance).

I don't want to add a script if I don't have to for ease of use.


texts:
INTRODUCTION: I have 5000 mammal species and I'm studying Extreme Climate Event impacts on them. I'm now analysing warm spell duration index (WSDI) which is defined as number of days in a year where Temp>90th percentile of the historical (1961-1990 in my case). I'm modifying the index calculating the 90th percentile of the historical within the range of the species, and not limited to the "temporal evolution" of the single cell, so that I have for every mammal a WSDI which is proper to that mammal.
PROBLEM : The whole problem of WSDI modification works, and it is not the topic of this question. The problem is that it takes really long time (1 month). I was suggested to use parallel, running at the same time 10x500 (or even less) instead of 5000 species all togheter. The problem is that for the process of modification I have to change region at every step, at every species, using also r.mask to calculate 90th percentile within range. Changing region during parallel seems not possible from the guide here: https://grasswiki.osgeo.org/wiki/Parallel_GRASS_jobs.
QUESTION: Is there a way to avoid that? I read on  the guide that something similar is possible saving the region at every time using WIND_OVERRIDE, is that exact? Otherwise, I maybe can use only r.mask without using g.region, so that the problem does not appear?


texts: Finding a certain length of a line from a central point

texts:
I am not exactly sure on how to describe what I am trying to do so hopefully this picture will help. What I am trying to figure out: I have a point and from that point I have, walking paths from that point. Say I want to know how far a half mile is on the paths, Is there a a method or tool that will aid me in this? By just doing a 1/2 mile buffer on the point does not achieve what I want since the paths are not straight lines.
I hope what I trying to do is clear enough


texts:
I'm stuck trying to get metadata for vector-layers into a PostgreSQL database with PostGIS-Extension. The ESRI-Environment allows to define metadata (ISO, FGDC, INSPIRE,...) and saves this information together with the features e.g. in Oracle, FGDB, Shape/XML,...). QGIS provides metadata-properties for layers, but these don't find their way into PG.

Does someone know a solution to save metadata together with the features in PG (based on the common standards)?


texts: How to write to a single table via Python multiprocessing

texts:
A comment in the post Can multiprocessing with arcpy be run in a script tool? got me thinking, as I often need to do exactly this:

  Just beware of deadlocking situations (two Insert cursors in the same
  table for instance)

My question is, how can you write to a single table when using multiprocessing?

Here's an example script, which iterates through the sample City layer, and uses multiprocessing to copy each city's values to an output table via an insertCursor (this simplifies the more complicated scenario I have in mind).

# Testing how to write to an output fGDB table via multiple threads
import os, sys, arcpy, multiprocessing
from multiprocessing import Process, Queue

def Worker(input, output):
    for func in iter(input.get, 'STOP'):
        inputs = func
        doProcess(inputs)

def doProcess(inputs):
    outTblName = inputs[0]
    city = inputs[1]
    pop = inputs[2]
    with arcpy.da.InsertCursor(os.path.join(arcpy.env.scratchGDB, outTblName), ["NAME", "POPULATION"]) as iCursor:
        try:
            iCursor.insertRow([city, pop])
        except:
            print("Problem inserting " + city + " : " + str(pop) + " : trying again" )
            doProcess(inputs)

if __name__ == '__main__':
    NUMBER_OF_PROCESSES = 8
    task_queue = Queue()
    done_queue = Queue()

    inFC = "C:\Program Files (x86)\ArcGIS\Desktop10.2\TemplateData\TemplateData.gdb\World\City"
    outTblName = "testTable"

    #Create the empty table
    outTable = os.path.join(arcpy.env.scratchGDB, outTblName)
    if(arcpy.Exists(outTable)):
       arcpy.Delete_management(outTable)
    arcpy.CreateTable_management(arcpy.env.scratchGDB, outTblName)
    arcpy.AddField_management(outTable, "Name", "TEXT")
    arcpy.AddField_management(outTable, "Population", "DOUBLE")

    #Iterate through the cities. Send each one to the multiprocessor
    with arcpy.da.SearchCursor(inFC, ["NAME", "POPULATION"]) as sCursor:
        for city in sCursor:
            cityName = city[0]
            pop = city[1]
            task_queue.put([outTblName, cityName, pop])

    for i in range(NUMBER_OF_PROCESSES):
        Process(target=Worker, args=(task_queue, done_queue)).start()

    for i in range(NUMBER_OF_PROCESSES):
        task_queue.put('STOP')

As expected it runs into problems when multiple threads try to access the output table simultaneously. Even when calling the doProcess function again recursively after errors are detected, the output table contains fewer rows than the input table.

An idea is for each thread to create a new table, and to append them all at the end. Are there any best-practise suggestions?


texts:
By what I understand, a "standard" slc from a SAR has complex values for each pixel. Here the phasor consists of the magnitude and the relative phase. What I don't understand is the complex nature of an interferogram and how this relates to the coherence.
Am I right in that in an interferogram the complex value consists of the phase-difference between the two SLCs and the coherence or this working completely different?


texts: Carto animation only showing data points in 1 min intervals

texts:
I'm trying to animate location data points, and instead of showing each point at its correct time, Carto seems to show a group of points every 1 minute (1 minute in the visualization, not real life).

My duration is set at 10,000 since this animation takes place over a 24 hour period and I don't want the objects in the animation to move too fast to be able to understand what's happening.
I've tested with the duration set at 5,000 and it's still a 1 minute difference.

Does anyone know if there's anything in Carto that would cause this issue?

Video of glitch

Here's my viz code:

        @duration: 10000
        @animation: animation(@timeSteps, @duration, fade(0, 2))
        @timeSteps: linear($time, time('2020-01-01T00:00:00Z'), time('2020-01-02T00:15:16Z'))
        color: ramp(@timeSteps, SunsetDark)
        filter: @animation
        width: 5
        strokeWidth: 1


texts: Concatenating fields using VBscript in ArcGIS Field Calculator?

texts:
I'm trying to get one column of text to move into another column that already has text.  I DON'T want to replace it, I want the transferring text to be added to the already existing text in a shp file attribute table.
Example
permit   acct id
123      567
324      456

I want to take the "acct id" numbers and bring them into the permit column and have it look like this...
permit
123567
324456


texts: Add the lengths of columns with the same id in arcmap

texts:
I have a polyline dataset that contains a unique feature (manhole) and the lines segmented by outflows from it. I need to get the distance from the manhole to each break on the line.

So in the MH Column I want

1 = 2.029
2 = 1+2 = 37.19~
3 = 1+2+3 = 58.3~
4 = 1+2+3+4 = 88.6~
a =20.31.62
b = a+b = 56

It doesn't seem as if I can use summary stats.

I can get the MH_Order to be MH_ID+SequentialNo. So it's SSG015-06_1, SSG015-06_2 etc instead of SSG015-06_676, SSG015-06_679 etc.

But the critical factor is how to I get it to do 1=1, 2=1+2, 3=1+2+3 etc? should be in field calculator python


texts:
I'm using QGIS 3.12's Raster -> Analysis -> Hillshade tool on a fairly coarse DEM provided by the British Columbia provincial government.

The DEM has elevation values between 500m and 2650m, so a range of around 2150m. The image below shows the output in my area of interest. On the mountain (left) side the hillshading mostly works as expected. However, on some parts of the mountain, and in the neighbouring valley (right), there is noise and some striping. This occurs in areas where cells differ by 1 or 2 metres, some of which may be real-world variations while others may be anomalies / errors in the data.

The following image is zoomed to an area of striping. From left to right elevation drops by only 5 metres and the slope is consistent, yet the hillshading algorithm generates a lot of noise.

Is there any way I can limit the hillshading output to only 'significant' elevation change? E.g. if elevation across the DEM changes by 2000m then only show changes where elevation changes by at least 20m?

I tried changing the Z factor up and down and switching between ZevenbergenThorne and Horn's formula but neither approach made a significant difference.


texts: Convert osm.bz2 file to shapefile for specific extent

texts:
I'm trying to get a historic shapefile (from 2010) of a road network for a specific country from OSM. I downloaded a planet osm.bz2 file from here
I'm wondering how to create a shapefile from this for a specific country? For example, is there a way to: (1) crop the osm.bz2 file to an extent then (2) convert to shapefile?
Just mentioning cropping the osm.bz2 before converting to shapefile, as I imagine converting the entire osm.bz2 file to shapefile would create a fairly large file. But good to know if this can't be avoided. Or if there's a better way to do this than download the whole osm.bz2 file good to know as well!


texts: How can I have the ground sampling distance from the images?

texts: I have some aerial images and I need to know how if I can have the ground sampling distance value from the image without calculate it from the plane fly's parameters.


texts: Using Leaflet.VectorGrid with local GeoJSON and addLayer

texts:
I am trying to test out Leaflet.VectorGrid from https://github.com/Leaflet/Leaflet.VectorGrid. My overall goal is to load a large 20MB or larger geojson quickly.

I am getting this error...

"Uncaught (in promise) TypeError: Cannot read property 'addLayer' of undefined
    at NewClass.addTo (leaflet-src.js?e11e:6559)

I am following this example https://github.com/Leaflet/Leaflet.VectorGrid/blob/master/docs/demo-geojson.html

import VectorTile from "leaflet.vectorgrid/dist/Leaflet.VectorGrid.js";
import Pbf from "leaflet.vectorgrid/dist/Leaflet.VectorGrid.js";
import euCountries from "../eu-countries.js";
 var highlight;
  var clearHighlight = function() {
    if (highlight) {
      vectorGrid.resetFeatureStyle(highlight);
    }
    highlight = null;
  };
  var vectorGrid = L.vectorGrid
    .slicer(euCountries.euCountries, {
      rendererFactory: L.svg.tile,
      vectorTileLayerStyles: {
        sliced: function(properties, zoom) {
          var p = properties.mapcolor7 % 5;
          return {
            fillColor:
              p === 0
                ? "#800026"
                : p === 1
                ? "#E31A1C"
                : p === 2
                ? "#FEB24C"
                : p === 3
                ? "#B2FE4C"
                : "#FFEDA0",
            fillOpacity: 0.5,
            //fillOpacity: 1,
            stroke: true,
            fill: true,
            color: "black",
            //opacity: 0.2,
            weight: 0
          };
        }
      },
      interactive: true,
      getFeatureId: function(f) {
        return f.properties.wb_a3;
      }
    })
    .on("mouseover", function(e) {
      var properties = e.layer.properties;
      L.popup()
        .setContent(properties.name || properties.type)
        .setLatLng(e.latlng)
        .openOn(myMap);

      clearHighlight();
      highlight = properties.wb_a3;

      var p = properties.mapcolor7 % 5;
      var style = {
        fillColor:
          p === 0
            ? "#800026"
            : p === 1
            ? "#E31A1C"
            : p === 2
            ? "#FEB24C"
            : p === 3
            ? "#B2FE4C"
            : "#FFEDA0",
        fillOpacity: 0.5,
        stroke: true,
        fill: true,
        color: "red",
        opacity: 1,
        weight: 2
      };

      vectorGrid.setFeatureStyle(properties.wb_a3, style);
    })
    .addTo(myMap);

I know I'm importing the geojson properly because I can console.log (Object.keys(euCountries.euCountries) and get the expected result {type, features} prints in the console.

My eu-countries.js looks like this...

    export default {
        euCountries: {
            "type": "FeatureCollection", "features": [
                {..all the features are here..}
]}}

What am I missing here?


texts: Polygon zone fell outside bounds of raster using ArcPy

texts:
I'm trying to perform zonal statistics, and get a message that:

  A polygon zone fell outside the bounds of the raster.

It is not possible because the polygon is just a 50-km buffer zone around the city, and the input raster is global.

Here is the part of the code:

 temp_table = (self.temp_dir + os.sep + 'zone_{0}.dbf'
                              .format(str(row[0]).zfill(2)))
 sa.ZonalStatisticsAsTable(out_layer, 'name', global_ras,
                                          temp_table, 'DATA', 'MEAN')


texts: OpenStreetMap/GeoServer fill polygon color on user interaction

texts:
I'm trying to build a system to show census tract on map. I used this example to import the census tract TIGER shapefile into DB, flatten the data on the column I need and managed to show a map with census tract using GeoServer, OpenLayers and OSM.

Now I need to add interaction where a user can select a census tract and be able to change the fill color. The example I followed creates a ol.Layer.Image layer which, I think, won't work here. I tried using a vector layer with a GeoJSON source from GeoServer but the browser is crashing trying to load 74K rows through JSON. I'm looking at vector tile example here but it doesn't cover dynamic style change from user interaction.

Is it possible to use OpenLayers/GeoServer with what I'm trying to do?

I'm new to these technologies and this is my first experiment with GIS.


texts: How to cast geom to geog in PostGis to find points within a given distance

texts:
I have a Geo table with geometry column storing point data SRD=4326

I need to find all points within a given distance in meters of a  Longitude, Latitude

Do I create a geography column from the geometry column or the Longitude, Latitude stored in the table

or

Do the conversion on the fly within the query

Could someone suggest how to perform both options and there pro & cons?


texts: QGIS: Joining data from multiple road lines onto polygons by location

texts:
I'm looking to add traffic data linked to a polyline shapefile of road centerlines to a multipolygon shapefile of land parcels. I want to join this data using proximity to roads instead of by matching addresses to create the most accurate traffic data (i.e. if a parcel is addressed on a side road but is cornered on a major thoroughfare). Inevitably, many parcels will have multiple traffic data inputs since one parcel can be near multiple roads. Many road lines will also be sharing information with multiple parcels.

I created a buffer for the road centerline to overlap with parcels nearby.  I then used the "Join attributes by location" to input traffic data onto the parcel shapefile.  I selected the "one to many" feature to get all relevant traffic info onto each parcel polygon.

The join created several superimposed polygons for each parcel, each with different traffic data from a different individual neighboring road (which I feared would happen). This functionally works for my purposes but is very inefficient. The file is multiple times the size it needs to be and is delayed in loading with each pan/zoom action.

How can I converge multiple traffic data points into one polygon, instead of assigning one polygon per traffic data point?

One acceptable, but less than ideal, solution will be to only accept the polygon with the greatest traffic flow, but I'm unsure of how to filter and eliminate polygons in that way.


texts: postgis problem with shortest distance calculation

texts:
while working with POSTGIS pgrouting, for calculateing the distance between two roads(lines) i got the shortest_path function.

But the logic is based on Start_point(Start_id) and end_point(end_id) but in my data the linestring contains so many internal points like ('linestring(1 1,2 2,3 3,4 4,5 5)' just for example..)

it is taking start point (1 1) endpoint(5 5)

if other line starting with (5 5) it is showing as route...like ('linestring(5 5,6 6)')

But line which crossing the point inside the linestring like(2 2,3 3,4 4) which is not telling as connected.. example

table roads: id name way 1 A linestring(1 1,2 2,3 3,4 4,5 5) 2 B linestring(5 5,6 6) 3 c linestring(2 1,2 2,2 3)

if i am applying shortest_path function from point(1 1) to (6 6) its showing the way but for (1 1) to (2 3) it is not showing anything...but there is a route for this (1 1,2 2,2 3)

can anyone please help me out for finding the solution..

Regards
Deepak M


texts: I'm new to GIS development, and am looking for a mapping solution for a web app. The map will be fairly complex with 4-5 different layers. I'm looking to convert shapefiles into a usable geo-coded format, to be served from a few different layer specific API's.


texts: Finding point in country furthest from boundary?

texts:
i'm trying to find the point in south africa furthest from the sea or a landward boundary.

my method thus far:

i obtained a first guess at the centre of the inscribed circle and the three segments of boundary this circle might touch using a compass and 1:2 500 000 map.
i then clipped these three segments of boundary from a data set of boundary coordinates supplied by national geo-spatial information (WGS84; decimal degrees; 9 decimal places; roughly one point for every 40 metres of boundary).
i converted the boundary coordinates to radians and used cartesian equations (distance between two points and centre of a circle through three points) to iteratively improve the centre of my circle and the three boundary points it touches. (obviously, i'm ignoring curvature; i'm effectively finding the centre of the base of the spherical cap defined by the inscribed circle -- i.e. a point approximately 16km underground.)

my questions:

is it reasonable to plug lat-long coordinates expressed in radians into cartesian equations to find the centre of the base of the spherical cap defined by the inscribed circle (given that the surface distance represented by one radian varies with latitude, and that the base of the spherical cap in question is roughly 800km in diameter and spans latitude 30S)? what order of error might be involved here?
is it reasonable to assume that the three boundary points tangential to the base of the spherical cap will likewise be the three boundary points closest to the actual point in the country furthest from a boundary?
i have the feeling i'm not even asking the right questions here. i guess there must be a better way to find the point in a country furthest from a boundary.


texts:
I've been looking for a way to enter a unique field name into an attribute table for up to 40 feature classes. I need the field names to be meaningful, not just 'field1', 'field2' etc. For example, if the feature class represents wetlands I would like a field called 'Wetland'.

I thought the easiest way of doing this would be to use TKinter. I can print the feature class name and then dynamically enter a field name with the TKinter data entry msgbox. Then keep looping through the feature classes until they all have a field with a unique name.

I've found some examples but am unsure about how to use data entry text as a variable in Python. If anyone out there had an example script or can make suggestions on this script I'd appreciate it:

import arcpy
from arcpy import env
from Tkinter import *

# Environment Workspace
env.workspace = r'G:\TEST\Python_Test\Test.gdb'

def receive():
    text = E1.get()
    print (text)
top = Tk()
L1 = Label(top, text="User Name")
L1.pack( side = LEFT)
E1 = Entry(top, bd =5)
E1.pack(side = RIGHT)
b = Button(top, text="OK", width=10, command=receive)
b.pack()
top.mainloop()


texts: Using Openlayers for GeoRSS from remote website

texts:
I am trying to get features from remote GeoRSS file with Openlayers... But I am stuck. Openlayers can not read remote file. I have searched on the net but I could not find any solution... Here is the code:

map = new OpenLayers.Map('map');

var wms = new OpenLayers.Layer.WMS("Osgeo WMS",
    "http://vmap0.tiles.osgeo.org/wms/vmap0", {
    layers: "basic"
});
map.addLayer(wms);

var earthquakes = new OpenLayers.Layer.Vector("T", {
    strategies: [new OpenLayers.Strategy.Fixed()],
    protocol: new OpenLayers.Protocol.HTTP({
        url: "http://brainoff.com/test/rss.xml",
        format: new OpenLayers.Format.GeoRSS()
    }),
});
map.addLayer(earthquakes);

var select = new OpenLayers.Control.SelectFeature([earthquakes]);
map.addControl(select);
select.activate();

var draw = new OpenLayers.Control.DrawFeature(
    earthquakes, OpenLayers.Handler.Point);
map.addControl(draw);

function toggle() {
    if (document.getElementById("toggle-id").checked) {
        draw.activate();
    } else {
        draw.deactivate();
    }
}
document.getElementById("toggle-id").onclick = toggle;
toggle();
map.addControl(new OpenLayers.Control.LayerSwitcher());
map.zoomToMaxExtent();


texts:
I got a tiny script:

def calc(value):
 if value is '123':
    return 'ok'
 elif value is '12 3':
    return 'not_ok'
 else:
    return value

if the value is without spaces it works perfectly, but if there is a space in it, it just returns the value.
Where is my error? :)


texts: Clarification on a previous Utilization Distribution Overlap Index (UDOI) answer

texts:
I'm really struggling to calculate the Utilization Distribution Overlap Index (UDOI) and could really use some clarification on a previous answer.

There is a comprehensive guide to what i am trying to achieve here:

How to Calculate Utilization Distribution Overlap Index (UDOI) in ArcGIS Desktop?

I have followed the instructions through to the calculation of the integrand but I am really stuck trying to calculate the zonal sums. This is where i have got to :

"This is the integrand. To integrate it, notice that because the grid you are using has rectangles all the same size, dx dy never changes. That lets us factor it out, leaving us to just do the sum. That is,compute the zonal sum of the product grid using the polygon of intersection as the zone.This gets multiplied by dx (the x cellsize) and by dy (the y cellsize, usually equal to the x cellsize). For the example, I obtain a value of 0.183814."

First can i clarify that the link in point 5 of Whuber's summary at the bottom of the answer is for the correct tool? The link takes you to the Focal Stats help when i would have thought that it should go to the Zonal Stats in the Zonal toolbox?
In Whuber's example his first zonal sum comes to 0.183814. Mine is coming to over 9000 so I am assuming that I am doing something wrong. Could someone confirm that i am looking for a figure in the 0.something order of mag?
Where it says "This gets multiplied by dx (the x cellsize) and by dy (the y cellsize, usually equal to the x cellsize". Is this explaining what the calc is doing or giving me an instruction?

I want to clarify that I am not criticizing Whuber's post as I have found it very useful.I just don't understand it enough to complete the calculation.


texts: How to animate a single raster layer with variable substitution SLD Geoserver?

texts:
I was trying to animate a single raster layer with a variable substitution in SLD. I've followed the Geoserver documentation http://docs.geoserver.org/stable/en/user/styling/sld/extensions/substitution.html#sld-variable-substitution
How can I animate a single layer by using a variable substitution on Geoserver Style? This is part of the style I was using:

<?xml version="1.0" encoding="UTF-8"?>
<StyledLayerDescriptorxmlns="http://www.opengis.net/sld" xmlns:ogc="http://www.opengis.net/ogc"mlns:se="http://www.opengis.net/se"xmlns:xlink="http://www.w3.org/1999/xlink"xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" version="1.0.0"si:schemaLocation="http://www.opengis.net/sldhttp://schemas.opengis.net/sld/1.0.0/StyledLayerDescriptor.xsd">
<NamedLayer>
  <Name>Umidade Relativa Mínima</Name>
  <UserStyle>
     <Name>Style</Name>
     <Title>Style for a Coverage layer.</Title>
     <FeatureTypeStyle>
        <Rule>
           <RasterSymbolizer>
              <ChannelSelection>
                 <GrayChannel>
                    <SourceChannelName>1</SourceChannelName>
                 </GrayChannel>
              </ChannelSelection>
              <ColorMap type="ramp">
                 <ColorMapEntry color="#ff0000" quantity="0.60" opacity="${env('op',0)}" label="20" />
                 <ColorMapEntry color="#00ff00" quantity="0.70" opacity="1.0" label="30" />
              </ColorMap>
           </RasterSymbolizer>
        </Rule>
     </FeatureTypeStyle>
  </UserStyle>

and this is the query I was using to try to animate:
http://localhost:8080/geoserver/wms/animate?&format=image/gif;subtype=animated&format_options=gif_loop_continuosly:true&layers=terrama2_4:view4&aparam=op&avalues=0;0.2;0.5;1


texts: What UTM like coordinate system is this?

texts:
I have this set of coordinates: 4622580.00000, 503749.00000

They should roughly correlate with 51°01'39.4"N 13°44'44.2"E in Dresden, Germany.

At first I thought they might be UTM, but that puts me off the coast of Somalia. Does anybody have an idea? Or maybe they are UTM but with special options? I am unfortunately not very knowledgeable with coordinate systems.


texts:
I'm struggling to get the layout to display normally. In the past, the layout has displayed the paper with an adjustable data frame which I can resize freely. Currently, when I switch to layout, the base map continues to occupy the whole screen, and there's no adjustable border and no sign of the actual paper (see image). Does anyone know what's going on here?

Edit: I'm using ArcMap 10.5.1 I've explained more in the comments. I can still export the layout map, and the exported image appears correctly, but the actual layout view is unchanged from the image.


texts: Getting XY screen coords from ArcScene 3D point using ArcObjects?

texts:
I already know how to convert screen point (X,Y in pixels) to 3D point:

ISceneGraph pSceneGraph = m_pSceneHookHelper.SceneGraph;
IPoint pPoint;
pSceneGraph.Locate(pSceneGraph.ActiveViewer, X, Y, esriScenePickMode.esriScenePickFastLocation, true, out pPoint, out pOwner, out pObject);

But how may I get X,Y screen point from 3D point?
I just use snapping to change 3D point coordinates, so I need to know new screen coordinates to use Locate function (because I need to find feature object).


texts: Geojson multipolygon doesn't work with Leaflet bindPopup?

texts:
I have polygon and multipolygon data in Geojson format and when I try to bind popup for polygon, it works perfectly, and when I try to bind it for multipolygon, it doesn't work.

Is that some Leaflet bug and how can I solve this?

Thanks


texts: Automating georeferencing of aerial photos?

texts:
I have tried using the "automatically place Passpoints" tool on the georeferencing toolbar in ArcGIS 10.3.1 but this only seems to work with Aerial photos which are near to identical. I have many colour infrared aerial photos and normal RGB aerial photos for the same areas which were taken 5 years apart from each other. This tools hardly works at all for these - either because of the colour difference or the fact that there are many discrepencies due to the time interval.

Does anyone know of other software which automates georeferencing? Or is there are way in which I could make it work better in ArcGIS 10.3.1?


texts: Reading in local GeoJSON file for OpenLayers

texts:
I am kinda new to OpenLayers and not (yet) very experienced in JavaScript.

I have a WordPress blog and I want to replicate and adjust the following script to my needs: http://jsfiddle.net/eqsrnjos/1/

With one big exception: I need to load the GeoJSON file from the server instead of defining it in the top (so replace the var geoJson).

However the most straight forward way does not seem to work for me:

function getData() {
return jQuery.getJSON("/json/studyCentroids.json");
}

var features = new ol.format.GeoJSON().readFeatures( getData() );

In Firebug I can see that the GeoJSON file is properly loaded, but it does not seem to appear for some reason?


texts: Simplifying a SpatialPointsDataFrame in ggplot world

texts:
I've been working on a map showing lines running across Europe. Each line is coming from a GPX file, parsed in like so:

spdf <- readOGR("data/file.gpx", layer = "track_points")

Now, if I plot the map using plot(), I can simplify each of these lines without any issue:

spdf_simple <- gSimplify(spdf, tol = 0.1)

However, in this particular instance I am using ggplot2 and need to keep a dataframe structure, as I'm doing some filtering (see below):

# load each file
for (file in list.files('data/') {
    # parse each file
    track_points <- readOGR(paste("data/", file, sep = ""), layer = "track_points")
    track_points@data$time_clean <- ymd_hms(track_points@data$time)

    # we need to preserve @data to filter the time value later on
    df1 <- data.frame(track_points@data)
    df2 <- data.frame(track_points@coords)
    df <- bind_cols(df1, df2)

    # filter out all records that happened after FINISH time
    df.filter <- df %>%
      filter(time_clean < FINISH) %>%
      select(time, coords.x1, coords.x2) %>%
      mutate(group = file)
    colnames(df.filter) <- c('time', 'lon', 'lat', 'group')

    paths <- rbind(paths, df.filter)
}

Subsequently I can plot my paths like so:

ggplot() +
    geom_path(data = paths, aes(lon, lat, group = group))

Now the result of gSimplify is a sp::SpatialPoints object which doesn't preserve what's in @data.

I'd like to simplify these paths, although I'm at a loss as to where to do it in ggplot world.


texts: Performing tag change detection in OpenStreetMap?

texts:
I'd like to have a tool for monitoring changes in OpenStreetMap (OSM) for tags only, not geometry. For instance, I'm interested in detecting any tag changes regarding one-ways in streets. This can be the 1) suppression of the tag oneway=*, 2) the modification of the value of the key (oneway:yes -> oneway:no), 3) the addition of a tag oneway=* to a highway, or 4) the creation of a new object tagged with oneway=* (the latter including the creation of a new geometry).

For now, I'm using a overpass query using the adiffargument, as e.g.:

[out:json][timeout:250][adiff:"2017-01-01T00:00:00Z","2017-09-01T00:00:00Z"];
{{geocodeArea:Ixelles}}->.searchArea;
(
    way["oneway"](area.searchArea);
);
out meta;
>;
out skel qt;

This query shows me all the changes associated to any objects tagged with oneway=* in my search area between the 2 dates. The problem is that I have a lot of false positives which are not strictly related to the oneway tags: addition of other tags (eg; surface=*), slight changes in the geometry like node alignments, etc.

How can I detect the changes that are only related to the tag "oneway"?


texts: Layers with same coordinate system should align/overlap in ArcMap but do not?

texts:
I put 4 different layers into my ArcMap document, but only one of them will show up on the map! I should be seeing the other 3 layers as well. When I do a Zoom to layer, I can see one layer, but I cannot figure out how to make the other 3 show up at once on top of each other.

Before I opened ArcMap, I used ArcCatalog where I right-clicked 3/4 of the shapefiles and  clicked Properties → Coordinate system → Import and then clicked on my Trees layer, which has the coordinate system I want the other 3 to have. But when I opened ArcMap the 4 layers did not overlap, nor could I move them above or below each other.


texts: Attaching vertex or polygon to line in QGIS

texts:
I want to connect the green polygon and the blue line. The buffer tool is not suitable because it makes the whole shape bigger. I'd like to make longer the length of the green one to attach to the line.
Is there any way to attach the vertex or polygon to the nearest line? I thought snap tool could help that, but the output was the same as the input geometries. I also tried the v.clean tool with dissolved layer, it was the same as the result of the snap tool.


texts: Removing No Data From Raster (*.sid) using ArcGIS Desktop?

texts:
I have an Image.sid that has a white border which is No Data.
I've tried using Clip Raster Tool with no success.
I've also changed the symbology value for No Data to be transparent.

I'd like to remove No Data if possible or simply trim / erase the white border.

How can I do this?


texts:
I have a project based on data from a water company, but the arrow heads are all pointing in the wrong direction to indicate the direction of flow.
Basically I need to turn all the SYMBOLROTA values to negatives. How can I change all the values of a set of features to become minus - it might be a search function of a script or something?
Better than that would be to be able to change the type of rotation (like you can do in AutoCAD) to be anti-clockwise rather than clockwise. I haven't been able to find out how to do that but that would be the easiest?


texts: Creating installer package from GeoServer source?

texts:
I have a project which required to develop by Geoserver Source (clone from GitHub).

Now, i have done my project and i want to create the installer from my project (such as in Visual Studio you can create setup package from your solution by adding new install wizard to create setup file...).

I do not know where to start.

Can anyone help me, or show me any thread, forum that discuss about it.


texts: CartoDB - How would I customize my legend to not be white and have blocks instead of bullets?

texts:
First off, CartoDB is amazing.

I am trying to modify the look of my legend to look something like the one as seen on the "Elephant Movement in a Human Landscape" example on the CartoDB website. (http://spaceforgiants.org/?q=blog/2013/02/28/elephant-movement-in-a-human-landscape)

Similar to that example, my map has both point data and polygon data. I want to represent the point data with the standard bullets provided with the default category legend, but for the polygon layers I want to create colored (and also preferably semi-transparent) blocks.

Is there any way I can make these stylistic changes within the UI/legends/custom html window?

One last question, how would I go about changing the color of legend background?

My code looks like this:

EcoRegions

         Cascade Crest Montane Forest

         Cascade Subalpine/Alpine

         Pluvial Lake Basins

         Cold Wet Pumice Plateau Basins

         Deschutes River Valley

         Ponderosa Pine/Bitterbrush Woodland

         Pumice Plateau Forest

         John Day/Clarno Uplands

         High Lava Plains

I thought maybe I could replace the word "bullet" with something like "box" or "block" but it didn't work. I am not a programmer so don't laugh - I thought it was a good try.

Anyway, thanks for your help!

El_Tigre


texts: Why the scale of an image generated from the linear regression is not same as its inputs in Earth Engine?

texts:
Here is a sample example from developers guide,
https://developers.google.com/earth-engine/guides/reducers_regression
// This function adds a time band to the image.
var createTimeBand = function(image) {
  // Scale milliseconds by a large constant.
  return image.addBands(image.metadata('system:time_start').divide(1e18));
};

// This function adds a constant band to the image.
var createConstantBand = function(image) {
  return ee.Image(1).addBands(image);
};

// Load the input image collection: projected climate data.
var collection = ee.ImageCollection('NASA/NEX-DCP30_ENSEMBLE_STATS')
  .filterDate(ee.Date('2006-01-01'), ee.Date('2099-01-01'))
  .filter(ee.Filter.eq('scenario', 'rcp85'))
  // Map the functions over the collection, to get constant and time bands.
  .map(createTimeBand)
  .map(createConstantBand)
  // Select the predictors and the responses.
  .select(['constant', 'system:time_start', 'pr_mean', 'tasmax_mean']);

// Compute ordinary least squares regression coefficients.
var linearRegression = collection.reduce(
  ee.Reducer.linearRegression({
    numX: 2,
    numY: 2
}));

// Compute robust linear regression coefficients.
var robustLinearRegression = collection.reduce(
  ee.Reducer.robustLinearRegression({
    numX: 2,
    numY: 2
}));

// The results are array images that must be flattened for display.
// These lists label the information along each axis of the arrays.
var bandNames = [['constant', 'time'], // 0-axis variation.
                 ['precip', 'temp']]; // 1-axis variation.

// Flatten the array images to get multi-band images according to the labels.
var lrImage = linearRegression.select(['coefficients']).arrayFlatten(bandNames);
var rlrImage = robustLinearRegression.select(['coefficients']).arrayFlatten(bandNames);

print(lrImage.projection().nominalScale(), "slope_regression image") //111319.49079327357
print(collection.first().select('pr_mean').projection().nominalScale(), "slope_original image") //927.6624229062147


texts: How to install GDAL 2.1.3 using virtualenv (Python 2.7) under Windows 10?

texts:
I use Windows 10, have Python 2.7 installed as part of ArcGIS 10.5. I have installed GDAL using a package manager that runs fine in Python, meaning I can import both gdal and ogr in a Python editor. But when I try to install it using virtualenv, things don´t work - I can´t install GDAL properly. I have tried various other solutions, such as specifying the GDAL version, but I cannot install GDAL using virtualenv. How can I install it using virtualenv?

(venv) C:\Home\venv\Scripts>pip install gdal
Collecting gdal
  Using cached GDAL-2.1.3.tar.gz
Building wheels for collected packages: gdal
  Running setup.py bdist_wheel for gdal ... error
  Complete output from command c:\home\venv\scripts\python.exe -u -c "import setuptools, tokenize;__file__='c:\\users\\eric\\appdata\\local\\temp\\pip-build-ox37ve\\gdal\\setup.py';f=getattr(tokenize, 'open', open)(__file__);code=f.read().replace('\r\n', '\n');f.close();exec(compile(code, __file__, 'exec'))" bdist_wheel -d c:\users\eric\appdata\local\temp\tmpfydmk9pip-wheel- --python-tag cp27:
  running bdist_wheel
  running build
  running build_py
  creating build
  creating build\lib.win32-2.7
  copying gdal.py -> build\lib.win32-2.7
  copying ogr.py -> build\lib.win32-2.7
  copying osr.py -> build\lib.win32-2.7
  copying gdalconst.py -> build\lib.win32-2.7
  copying gdalnumeric.py -> build\lib.win32-2.7
  creating build\lib.win32-2.7\osgeo
  copying osgeo\gdal.py -> build\lib.win32-2.7\osgeo
  copying osgeo\gdalconst.py -> build\lib.win32-2.7\osgeo
  copying osgeo\gdalnumeric.py -> build\lib.win32-2.7\osgeo
  copying osgeo\gdal_array.py -> build\lib.win32-2.7\osgeo
  copying osgeo\gnm.py -> build\lib.win32-2.7\osgeo
  copying osgeo\ogr.py -> build\lib.win32-2.7\osgeo
  copying osgeo\osr.py -> build\lib.win32-2.7\osgeo
  copying osgeo\__init__.py -> build\lib.win32-2.7\osgeo
  running build_ext
  building 'osgeo._gdal' extension
  creating build\temp.win32-2.7
  creating build\temp.win32-2.7\Release
  creating build\temp.win32-2.7\Release\extensions
  C:\Users\Eric\AppData\Local\Programs\Common\Microsoft\Visual C++ for Python\9.0\VC\Bin\cl.exe /c /nologo /Ox /MD /W3 /GS- /DNDEBUG -I../../port -I../../gcore -I../../alg -I../../ogr/ -I../../ogr/ogrsf_frmts -I../../gnm -I../../apps -IC:\Python27\ArcGIS10.5\include -Ic:\home\venv\PC -Ic:\home\venv\lib\site-packages\numpy\core\include /Tpextensions/gdal_wrap.cpp /Fobuild\temp.win32-2.7\Release\extensions/gdal_wrap.obj
  gdal_wrap.cpp
  C:\Users\Eric\AppData\Local\Programs\Common\Microsoft\Visual C++ for Python\9.0\VC\Include\xlocale(342) : warning C4530: C++ exception handler used, but unwind semantics are not enabled. Specify /EHsc
  extensions/gdal_wrap.cpp(3085) : fatal error C1083: Cannot open include file: 'cpl_port.h': No such file or directory
  error: command 'C:\\Users\\Eric\\AppData\\Local\\Programs\\Common\\Microsoft\\Visual C++ for Python\\9.0\\VC\\Bin\\cl.exe' failed with exit status 2


texts: What tool was used to produce these Cross Sections?

texts:
To start with the actual question: What software was used to create these cross sections? (PDF, aprox. 10MB)

Background:

I have been really struggling to find a good software package that can generate professional cross sections. Most software focuses on 3D visualization but analysts want to see actual cross sections.
I've already asked a similar question here but this time my requirements are more complex. Where before I only cared about visualization of straight-line boreholes now I need to visualize mining drillhole survey data within section slices of variable thickness.  These drillholes are mostly inclined and curved, often slightly spiraled in  attempts to intersect subsurface mineralized zones.

Basic Parameters:

Import drillhole and visualize drillhole surveys & data; collars, surveys, assays, lithology, etc.
Must support deviated and inclined drill holes
Create subsurface geology models from drillhole data and display on sections
Import and visualize 3D DXF and other 3D (Leapfrog) subsurface models in sectional view
Define section slices (variable thickness) in plan view
able to view sections from various directions (Looking North, South, any azimuth)
Edit completed section output (adjust labels, colors, layout)
3D view and navigation is NOT necessary

It seems that the above linked sections satisfy all or most of these parameters. I would really appreciate If someone can suggest software that has these capabilities.

I've already tried and dismissed: Gemcom GEMS, Geosoft Target, CrossView for ArcGIS

Here is an example from the linked document. (Detour Gold Corporation):


texts:
I have a vector layer that I need to filter and don't to view the other registers. That filter is like the picture, I find a code Python and I executing in console but instruction return other element.

python' expression = "\"id_parroqu\"= 1" '
python' QgsFeatureRequest().setFilterExpression(expression)'

<qgis._core.QgsFeatureRequest object at 0x000001A24D00DB88>

'


texts:
Before I sit down and make the mistake of start writing my own, I wanted to know my options for good, modern, fleet tracking options that are not too expensive. What I am looking for is something that doesn't look like it was built in 1995 and that will not cost me an arm and a leg as I scale up. Brownie points if the project is Open Source, but using a hosted solution will do just fine provided that it is not too expensive.

What are people using to track on the field assets in realtime?
If price was not an issue, what would you use?

Update: Hardware is not much of an issue for this project. There are several spare servers  with plenty of bandwidth and redundant drives connected to the Internet 24/7. As far as client hardware, there are several of-the-shelf devices with GPS and connection to the internet, not to mention that some of these devices will be modified smart phones. My concentration is now more related to the server and web front end software I will use for this.

Update 2: I was asked for a use case - so I will give a simplified version. A company that does deliveries wants to know where their fleet is. There are several managers and each of them manage a chunk of delivery trucks. Each of the managers (there are about 100) is only interested in their delivery trucks (each has around 20) and wants to see them in realtime at - the same time -on a web map. Each driver wants to see where their partners are, too - but this is secondary. Privacy issues are just that each manager will be able to see their own drivers only. I have any other data required (streets, geofences, etc - but that is a separate discussion - I am looking only for this basic functionality). Yes, some time-based reports about where the trucks are/were will be necessary.


texts:
I am a couple of rungs on the ladder below junior novice in GIS. That said, I have a map, have about 100 lat/long points plotted on the map, and for those points for which I have pictures (typically 2 pictures per point).

When I  mouse roll over a point, the pictures and a brief text message pops up just as I want. However, to make this work sensibly, I had to re-save all of the pictures at a very low resolution, 100 - 200KB. If I left the pictures at their original 4 - 6 MB resolution, the screen would flood with a very zoomed in tiny portion of a picture.

Obviously re-saving everything is a solution, just not very elegant, so, how do I control the size of the popup window and make sure everything fits inside of it? And naturally some of the pictures are landscape, some portrait and the popup may have 2 landscape, 2 portrait or one of each.


texts: Summing all values with same ID in ArcGIS Desktop?

texts:
I did a one to many spatial join and ended up with an attribute table as follows:
Lets say we have 5 different IDs that have numerical values in 500 rows. I need to sum the values for each ID, and end up with 5 rows of data containing the ID and the total of all values that had that ID. Is group by the right tool for this? If so, how exactly do I use it?

Please note that these numbers are simplified and the actual number of rows is much larger.


texts:
I am using arcmap 10.2.

How do i insert survey points that were collected in PNEZD format?


texts:
When I add a WMS layer with a tile size restriction (2048 2048) the resulting WMS datasource in the QGIS project file is

<datasource>tiled=2048;2048,url=http://openmaps.gov.bc.ca/mapserver/base3</datasource>

This does not look to be a properly formed GetMap request, which I believe is resulting in the WMS server sending out a tiled WMS feed while ignoring the fact that the feed is transparent.  Using the transparency tools in layer properties results in the entire layer (features and background) becoming transparent.  If the tiled portion of the datasource is removed, the WMS feed displays correctly but does not print.   If you attempt to demo this, use the Layer group called "BASE_MAPPING_20K"

So is there a way to edit the GetMap request to be properly formatted with Width & Height & Transparent?


texts: Merge Zip Codes to Create new sales Markets

texts:
Let me preface this question by saying I am very novice with respect to GIS.  With this in mind, I am looking to merge zip codes into sales territories. Instead of plotting data by zip code, I want to have a new file that lets me plot data for the respective territories.

Currently, I get around this by plotting the same value for all zip codes in the territory.

What is the best/easiest way to do accomplish this task?  I am quite familiar with the programming language R, but usually can work through examples if this requires programming in other languages.

Many thanks in advance


texts: Getting feature coordinates from dynamically loaded layer using OpenLayers?

texts:
i want to load the coordinates of a point that is stored inside a vector layer and the vector layer is stored inside a layergroup.
I do load the vector layers dynamically to the layergroup "Runs". So i can't just do:

var source = layer.getSource();
source.forEachFeature(function(feature){
  var coord = feature.getGeometry().getCoordinates();
  // ...
});

At the moment it is totally unclear how to filter for the layer and then for the point by name to get the coordinates.

So what i have available is the layer name, layergroup name and the feature name. I tried different approaches like:

My layergroup is called 'runs' so to get the layers inside runs i did:

i tried also:

runs.getLayers().getSource();

But i get:

  Uncaught TypeError: runs.getLayers(...).getSource is not a function.

I'm having a hard time to get my head around this.

So i made it a little furhter:

      function getActiveLayer(layerName,picID) {
            // get layerGroups
            var layers = map.getLayers();
            var length = layers.getLength(), l;
            for (var i = 0; i < length; i++) {
                l = layers.item(i);
                var lt = l.get('title');
                // check for layers within groups
                if (lt === 'Runs') { // Title of Group
                  // get layers from Runs
                  var layers = l.getLayers();
                  //get length
                  var length = layers.getLength(), l;
                  if (length > 0) {
                    for (var i = 0; i < length; i++) {

                      l = layers.item(i);
                      var lt = l.get('title');
                      // check for Layer Title
                      if (lt === layerName) { // Title if Layer
                        var innerLayers = l.getSource().getFeatureById(picID).getGeometry().getCoordinates();
                        return innerLayers;
                      }
                    }
            }
                }
            }

      };

With this function i get the LayerGroup list that i can parse and identify the Group "Runs" i want. Until this point it is working.

The Problem starts at the second Point. The layers.getLength() is always zero for the layers that are inside the "Runs" Group.

The Upper one is the getLayers from map.
The lower one is the getLauers from "Runs".

Why is length 0 if there are Layers inside the array?

EDIT

      function getActiveLayer(layerName,picID) {
        var layers = map.getLayers();
        console.log(layers);
            var length = layers.getLength(), l;
            console.log(length);
            for (var i = 0; i < length; i++) {
                l = layers.item(i);
                var lt = l.get('title');
                // check for layers within groups
                if (lt === 'Runs') { // Title of Group
                  var layers = l.getLayers().getArray();
                  console.log(layers);
                  //if (length > 0) {
                    for (var i = 0; i < 99; i++) {

                      l = layers[i];
                      console.log(l);
                      console.log(layerName);
                      var lt = l.get('title');
                      console.log(lt);
                      // check for Layer Title
                      if (lt === layerName) { // Title if Layer
                        var innerLayers = l.getSource().getFeatureById(picID).getGeometry().getCoordinates();
                        return innerLayers;
                      }
                    }
                }
            }

      };

This function works. The only Problem i have that i can not figure out how the get the length of the object to config the second for loop.


texts: Create a polygon from coordinates in GeoPandas with Python

texts:
I was trying to create a simple Polygon from a list of coordinates and export it to multiple formats (including shapefile and GeoJSON) and wanted to automate it using Python. There were some good answers on creating polygons from coordinates in Python:pyshp and in Python:gdal/ogr, but I prefer using GeoPandas.

I came up with some simple code to create a simple polygon from a list of coordinates, but other users on GIS StackExchange helped to improve the code. The code is working perfectly in Python 3.6 with GeoPandas 0.4.0:

import geopandas as gpd
from shapely.geometry import Polygon

lat_point_list = [50.854457, 52.518172, 50.072651, 48.853033, 50.854457]
lon_point_list = [4.377184, 13.407759, 14.435935, 2.349553, 4.377184]

polygon_geom = Polygon(zip(lon_point_list, lat_point_list))
crs = {'init': 'epsg:4326'}
polygon = gpd.GeoDataFrame(index=[0], crs=crs, geometry=[polygon_geom])
print(polygon.geometry)

polygon.to_file(filename='polygon.geojson', driver='GeoJSON')
polygon.to_file(filename='polygon.shp', driver="ESRI Shapefile")

Possibly this code helps others and it would be nice to have suggestions how to create a better Python script that creates a polygon from coordinates in GeoPandas.

I added improvements by ImanolUr, gene and reevesii to the code! I made a stupid mistake to reverse latitude and longitude order. The location of the polygon can be checked by using folium in Jupyter Notebook:

import folium
m = folium.Map([50.854457, 4.377184], zoom_start=5, tiles='cartodbpositron')
folium.GeoJson(polygon).add_to(m)
folium.LatLngPopup().add_to(m)
m


texts: Convert from EPSG:4326 to UTM in PostGIS

texts:
I would like to convert from EPSG:4326 to UTM (30N/EPSG:32630 or 29N/EPSG:32629) in PostGIS. I do the following query but I get wrong results:

SELECT ST_AsText(ST_Transform(ST_GeomFromText('POINT(36.5277099609375 -5.86424016952515)',4326),32630)) As check;

I get "POINT(5262418.33128724 -839958.963432011)" when it should be something approximate to 243625.00,4046330.00 in UTM 30N. If I do the conversion from 4326 to UTM I get the right result but not from UTM to 4326.

What's wrong with the query?
And are there anyway to get the UTM timezone from the coordinates in
EPSG:4326 because I don't know if they belong to 30N or 29N?


texts: Different output for Fiona Shapefile CRS

texts:
I'm trying to read the crs property on a shapefile to validate projection. Below I contrast fiona documentation with what is output on my console.

Either output would be fine for my use case, but any idea why they are different?

Documenation states:
http://toblerity.org/fiona/manual.html#format-drivers-crs-bounds-and-schema

>>> c = fiona.open('docs/data/test_uk.shp')
>>> c.crs
>>> {'no_defs': True, 'ellps': 'WGS84', 'datum': 'WGS84', 'proj': 'longlat'}

I get, as if I were using from_epsg:

>>> {'init': u'epsg:4326'}

Are there scenarios where different shapefiles / prj formats yield different crs objects?

My .prj file:
  GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]]

In looking at an additional shapefile with prj:
    PROJCS["USA_Contiguous_Albers_Equal_Area_Conic",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Albers"],PARAMETER["False_Easting",0.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-96.0],PARAMETER["Standard_Parallel_1",29.5],PARAMETER["Standard_Parallel_2",45.5],PARAMETER["Latitude_Of_Origin",37.5],UNIT["Meter",1.0]]

I get the following output from c.crs:

{u'lon_0': -96, u'datum': u'NAD83', u'y_0': 0, u'no_defs': True, u'proj': u'aea', u'x_0': 0, u'units': u'm', u'lat_2': 45.5, u'lat_1': 29.5, u'lat_0': 37.5}

If I am just validating projection, is there a uniform interface besides .crs to get that info?


texts: pair value comparison to generate new raster using arcpy/python

texts:
I want to write a function which takes in two arguments (from two rasters) and compare those 2 values against a set of hard coded pair values inside function. Then return a value which i want to use to generate a new grid. I am not quite there to complete the function. Any suggestions?

def grid(val1, val2):
    table = {'10','1':'75','10,2':'75','10,3':'83','10,4':'87'}
    for i, j, k in table.iteritems():
        if val1 == i:
            if val2 == j:
                val = k
    return val

t = grid(10,1)


texts: Sampling land cover map using ArcGIS Desktop?

texts:
I have a little bit of a sampling problem in ArcGIS Desktop. I have a classified map with 7 Landcover classes, and this classes are as follows with their percentage:  Forest=20.9%, agriculture=6.3%, grassland=67.8%,water=1.50%, wetland=3.4% etc. I am generating 10000 random points and I want to distribute it according to the % or area of each class to total class (990615ha). e.g for forest 20.9/100*990615= 10000, agriculture 5.9/100*10000=622, grassland 67/100*10000=6651,

How do I generate the sample points assigning specific samples to each class ?

If I generate random points with the landcover map, it assigns randomly and disproportional to the class area and total area.

I have an X tool, I read it is possible with x tools but I cannot figure it out.


texts: How to check if parcel is located by golf course

texts: Is there any simple way for me to run millions of parcels through some type of script to see if a shape file parcel is on a golf course? For example if I have a layer that has all the golf courses and parcel layers and let's say the max distance would be 100 feet then how would i be able to extract just those parcels?


texts: How do I convert points to line in QGIS 3?

texts:
I have a dataset of gps points with DateTime and id field and I would like to create a line vector so that I can use this vector for interpolation on a raster. My goal is to find the mean of the altitude between a location and the next one. Is it possible in QGIS?
I manage to create the trajectories with the Points to path function, but now I don't know how to use this "linestring" file. Can someone help me?

EDIT: The altitude is stored in a different file, a raster. What I would like to achive is to create a line with the gps locations and then use this line to find on my altitude raster the mean of the altitude between two successive locations. My output should be something like this:

x          y        DateTime          id   mean.altitude
716310.48  7343903  12/07/2010 00:00  05   356
718398.12  7341580  12/07/2010 02:00  05   301

where mean.altitude is the mean altitude of the imaginary elements (where the line lies) between the first point and the second one.


texts: Classify point cloud from land cover raster values with lidR in R

texts:
I'm working with the R package lidR to process point clouds.    I have a raster with 4 coarse land cover classes and the corresponding tile with the LiDAR data in las format (both in the same coordinate system). Now I want to attach the land cover classes to the point cloud with the  lasmergespatial tool from the lidR package.

#import las
mylas<-readLAS(files="xyz.las", select = "xyzia", filter = "")

#import raster as rasster (not as stack)
landcover<-raster("LC_raster.tif")

#classify las with raster
class_las<-lasmergespatial(mylas, landcover, attribute = "class")

So basically I don't get errors but the attribute class in the new created las file is empty with NA. What might be a problem?


texts: PostGIS pgr_dijkstra function not finding node that appears in ways table?

texts:
I am having a problem with pgr_dijkstra algorithm not finding some of the node vertices in my table.  I have a large dataset of survey travel diaries and am trying to determine their shortest paths between each place they visit.  I linked each place visited to the closest node in my OSM network using the ways_vertices_pgr table.  I can find the nodes that are returning a blank shortest route in either the source or target column of the ways table (my edge table) but for some reason pgr_dijkstra is not finding the node when I call the function for about half of my routes.

Do I need to do something to tell the algorithm to look in both source and target columns of the ways table?

Here is my query that fills in the routes, which worked for about half of my routes.  I should mention that the prg_dijkstra function I use here is the prg_dijkstra wrapper function given on the workshop website on the pgrouting homepage.  I have tried with the full function as well and it doesn't seem to make a difference.  I have also done pgr_analyzegraph and the network topology returns OK.

update analysis.routing
set
routeline2 = (select st_multi(st_union(pgr.geom)) from (select * from
pgr_dijkstra('ways',analysis.routing.source,analysis.routing.target)) as pgr) ;

Do I need to do ST_Reverse or some other transformation to get the function to find the node?


texts: How to create dynamic text in ArcPy that recursively reflects the exported layer name from a list of layers

texts:
I have a script shown below with the help from @Midavalo that exports irrig_0, citylyr, boundlyr as a GIF. Once this is done the script then exports irrig_1, citylyr, boundlyr as another GIF and so on. I am now trying to include a text element that reflects the 'Irrig' layer name (i.e. irrig_2) and that changes for each corresponding exported GIF.

import arcpy

projFolder = r"E:\IrrigatedLands\FC_test"

mxd = arcpy.mapping.MapDocument(r"{}\irrig2.mxd".format(projFolder))
df = arcpy.mapping.ListDataFrames(mxd, '')[0]

GIFPath = r"E:\IrrigatedLands\FC_test"

layers = ['irrig_0','irrig_1','irrig_2'] # List of your output layers

lyrlist = arcpy.mapping.ListLayers(mxd, '', df)

for layer in layers:
    layerOn = "" # When required layer is turned on this will be the layer name to output
    for lyr in lyrlist:
        lyr.visible = False
        #if lyr.name in layers and lyr.name <> layer:
        #   lyr.visible = False
        if lyr.name == "citylyr":
            lyr.visible = True
        if lyr.name == "boundlyr":
            lyr.visible = True
        if lyr.name == layer:
            lyr.visible = True
            layerOn = lyr.name

        for elm in arcpy.mapping.ListLayoutElements(mxd, "TEXT_ELEMENT"):
            if elm.name == 'TestName':
                elm.text = lyr.name

            else:
                elm.text = ""
    if layerOn:
        arcpy.mapping.ExportToGIF(mxd,r"{}\{}.gif".format(GIFPath,layerOn))

del mxd


texts: Assigning existing geometry to another layer in QGIS?

texts:
I have several line layers that run parallel to each other and would like to assign the geometry of the central line to all the others without losing information in the attribute tables.

Is there a way to automatically move all lines towards the position of my central line?


texts: arcpy.CalculateField_management doesn't recognize variable in code block

texts:
I found a useful field calculator code block that I'm trying to add to a custom python tool. It works great in the manual field calculator; however, when inserted into the code block text to be used as a script, I get "NameError: name 'fc' is not defined". Apparently the script isn't seeing my variable. My question then is: How do I get the variable inside the code block?

I've tried everything I can think of including but not limited to setting 'fc' as global, escaping 'fc' from the code block string, and defining 'rows' outside of the code block, but I seem to be getting nowhere. I feel like I'm missing something obvious, but I'm out of ideas for today.  Below is my code and the error:

import arcpy

fc=arcpy.GetParameterAsText(0)

DupNumCode="""uniquelist = {}
rows = arcpy.SearchCursor(fc)
for row in rows:
  value = row.getValue('SpatialID')
  if value not in uniquelist:
    uniquelist[value] = 1
  else:
    uniquelist[value] = uniquelist[value] + 1

def findIncidence(inValue):
  return uniquelist[inValue]"""

DupNumExpression="findIncidence(!SpatialID!)"

arcpy.CalculateField_management(fc, "DupNum", DupNumExpression, "PYTHON", DupNumCode)

I get the following:

ExecuteError: ERROR 000539: Runtime error
Traceback (most recent call last):
  File "<string>", line 2, in <module>
NameError: name 'fc' is not defined

Failed to execute (CalculateField).

Here is my parameter setup:


texts: Why does my Validation Script continue to run when adding the result to the map?

texts:
I am in the process of developing a python script that creates a list available imagery for our project areas that would then allow the user to choose which imagery source they would like to add to a map.  My idea to accomplish this was to use a validation script to add an index layer to the map and create a polygon from an array and then do a select by location. A list of available imagery would be created.  I have gotten the validation script to create the list, but the problem comes in after selecting the imagery to add.  As far as I can tell the validation script runs again to every selection and after the initial run.  Below is the validation script:

class ToolValidator(object):
  """Class for validating a tool's parameter values and controlling
  the behavior of the tool's dialog."""

  def __init__(self):
    """Setup arcpy and the list of tool parameters."""
    import arcpy
    self.params = arcpy.GetParameterInfo()

  def initializeParameters(self):
    """Refine the properties of a tool's parameters.  This method is
    called when the tool is opened."""
    # Define the map document and data frame
    mxd = arcpy.mapping.MapDocument("CURRENT")
    df = arcpy.mapping.ListDataFrames(mxd, "*")[0]
    frameExtent = df.extent

    #Data Frame Array
    XMAX = frameExtent.XMax
    XMIN = frameExtent.XMin
    YMAX = frameExtent.YMax
    YMIN = frameExtent.YMin
    pnt1 = arcpy.Point(XMIN, YMIN)
    pnt2 = arcpy.Point(XMIN, YMAX)
    pnt3 = arcpy.Point(XMAX, YMAX)
    pnt4 = arcpy.Point(XMAX, YMIN)
    array = arcpy.Array()
    array.add(pnt1)
    array.add(pnt2)
    array.add(pnt3)
    array.add(pnt4)
    array.add(pnt1)
    polygon = arcpy.Polygon(array)

    #Add ImageryIndex Layer
    imgIndexLayer = arcpy.mapping.Layer(r'\\mdcfcs2\gisdept\GIS_Department\GIS_Users\HartmanB\Projects\2015\!ImageryProject\ImageryIndex.gdb\ImageryIndex')

    #Select ImageryIndex features by data frame corners
    arcpy.SelectLayerByLocation_management(imgIndexLayer, 'INTERSECT', polygon)

    field = "ImagePath"
    cursor = arcpy.da.SearchCursor("ImageryIndex", field)
    imageryList = [str(val) for val in
                     sorted(
                       set(
                         row[0] for row in cursor
                           )
                         )
                   ]
    self.params[0].filter.list = imageryList

    #Remove Layers
    for lyr in arcpy.mapping.ListLayers(mxd, "*", df):
        if lyr.name == "Polygon_Extent" or lyr.name == "ImageryIndex":
            arcpy.mapping.RemoveLayer(df,lyr)
    return

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

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


texts: ol-cesium throws exception when map.addLayer() is called with a vector layer

texts:
I have a strange behavior here using OpenLayers 5.2.0, Cesium, and the compatibility library ol-cesium. When I try to add a layer to the map, it throws an exception:

TypeError: invalid 'instanceof' operand type[En savoir plus] asserts.js:6
    assertInstanceof asserts.js:6
    createSingleLayerCounterparts VectorSynchronizer.js:106:4
    _loop AbstractSynchronizer.js:127:24
    addLayers_ AbstractSynchronizer.js:100:28
    listenForGroupChanges_/listenAddRemove</contentKeys< AbstractSynchronizer.js:253:12
    h events.js:41:11
    i</e.prototype.dispatchEvent Target.js:95:12
    M</t.prototype.insertAt Collection.js:178:4
    M</t.prototype.push Collection.js:203:4
    Rs</t.prototype.addLayer PluggableMap.js:498:11
    draw_graphMap http://127.0.0.1:8000/static/carto/js/map_graph.js:166:5
    success http://127.0.0.1:8000/static/graph/js/completion.js:135:17
    u http://127.0.0.1:8000/static/js/jquery-3.3.1.min.js:2:27452
    fireWith http://127.0.0.1:8000/static/js/jquery-3.3.1.min.js:2:28202
    k http://127.0.0.1:8000/static/js/jquery-3.3.1.min.js:2:77649
    n/< http://127.0.0.1:8000/static/js/jquery-3.3.1.min.js:2:79907

And the weirder thing is that it happens also when I try to add an empty layer to the map.

var GraphGeoSource = new ol.source.Vector({});

var GraphLayer = new ol.layer.Vector({
    source: GraphGeoSource
});

map.addLayer(GraphLayer);

I noticed that it happened only when I try to add a vector layer to the map. WMS tiles for instance don't cause any problems and displays fine in Cesium's 3D view.


texts: Using AttributeInspector to Create New Record

texts:
This is a follow up question to another questions I posted earlier:
Adding Record to Related Table using JavaScript API

Without repeating too much information from that previous question, I was able to figure out how to use ApplyEdits to add a (hard-coded) record to the related table. Now, when the user clicks on the button to "Add Record", I want a popup/info window to appear that has all of the attributes from the table, and that the user can enter the new information into for each attribute. Doing some research, I saw that the AttributeInspector will bring up the fields/values for a feature, but it seems that it must already be created.

Does anyone know a way to use AttributeInspector or another method to retrieve all of the fields (and their domains) to create a new record in a table and not just edit an already existing one? Thanks.


texts: QGIS QgsLayerTreeGroup signal "currentLayerChanged"?

texts:
In a pyqgis context, I use the "currentLayerChanged(QgsMapLayer *)" signal of the QgisInterface class to do something when the current layer is changed, but I need now to do the same when I select a group, not a layer. The class QgsLayerTreeGroup derives from QgsLayerTreeNode who has signals reflecting a change in structure, or a change in active/inactive and expanded/collapsed, but no selection signal.

Is there a way to do that ?


texts: Google maps : Is it possible to draw a polygon and retrieve all street names within it?

texts: Using Google Maps' API v3, is it possible to draw a polygon on a Google maps and get a list of all the street names that fall within that polygon?


texts: Polyline length within corresponding polygon

texts:
ArcGIS Pro 2.3: I have a feature class of polylines and for each polyline I want to analyze against its corresponding polygon. The link is the FLIGHTID field and it is a 1-to-1 relationship. I am attempting to use the Summarize Within tool to determine the length of polyline within its corresponding polygon but cant seem to figure out how to setup the inputs. Here is what my data looks like:

I am confused as to what to use for inputs to get the 1-to-1 analysis of segment length within a polygon


texts:
I have a tool that I created from a Python script which I update from time to time. Usually I will right-click the tool and choose Import Script, but lately this has been unavailable, as shown below:

This has been happening to almost all the tools I create, usually after a few imports. The problem is that if the script is cached the only workaround is to create a whole new interface, which can be pretty time-consuming. Does anyone know why this happens?

I found this post on the Esri forums describing the same problem:
http://forums.arcgis.com/threads/70829-Python-Import-Script-grayed-out

Thanks!


texts:
I have two layers, a set of households and burned area (vector data). I created buffer zones around each household and now want to estimate the burned area within a buffer. Unfortunately "intersect" doesn't help and always gives me the error code 506 - Object (506) has invalid geometry.

I am working in QGIS 3.4.


texts: How can I make my data efficient and "light" to do 3D Visualizations in Google Earth?

texts:
As a test run, I'd like to visualize population on census blocks for a few states (New England for example). I assume that the kml to do this will be huge, when you factor in the census boundaries and population data.

How can I make my data efficient to be able to do such visualizations keeping in mind that I might extend this to having multiple census periods, and doing a time lapse combined with a flythrough?

Also note that usage of Google Earth isn't necessary. I'm open to using QGIS


texts: How to make a local copy of a large WFS dataset?

texts:
We are writing a application in which a part of the job is to fetch a relatively large (>10000 features) and often updated dataset from a WFS server, store a local copy, and do geospatial analyses on it. We are stuck on the downloading part, as we hit the default 1000 features limit on the server side, which we do not control.

We see that QGIS for example can access this dataset completely, from which we can in principle manually save it in a suitable format for our application, but this is not the desirable way to go, as we want a completely automatic procedure.

We are using OGR (ogr2ogr) and can get local copies of up to 1000 features. We wonder if there is an easy way to get the entire dataset using this approach. Hints to other approaches are also welcome.


texts: Extracting pixels from a raster based on defined zones in another raster - ArcGIS

texts:
This seems like it should be a very simple task but for some reason is more difficult than I thought.

Basically I have rasters with water table depth information, with each pixel representing a depth in meters.

I have a global land-use raster, with each pixel representing a number from 1-6 (different land-use types).

I basically want to know what is the water depth in each of these 6 land-use types. Zonal statistics gives me pretty much all the information I need, EXCEPT the quartiles. And Extract by mask doesn't summarize the data by zone.

I want to know the mean, median, quartiles, and standard deviation of water table depth within each zone (land-use type).

What should I be doing?


texts: Downloading LASTools on Windows 10, QGIS 2.18.24 - Error Message

texts:
Following the instructions on this page;
How to install LAStools toolbox in QGIS

I came across this error when I tried to define the path to my LASTools folder in the QGIS Processing Toolbox Menu [C:\Lastools]

Please could anyone explain to me how I can fix this error?


texts:
I am new to Leaflet but so far so good. I managed to add the plugins I needed and they work. Except the shapefile editing.

I am using leaflet.pm to draw, snap, cut polygons (only polygons). I also use
the Leaflet.shapefile to upload a shapefile in my map and convert it as a GeoJson. All of this works. My shapefile uploads perfectly and I can see my polygons. However, the editing does not. I tried to follow this post and comments associated without success.

I also followed this post but I am clearly stucked. I turned this in every directions I could tink off.

I saw in leaflet.pm page the part about editing but I am not sure how to use it / adapt it to my uploaded file.

This is my code: https://codepen.io/mscfourn/pen/YRZebJ

Can someone put me in the right direction using the edit mode of leaflet.pm?

PS: I checked the similar questions, they concern the plugin leaflet.draw, so I apologize in advance if you consider this as a duplicate.


texts: Identifying an unknown coordinate system

texts:
I'm currently doing some software work on a set of exports of a maritime database. But I have no idea how to approach this, as I do not know the coordinate system in use. I am hoping that it is a known (although somewhat obscure) system, and not something homebrewed.

Here's an extract:

4111.099998 898.500000
4110.559998 899.800002
4111.240002 902.910000
4111.840002 901.500000

I have reason to believe that these correspond to:

68.5306667N, 15.0250000E
68.5183333N, 14.9750000E
68.5206667N, 15.0485000E
68.5306667N, 15.0250000E

(somewhere in the Norwegian or Barents Sea)

Anyone able to identfy the system, so that I can work with these coordinates?


texts: Importing Survey123 SQLite Database from Device

texts:
I need to import an SQLite database generated by ESRI Survey123 to a csv or some other useful format.  I modified this script so that I can use arcToolbox to acquire input and output parameters.  I am not familiar enough with the SQLite Survey123 format to understand how to address the error below the script:

import csv, sqlite3, json, os, sys, arcpy

def readS123db(inDB):
    conn = sqlite3.connect(inDB)
    cur = conn.cursor()
    surveys = {}

    for row in cur.execute('SELECT name, feature, status from Surveys where status = 1 or status = 3 or status = 2'):
        arcpy.AddMessage(row)
        arcpy.AddMessage ('-----------------')
        surveyName = row[0]
        if surveyName not in surveys.keys():
            surveys[surveyName] = {"adds":[], "updates":[]}
        jTransaction = json.loads(json.loads(row[1]))[0]
        arcpy.AddMessage(jTransaction)
        #Case for Adds
        if "adds" in jTransaction.keys():
            jRow = jTransaction["adds"][0]
            outRow = jRow["attributes"]
            #Add Geometry
            outRow[u"x_geometry"] = jRow["geometry"]["x"]
            outRow[u"y_geometry"] = jRow["geometry"]["y"]
            outRow[u"z_geometry"] = jRow["geometry"]["z"]
            #Add Attachments
            if "attachments" in jTransaction:
                jAttach = jTransaction["attachments"]
                for jAttRow in jAttach:
                    if jAttRow != None:
                        for jAtt in jAttRow:
                            outRow[jAtt["fieldName"]] = jAtt["fileName"]
            surveys[surveyName]["adds"].append(outRow)
        if "updates" in jTransaction.keys():
            jRow = jTransaction["updates"][0]
            outRow = jRow["attributes"]
            #Add Geometry
            outRow[u"x_geometry"] = jRow["geometry"]["x"]
            outRow[u"y_geometry"] = jRow["geometry"]["y"]
            outRow[u"z_geometry"] = jRow["geometry"]["z"]
            if "attachments" in jTransaction:
                jAttach = jTransaction["attachments"]
                for jAttRow in jAttach:
                    if jAttRow != None:
                        for jAtt in jAttRow:
                            outRow[jAtt["fieldName"]] = jAtt["fileName"]
            surveys[surveyName]["adds"].append(outRow)
        #arcpy.AddMessage(outRow)
    return surveys

if __name__ == '__main__':

    inDB = arcpy.GetParameterAsText(0)
    outDir = arcpy.GetParameterAsText(1)

    surveys = readS123db(inDB)

This gives me the following error:

Traceback (most recent call last):
  File "T:\Projects\Work_in_Progress\Guphy\2017\GRPAppData\GRPAppData\ReadDB_tool\ReadDB\readDb.py", line 100, in <module>
    surveys = readS123db(inDB)
  File "T:\Projects\Work_in_Progress\Guphy\2017\GRPAppData\GRPAppData\ReadDB_tool\ReadDB\readDb.py", line 21, in readS123db
    jTransaction = json.loads(json.loads(row[1]))[0]
KeyError: 0

Failed to execute (ImportSurvey123).

Is this a problem with my code, the SQLite database, or both?