texts: Specifying attributes for network dataset?

texts:
Using ArcGIS 9.3 to try and create some routes with network analyst... currently following this tutorial: Network analyst tutorial for 9.3

I have a street dataset from NYS clearinghouse. I can't seem to figure out how to add the FT_minutes attribute. Not sure if I can build routes based on time without it.

Any ideas?


texts: GDAL rastercalc giving unexpected results

texts:
I have a raster file containing the following values:

data values from 0 to 8488
nodata values
values equalling -3.40282347e+38

The values equalling -3.40282347e+38 are presumably supposed to be nodata values. I would like to generate a new raster where they are all equal to "nodata". For this I first tried to isolate them using the gdal raster calculator. I used the following expression containing functionality borrowed from numpy

where(isnan(A), 1,0)

I expect every value equalling nodata to become 1 in the new raster and every value not equalling nodata to become 0.

However what occurs is this:

nodata values are converted to nodata
-3.40282347e+38  values are converted to 0
data values from 0 to 8488 are all converted to 9.

I'm bewildered as to why. What am I doing wrong here?


texts: Creating triangles from points with given bearing and leg length using ArcGIS Desktop

texts:
I have a spreadsheet that contains ID, LAT_DD, LONG_DD, Bearing and LEN fields. I need to automate a way to create a triangle for every point identified in the table and I have no idea how. I know the steps should be as follows:

Create points from X,Y (Creating the known vertex, beginning, of my triangle)
Create triangles

Use bearing for triangle bisector
Use LEN for length of legs radiating from known vertex.

I'm fairly new to GIS programming and don't even know where to start. I can create the points, but that's it.


texts: zonal histogram fails to create output table when looping

texts:
Using python I'm having problems when trying to attempt a FOR LOOP with the ZonalHistogram function of ArcGIS 10.1. It works fine if I use it only for one image but it gives me this error while for looping:

ExecuteError: ERROR 999999: Error executing function.
Create output table failed
Failed to execute (ZonalHistogram)

Here is my very simple script:

import arcpy
from arcpy.sa import *

arcpy.env.workspace = 'G:\\Chile_2012\\MODIS10A2.V005\\prova'

elev_ranges = 'G:\\Chile_2012\\MODIS10A2.V005\\DEM\\elev_ranges'
list = arcpy.ListRasters ()
out_dir = 'G:\\Chile_2012\\MODIS10A2.V005\\prova\\hypso\\'

arcpy.CheckOutExtension ('Spatial')

for raster in list:
    ZonalHistogram (raster, 'Value', elev_ranges, out_dir+ raster+ '.dbf')


texts: Spatial databases learning resources for newbies?

texts:
Would be taking up a course on spatial databases next semester.

Would anybody know a newbie friendly introductory resource?

A vendor neutral one would be preferable.


texts: Seeking shapefile for Pacific Ring of Fire?

texts: I have been googling around for the shapefile that can be used to produce a map of the Pacific ring of fire, e.g. the pink region in the following map. All I could find were the map images, but I need the coordinates of the polygon that defines the boundary of the region. Where can I get hold of such data?


texts: Making Scale Bar end in zeros for Data Driven Pages?

texts:
This is related to Making Scale Bar in ArcMap end in zeros?

I want the scale bar to end in round numbers continually in data driven pages without overrunning the data frame. With the Adjust width setting, I get the round numbers and with the adjust division value, I get crazy numbers. I did see on the data driven pages setup page (extent), best fit option the ability to round scale to the nearest but its rather finicky and I can't find any documentation on how to get it to relate to the map scale bar.

See picture below, but I would prefer the miles to end in a whole number and the feet to end in 1000. The polygons that are driving the data driven pages are municipalities and range in size from small to large.

I'm using ArcMap 10.5


texts: ArcMap Form/Button VB.NET MxDocument hook

texts:
---UPDATE--- see below

I have a scaling raster tool that I'm attempting to migrate from ArcGIS 9.3.1 VBA to VB.NET (Visual Studio 2010).  I have a BaseTool that is used to put two graphic points defining the left and right extent of the image that I want to scale.  After the second point is added to the view a form opens to input the scaling parameters.  My issues is when I go to click the scale command button (within the form) I get the "Null Reference Exception was unhandled by user code, Object reference not set to an instance of an object", error.  The error occurs on this line:

pMap = m_pMxDoc.FocusMap

I have my ESRI references imported into the project for both pMap and m_pMxDoc.  It seems like the form is not hooked to the ArcMap document?  I tried to hook the form on load via the code below, but I get the same error.

Private Sub frmScaling_Load(ByVal sender As Object, ByVal e As System.EventArgs) Handles Me.Load
            Dim hook As Object = Nothing
            m_application = hook
            m_pMxDoc = m_application.Document
End Sub

Any suggestions of what I'm missing would be greatly appreciated.

Thanks

---UPDATE---

I created the docable window class and memebers.  The On Create sub now is:

Public Sub OnCreate(ByVal hook As Object) Implements ESRI.ArcGIS.Framework.IDockableWindowDef.OnCreate

        m_application = hook
        m_pMxDoc = m_application.Document

End Sub

Now when I click the command button within my form to scale the image I get a LoaderLock error (Attempting managed execution inside OS Loader lock. Do not attempt to run managed code inside a DllMain or image initialization function since doing so can cause the application to hang).  The error is still on this line:

pMap = m_pMxDoc.FocusMap

Is this error present because my form is open while the scaling action is trying to take place?  So, the hook is now on the member "On Create", do I also need to hook within the button click on the form?

Are there any other good VB.NET form/code examples within ArcMap that contain some sort of button that interacts with the ArcMap view?  I have only found this one,


texts:
I know that this is definitely an ambiguous question. However, I initially had my Python toolbox working great. I then wanted to have additional functionality with it and added some more coding. Unfortunately, I ran into some indentation errors and ultimately redid the indents from scratch (tried multiple tools).

Now when I run a syntax check in ArcMap, it passes fine. However, when I double click on my tool, called Apply_LayerFeatures, then nothing happens. I get absolutely no feedback from ArcMap.

Does anyone have any idea why this may be happening? The purpose of the tool is to be able to apply the same symbology to multiple layers. The additional feature I was trying to implement was to check if the input is a set of points, then to enable labels using the same font format.

The code I have been working on is shown below:

import arcpy

class Toolbox(object):
def __init__(self):
    """Define the toolbox (the name of the toolbox is the name of
    the.pyt file)."""
    self.label = "Toolbox"
    self.alias = ""

    # List of tool classes associated with this toolbox
    self.tools = [Apply_LayerFeatures]

class Apply_LayerFeatures(object):
def __init__(self):
    """Define the tool (tool name is the name of the class)."""
    self.label = "Apply_LayerFeatures"
    self.description = "Applies feautures to multiple layers. Customize a single layer and use it as a source layer in the tool then select layers to be modified by the tool"
    self.canRunInBackground = False

def getParameterInfo(self):
    """Define parameter definitions"""

    # First parameter = layer with customized features
    source = arcpy.Parameter(
        displayName = "Source layer with customized features",
        name = "source",
        datatype = "Feature Layer",
        parameterType = "Required",
        direction = "Input")

    # Second parameter = layers to which customized features will be applied
    updlayers = arcpy.Parameter(
        displayName = "Layers with features to be updated",
        name = "updlayers",
        datatype = "GPValueTable",
        parameterType = "Required",
        direction = "Input")

    labelFontSize = arcpy.Parameter(
        displayName = "Label Font Size",
        name = "labelFontSize",
        datatype = "GPDouble",
        parameterType = "Required",
        direction = "Input")

    updlayers.columns = [['Feature Layer', 'Features']]

    params = [source,updlayers,labelFontSize]
return params

def isLicensed(self):
    """Set whether tool is licensed to execute."""
return True

def execute(self, parameters, messages):
    """The source code of the tool."""

    #Get Inputs
    updlayers = parameters[1].values
    labelFontSize = parameters[2].values

    #labelFontSize = parameters[2].values

    for layers in updlayers:
        for labels in layers:
            labels.showLabels = True
            for a in labels.labelClasses:
                    a.expression = '"%s" + [Name] + "%s"' %  ("<BOL><FNT name='Arial' size=labelFontSize>","</FNT></BOL>")

    source = parameters[0].valueAsText
    updlayers = parameters[1].valueAsText.split(";")

    for app in updlayers:
        arcpy.ApplySymbologyFromLayer_management(app.strip("'"),source)

    arcpy.RefreshTOC
    arcpy.RefreshActiveView
return


texts: I am trying to join traffic counts to a roads layer. I have done this with a spatil Join but I have to manualy cut the road into segments before, so its 1 road segment joined to 1 traffic count.That way each road segment has a traffic count. I want to continue this but automate the cutting of road segments as I have over 100 traffic count locations. Qgis or ArcMap.


texts: Importing CSV matrix of z values to raster with no x y coordinates using ArcGIS Desktop?

texts:
I have a CSV file of depth values that have no spatial information. I know the "cell size" of the data but not where it should be put in space.

How do I convert this into a data format, such as tif, that I can import to ArcGIS Desktop to georeference as though it is an image?


texts: Problem with direction calculation in ArcGIS

texts:
I want to calculate some values with field calculator using Python in ArcGIS table. I have used this link: http://support.esri.com/technical-article/000012478 but with different values of angles. I used this:

if (angle >= 67.5 and angle <=90):

    return "N"
 if (angle > 90 and angle <=112.5):

    return "N"
 if angle >= 22.5 and angle <67.5):

    return "NE"

 if (angle >=0 and angle <22.5):

    return "E"

 if (angle >= -22.5 and angle <0):

    return "E"
 if (angle >= -67.5 and angle <-22.5):

    return "SE"
 if (angle >= -90 and angle <-67.5):

    return "S"
 if (angle >= -112.5 and angle <-90):

    return "S"
 if (angle >= -157.5 and angle <-112.5):

    return "SW"
 if (angle > -157.5):

    return "W"
 if (angle >112.5 and angle <=157.5):

    return "NW"
 if (angle >157.5):

    return "W"

but that doesnt work.I have set all parameters like it is shown in http://support.esri.com/technical-article/000012478
The error message is:
,,There was a failure during processing, check the Geoprocessing Results window for details".
What is the problem?


texts: GDAL does not recognize EPSG code of a shapefilefile reprojected by ogr2ogr

texts:
I am using ogr2ogr to reproject shapefiles. For example:

ogr2ogr -s_srs EPSG:1234 -t_srs EPSG:3857 dest.shp srs.shp

ogr2ogr then creates dest.shp, dest.shx, dest.dbf and dest.prj as expected.

The problem occurs when I try to find the EPSG code of dest.shp using GDAL. GDAL cannot identify that the file is in EPSG 3857.

Why can't GDAL identify a prj's EPSG when it created the prj from a EPSG?

Please note that I am using GDAL C# library.


texts: Adding 'marks' to point file in spatstat?

texts:
Using help on this site I have already got the following code from where I can perform certain analyses.  However, I'd now like to perform analyses on a certain variable (column z) within the point shapefile.  To do this I need to establish 'marks' on column z (the 7th column).  This column provides either a '4' or a '10', defining two types of point.

Can anyone assist?

##Run libraries:
library(“spatstat”)
library(“maptools”)
library(“sp”)
library(“rgdal”)

##set the working directory
setwd(“Z:/CP/spatstat/data/”)

##read in the shapefiles
Points <- readOGR(getwd(), "MyPointShapefile")
StudyArea <- readOGR(getwd(), "MyPolygonShapefile")

##create the window from polygon shp
w <- as.owin(StudyArea)

##coerce the points to window
pts.ppp <- as.ppp(coordinates(Points), w)
plot(pts.ppp)


texts: Iterate Field Value and Select by Attribute in ArcMap does not select rows after first iteration

texts:
I'm running ArcGIS 10.3 on Windows 10.
I need to find all the Block names and identify the record with the largest area. Then I need to calculate the Species across to the last field. However, it doesn't seem to do anything after the first iteration.
To do this I set up a basic iteration on field values, using unique values in the Blocks feature class below (which has the same Block names as the input feature class for the "Make Feature Layer" tool), selecting the first BlockId in the Feature Layer and then executing a nested subset selection (SUM_SourceHa=(SELECT MAX( SUM_SourceHa) FROM TestSource_Dissolve)) on the records for that Block to find the instance with the maximum area. Then I simply use that selection to calculate the Species value across to the SpeciesSelect field. The problem is that it calculates this for the first iteration, but not thereafter.
I suspect it's holding on to the first selection somehow and when it runs the second time, nothing is selected when it's looking for Block A2 inside the initial selection for A1 and then it gives me the green error saying that it's empty when it tries to select the maximum. As you can see, I've even tried to add the "Clear Selection" operation, but it's got no affect. Strictly speaking I shouldn't have to do that as it's a new selection on each iteration, but I thought it won't hurt.
What am I missing?


texts: Rebuild PostgreSQL index during import osm2pgsql

texts:
I started the full planet import, and after a week of work the job was ALMOST done, failed during building indexes after the import as the harddrive run out of space. Is there a way how I can just rebuild the indexes without reimporting the whole set?
Here is the log.
sudo -u renderer osm2pgsql -d gis --create --slim -G --hstore --tag-transform-script /home/renderer/src/openstreetm                                                                                ap-carto/openstreetmap-carto.lua --number-processes 8 -S /home/renderer/src/openstreetmap-carto/openstreetmap-carto.s                                                                                tyle /data.osm.pbf -C 32000 --flat-nodes /nodes/flat_nodes.bin
2021-01-30 17:39:17  osm2pgsql version 1.4.0
2021-01-30 17:39:17  Database version: 12.5 (Ubuntu 12.5-1.pgdg18.04+1)
2021-01-30 17:39:17  PostGIS version: 3.0
2021-01-30 17:39:17  Node-cache: cache=32000MB, maxblocks=512000*65536, allocation method=11
2021-01-30 17:39:17  Setting up table 'planet_osm_point'
2021-01-30 17:39:17  Setting up table 'planet_osm_line'
2021-01-30 17:39:17  Setting up table 'planet_osm_polygon'
2021-01-30 17:39:17  Setting up table 'planet_osm_roads'
2021-02-07 22:53:59  Reading input files done in 710082s (197h 14m 42s).      4/s)
2021-02-07 22:53:59    Processed 6634698177 nodes in 5524s (1h 32m 4s) - 1201k/s
2021-02-07 22:53:59    Processed 733743019 ways in 303789s (84h 23m 9s) - 2k/s
2021-02-07 22:53:59    Processed 8593984 relations in 400769s (111h 19m 29s) - 21/s
2021-02-07 22:54:00  Clustering table 'planet_osm_polygon' by geometry...
2021-02-07 22:54:00  Clustering table 'planet_osm_point' by geometry...
2021-02-07 22:54:00  Clustering table 'planet_osm_line' by geometry...
2021-02-07 22:54:00  Clustering table 'planet_osm_roads' by geometry...
2021-02-07 22:54:01  Building index on table 'planet_osm_ways'
2021-02-07 22:54:01  Building index on table 'planet_osm_rels'
2021-02-08 05:14:21  node cache: stored: 3812666007(57.47%), storage efficiency: 90.90% (dense blocks: 472519, sparse nodes: 161717496), hit rate: 56.10%
2021-02-08 05:14:24  ERROR: Database error: ERROR:  could not extend file "base/16385/5832519.8": No space left on device
HINT:  Check free disk space.


texts: How to visualize and validate public transportation data in QGIS?

texts:
I am building a routing solution app for my city using public transit. Currently I have data regarding bus stops and routes.

Currently its an sqlite table with one table containig busstops with their latitude and longitude, and other containing all the bus routes, linking bus stops and having the bus time table.

How should I go about visualising the data using qgis?  I need to visualise the data to validate my data set and to extract more information to develop my routing solution.


texts: How would one go about displaying the time that a feature was edited or digitized in a date field so as to display MM/DD/YYYY 12:00:00 AM in the attribute table?


texts: QGIS 2.18 launch Python console at startup

texts:
I'm using QGIS 2.18 on both Windows and OSX.  I have created a startup.py script in my user directory within the .qgis2/python directory.  QGIS is reading the startup script correctly because it will execute print statements and other simple tests.

I'm trying to automatically launch the Python Console when QGIS starts using the following code:

from console import console
console.show_console()

I don't get any errors.  Nothing happens.  I can however trigger exceptions with intentionally incorrect Python code.

I know QGIS will keep the console if it was open at shutdown, but I want to guarantee it's available at startup.

Launching the console at startup is fairly well documented around the Net but doesn't seem to be working in QGIS 2.18.  Am I doing something wrong or can somebody verify similar behavior in this version?


texts: Software needed to cross-reference geocodes to a Shapefile

texts:
I do not have any GIS software except Microsoft MapPoint.

I  have State district KML and Shapefiles of the districts.

I am assuming these are for ArcGIS.

I do not want to draw maps.

Instead, I am looking for a program that, if I give it one of 500,000 geocodes, it will tell me from a database perspective which District it is in.

I am looking for analytical information not visual.


texts: Leaflet : Differentiating colors dynamically

texts:
Here's what I want to do:

Show all features of the GeoJSON file using Leaflets. Connect to the database, in a database there is a table that has an ID column that matches the GeoJSON file feature, and there is a column containing the range (1-4).

IF
range = 1: red;
range = 2: yellow;
range = 3: green;
range = 4: blue;

Does anyone know this example or the idea how this should be done?
I'm following this method Choropleth but I can only do this by adding a property and entering the value into the GeoJson file. But now I want to do it dynamically, by storing the range number values in the database.


texts: Datum transformations, do we use wgs84 as an intermediate common reference?

texts: I was transforming some data the other day using the 7 param helmert transformation and it got me thinking. I have a set of transformation parameters, but these only describe the transformation from WGS84 datum, to Bessel 1841 datum. In practice... using my software, i can technically transform this data from/to any other datum i want, but does this mean: that all datums that exist actually use WGS84 (or some other) as a common denominator to which they all have known transformation parameters, and that the transformation actually happens (in the software) in the form of Datum A -> WGS84 -> Datum B, because there are no known direct parameters from A to B. In my case, since A is already WGS84 it isnt noticable of course (WGS84 -> System B)


texts: How to create an inlet drainage basin layer

texts: Using a point layer of 8,375 inlets and a DEM, how do I create an inlet drainage basin layer with contiguous basins, each delineating water drainage for and into each inlet?  I used an integer DEM, from which I used the Flow Direction tool to make an input flow direction raster.  Then I used the Watershed tool which took the input flow direction raster and the feature pour points (from the inlets layer) to create an output raster which I thought would yield contiguous inlet drainage basins but there's a lot of blank space between the basins (see screen shot), much more than make up the basins themselves.  I would think it would look like patchwork adjacent polygons (leaving no areas of the DEM blank) showing surface drainage into one of inlets at any given point.


texts: Are licensed file geodatabase restrictions/permissions inherited when converting in ArcMap?

texts: If I create a licensed file geodatabase and I convert a feature class within in it as KMZ/KML file through the conversion tool, will it inherit the restrictions that I created within the licensed filegeodatabase?


texts: How to create regular soil sampling grids?

texts: I need to be able to set up a soil sampling grid at regular spacings of around 50 - 100 meters. Then export the points to a GPS so I can be guided to the points. I am new to GIS and have tried to do this with QGIS and ArcMap 10 - so far with no luck.


texts:
How can I round the decimal number in a table in Arcmap:

e.g.

how can I set it to 2.60?


texts: Find which line has the highest Z value of a set of lines intersecting a polygon

texts:
I've a dataset of polygons representing buildings and a dataset of polylines representing buildings. The dataset  with lines has a Z value for each line, the polygon dataset doesn't. There are several lines of different object types that make up a building, each with different z values at they're a different part of the building.

I want to find the one line with the highest Z value for each building so I can then transfer the Z attribute of the highest line to the polygon file and visualise the max height of all buildings in 2D. So far I've tried just doing a one to many spatial join of lines within a building polygon to each building polygon, choosing max on the Z attribute as a merge rule without luck. Any ideas for how this can be done?


texts:
I am a high school student trying to decode a pycom code I was given, and within it I found this:

  $PMTK001,353,3,1,1,1,1,0'
  I have been searching a lot about this online, and I know that for PMTK001

  0 = Invalid command / packet.

  1 = Unsupported command / packet type

  2 = Valid command / packet, but action failed

  3 = Valid command / packet, and action succeeded

while for PMTK353,

  ‘0’: disable (DO NOT search GPS satellites)

  ‘1’ or non-Zero: search GPS satellites

Unfortunately, I've tried to put these messages together, but they don't seem to make much sense to me as of now. Could someone help me understand what the "3,1,1,1,1,0'" part says?


texts: How to emulate Google Maps driving directions using pgRouting?

texts:
Would like to emulate google maps' driving directions functionality but use pgRouting instead. google's data over the project area & scale is quite a bit off so we want to use our own transport network data in postGIS.

any hints?

thanks.


texts:
Ok. I've a server which send me a geojson data, according to the specifications. I'm using mapbox gl js, and not sure what type of method to choose: 1) creating a layer with geojson properties of geometry + adding an image to it and creating label; 2) creating a marker for each of this geojson feature.

Which method influence performance and which variant is the best, taking into account that the geojson has 1000 features inside?


texts: Calculating point density of distinct points using QGIS?

texts:
I want to calculate the density of a point layer with QGIS, but these points are dependent on one attribute. (I have different kind of points)

I do not only want the normal heatmap; I need the density of distinct points, because I want to discover hotspots of mixed points vs hotspots of homogeneous areas.

Do you know, if this is possible with normal QGIS?

I have not tried to apply GDAL/OGR in a Python script in order to achieve the objective.


texts: I have raw data with x, y, z coordinates before work and after work. I have to calculate mass volume between two of them. I have converted two files into tin. I have used surface difference (3D Analyst) tool in Arcscene. Output file contain column Volume in m3. Is this correct steps of calculation ?


texts:
Are the datum and the ellipsoid the same? What is a Datum?

  The Earth is shaped like a flattened sphere. This shape is called an
  ellipsoid. A datum is a model of the earth that is used in mapping.
  The datum consists of a series of numbers that define the shape and
  size of the ellipsoid and it's orientation in space. A datum is chosen
  to give the best possible fit to the true shape of the Earth.  in this
  article....


texts: Unable to calculate field correctly using ArcPy

texts:
I am attempted to calculate a field in a feature class ("Feat_Class_1", by using the count of selected features within another feature class ("Feat_Class_2") as a numerator, while using the value of a field in Feat_Class_1 as the denominator. I am iteratively selecting each row in Feat_Class_1 to calculate this. The script executes successfully right up until the end, when the selected row retains 0 as its value (I am expecting a proportional value). My script is as follows:

## Select Routes which are less than 800m, count the values and retain as denominator
feat_class_2_subset = str(feat_class_2+"_short")
print feat_class_2_subset
arcpy.MakeFeatureLayer_management(feat_class_2, feat_class_2_subset)
arcpy.SelectLayerByAttribute_management(feat_class_2_subset, "NEW_SELECTION", "Total_Kilometers < 0.8")
result = arcpy.GetCount_management(feat_class_2_subset)
count = str(result.getOutput(0))
print count ##<-to confirm that the subset has been successfully selected
feat_class_1_subset = str(feat_class_1+"_individual_record")
print feat_class_1_subset
arcpy.MakeFeatureLayer_management(feat_class_1, feat_class_1_subset)
selection_sub_expression = str("'"+SA2_String+"'")
print selection_sub_expression
selection_expression = str("SA2_String = "+selection_sub_expression)
print selection_expression
arcpy.SelectLayerByAttribute_management(feat_class_1_subset, "NEW_SELECTION", selection_expression2)
test_results = arcpy.GetCount_management(feat_class_1_subset)
test_count = str(test_results.getOutput(0))
print test_count ##<-to determine if the subset has been selected properly
exp = str(count+"/!Total_Properties!")
print exp
arcpy.CalculateField_management(feat_class_1_subset, "Train_800m_Perc", exp, "PYTHON", "")

This works right up until the "test_count" component but then fails to actually execute the calculation.

Any ideas on where I'm going wrong? I've been careful to ensure my expression is a string for inclusion in the function.


texts: Adding Basemaps in Map Control using ArcGIS Engine

texts:
I am developing a .Net application that is using ArcGIS Engine 10.8. I would like to add programmatically the "Add Basemap" functionality into the application's Map Control.
How can I add this command to my application?
Command: ESRI_BaseMap_BasemapCommand
GUID: {76B26001-63CF-4CAC-A7E0-91CEE8686FAC}
BaseMap.BaseMapCommand


texts: How to add the current project path in Print Composer?

texts:
I'm trying to demonstrate to my company that we can begin using QGIS in conjunction with (and sometimes even instead of) ArcGIS and MapInfo.  I've started setting up some layout templates but have run into a.  In Print Composer, is there a way to add the current project path?  In ArcGIS, the dynamic text is:

<dyn type="document" property="path"/>

I'd prefer not to have to copy/paste the document path every time I make a map.

EDIT:

To address answer and comments by artwork21 and Nathan W:

I'm using QGIS 2.0.1.  I would also like this path to appear to anyone who may need to open the project for editing without modifying their startup.py file.


texts: Adding ArcGIS Rest Service Tools/Special Features to ArcGIS Desktop?

texts:
Adding from the ArcGIS Services Directory shows how to add ArcGIS Rest Service Map Layers to ArcGIS Desktop.  However, it does not mention adding tools/Special Services.  I am attempting to add the GetLatLon service from BLM to my ArcGIS Desktop program but cannot figure it out.  I've attempted to add a Tool, and navigate to the server https://gis.blm.gov/arcgis/rest/services/ but no tools are available.

How do I add ArcGIS Rest Service Tools/Special Features to ArcGIS Desktop?


texts: How to construct a grid over a given position from polygons?

texts:
I have the following situation: I have a MultiPolygon Layer that shows some squares with the size of 100mx100m. The data looks like a raster:

I want to construct a grid with the size of 100mx100m (rectangle(polygon)), so that it fits over the MultiPolygons. I have tried it with the grid-tool in QGIS. In the settings for the expansion parameter I have tried all options, but I always get a result like in the following picture:

Does anyone have an idea how to fix the grid over the squares so that the grid fits perfectly over them?


texts: Grass is very useful tool for pre processing of imagery. I want to use the i.atcorr module separately. Will it be possible? I have compiled using .cpp files (AerosolConcentration.cpp, AerosolModel.cpp, Altitude.cpp, AtmosModel.cpp, etc) given but its asking grass functions to be declared. I dont want to use entire grass, only i.atcorr is required. How can i build that single module? Please suggest.


texts: Open Source Tool capable of reading Esri File Geodatabase (.gdb)?

texts: I have an ESRI File Geodatabase (.gdb) file that I would like to convert to shapefile (or any other non-proprietary format).  Are there any Open Source tools to do so?


texts: How do I get the "accordion" effect for input parameters in a Custom Tool?

texts:
I want to get the "Accordion" effect you find in the some tools dialog boxes in my own Python Script tool.  I can't find anything on the ArcGIS Resources website about this.

An example is shown below:

I am referring to the feature shown in "Additional Model Parameters" and the other additional toolsets.


texts: Convert WKT (or geojson) geometry to numpy array of coordinates

texts:
I am looking for the fastest way to convert WKT or geojson geometry from string to numpy array of coordinates.

Is it possible to do it strictly in Python and without creating geometry object such Shapely shape, OGRGeometry or GEOSGeometry?


texts: Implementing numpy.mean in Expression using Calculate Field Tool in ModelBuilder

texts:
I am using ArcGIS 10.4 ModelBuilder.
I developed a model that iterates through a bunch of feature datasets.  I need to add a field to each one of these features and apply one formula on it. My formula consists on dividing one column by its mean.
I have implemented a formula on the "expression" part of the "Calculate Field" tool.  As it is shown in the image:

The result of this is always 1. It seems that the mean is calculated per row and not for the whole column.  So it is dividing the value by itself giving 1 as a result.

I have also tried to define my own mean, but so far with out success
def my_mean(sample):
 return sum(sample) / len(sample)

How should I proceed?


texts: Finding number of links where buffer and link intersect using PyQGIS?

texts:
Created a buffer of 0.001 size at the beginning and end of each link. I want to extract the starting point and the ending point of which the number of the buffer and the link intersect is one.

layer = iface.activeLayer()
for f in layer.getFeatures():
  geom = f.geometry().asPolyline()
  start_point = QgsPoint(geom[0])
  end_point = QgsPoint(geom[-1])


texts: Alternative to using ArcGIS Online Geocoding

texts: I have 1.4 million addresses that I need to batch geocode. I can't find a free service to do so. Does anyone know of something out there that I can use? Oh, the service needs to be were I can upload a file in a .txt file formate, not a copy and paste. Copy and paste services crash the 661 mb size file.


texts: Copying features from a shapefile into feature class using ArcPy with ArcMap gives ERROR 000725

texts:
I have a shapefile which contains some features and I need to copy these features and paste them into an EXISTING feature class using ArcPy. I have tried the (arcpy.CopyFeatures_management(..)) but it seems it is used only to copy data from a shapefile to a NEW feature class.

import arcpy
input_ShapeFile ="C:\Users\Moh\Desktop\file.shp"
output_FeatureClass ="C:\Users\Moh\Desktop\Data.mdb\fc"
arcpy.CopyFeatures_management(input_ShapeFile, output_FeatureClass)

...ERROR 000725: Output Feature Class: Dataset "C:\Users\Moh\Desktop\Data.mdb\fc"
 already exists. Failed to execute (CopyFeatures).

I am using ArcMap 10.5.


texts: How to match two feature classes with same geography and use part of the name of one to name the other

texts:
So I have a polygon feature classes with unique identifying names to them that I want to use to rename some other feature classes that are of the point feature type and have the same geography of one of the matching polygon feature classes.  Basically the points from the second feature class fit within the polygon feature class.  The point feature class has no name identify to it so I would like to change that and add at least part of the name of the polygon feature class to the name based on that fact they share the same geography.

I was thinking I could use the Select by Location tool but the problem is I basically just have a bunch of unnamed point feature classes in a folder so it's not like I can select individual feature classes.  So some how I need to iterate thru all point feature classes and rename them based on the polygon feature class they match up to.  I'm guessing I may use Arcpy to do this.  Any help or direction would be appreciated.


texts: Setting GPFeatureLayer parameter programmatically in Python Toolbox (.pyt)?

texts:
I've created a python toolbox (.pyt file) with 2 tools in it: Prep() creates various derived datasets; Track() does processing using the prepped data. The idea is for Prep to be run just once, interactively (from the user dialog), while Track will be run many times, interactively. Both tools are working perfectly when run interactively.

I would like, however, for Prep to call Track a single time after the prep work is done, with a default selection so that there is some processed output as well as the prepped data when it has completed -- this entails an internal call to Track()  rather than an interactive call from its user dialog.

At the end of Prep's execute() method, I create an instance of Track, get its parameter info, assign parameter values, and invoke its execute() method with those parameters. The GPString and GPBoolean parameters work fine, but the GPFeatureLayer parameter does not: when I assign a Feature Layer object to the parameter value, nothing changes: the parameter value is still None after the assignment (in contrast, the other parameters show their assigned values after the assignments).

Here's how I'm assigning the parameter values and calling Track:

    track = Track()
    trackParams = track.getParameterInfo()
    trackParams[0].value = 'default'           # GPString
    # targetpoints_layer holds a point FeatureLayer returned from arcpy.mapping.ListLayers()
    trackParams[1].value = targetpoints_layer  # GPFeatureLayer
    # REPLACE_EXISTING holds a boolean value
    trackParams[2].value = REPLACE_EXISTING    # GPBoolean
    track.execute(trackParams, None)

Here's getParameterInfo() for the Track tool:

def getParameterInfo(self):
    params = []

    p = arcpy.Parameter(
        displayName='Track name',
        name='TrackName',
        datatype='GPString',
        parameterType='Required',
        direction='Input')
    params.append(p)

    p = arcpy.Parameter(
        displayName='Target point layer',
        name='TargetPointLayer',
        datatype='GPFeatureLayer',
        parameterType='Required',
        direction='Input')
    p.filter.list = ['Point']
    params.append(p)

    p = arcpy.Parameter(
        displayName='Replace existing outputs?',
        name='ReplaceExisting',
        datatype='GPBoolean',
        parameterType='Required',
        direction='Input')
    p.value = False
    params.append(p)

    return params

The first sign of trouble was when I attempted to get the name of the target point layer in Track.execute(): it threw an exception "AttributeError: 'NoneType' object has no attribute 'name'". I followed this back to what seems to be the problem: the assignment of the parameter value in Prep before the call to track.execute().

Inspecting trackParams[1].value before and after assigning it targetpoints_layer shows no change -- it's value is None both before and after. In contrast, inspecting the other parameter values before and after assignment shows them changing from None to 'longest' and False, respectively. Inspecting targetpoints_layer itself shows it is indeed a FeatureLayer object (for example, targetpoints_layer.isFeatureLayer is True). I can successfully assign targetpoints_layer to a new variable, but not to a parameter's value attribute.

Any ideas why the assignment is having no effect (nor throwing an exception)?

Do you see anything I'm doing wrong or failing to do?

Why do the other assignments work but not the feature layer one?


texts: Preserving Geometry of Shape file after merging edges using Networkx

texts:
I have a stream network and I want it to convert in a graph like data structure having vertices with only degree 1 and 3. The degree of the vertex is the number of edges incident to that vertex. The current stream network also contains degree 2 vertices which are to be eliminated. The elimination is expected as the merging of adjacent edges of the vertex with degree 2 node as shown below. The label in the figure represents the stream vertex id and degree of the vertex (ID-Degree)

Input shapefile--Stream Network and Point layer

Expected Output:-

It can be seen from the above figure that initially the shapefile has stream segments as 60-> 61, 61->64, 64->65 and the output stream network has 60-> 65 as the stream segment with edges corresponding to degree two nodes are merged.

Below is my code using Networkx which reads a shapefile and eliminates the degree two nodes from the network.

import networkx as nx
import matplotlib.pyplot as plt
import csv
import warnings

def fxn():
    warnings.warn("deprecated", DeprecationWarning)

with warnings.catch_warnings():
    warnings.simplefilter("ignore")
    fxn()

path='/home/swapnil/Desktop/code/Final_files'
G=nx.read_shp(path + '/Stream_Segment_Vect.shp' , geom_attrs =  True , strict = True)
pos = {k: v for k,v in enumerate(G.nodes())}
X=nx.Graph() #Empty graph
X.add_nodes_from(pos.keys()) #Add nodes preserving coordinates
# print(pos.keys())
l=[set(x) for x in G.edges()] #To speed things up in case of large objects
edg=[tuple(k for k,v in pos.items() if v in sl) for sl in l] #Map the G.edges start and endpoints onto pos
# print(edg)
nx.draw_networkx_nodes(X,pos,node_size=100,node_color='r')
X.add_edges_from(edg)

remove = [node for node,degree in dict(X.degree()).items() if degree == 2]

# For each of those nodes
for node in remove:
    # We add an edge between neighbors (len == 2 so it is correct)
    # print(node)
    X.add_edge(*X.neighbors(node))
    # And delete the node
    X.remove_node(node)

X.remove_nodes_from(remove)

nx.draw_networkx_edges(X,pos)

# print(X.degree[1])

plt.show()

But this code creates a straight line between merged edges as shown below instead of preserving the geometry which stream network had.

Doubts:-
1)Does networkx preserve shapefile geometry when imported to QGIS, irrespective of how the image output is displayed, if yes, how?

2)How can I import geometry attributes in QGIS of the image output of networkx?


texts: Creating labels based on attribute field in shapefile with Arcpy?

texts:
I want to be able to label point features on my map in a python loop.

How do I access the ArcGIS label engine with ArcPy?


texts: Aligning shapes to terrain in CityEngine?

texts:
I have a problem. As soon as I click "Align shapes to Terrain", my shapefile rectifies and falls down the edge.

In the first picture(Border_Issue1) you can see that all shapefiles have the right size, in the second picture(Border_Issue) you can see the mentioned problem.

Does anyone have the same problem and can help me?


texts: Clip one folder of shapefiles by another with same name endings

texts:
I have two folders with shapefiles:
r'E:\SO\Polygons' and
r'E:\SO\Districts'

Within each folder I have 215 shapefiles named with similar naming conventions:
poly_100 or dist_100. I need to clip the polygons by the districts if their ending digits are equal. I am new to python and at a loss:

    #folder for clipped Voronoi polygons
Voronoi_file = str(maindir) + "/" + "Voronoi_polygons"
if not os.path.exists(Voronio_file): os.makedirs(Voronoi_file)

env.workspace = r'E:\SO\Polygons'
outputworkspace= Voronoi_file
count=0
for fc in arcpy.ListFeatureClasses():
    end1 = str(fc[4:])
    env.workspace = r'E:\SO\Districs'
    for fc in arcpy.ListFeatureClasses():
            end2 = str(fc[4:])
            if end1==end2:
                out_name = str((outputworkspace) + "/" + "dis_poly" + str(count))  # Assemble the output point name and path
                arcpy.Clip_analysis(end1, end2, out_name)
                count = count + 1


texts: Writing collector cable optimizer for ArcMap?

texts:
We're going to create an extension within ArcMap that uses python to process information within a shapefile that the user drags into a toolbox. The program is supposed to determine as final output the best routes for laying collector cables.

The model will need to account for, among other things, road allowances, properties signed onto the project, the maximum number of Wind Turbines per electrical circuit, environmental features, cost for running x number of collector cable circuits overhead, cost for running x number of collector cable circuits underground etc

It's a pain because i've got to learn a new language from scratch. Anyway, i understand that such models exist already and i was wondering if someone could shoot me off in the right direction by pointing one or more out and possibly any other useful resources as well.


texts:
I've got a ".shp" file (plus the additional .dbf, .fix, .prj, .shx files) that was created using the geotools library and want to view its data using QGIS 3.8.2.

QGIS properly displays all of the points in the .shp file and clicking on the "i" button and then on one of the points displays its information at the right-hand side:

Name
ID
X
Y
Z

The only thing missing is the "M" value. I know this is probably a very basic question but I can't find any information about it anywhere: How do I get access to the M value?

I don't want to edit it or create new ones, I just want to see the existing one. Is there maybe any type of table that has a list of every point in the .shp file and all of its values (XYZ, M, name,...)?


texts: Using xarray to resample and merge two datasets

texts:
I have one dataset of satellite based solar induced fluorescence (SIF) and one of modeled precipitation. I want to compare precipitation to SIF on a per pixel basis in my study area. My two datasets are of the same area, but at slightly different spatial resolutions. The SIF is a little lower resolution than the rainfall. I can successfully plot these values across time and compare against each other when I take the mean for the whole area, but I'm struggling to create a scatter plot of this on a per pixel basis.

I'm not sure if this is the best way to compare these two values when looking for the impact of precip on SIF so I'm open to ideas of different approaches. As for merging the data currently I'm using xr.combine_by_coords but as described below it doesn't seem to be able to resample on the fly or I'm not using it correctly. I could also do this by converting the netcdfs into geotiffs and then using rasterio to warp them, but that seems like an inefficient way to do this comparison. Here is what I have thus far:

import netCDF4
import numpy as np
import dask
import xarray as xr

rainy_bbox = np.array([
    [-69.29519955115512,-13.861261028444734],
    [-69.29519955115512,-12.384786628185896],
    [-71.19583431678012,-12.384786628185896],
    [-71.19583431678012,-13.861261028444734]])
max_lon_lat = np.max(rainy_bbox, axis=0)
min_lon_lat = np.min(rainy_bbox, axis=0)

# this dataset is available here: ftp://fluo.gps.caltech.edu/data/tropomi/gridded/
sif = xr.open_dataset('../data/TROPO_SIF_03-2018.nc')

# the dataset is global so subset to my study area in the Amazon
rainy_sif_xds = sif.sel(lon=slice(min_lon_lat[0], max_lon_lat[0]), lat=slice(min_lon_lat[1], max_lon_lat[1]))

# this data can all be downloaded from NASA Goddard here either manually or with wget but you'll need an account on https://disc.gsfc.nasa.gov/: https://pastebin.com/viZckVdn
imerg_xds = xr.open_mfdataset('../data/3B-DAY.MS.MRG.3IMERG.201803*.nc4')

# spatial subset
rainy_imerg_xds = imerg_xds.sel(lon=slice(min_lon_lat[0], max_lon_lat[0]), lat=slice(min_lon_lat[1], max_lon_lat[1]))

# I'm not sure the best way to combine these datasets but am trying this
combo_xds = xr.combine_by_coords([rainy_imerg_xds, rainy_xds])

Currently I'm getting a seemingly unhelpful RecursionError: maximum recursion depth exceeded in comparison on that final line. When I add the argument join='left' then the data from the rainy_imerg_xds dataset is in combo_xds and when I do join='right' the rainy_xds data is present, and if I do join='inner' no data is present. I assumed there was some internal interpolation with this function but it appears not.


texts: Trimble Access vs. TerraFlex: Which is simpler for creating data entry forms?

texts:
I am working on an archaeological project and need to create forms for users in the field to use and take photos. I am renting a Trimble GNSS unit and can get Trimble Access for free with the unit, or I can pay for Terraflex separately. Terraflex seems simpler to use, but I want to keep costs low.

Does anyone have experience with both Trimble Access and TerraFlex? How hard is it to create forms with Trimble Access compared to TerraFlex? I have created a TerraFlex trial account and successfully tested creating forms and tasks on the office/web side which can be pushed out to mobile users.


texts: GDAL and Python: how to create new files for (e.g.) clipping with gdalwarp

texts:
I'd like to learn using GDAL, and I want to use it in Python.
Following the tutorial (http://www.gdal.org/gdal_tutorial.html), I succeded in opening vector and raster files using ogr.Open and gdal.Open () calls respectively.
I then tried to understand how to create new files, but I didn't figure it out so far.
For instance: if I want to use "gdalwarp" function to (let's say) clip a raster with a polygon, it does its job perfectly and clearly within QGIS, but using GDAL in Python how do I:

1) specify the output folder?

2) create the new file (I think it's done using the Create or CreateCopy () methods, is it?)?

I hope someone could help me, thanks in advance to the community.

UPDATE

I came up with creating a file in a specific folder. Anyway, it is just a file to be filled, namely, now I want to know how can I do my calculation and output to this "empty" file. Here is my code to let you understand me:

from subprocess import call
call(["ls", "-l"])
import gdal, ogr
from gdalconst import *

input= 'C:/Users/unimi/Desktop/sept2000_B5_refl.tif'
clip_area= 'C:/Users/unimi/Desktop/study_area_foscagno.shp'
output= 'C:/Users/unimi/Desktop/clip_sept2000_B5_refl.tif'

shp= ogr.Open (clip_area)
if shp is None:
    raise SystemExit("The shapefile could not opened")
raster= gdal.Open (input, GA_ReadOnly)
if raster is None:
    raise SystemExit("The raster could not opened")
print raster.GetProjection ()
format= 'GTiff'
driver= gdal.GetDriverByName (format)
metadata= driver.GetMetadata ()
dst_ds= driver.Create (output, 2133, 1273, 1, GDT_Float32)

ClipByPolygon= 'gdalwarp -dstnodata 0 -q -cutline %s -crop_to_cutline -of GTiff %s %s' % (shp, raster, dst_ds)
call (ClipByPolygon)

dst_ds= None

The code doesn't raise any error. However, the "empty" file from the "Create" call is produced, but gdalwarp has no effects afterwards.


texts: Multiple Ring Buffer for Multiple Points

texts: I am creating multiple ring buffer around points. I am using Multi Ring Buffer plugin of QGIS. When I am creating 5 ring buffers of a single point its returning a shape I desire but in Attribute Table there is no trace of point information from which the buffer was created. Moreover when I am creating 5 ring buffers of multiple points its returning similar single buffer shape using set of points. But I am trying to get separated 5 ring for each point. Is there any way to get separated multiple ring buffer for each point and having the referenced point information in QGIS?


texts: ArcGIS 10.3.1: Spatial Join & Count of points on a regular grid

texts:
Following the indication given in Counting number of points within a grid I was able to count how many points of a shapefile fall within each cell of a fishnet in my case study.

After loading the shapefile, I have created the fishnet and have used Analysis Tools->Overlay->Spatial Join. The target feature was the fishnet, and the join feature was the point shapefile. I have opted for the Join-One-to-One option, and have added a count field to know how many points are within each cell of the fishnet.

My problem: when I select one cell with multiple points (e.g., the red cell in the image below) I get the correct total number of points (3), but I do not see their IDs.

So I only see the cell ID (117), but not the IDs of the points inside of it. Is there any way I can get these IDs? I need to derive some sort of phonebook where I know that cell i is associated with points p1, p2, ....


texts:
I have a point map of trees along with their size (diameter) as an attribute and I want to create a raster that with the output of a neighborhood analysis

.

Basically for every point in a raster of specified resolution, I'd like to make the following calculation

For every tree within radius r, sum up for all trees.

(tree size)^2 / (distance from center)^2

Note that this is getting the quotient for each tree individually, then summing the up and assigning the middle pixel that value.

I can't picture this working with with the raster calculator or focal analysis tools, so I think I'd have to program something, but I'm not very experienced using Python.

Any suggestions?

I am using ArcGIS 10. And have tried various ways of looking at both Focal Statistics and Point Statistics. The problem with both of these is that they summarize on a single attribute. That is, I can get point densities for a given area, or I can sum up tree sizes for a given area... Dividing those two layers dones't quite give me what I need.

If I were to do this by hand, I would repeat these steps in a grid every 5 meters

Identify all trees within radius
Get size of each tree
Divide each tree by the distance^2 from the center
Add up the quotients for all trees within radius
Assign that sum to a raster pixel at the point center

I've been able to do this analysis in R using a loop that goes through each point, but I'm not familiar with Python, so I hoped I could do this out of the box


texts:
I have a Trimble GeoExplorer Series GeoXM that runs Windows PocketPC version 4.2 (Built March 15 2004 for a bit more resolution there).

Lately when the unit has been out on a job it will periodically disconnect from the internal GIS receiver (TISP on COM3) and give an error message that says:

"Failed to connect to the GPS receiver because the COM port in in use by another application"

My field techs tell me that this is happening randomly, sometimes several times a day.  So far the only way I have to deal with it is to tell them to just keep trying to connect and after a few tries (generally around 20 minutes of waiting) it comes back online and they carry on with their task.

This unit is used exclusively to run Terrasync (v3.01) - it doesn't do anything else so I don't understand what other program could be occupying the COM port when Terrasync needs it.

I have tried to uninstall and reinstall Terrasync and the problem has persisted.


texts: Which SQL Server spatial data type to use for UK data in SQL Server?

texts:
Which spatial data type should be used for UK data in SQL Server?, i.e. GEOMETRY or GEOGRAPHY?

I'm aware of OS National Grid and my understanding is it's an X/Y coordinate system, covering from Sicily Isle to Shetlands.

So it looks like best option would be GEOMETRY

But if we store our data in GEOMETRY, will it work with things like Open Street Maps, e.g. using LeafletJS with OSM base, and our data as overlays


texts:
Can Nominatim search only boundary:administrative, because I get for every place two or three result, but I only need one.
Examples:
Bobovo
Stara Pazova


texts: How to calculate shortest path dynamically considering other shortest paths calculated in same run?

texts:
My topic is to minimize the cost of construction line digging. Starting from a source point I have several households that need to be connected to the that line. My plan was to calculate the shortest path for each source -> house and then cutout the overlapping lines.
In the example below I would then have a route on the Eastern street that goes from house 1-4, where the line for house 4, for instance, only reaches to house 3 and not the whole way to the source point (red) as I cutout its shortest path with the shortest pathes from houses 1 to 3. However, for house 5 the shortest path would have been the Western road, which is why it doesn't intersect with any other shortest path. But practically a connection from house 4 would have been more cost efficient. How can I include this consideration into my algorithm using PostGIS and pgrouting?


texts: Looking for software to interface GPS data with PC

texts:
I have GPS Unit from Garmin (GPS 18xPC).

I'm looking for software I can install on my PC that will allow me to interface with the data I collect on my GPS unit.

Any help will be great


texts: ArcGIS Length and Area Calculation Scenarios

texts:
Trying to get a handle on how length and area are calculated in different scenarios in ArcGIS. I don't know why I can't find a answer on the feature class fields, but I can't find a precise answer, unless I'm not understanding something and I know there is some history. Can you help me fill in the question marks? Or tell me why I'm going about this all wrong ; )

GCS = Geographic coordinate system
PCS = Projected coordinate system
All links are to 10.1 help docs --

Feature Class shape_length and shape_area fields
a. GCS - ?
b. PCS - Using simple planar
c. Is it always auto-updated, except for shapefiles? yes
ArcMap Measurment Tool
a. GCS - default geodesic, alternatives are Loxodrome and Great Elliptic, but not planar. Area calculation not available!
b. PCS - default planar, alternatives Geodesic, Loxodrome, and Great Elliptic
http://resources.arcgis.com/en/help/main/10.1/index.html#//00s500000022000000
Attribute Table Calculator
a. GCS - not available
b. PCS - planar
http://resources.arcgis.com/en/help/main/10.1/index.html#//005s00000027000000
Calculate Field Tool (Data Management toolbox)
a. GCS - geodesic linear, area available but questionable
b. PCS - planar
http://resources.arcgis.com/en/help/main/10.1/index.html#//00170000004m000000
Buffering Tool (and other tools coming)
a. GCS - geodesic
b. PCS - planar or specify GCS output
http://resources.arcgis.com/en/help/main/10.1/index.html#//000800000019000000
Javascript API Clientside
a. GCS - geodesic area and length functions
b. PCS - can convert from web mercator to geographic (or use geometry service)
http://help.arcgis.com/en/webapi/javascript/arcgis/help/jsapi/namespace_geometry.htm
Flex API Clientside
a. GCS - geodesic area and length functions, "The length [or area] will be calculated using a custom cylindrical equal-area projection". This is not mentioned in the javascript api!!
b. PCS - can convert from web mercator to geographic
http://resources.arcgis.com/en/help/flex-api/apiref/com/esri/ags/utils/GeometryUtil.html
ArcGIS Server REST API - Geometry Service
a. GCS - geodesic
b. PCS - planar
http://help.arcgis.com/en/webapi/javascript/arcgis/help/jsapi/geometryservice.htm

Another question, what exactly is a geodesic measurement? I thought it meant a 3D trig formula on a spheroid (haversine?). And is it too slow to use in calculating an area and that is why equal area projections are used?

Another question, when determining length and area -- is an equal area projection more accurate than a geodesic calculation using the same datum, spheroid? And briefly why?


texts: Isochrone along direction map of watershed

texts:
I am interested in producing a isochrone map given a raster speed map and a direction map (like the output of r.watershed), but I do not see how to produce it.
I have studied r.walk and r.cost wich calculate accumulated cost on a friction map, but it is not intended to work with a direction map.
I have been trying to see if there is a way to produce the inverse of an accumulative map like the result from r.drain but backword, so given the outlet and the direction maps, assign a cumulative cost for each cell of the catchment area from the sarting point (downstream).
I am familiar with GRASS and QGIS, but any solution with a free software is welcome.


texts: Placing a cost raster in Arcmap Table of Contents in 10.2.2

texts:
We have a set of routines that have been used for writing rasters to the ArcMap table of contents from our tools.  The code that always worked is this:

    /// <summary>
    /// Used for testing only.  Do not execute this in published server code.
    /// </summary>
    /// <param name="rasterDataset"></param>
    /// <param name="layerName"></param>
    internal static void PutToMap(IGeoDataset rasterDataset, string layerName)
    {
        if (rasterDataset == null)
        {
            throw new ArgumentNullException("rasterDataset");
        }
        RasterLayer layer = new RasterLayer();
        IApplication baseApp = (IApplication)Activator.CreateInstance(Type.GetTypeFromProgID("esriFramework.AppRef"));
        IMxDocument doc = (IMxDocument)(baseApp.Document);

        IMap activeMap = doc.ActiveView.FocusMap;
        if (rasterDataset is Raster)
        {
            layer.CreateFromDataset(((IRaster2)(Raster)rasterDataset).RasterDataset);
        }
        else
        {
            layer.CreateFromDataset((RasterDataset) rasterDataset);
        }
        layer.Name = layerName;
        layer.Visible = true;
        activeMap.AddLayer(layer);
    }

    /// <summary>
    /// Used for testing only.  Do not execute this in published server code.
    /// </summary>
    /// <param name="layer">The raster layer or subclass that will be displayed</param>
    /// <param name="layerName">The name to be put in the TOC</param>
    internal static void PutToMap(IRasterLayer layer, string layerName)
    {
        IApplication baseApp = (IApplication)Activator.CreateInstance(Type.GetTypeFromProgID("esriFramework.AppRef"));
        IMxDocument doc = (IMxDocument)(baseApp.Document);

        IMap activeMap = doc.ActiveView.FocusMap;

        layer.Name = layerName;
        layer.Visible = true;
        activeMap.AddLayer(layer);
    }

    internal static void PutToMap(IRaster raster, string layername)
    {
        IRasterLayer layer = new RasterLayerClass();
        layer.CreateFromRaster(raster);
        PutToMap(layer, layername);
    }

We are now trying to debug a problem but when I try to write the cost raster out, these methods return an error when I try getting the singleton for ArcRef.  The Activator.CreateInstance throws an exception in trying to load the module.  Has anyone done this in 10.2.2 and is there a way to get things running that I'm missing???

I have tried changing the CreateInstance from Type.GetTypeFromProgID to Type.GetTypeFromCLSID and both approaches throws an exception.  It returns an error of 8000ffff or something like that.


texts: Get the begin and end coordinates of city municipality

texts:
I am not sure if this is the right place to ask this question. I would like to have an idea of precise coverage (extend) for the latitude/longitude of the city of Porto, Portugal.
This is in connection with my research that I need to filter out geolocation data for the cities of Porto (Portugal) and Beijing (China). Luckily, I found one for Beijing city here  as "ranges from 39° 27' to 41° 03' N and in longitude from 115° 25' to 117° 30' E" which is what I want.
For the city of Porto, all I get is a one-point coordinates (41.1579° N, 8.6291° W).
Any pointer to some resources?


texts: How do I configure osm2po.config to parse both car and foot routes?

texts:
I have tried to play around with the osm2po.config file, to parse both car and foot routes but have had little success. Could someone please tell me how to do it?

Thanks


texts: ESRI Workspace Domains sometimes not implementing ICodedValueDomain2

texts:
I am attempting to use the .SortByName() and .SortByValue() methods of the interface ICodedValueDomain2.

I use the same extension methods from an IWorkspace, however if I open the IWorkspace from an .sde connection file and try to cast the domain to an ICodedValueDomain2, I get the exception:

Unable to cast COM object of type 'System.__ComObject' to interface type 'ESRI.ArcGIS.Geodatabase.ICodedValueDomain2'. This operation failed because the QueryInterface call on the COM component for the interface with IID '{057E88C6-0E61-426E-B4BB-D192C481CF1F}' failed due to the following error: No such interface supported (Exception from HRESULT: 0x80004002 (E_NOINTERFACE)).

// *** Extension Methods ***

//This Method Works With .sde files and through ArcMap
public static ICodedValueDomain GetCodedValueDomain(this IWorkspace pWorkspace, string pName)
{
    IWorkspaceDomains workspaceDomains = (IWorkspaceDomains)pWorkspace;
    IDomain domain = workspaceDomains.get_DomainByName(pName);
    ICodedValueDomain codedValueDomain = (ICodedValueDomain)domain;
    return codedValueDomain;
}

//This Method does not work with .sde files, only through ArcMap
public static ICodedValueDomain2 GetCodedValueDomain2(this IWorkspace pWorkspace, string pName)
{
    IWorkspaceDomains2 workspaceDomains = (IWorkspaceDomains2)pWorkspace;
    IDomain domain = workspaceDomains.DomainByName[pName];
    ICodedValueDomain2 codedValueDomain = (ICodedValueDomain2)domain;
    return codedValueDomain;
}

// *** Workspace Retrieval ***

//Cannot use .GetCodedValueDomain2()
private static IWorkspace WorkspaceFromSDEFile(string pPath)
{
    SdeWorkspaceFactory wsFactory = new SdeWorkspaceFactory();
    IWorkspace ws = wsFactory.OpenFromFile(pPath, 0);
    return ws;
}

//Can use .GetCodedValueDomain2()
private static IWorkspace WorkspaceFromArcMap()
{
    IWorkspace ws = ((IDataset)((IFeatureLayer)ArcMap.Document.FocusMap.get_Layer(0)).FeatureClass).Workspace;
    return ws;
}

Why would the IWorkspace returned from wsFactory.OpenFromFile() not implement this interface, but the IWorkspace retrieved from the map does?

Here is the ArcObjects documentation for the ICodedValueDomain2 interface, present since 10.1: http://resources.arcgis.com/en/help/arcobjects-net/componenthelp/#/ICodedValueDomain2_Interface/002500000nrt000000/


texts: OpenLayers 2 cluster map convert from WMS to OSM

texts:
I have implemented a simple OpenLayers map to plot stores on the map. I've done this using a WMS.

Fiddle: https://jsfiddle.net/7n42uaux/3/

Now I want to use OpenStreetMap(OSM) maps.

new OpenLayers.Layer.OSM();

When I do change it, the map gets drawn but the clusters don't get plotted properly. They are all shown as a summation and in the wrong location in the map.


texts: Getting full list of cities and regions from OSM

texts:
I need to match cities/postal codes from a large dataset with their regions (NUTS3 if available, otherwise administrative divisions as can be found here https://en.wikipedia.org/wiki/List_of_administrative_divisions_by_country).
Cities can be all around the world and have all sizes.

I've come across that OSM provides regions for all(?) cities.
Searching for a small village in Germany ('Unterkirnach'), I get the following result:
    Village Unterkirnach, Verwaltungsgemeinschaft Villingen-Schwenningen, Schwarzwald-Baar-Kreis, Regierungsbezirk Freiburg, Baden-Württemberg, Germany

That's all the information I need, namely administrative divisions on several levels, and it also seems to work for small towns, e.g., in China.

However, I would need to extract a list/table of towns and their assignment to regions worldwide in order to upload it to my database and to look up the regions of the cities in my dataset.
At this point, I do not want to do any 'fancy' stuff I do not understand, I simply want to create a table with as many cities as possible and columns for their administrative regions on different levels.

Does anyone have an idea how to do this?


texts: Why doesn't QGIS automatically find the new plugin update during startup?

texts: I have done a plugin for QGIS and I made a repository (xml) for the end users to install it. Then, I made a new version and update also my repository file with the new one but Qgis doesn't automaticaly find the new update during startup. Where should I flag the new version? I flag the new version in xml file like that:from <version>1.0.0</version> changed to <version>1.0.1</version> but Qgis doesn't recognize any new version.


texts: Failure of the code calling the function “rasterize” in R: NAs or no slots

texts:
I have some problems to convert a shapefile into a raster. I used the function rasterize (package raster) but all raster values are equal to NA.

Following the answers to this question rasterize produces a raster full of NA values. The raster aligns with the SpatialPolygons:

r.ext <- raster(shp.prj, res=30)
rfn <- rasterize(shp.prj, r.ext, field="land_cover_type")

class(shp.prj$land_cover_type)
# [1] "character"

But I have always NAs in the output raster (rfn).

> rfn
class       : RasterLayer
dimensions  : 1, 1, 1  (nrow, ncol, ncell)
resolution  : 30, 30  (x, y)
extent      : -74, -44, 15.875, 45.875  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=NAD83 +no_defs +ellps=GRS80 +towgs84=0,0,0
data source : in memory
names       : layer
values      : NA, NA  (min, max)

Here are some information on the input shapefile (shp.prj) and the input raster (r.ext):

> shp.prj
class       : SpatialPolygonsDataFrame
features    : 78152
extent      : -74, -71.5, 45, 45.875  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=NAD83 +no_defs +ellps=GRS80 +towgs84=0,0,0
variables   : 3

> r.ext
class       : RasterLayer
dimensions  : 1, 1, 1  (nrow, ncol, ncell)
resolution  : 30, 30  (x, y)
extent      : -74, -44, 15.875, 45.875  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=NAD83 +no_defs +ellps=GRS80 +towgs84=0,0,0

I also tried to reduce the raster resolution (res(r.ext) <- 0.05) but I obtained this error message:

> rfn <- rasterize(shp.prj, r.ext, field="land_cover_type")
Error in cbind(poly, rbind(poly[-1, ], poly[1, ])) :
  trying to get slot "coords" from an object of a basic class ("NULL") with no slots

Following the answers to this question r raster package rasterize failing:no slots, I also tested:

rfn <- rasterize(shp.prj, r.ext, shp.prj$land_cover_type)

but I have always the error message:

Error in cbind(poly, rbind(poly[-1, ], poly[1, ])) :
  trying to get slot "coords" from an object of a basic class ("NULL") with no slots


texts: Using parentheses in query manager of ArcGIS to define query order?

texts:
I would like to make a nested query in a single table by means of parentheses in the Query Builder of ArcGIS.

More specific: I only want to select those Columns which start with TR- and in those I want to select those columns that have the attribute 0 (FieldD) and V and E (FieldB) or own (FieldC).

I used the following Syntax:

SELECT * FROM Table WHERE:
FieldA LIKE 'TR_%' AND FieldD = 0 AND FieldB in ('V','B') OR FieldC = 'own'

the problem is that when I use the Operator OR for fieldC ArcGIS also considers KS-


texts: How to generate raster from vector in QGIS with resolution matching existing raster?

texts:
Is there any way to use open source tools (preferably QGIS) to rasterize a vector to create a raster with the same coordinates and resolution as an existing raster?

I know that the map size can be set using gdal_rasterize, but I'm hoping for something more precise...


texts: How to rectify the walls of buildings in PostGIS?

texts:
I need an operations in PostGIS to eliminate intermediate vertexes in a polygon to rectify the walls of a building.

I was revising http://trac.osgeo.org/postgis/wiki/UsersWikiExamplesSpikeRemover but it only works to eliminate spikes.

example:


texts: Edit array type attributes in attribute table

texts:
I loaded a PostgreSQL layer in QGIS project. The layer has some array type attributes like integer[].
I'd like to edit these attributes in Attribute Table. But it doesn't work, I could not type any value in the table cell.
I also tried changing the Fields Widget of the attribute but it seem doesn't support multidimensional array value.


texts: Running non blocking python script in QGIS?

texts:
I am new to python programming with QGIS.

I'm looking for a solution to start in QGIS a socket client to receive continuously data from a server which sends telegrams with point data if an event occurs. The points should be then visualized and updated on a temporary layer.
In a standalone python script the socket client runs in an endless loop and processes the data until it's interrupted with + or the sever stops.
If I start the script in the QGIS python console, it blocks the application until the server connection breaks.

How do I run such a script in QGIS, which doesn't block the program itself?


texts: ogr2ogr "Error: pointer being freed was not allocated" on mac mavericks

texts:
I've successfully installed GDAL via homebrew on mac

But while I'm using command ogr2ogr, keep getting error below

$ ogr2ogr -f geoJSON tl_2010_20_tract10.json tl_2010_20_tract10.shp

CPLAcquireMutex: Error = 22CPLMutexHolder: Failed to acquire mutex!
ogr2ogr(2439,0x7fff7743a310) malloc: *** error for object 0x7fb7a0d00560: pointer being      freed was not allocated
*** set a breakpoint in malloc_error_break to debug
[1]    2439 abort      ogr2ogr -f geoJSON tl_2010_20_tract10.json tl_2010_20_tract10.shp

and it also happened while using ogrinfo.

Do anyone knows what this error means, and how could I fix it? Thanks.


texts: Querying for All OSM Street Intersections in PostGIS

texts:
I created a GIS database using an OSM datafile for a single city (from Geofabrik) and imported it into PostGIS using osm2pgsql. I have around 180K rows in osm_planet_line database, a manageable sie. Now I am trying to query all intersections in this city. This query

SELECT st_asText(st_transform(st_intersection(p1.way,p2.way), 4326))
FROM planet_osm_line as p1 inner join planet_osm_line as p2
     on st_intersects(p1.way,p2.way)
where 1=1
and p1.osm_id <> p2.osm_id
and p1.name <> p2.name

is taking very long (I waited 30 mins then stopped the query). I know the query works in general, because if I add the following

and st_dwithin(st_transform(p1.way, 4326),st_geomfromtext('POINT(13.4270146842 52.5408961579)',4326)::geography,1000)
and st_dwithin(st_transform(p2.way, 4326),st_geomfromtext('POINT(13.4270146842 52.5408961579)',4326)::geography,1000)

which only gives me intersections within a kilometer of a certain region, then after 4-5 mins I get a result, and they seem to be correct.

Is there a faster way of retrieving all intersections, using a different join perhaps, an additional filter, etc?


texts: OpenLayers 4.6.5: Can't display a Point feature

texts:
I'm starting to use OpenLayers, and I can't use version later than 4.6.5

I'm just trying to display a point with this piece of code, with no success so far.
The map is displayed, is centered where I want, but no point in sight.  No error in console.
Somebody could point me the problem?


texts: Flex popup from outside resource weather warnings

texts:
I am using the USGS service of NHSS NOAA weather warnings in an application.  I have the service displaying in my map using
http://rmgsc.cr.usgs.gov/ArcGIS/rest/services/nhss_weat/MapServer
The service shows that identify is supported.
The question is how do I get a popup from this service with the warning type and other info?
If I try to use the same service I get no popup, If I point it at my own server I get a popup with no information.
I have also tried using the atom described here but not sure of the exact usage in this context.

Any help would be appreciated.


texts: How does OpenLayers GetFeature filter strategy work?

texts:
I want to select on map that "res_id = 182" but that is selecting all objects on the map.

var filter = new OpenLayers.Filter.Comparison({
    type: OpenLayers.Filter.Comparison.EQUAL_TO,
    property: "res_id",
    value: "182"
});

var  filterStrategy = new OpenLayers.Strategy.Filter({filter: filter});

var control = new OpenLayers.Control.GetFeature({
    protocol: OpenLayers.Protocol.WFS.fromWMSLayer(mylayer),
    strategies: [new OpenLayers.Strategy.Fixed(),filterStrategy],
    box: true,
    multipleKey: "shiftKey",
    toggleKey: "ctrlKey"
});

I actually want to select an area on map and filter selected items by attribute from combobox.


texts: Bringing multiple non-projected CAD .dwg files over to ArcMap as a projected ESRI .shp

texts: Is there a streamlined method to bring over CAD .dwg files that have unknown coordinate systems into ArcMap as ESRI .shp files and project the coordinate system of each new .shp accordingly? I have only been able to re-draw polygons or polylines in Auto CAD, project them there, save as .shp and then bring over that way with any luck. It just seems there would be a better method for this. There will be multiple layers and multiple .dwgs to bring over and preserve integrity of.


texts: Space time analysis. Finding spatial patterns (cluster&disperse) over time in point data sets

texts:
I am trying to describe/find statistically spatial patterns in a point data set of industrial sites. The points feature apart from

their location,
a date for their opening,
a date for their closure and
the length of operation, but no other value.

I would like to see at which moment in time certain patterns emerge and disappear, such as sudden dense cluster or dispersal. I investigated the data set 'manually' and found several moments in time, which show the kind of patterns I am looking for (see image 1). I would however like to explore and describe these statistically, also most likely there are relations through out time, which I wouldn't find else wise.

(http://i.stack.imgur.com/HJe7L.jpg)

I am currently working with ArcMap 10.2.2 and have explored some possibilities of space-time analysis. I tried several methods, such as Average Nearest Neighbour, Grouping or Multi-Distance Spatial Cluster Analysis (Ripleys K Function), the difficulty is that these method do not incorporate the time dimension. The best I could find is as described here (http://video.esri.com/watch/1681/spatial-pattern-analysis-mapping-trends-and-clusters) a time space analysis using a Spatial Weight Matrix (Although here I have to introduce a bias in the length of relation [how long and how far away events should be so that the analysis considers them as related]), in combination with a HotSpot analysis, followed by 3 dimensional representation.

The problem I have with this is, that my data set doesn't feature a value or magnitude field against which I could apply a HotSpots analysis. The only values I have are Dates, or respectively time length, which I can not apply a hot spot analysis on.

Does anyone of you maybe has a good idea to analyse such a point data set?


texts: ArcGIS Collector offline synchronous update

texts:
When I work offline in the ArcGIS Collector application, the changes made to the data I collect (including attachments to this) are not updated in the service when I say synchronize the data.
The Collector says that the data changes are being made successfully and the data is gone. However, I cannot see these changes in my map in the service.
How can I solve this problem?


texts: Making simple, web based map without uploading data to 3rd party websites?

texts:
I have a lot of point data (Lat, Long) that I want to display on a map. I want the map to be accessible to other users on my network, ideally via a web browser. However, I don't want to (read will not be allowed to) upload this data to ArcGIS online or any similar website.

I would also like to be able to customize the points, have them showing a little pop up window (giving information on each point - number of samples etc) and use clustering (because some points have been studied several times and will have different records at exactly the same coordinates).

I need a free (or very cheap) solution, and I have no programming skills. I've tried the followwing -

ArcGIS Viewer for flex - looked good, but the only way to add local data is via a widget, but then you can't customize the display or add pop ups etc?

ArcGIS Explorer - Again this looked really good, but I couldn't seem to get any control over the pop up information or how the points displayed

ArcGIS for SharePoint - this was the best one, but doesn't seem to be able to do the clustering.

I know I'm asking for quite a lot of free software here, but is there anything that will do it? The ArcGIS API for Silverlight looks like it would do the job, but I'm not sure I'm smart enough to use it. I have access to ArcGIS 10 Desktop, but not ArcGIS Server, is there anything I can export from here? Or have I missed something obvious somewhere?


texts: How to provide meanings (metadata) of attributes/fields with shapefiles?

texts: Usually when doing GIS projects, especially when working as a freelancer, you're just provided with a bunch of shapefiles and you have to figure out yourself what all the table attributes mean. Because shapefiles have a *.dbf limited character number for naming attribute fields, the names are usually abbreviations. Is there a way to add a description of the attribute field names when creating them, or how do you add attribute field metadata. What do you do when you don't understand the abbreviations of a shapefile attribute table metadata that has been created by someone else?


texts: Problems with installing the GeoServer Metadata extension

texts:
I have recently tried to install the GeoServer Metadata extension. In the installation instructions indicates that any required configuration must be made in the servlet container, but what changes must be made?

I downloaded the extension from the GeoServer Download Page, the file name is called geoserver-2.17-SNAPSHOT-metadata-plugin.zip, I extracted these files and placed the JARs in WEB-INF/lib, but when I start the service the following message appears

  HTTP 404 Status - Not Found
  Type Status report
  message The required resource [/ geoserver] is not available
  Description The required resource is not available.

I've already done several changes to the code in web.xml file

<servlet>
      <servlet-name>metadata</servlet-name>
      <servlet-class>org.geoserver.community.gs-metadata</servlet-class>
 </servlet>

<servlet-mapping>
     <servlet-name>metadata</servlet-name>
     <url-pattern>/webapps/geoserver/WEB-INF/lib/*</url-pattern>
 </servlet-mapping>

The installation I have is Windows 10, GeoServer 2.17-SNAPSHOT, Apache-Tomcat/8.5.55, Java jre1.8.0_231.


texts: How do I get the population by country when using the Esri Demographics Geoprocessing Service?

texts:
I am using the ArcGIS JavaScript API and Geoprocessing Services to get population. The service takes an input circle and returns the population inside.

Some times, a circle can cover cover more than one country (for example Spain and Portugal) but I need to get the population of a single country using a filter.

I am using this function to get the area using Google Maps:

function getPop()
myFunctions.ftn_findpop1(selectedShape.center.lat(),selectedShape.center.lng(),selectedShape.getRadius());
}
Then...
function findPop(ring) {
        var gp = new esri.tasks.Geoprocessor("https://sampleserver1.arcgisonline.com/ArcGIS/rest/services/Demographics/ESRI_Population_World/GPServer/PopulationSummary");
...


texts: Comparing several values within each cell of a DEM

texts:
I am working with a DEM and got a GPS line of my working area.
I converted the GPS line into single measurement points. Now I want to compare each of these GPS points with the corresponding elevation value of my DEM.

I tried to do this manually by hand but working with several thousand points makes it nearly impossible to accomplish that.
I thought about calculating average values of the GPS values within a DEM cell, but what if the elevation differs significantly within that specific cell?

Is there a tool either in ArcGIS or QGIS which can help me comparing my values?


texts:
I would like to add a line after the transformation of my shapefile in order to save it in a temporary file.

uploadShpfile <- reactive({
    if (!is.null(input$shp)) {
      inFile <- input$shp
      inFolder <- substr(inFile$datapath, 1, nchar(inFile$datapath) - 5)
      unzip(inFile$datapath, exdir = inFolder)
      shpDF <- input$shp
      pwd <- getwd()
      updir <- dirname(shpDF$datapath[1])
      setwd(updir)
      for (i in 1:nrow(shpDF)) {
        file.rename(shpDF$datapath[i], shpDF$name[i])
      }
      shpName <- shpDF$name[grep(x = shpDF$name, pattern = "*.shp")]
      shpPath <- paste(updir, shpName, sep = "/")
      setwd(pwd)
      shpFile <- readOGR(shpPath)
      shpFile <- spTransform(shpFile,CRS("+proj=longlat +datum=WGS84"))

Any ideas ?


texts: How to compute the outline of multiple tracks?

texts:
I have multiple gps tracks with depth data gathered from the lake. The tracks are overlapping/intersecting. I want to get an outline polygon.
I have a multi-line shp. I've tried to diffuse it in qgis and buffer. But buffering amounts of data I have simply stalls the PC, after eight hours of work buffer was still processing data.

Is there a simpler and more productive way to get outline from multiple tracks? (Preferably using qgis, gdal, or any other free software)

This is what I have:

This is what I need to get:

Another part of the goal is to preserve the precision of vector data. And hopefully the solution should work for large majority of cases automatically. And it probably shouldn't involve rasterizing image to keep memory use as low as possible.

It should be noted that simple approach to do some kind of approximation causes the solution to be sensible to data. I.e. each data set will have to be treated individually to remove errors related to approximation.

EDIT: Images updated


texts: ArcMap ModelBuilder: Output of Raster Iterator only runs one Output

texts:
I built a model in ArcMap 10.3.1 to quickly repeat a process with multiple raster images stored in the same folder.

The idea was to have the model iterate through the rasters , then have the pixels in all of the rasters be multiplied by a certain constant value using Times, after which the model would run those output rasters through the Zonal Statistics Function and output a Final Table. I included the Parse Path step in there because I was having trouble running the Zonal Statistics step with the .tif at the end of the filename when I was simply using %Name% for my inline variable substitution.

I was having a lot of issues with this model, and finally managed to get it to (sort of) work, except now, instead of iterating through all the rasters, it only creates an output table for the first file in the "Actual" folder where the input data is stored. Instead of saying "Raster" on the output bubble from the iterate rasters step, it just says the filename for that first file, "A010216.tif", and only generates its output at the end.

How do I solve this?

I managed to fix the model in two different ways: (1) By removing the iterator step and just having the model accept a single raster as an input, it worked just fine (2) for some reason adding a Set Null step to change zero values to NoData/Null worked really well, and it let me use the Iterator without any trouble.

Here's a picture of solution number 2 (I also added a couple more steps to the model and those worked fine as well!


texts:
I am using ArcMap 10.3.

I am trying to create a new raster from 12 pre-existing rasters using the Raster calculator and I am having some trouble. You see, I want to create a new raster that represents the Precipitation Concentration Index, which uses the following formula:

PCI = (Summation (p^2) / (Summation (p)^2)) * 100;
where p = monthly precipitatoin average

I have all the rasters for world precipitation means by month, and hence I am using the following formula on Raster Calculator:

((("prec_1" ^ 2) + ("prec_2" ^ 2) + ("prec_3" ^ 2) + ("prec_4" ^ 2) + ("prec_5" ^ 2) + ("prec_6" ^ 2) + ("prec_7" ^ 2) + ("prec_8" ^ 2) + ("prec_9" ^ 2) + ("prec_10" ^ 2) + ("prec_11" ^ 2) + ("prec_12" ^ 2)) / (("prec_1" +  "prec_2" + "prec_3" + "prec_4" + "prec_5" + "prec_6" + "prec_7" + "prec_8" + "prec_9" + "prec_10" + "prec_11" + "prec_12") ^ 2)) * 100

As a result, I was expecting something like this:

But instead I got this:

I have tried changing the symbology of the raster layer to "Streched" and adjusting the values to a suitable scale (they range between 0 and 50) following the example of the monthly precipitation layers - but then the map is just blank.

Anyway, since I do not know much about coding, I was hoping someone could  help doing this the simplest way possible.

This is the exact message I got from ArcMap:

Messages
Executing: RasterCalculator "((("prec_1" ^ 2) + ("prec_2" ^ 2) + ("prec_3" ^ 2) + ("prec_4" ^ 2) + ("prec_5" ^ 2) + ("prec_6" ^ 2) + ("prec_7" ^ 2) + ("prec_8" ^ 2) + ("prec_9" ^ 2) + ("prec_10" ^ 2) + ("prec_11" ^ 2) + ("prec_12" ^ 2)) / (("prec_1" + "prec_2" + "prec_3" + "prec_4" + "prec_5" + "prec_6" + "prec_7" + "prec_8" + "prec_9" + "prec_10" + "prec_11" + "prec_12") ^ 2)) * 100" C:\Users\Documents\ArcGIS\Default.gdb\rastercalc
Start Time: Tue Feb 11 10:56:15 2020
(((Raster(r"prec_1") ^ 2) + (Raster(r"prec_2") ^ 2) + (Raster(r"prec_3") ^ 2) + (Raster(r"prec_4") ^ 2) + (Raster(r"prec_5") ^ 2) + (Raster(r"prec_6") ^ 2) + (Raster(r"prec_7") ^ 2) + (Raster(r"prec_8") ^ 2) + (Raster(r"prec_9") ^ 2) + (Raster(r"prec_10") ^ 2) + (Raster(r"prec_11") ^ 2) + (Raster(r"prec_12") ^ 2)) / ((Raster(r"prec_1") + Raster(r"prec_2") + Raster(r"prec_3") + Raster(r"prec_4") + Raster(r"prec_5") + Raster(r"prec_6") + Raster(r"prec_7") + Raster(r"prec_8") + Raster(r"prec_9") + Raster(r"prec_10") + Raster(r"prec_11") + Raster(r"prec_12")) ^ 2)) * 100
Succeeded at Tue Feb 11 11:03:31 2020 (Elapsed Time: 7 minutes 16 seconds)


texts: QGIS Print Composer doesn't save layout of attribute box

texts:
I've added a square map to a landscape page and to the side of it I have a multi column index of street names. All is good with how it looks.
I use the "Extend to Next Page" option for the index attribute table, under the frames menu and when the table gets extended onto multiple blank pages I drag the tables up and onto the side of the map in the first page to get multiple columns of index text. This is great until I save and re-open.
When I open the layout again, all my attribute segments have expanded back to being on multiple pages and I have to reorganise them again.
I've tried locking the layers, using every other option in the Frames menu, tried all the other items in the Item Properties box that have anything to do with the attribute table. But nothing works. Can someone point me in the right direction?
I've had to devise a convoluted way of creating this index which is all fine right up until this last step. If there's a particular work around or fudge, that's fine.


texts: Turning shapefile into mask and calculating mean

texts:
I have been using Python to process huge arrays stored as NetCDF files. I would like to calculate the mean of an area defined by a shapefile.
I have just installed GDAL but if there are other tools I should use please let me know.
So my main question is how can I turn my shapefile into a mask?
Then use my mask to calculate the NumPy mean of my array within the area of the mask?
shp="country.shp"
array=mynumpyarray
mask=useGDALtochangeshapefiletomask
meanofarea=N.mean(array,limits=mask)

I have been able to convert my shapefile to a raster using gdal_RasterizeLayer link supplied and then a masked array using
maskarray=mask_ds.GetRasterBand(1).ReadAsArray()

Now it is just a question of matching my NetCDF extent with my raster/maskarray. Here are my extents as requested.
my shapefile/raster extent: x_min, x_max, y_min, y_max 140.962408758 149.974994992 -39.1366533667 -33.9813898583

my netcdf file extent: min longitude: 139.8 max longitude: 150.0 min latitude -39.2 max latitude: -33.6 LAT size 106 LON size 193


texts: Getting Annotation Text Geometry Only (without leader line) as Polygon Using ArcObjects?

texts:
I have an annotation feature that also has a line callout (see picture).  I need to use ArcObjects to get the polygon that surrounds ONLY the text (not the line callout).  The following link is very similar to the problem I am having.

Getting Annotation Text Geometry as Polygon Using ArcObjects?

However the solution is to use the QueryOutline.  The problem with the QueryOutline method is that I must pass a IDisplay parameter.  When I pass the IDisplay of the mxd, the polygon that is returned has a different size based on the zoom level of my mxd.

So if I'm zoomed in at 1:1000 the polygon that is returned is bigger than when my mxd is zoomed in at 1:2400.  How do I know what level to zoom in so I can get the exact size of the polygon?  Or if there is another to only get the text polygon that would work for me as well.


texts: Deploying localization files for Geoprocessing service in ArcGIS Server

texts:
I have a Python geoprocessing tool, it is split in several Python files, When I use  ArcGIS Pro to publish it as a service into the ArcGIS Server, all the dependencies are packed and deployed to the service.

Now the tool is internationalized and localized using gettext and have the following structure:

./Alignment.py
./locale/Spanish/LC_MESSAGES/Alignment.mo
./locale/Spanish/LC_MESSAGES/Alignment.po
./locale/German/LC_MESSAGES/Alignment.mo
./locale/German/LC_MESSAGES/Alignment.po

The Python files are deployed but neither *.mo nor *.po files are deployed.

How can I make that those files are included in the set of files deployed to the ArcGIS Server?

[UPDATE]
If In the file Alignment.py I make a horrible hack like this i:

locale_path = os.path.join(os.path.dirname(os.path.realpath(__file__)), "locale/Spanish/LC_MESSAGES/")

the folder and its content are added to the service, but the folder structure changes a lot:

MyService.sd
 |-p20
     |-Service
     |   |-Alignment.py
     |
     |-lc_messages
        |-Alignment.mo
        |-Alignment.po


texts: Mark latitude and longitude values on Google Earth

texts: I have a csv file of latitudes and longitudes. Is it somehow possible that all these lat-long values can be marked on Google Earth?


texts: How to fix Quickfinder OperationalError: table quickfinder_info already exists?

texts:
This is what I am getting when I try to search in project layers and create a .qfts file to use quickfinder on a local shapefile.

An error has occured while executing Python code:

  Traceback (most recent call last):
    File "C:/Users/xxxxxx/.qgis2/python/plugins\quickfinder\gui\configurationdialog.py", line 120, in createQFTSfile
      createFTSfile(filepath)
    File "C:/Users/xxxxxxxx/.qgis2/python/plugins\quickfinder\core\projectfinder.py", line 55, in createFTSfile
      cur.executescript(sql)
  OperationalError: table quickfinder_info already exists

  Python version:
  2.7.5 (default, May 15 2013, 22:44:16) [MSC v.1500 64 bit (AMD64)]

  QGIS version:
  2.8.2-Wien Wien, 1b929ef

  Python path: ['C:/Users/xxxxxxx/.qgis2/python/plugins\processing', 'C:/PROGRA~1/QGISWI~1/apps/qgis/./python', u'C:/Users/xxxxxxx/.qgis2/python', u'C:/Users/xxxxxxxx/.qgis2/python/plugins', 'C:/PROGRA~1/QGISWI~1/apps/qgis/./python/plugins', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\matplotlib-1.3.1-py2.7-win-amd64.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\nose-1.3.3-py2.7.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\tornado-4.0.1-py2.7-win-amd64.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\backports.ssl_match_hostname-3.4.0.2-py2.7.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\certifi-14.05.14-py2.7.egg', 'C:\PROGRA~1\QGISWI~1\bin\python27.zip', 'C:\PROGRA~1\QGISWI~1\apps\Python27\DLLs', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\plat-win', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\lib-tk', 'C:\PROGRA~1\QGISWI~1\bin', 'C:\PROGRA~1\QGISWI~1\apps\Python27', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\PIL', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\jinja2-2.7.2-py2.7.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\markupsafe-0.23-py2.7-win-amd64.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\pytz-2012j-py2.7.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\win32', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\win32\lib', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\Pythonwin', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\Shapely-1.2.18-py2.7-win-amd64.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\six-1.3.0-py2.7.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\wx-2.8-msw-unicode', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\xlrd-0.9.2-py2.7.egg', 'C:\PROGRA~1\QGISWI~1\apps\Python27\lib\site-packages\xlwt-0.7.5-py2.7.egg', u'C:/Users/xxxxxxxxx/.qgis2//python', 'C:\PROGRA~1\QGISWI~1\apps\qgis\python\plugins\fTools\tools']


texts: ArcPad 8.0 and Postgis Data - Viewing & Editing Problems!

texts:
I have data held on a postgis which I want to take out onto site to update.

I exported data from Postgis using dB manager in QGIS 3 into a file which I then transferred it onto my tablet for site use.  My tablet has ArcPad 8.0 which has a project set up to read and project the above mentioned data.

When viewing the data on arcpad however it is extremely small, when zooming in the size does not change. Also I can no longer change the symbology and no longer toggle editing as the tick box for this function has disappeared.

This exercise was never an issue as I used to store my data locally and used ArcInfo for all my editing/viewing/querying/exporting

We now have been asked to store everything on a Postgis and start using QGIS.  I tried exporting the data from postgis using Arc however this did not seem to work either. Only older shape files stored locally worked, so I know it's not ArcPad. It seems as though since uploading all my data onto Postgis is doing something to the data not just storing it!


texts: Python FindLabel Function for Label Expression with Formatting Font?

texts:
I originally had this FindLabel function (which works)...

def FindLabel ( [NAME], [Count_]  ):
  if [Count_] != None:
    return [NAME] +'\n' + "Well Log Count: " + [Count_]

I wanted to format the font of the second line only, after '\n', using formatting tags. This is what I tried, and it did not work...

def FindLabel ( [NAME], [Count_] ):
  if [Count_] != None:
   return"[NAME] +'\n' +  "<FNT name='Arial' size='13'>" & "Well Log Count: " + [Count_] & "</FNT>

I also tried to set a variable, and this did not work either...

def FindLabel ( [NAME], [Count_] ):
  myCount = "<FNT name='Arial' size='13'>" &  "Well Log Count: " + [Count_] & "</FNT>
  if [Count_] != None:
   return"[NAME] +'\n' +  myCount

Any thoughts?


texts: How can I determine what fonts are used by an mxd?

texts:
A client has sent me an mxd and a file gdb.

When I look at the symbols, they don't look good.  I think they need to send me some fonts.  Is there some way I can use arcmap to determine the name of the font file I need so that I can request it from them?


texts: Thematic map exported using QGIS2Web is extremely laggy

texts:
I have made a world thematic map and added layers of data for creating graduated color scheme on the country level. I need a interactive map as an output which shows values when user hover over any country.

So while I am exporting it in Leaflet format using QGIS2web it is not working properly with an extreme lag. Also my data is not that big (final output is only 2.5 MB's and only around 500 features in total).

So is there anyway I can remove this lag in the Leaflet format, as the other format (OpenLayers 3) cannot be used as it is not showing legends which are a must for this map?


texts:
I am looking to match both the projection and extent of two separate .tif files. The file I want all other .tif files introduced into raster calculations is as follows:

The other raster image covers the globe, however I am only interested in matching the extent of the first image, this second raster image looks like the following:

I attempted to reproject both rasters to the WGS 84 coordinate system (EPSG:4326) through the following GDAL code, however when I placed both rasters into QGIS to verify they were matched up, the second image was in the bottom left hand corner of the first, which verified my suspicion the code did not work as intended.

gdalwarp -t_srs EPSG:4326 -overwrite gfs.tif gfstest.tif

Other info I can provide from this is the layer extent of both images, which I hope can provide some clues:
Image 1:-1269.8514999999999873,-1269.8515000003390014 : 5446393.0834999997168779,3495901.1794999996200204

Image 2: -0.1250000000000000,-90.1250000000000000 : 359.8750000000000000,90.1250000000000000