texts: Editing ArcFM ArcSDE table in version with ArcPy?

texts:
I am trying to edit data that is in an ArcSDE Oracle database with updates in a file geodatabase. I want to be able to make these changes on a separate version to isolate them.

If I connect to a version created through ArcMap 10.2.1, I am able to edit and apply changes and they can be seen through the GUI. I want to be able to create a version through the python script and then apply the changes to that version.

However, when the code invokes the edit = arcpy.da.Editor(myconnect) command it returns the error:

edit = arcpy.da.Editor(myconnect)
RuntimeError: The database was not found.

There doesn't seem to be much I can find on the web about this error. The code below is how I create the version but is currently commented out to simplify the debugging process.

Initially I would connect without trying to connect to the version directly and as annoying as it seems, when editing each feature, I would create a layer for that, then change the version for that layer. However, I was unable to get that to work and then thought I could just connect to the version directly. So much for that idea! The objects are ArcFM features, I wouldn't think that has anything to do with it but as I've found out before, I could well be wrong!

Here is the abbreviated code I have:

app = Dispatch('Miner.Framework.Dispatch.MMAppInitializeDispatch')
au = Dispatch('Miner.Framework.Dispatch.MMAutoupdaterDispatch')
runtime = Dispatch('Miner.Framework.Dispatch.MMRuntimeEnvironmentDispatch')
runtime.RuntimeMode = mmRuntimeMode.mmRuntimeModeArcServer
app.Initialize(mmLicensedProductCode.mmLPArcFM)
au.AutoUpdaterMode = mmAutoUpdaterMode.mmAUMNoEvents
try:
    versionBase="script-test"
    versionName="%s.%s" % (user, versionBase)
##    myconnect = sdeConnection.connect("", "FMSQA", user, password,"")
##    found = False
##    for version in arcpy.da.ListVersions(myconnect):
##        if versionName in version.name:
##            found = True
##            break
##    if (found == False):
##        arcpy.CreateVersion_management(myconnect, "SDE.DEFAULT", versionBase, "PUBLIC")
    myconnect = sdeConnection.connect("", "FMSQA", user, password, versionBase)
    arcpy.env.workspace = myconnect
    edit = arcpy.da.Editor(myconnect)

My database connection code is basically taken from here https://gis.stackexchange.com/a/17028/14120


texts: I have some points overlayed on top of a raster.  Is there a tool in ArcGIS Desktop to add the values of the raster that the points fall on top of into the point attribute table?


texts:
I'm working on a Qgis plugin with python,

I'm creating a feature and I show a feature form to edit the attributes. But after the openFeatureForm function call the program crash and close

line = QgsFeature()
line.setGeometry(QgsGeometry.fromPolyline(transformedLine))
pr = parcoursLayer.dataProvider()
pr.addFeatures([line])
qgis.utils.iface.openFeatureForm(parcoursLayer,line)


texts:
I have a bunch of jp2 files with size going from 40MB to 1GB and I need to convert them all to GeoTiffs (.tif format). When using gdal_translate for the small files (having the jasper driver to support jp2 files), this works. But when applying this on files with size about 400MB, I get segmentation errors.

gdal_translate input.jp2 output.tif

In another attempt I've used opj_decompress (from the JP2OpenJPEG library) to convert jp2 to tif files, but again resulted in an error.

opj_decompress -i input.jp2 -o output.tif

Is there a way to handle large jp2 files with GDAL (or another library) to convert them to .tif?


texts:
i use the SQL-API as described in here:

And i find the following issues undocumented:

How many requests per second can i send to the SQL-API, before the user gets slowed down/blocked? I would like to throttle myself and also make calculations, how long those requests will take.
How many statements can be put in a batch? Or is this purely a limit on the POST request? Whatever fits in the POST request, will be executed?

Does anybody know or can clarify?


texts: Converting ASCII to Raster then building attribute table?

texts:
I recently obtained precipitation normals from the PRISM dataset online. I'm trying to import the ASCII files into ArcMap, but have stumbled across a few issues.

I converted the ASCII to raster. I did it first as float, then as integer.
The two rasters (float vs integer) had identical values, except the float raster values went to the hundredth decimal place, while the integer raster was just whole numbers. Not a problem; the ASCII is floating point anyway.

My problem then comes in with the attribute table. The float raster does not have an attribute table, while the integer raster does.

Can someone help me understand why the float raster didn't create an attribute table, and how I can build one?

My GIS knowledge is rather limited.  I'm using ArcMap 10.2.


texts: How to invoke Special Auto Updater via CustomAutoUpdater in ArcFM

texts:
I am trying to invoke SpecialAU(Auto-Updater) ArcFM with below mentioned code, but I am not able to trigger it.
I took help from
http://gisstudycenter.blogspot.com/2015/01/arcfm-auto-updaters-to-update-shape.html
using ESRI.ArcGIS.Carto;
using ESRI.ArcGIS.Geodatabase;
using Miner.ComCategories;
using Miner.Interop;
using System;
using System.Collections.Generic;
using System.Runtime.InteropServices;
using System.Text;

namespace DesktopClassLibrary2
{
    [ComVisible(true)]
    [Guid("fecebee8-67a5-47ad-a1a7-314ba50b2c9c")]
    [ClassInterface(ClassInterfaceType.None)]
    [ProgId("DesktopClassLibrary2.Class1")]
    [ComponentCategory(ComCategory.SpecialAutoUpdateStrategy)]
    public class Class1 : IMMSpecialAUStrategyEx
    {
        void IMMSpecialAUStrategyEx.Execute(IObject pObject, mmAutoUpdaterMode eAUMode, mmEditEvent eEvent)
        {
            int n = pObject.OID;
            IObjectClass objClass = pObject.Class;
            IFeatureLayer fl = new FeatureLayer
            {
                FeatureClass = objClass as IFeatureClass
            };
            eAUMode = mmAutoUpdaterMode.mmAUMArcCatalog;
        }

        string IMMSpecialAUStrategyEx.Name => "class1";

        public bool get_Enabled(IObjectClass pObjectClass, mmEditEvent eEvent)
        {
            return eEvent == mmEditEvent.mmEventFeatureCreate;
        }
    }
}


texts: Finding the area of polygons within other polygons in R

texts: Is there a function in R which will loop through a layer of zonal polygons (i.e. a polygon grid) and return the area of a second polygon layer that is contained in each grid cell? I have a layer of polygons representing areas of forest (in blue), and a second layer representing ordinance survey grid cells (in pink - each grid cell has a unique "name"). Both layers are polygon shapefiles. I would like to calculate how much forest is present in each grid cell but I cannot find a tool that will do this in R. I'm looking for something equivalent to tabulate area in ArcGIS. As you can see, the forest polygons are complex and in many cases there are several of them per grid cell, which is making it trickier.


texts: How to add a URL as a QGIS Favorite in browser

texts:
I have GeoJSON/TopoJSON data stored in the cloud (on my S3 bucket).

Currently I have to instruct users to use the entire URL to add a JSON file to a project.

I would like to save the main portion of the URL as a favorite, so I can instruct users to use only the actual JSON file name. The "favorites" in the Browser Panel only accept paths to local data layers.

Is there a way in QGIS to save a URL as a favorite in the Browser Panel?


texts: QGIS 2.18.2 - Digitizing Error geometry type is not compatible with the current layer

texts:
I have an issue with digitizing lines in QGIS 2.18.2. I was using cadtools and after digitizing lines in my Lines layer on saving i get the message:

--Could not commit changes to layer (insert layer name)

Errors: ERROR: 1 feature(s) not added - geometry type is not compatible with the current layer.--

I have digitized a line in a line layer - that is certain.

I have tested without using the Cadtools and it appears that even if i use the snapping tool to digitize to a feature it throws the above error but if i just digitize a random line not using snapping then it is able to save.


texts: Using osm2po to route using places (saved as multipolygon)

texts:
I use osm2po to create a routable network from OSM data and this works great.
Although I feel quite familiar with osm2po nowadays, there is one thing I can not solve using PgRoutingWriter's output:
routing using mulipolygon relations.

What I am trying to achieve is taking places/areas that are tagged with area=yes and foot=yes into account.

A location where this comes into play is the Waltherplatz in Bozen/Bolzano (but I am also sure that there are a ton of others as well)... ways end at the border of the place, one side is even bound by a street, but even when routing with type "foot" I can not get from one side to the other.

Can I somehow tell osm2po that it should use relations as well (not only ways) when creating my routable network (uncommenting the relationTagResolver.class in the config file does not seem to be enough)?

A very old question about a similar topic can be found here, but since it is 5 years old I guess it is a bit outdated (hopefully).

I know that this is still a problem that is under research, but I would be happy with a very basic solution (something like "walk along the outer line of the multipolygon"). Right now I do not need a better solutions like building all possible routes across a place...


texts: I started to use GeoServer, I uploaded shapefile and I want to publish it but I don't want to show all of the "Feature Type Details". For example, I want to publish only geometry. I know it's possible with SQL Views, but what about simple shapefile?


texts: Import ol-ext to a React - Open Layers 6 project

texts:
I'm trying to set up a React project where I want to use some of the tools provided by the ol-ext. I looked around and found some Codesandbox projects from the creator of the library but I was unable to make them work properly because of the way the library is imported in the projects causes an error.
Is there a problem with the syntax or incompatibility with ol-ext and the other Open Layers versions?


texts: Reading default qgis-auth.db from a script

texts:
I try to read username and password from the qgis-auth.db, the latter having been configured under QGIS (3.6.1, Ubuntu).

Here is my code :

import os
import sys
from qgis.core import QgsApplication

argvb = list(map(os.fsencode, sys.argv))
qgisapp = QgsApplication(argvb, True)
qgisapp.initQgis()

authM = qgisapp.authManager()
authM.setMasterPassword("myMasterPassword", True)  # returns True
print(authM.authenticationDatabasePath())
# returns '/home/me/.local/share/profiles/default/qgis-auth.db'

it's not the QGIS3 default profile location, which is :

/home/me/.local/share/QGIS/QGIS3/profiles/default/qgis-auth.db

I'd like to avoid copying the auth database and read it directly from the default profile.


texts:
I have two Lands Surface Temperature images, one from Landsat 7 and the other one is from MODIS, both for the same day. After projecting MODIS to UTM WGS84 the two images are slightly shifted, i.e. they do not well overlap.

Is there any suggestion to fix the matter?


texts:
Using QGIS, after run R script like Histogram it returns an empty result in "Results Viewer" window and no image are created inside html file.
Any ideia of missing configuration or package install? By the way, I'm unable to update R scripts.
(QGIS 2.18.21, Linux Mint Cinnamon)


texts: White line appears between the tiles using ImageMosaic / ImagePyramid technique on GeoServer

texts:
I have created image Pyramid of 4 levels using gdal retile and published it on GeoServer using ImagePyramid technique. When I am accessing it in QGIS, it seems to work fine but at the zoom level from 1:2000 to 1:900. I am getting some white line between the tiles.

also, I have published each folder of pyramid separately on Geoserver using image mosaic technique. when I am accessing these layers in Qgis. in folder layer 0, 2 and 4 I haven't found the white lines between tiles at any zoom level.
 but when I am loading folder layer 1, whites lines are appearing between the tiles at the same zoom level.
I have also tried the "Input/output Transparent Color" set to white but it's not helpful.

Can anybody suggest why white lines appear and how to remove it?


texts: Does GeoServer and/or MapServer have a services directory similar to ArcGIS's services directory?

texts:
Does GeoServer or MapServer have anything analogous to the ArcGIS Server services directory? Is there a convenient endpoint where you can see all the services hosted on a particular machine?

Edit:  I'm looking for the GeoServer or MapServer equivalent of this:  http://sampleserver3.arcgisonline.com/ArcGIS/rest/services


texts: Identifying layers from multiple ArcGIS Server instances using ArcGIS API for JavaScript 2?

texts:
This applies to ArcGIS Server 9.3 with the JavaScript API 2.4.

My goal is to click on the map and have the Identify task run on all visible Dynamic layers.

IdentifyTask requires a URL to the REST endpoint, with IdentifyParameters specifying which layerIds the identify operation should be performed on.

It seems the Identify tool is expecting all layers to be accessible from the same REST endpoint (ie, the same ArcGIS Server).

In my case, layers are being served from multiple ArcGIS Server instances - how can the Identify tool support this? For example, if the Buildings and Parcels layers in this map were coming from separate ArcGIS Servers.

(This is related to my previous question, but I realise now I need to answer this question first)


texts: Is it possible to change the default raster style of newly created rasters in QGIS?

texts:
If I do calculations on rasters with the raster calculator or other tools, the resulting rasters displayed in my canvas are always in the same default style (singleband gray, black to white...) so that I have to change it afterwards to my preferred style. Being able to change that default raster style to my own choice would save me a lot of clicks.

@Joseph:

    Traceback (most recent call last):
      File "C:/Users/Gumbo/.qgis2/python/plugins\multiqml\multiqml.py", line 98, in on_pbnApplyStyle_clicked
        self.myPluginMessage( QCoreApplication.translate("MultiQmlDlg", "Unable to apply vector qml style \"%1\" to raster layer \"%2\".")\
    AttributeError: 'unicode' object has no attribute 'arg'


texts: Using line feature as stops input in network analyst?

texts:
I'm trying to build a certain route with network analyst on a Network Dataset. As an input to create the route I want to use line features. I am trying to do this with network analyst in ArcGIS 10.0 to connect these lines into the shortest route. But this seems to give me only the opportunity to insert stops as the network analysis objects.

I tried to use points for the stops, by making a point on both ends of those separate line features and add them as the input locations. But then I get in some cases an insufficient output because one of the lines which I want in my route is not included. This is for example because the points are in a square and the longest line is not included in the route then because I search for the shortest route, but actually I want that line in. I also tried to find a tool to make instead of points on both ends of the lines a lot of points on the line with a very small interval to be shore the lines that I want to be in the route will be included in the route by the network analyst tool. But I cannot find such a tool for ArcGIS 10.0.

Can anybody give me some solutions for this problem?

Regarding the last comment of Chris W, I made a new question: Is it possible to make stops not in a certain sequence in order to find the shortest route with the network analyst?


texts: Google Earth Engine: precipitation data, summing monthly data is not matching annual sum

texts:
In essence, I used GRIDMET to pull both monthly and annual precip over the same area, same time, and same database, but summing each month for the same water year is not matching:

NEW NOTE: As pointed out in Answer#1, running the code using a point is fine; the issue is running the code over an area. When troubleshooting, create a geometry poly and then run function over the poly instead of point.

Approach: I created a sequence of years and then used ee.Date.fromYMD and the built-in .advance function to pull GRIDMET 'pt' data as an annual mean (for water years [Oct-Sept]) over user-defined areas (i.e., HUC6):

This works; however, I also needed to pull the same data as monthly means over the same areas; i.e., I changed .advance(1, "year") to .advance(1,"month"). The sum of the monthly data over the same water year should, theoretically, be identical to the annual mean (same dataset, same area, same timeframe). However, they are different. The pattern is very similar and co-varies; but why would the 'pt' be different? I need to know if the difference is due to the built-in code (e.g., does "month" correctly advance each month 30, 31, 30, 31 days?) or whether it is on the database side (TERRACLIMATE/GRIDMET).

Code: (create a geometric area in GEE instead of using a point)

var point = ee.Geometry.Point (-110.93, 32.28);

var climateYEAR =  ee.ImageCollection(ee.List.sequence(1994, 2000).map(function (year){
  var date_start = ee.Date.fromYMD(year, 1, 1);
  var date_end = date_start.advance(1, "year");
  return  ee.ImageCollection('IDAHO_EPSCOR/GRIDMET')    //changed to/from IDAHO_EPSCOR/TERRACLIMATE
                        .select(['pr'])
                        .filterDate(date_start, date_end)
                        .filterBounds(point)
                        .sum()
                        .set({year: year, 'system:time_start':date_start});
}));

var climateYEARmean = climateYEAR.map(function(image) {
  return image.reduceRegions({
    collection: point,
    reducer: ee.Reducer.mean(),
    scale: 4000,
    crs: 'EPSG:4326',
    tileScale: 1
  }).map(function(f){return f.set({'year': image.get('year')})});
  });
print(climateYEARmean, 'climateYEARmean');
var climateMEAN = (ee.FeatureCollection(climateYEARmean).flatten());
print(climateMEAN, 'climateMEAN');

Also, I likely should post as a different question, but while troubleshooting, I pulled the same data from TerraClimate and found that TC is computing much lower annual means for the same area, same timeframe (and they also covary very well)??


texts: Iterate over multiple shapefiles, changing year from "2014" to "2015" in file name

texts:
I have about 50 shapefiles  I would like to change each of their file names which contains "2014" to "2015". I would like to use a python script but have not been able to figure out how to do it. I have tried using the sample script from Esri but to no avail. I've tried to use other examples I've found Looping through list of shapefiles to rename them using ArcPy? but unable to work it into my problem.

import os
import arcpy
arcpy.env.workspace =
fcs = arcpy.ListFeatureClasses(fc)
if fc.contains("4"):
  arcpy.Rename_management("5")


texts: Using GDAL ogr2ogr correctly to generate vector MVT tiles

texts:
I am a bit lost here. I am trying to convert JSON to MVT tiles for Mapbox (directory format) using ogr2ogr. My problem is that no matter what I do on the Mapbox/browser side, the tiles do not show up. I am thinking I am generating the tiles wrong.

What I have is weather radar data, which looks like this in the JSON file :

   {
    "type": "Feature",
    "geometry": {      "type": "Polygon",
      "coordinates": [
        [
          [-86.78203526482261, 33.14694835752057],
          [-86.7860454984611, 33.13859751774804],
          [-86.78673620526779, 33.138838036230396],
          [-86.78255333547662, 33.1471287631669],
          [-86.78203526482261, 33.14694835752057]
        ]
      ]
     },
     "properties": {
        "value": 18.0,
        "colorIndex": 6
     }
   }

I see documentation for doing this at : https://gdal.org/drivers/vector/mvt.html. Here is the command I have used :

ogr2ogr -f MVT -dsco FORMAT=DIRECTORY -dsco MAXZOOM=10 target.mbtiles test.geojson

I am not sure what -dsco means and I can't find it in the documentation, but the script runs with it.

The output generates tiles in directories like I would expect it to look, but it will never load on a map.

Is there something wrong with the way I am converting the data? Something I am missing? I am not an expert here so I think that might be the problem. I am sure it's not a projection issue as the default settings are in EPSG:3857. Any help is greatly appreciated.


texts: How to use gdalwarp to reproject a raster without changing band labels?

texts:
Procedure
I used QGIS to reproject a raster data. QGIS does this via gdalwarp like so:
gdalwarp -overwrite -t_srs EPSG:4326 -of GTiff /path/to/input.tif /path/to/output.tif

Comparison of input and output files
input.tif

output.tif

Is it possible to use gdalwarp to reproject a raster without changing band labels? If so, how?


texts: Polygons will not buffer in ArcMap and instead give WARNING 000117: Warning empty output generated?

texts:
I have been trying to buffer 15 polygons (-15m) in ArcMap.  I did a batch buffer and 13 polygons worked fine, the other 2 have not buffered, an error occurs:
WARNING 000117: Warning empty output generated.
Executed (Buffer) successfully.
End Time: Fri Jul 06 18:58:32 2012 (Elapsed Time: 1.00 seconds)

All the input polygons are similar with full attribute tables so the output should not be empty for these 2!


texts: Looking for a template project for OSM data

texts:
I have a PostGIS 2 database filled with OSM data for my country via Osmosis. Is there perhaps a template/project (!styling) for QGIS i could download, where i would simply change the provider of data to my own database (different database name and credentials) so i wouldnt have to deal with the tedious styling.

Thanks


texts:
I am trying to import processing toolbox in my standalone PyQGIS application. But I am in trouble. I follow this instruction:
Importing QGIS Processing in stand-alone python script? and this Using QGIS3 Processing algorithms from standalone PyQGIS scripts (outside of GUI)

sys.path.append('/home/' + getuser() + '/.qgis3/python/plugins')
from processing.core.Processing import Processing
Processing.initialize()

But PyCharm doesn't recognize processing still. :/ I suppose it is because of the path, but I don't know where to point anymore.
I am using QGIS 3.0 and Linux Ubuntu.

UPDATE: I tried now this path: /home/usrer/.qgis3/python/plugins
After that, I got the option to install package, which I run. But after running it, I got the error:

Collecting processing
  Using cached https://files.pythonhosted.org/packages/3b/2d/a6f17cc99d9c45c33eb3eccd6999505d9197b31f0845a845919032262a01/processing-0.52.zip
Complete output from command python setup.py egg_info:
Traceback (most recent call last):
  File "<string>", line 1, in <module>
  File "/tmp/pycharm-packaging401/processing/setup.py", line 12
    raise ValueError, 'Versions of Python before 2.4 are not supported'
                    ^
SyntaxError: invalid syntax

----------------------------------------

Command "python setup.py egg_info" failed with error code 1 in /tmp/pycharm-packaging401/processing/

And installation failed.
And when I try it to install from terminal with:
    sudo pip3 install processing
I get the same error.


texts: Difference between arcpy.analysis.Buffer and arcpy.Buffer_analysis?

texts:
I'm running Arc 10.2 for Desktop with Python 2.7 installed at C:\Python27\ArcGIS10.2 (I also have Py Ver 3 as well).

Can someone enlighten me in a few words on difference between, say, arcpy.analysis.Buffer and arcpy.Buffer_analysis and where the applicable scripts reside and how they work.

I've found the analysis.py script at Program Files\ArcGIS\Desktop10.2\arcpy\arcpy\Geoprocessing which contains the Buffer class within the script and this would run. I assume this essentially sits behind the ArcToolbox\Analysis\Proximity\Buffer tool.

However the Buffer_analysis also comes up as an option when I type but cannot find where this sits. I know sometimes different versions of the programs have changes.


texts: ArcGIS Pro: How can I add layer effects to all Values with 'Unique Value' Symbology

texts:
I have a polygon feature class I want symbolized by a unique value attribute, in ArcGIS Pro.
Here is my example of this basic operation:

When I select one of these values, go to the 'More'-> 'format symbol(s)', I have the option to add a 'Donut effect' to the solid fill of my layer as shown here:

As you can see it edits the symbol for the feature I had selected.
I want to do this to ALL of the unique values of this feature class. But when I go to the 'More' -> 'Symbols' -> 'Format all symbols' for example, the option to add a layer effect is not there. Same for when I select more than one value and a time and do 'More'-> 'format symbol(s), the option to add a layer effect is not there.

Is it possible that the only solution to this is to edit each value's symbol one by one? I have hundreds in this feature class I need to edit.


texts: grey lines around map cells in pdf export from qgis 2.18.5

texts:
I am working on a map from raster files in the map composer of QGIS 2.18.5 from a 5 min resolution raster - export resolution set at 600dpi:

.

I try to save the map as pdf, but for some reason, a gray line appears around -most- cells:

.

Any ideas? I have created a new project with just one raster but I get the same problem again.


texts: Request values within an area from geoserver

texts:
I host some raster data on geoserver and want to access all pixel values within a given area (e.g. a geojson polygon).
I can access single pixel values via WMS GetFeatureInfo but this approach is limited to a single point, not areas.

WCS seems to work for area requests but the available operations don't serve the purpose of retrieving values.

Is there a way to run a WMS/WFS/WCS like request for raster data using an area as request?

Are cql filters a way to go?


texts: How to keep the image of a collection that are integrally contained in a geometry boundary in Google Earth Engine?

texts:
I have an image collection and a geometry ROI

var ROI = ee.Geometry.Rectangle(-78.20688664502586, 34.67650023609108, -77.65207707471336, 35.216816430066366);

var l5 = ee.ImageCollection("LANDSAT/LT05/C01/T1_TOA")
      .filter(ee.Filter.calendarRange(1984,1985,'year'))
      .filter(ee.Filter.calendarRange(1,12,'month'))
      .filterBounds(ROI)
      .select('B1')

When I print my image collection I get 190 images. However, some of the images cover only 10% of the ROI area. I would like to keep only the images that are integral contained in ROI. I have no idea how to do it.


texts: Setting proj parameters for Second World War military grid systems

texts:
I'm trying to georeference some WWII maps available at the University of Texas repository (e.g., example map).

The problem is that these maps refer to an unknown proj system (at least for me) and I can not figure out how to set up a custom reference system.

These maps refer to a "Mediterranean Grid" and "Caucasus Grid" as follow:

Is anybody aware of how to manage these kinds of grid systems and how to create a custom projection (and a graticule) on top of them?

In one of the adjoining map sheets, I found the following information:

It introduces new clues and allows the benefiting of this discussion: Setting proj parameters of old map.

EDIT April 5th:

On the basis of the different answers (thanks to all), I was able to generate the following custom projection:

+proj=lcc +ellps=clrk80 +lat_0=39.50 +lon_0=45 +lat_1=39.50 +x_0=2155500 +y_0=675000 +units=m +no_defs

Then I georeferenced this sheet. I tried either with Helmert or Polynomial or Thin Plate Spline transformation algorithm (with just few differences). Yet, the map appears inaccurately georeferenced. See the following example:

curl -L "http://legacy.lib.utexas.edu/maps/ams/turkey/txu-pclmaps-oclc-10127923-palu-e14.jpg" > ./output.jpg
# add control points
gdal_translate -of GTiff -gcp 751.581 2077.58 1.72e+06 600000 -gcp 6302.39 1037.18 1.81e+06 600000 -gcp 6510.44 5478.53 1.8e+06 530000 -gcp 1461.03 5785.78 1.72e+06 540000 ./output.jpg ./output_gcp.tif
# apply the proj parameters
gdalwarp -r near -tps -co COMPRESS=PACKBITS  -t_srs "+proj=lcc +ellps=clrk80 +lat_0=39.50 +lon_0=45 +lat_1=39.50 +x_0=2155500 +y_0=675000 +units=m +no_defs" ./output_gcp.tif ./output_mm_warped.tif

I'm afraid that, as suggested by @GabrielDeLuca, the missing of datum transformation is a big issue. But I have no idea how to collect this missing piece of the puzzle.

In addition, whereas other sheets of the same series report also the scale factor (see for instance the Levant Grid parameters below), to my knowledge, there are no maps of Caucasus Grid that offer this information.


texts:
I am using ArcPy to clip from a huge raster (about 100.000km²) a smaller part (about 20.000m²) based on the extent of a input polygon.
Running the analysis takes forever and although I let it run for hours, the clipping operation never seems to come to an end.
What settings can I use to improve the performance?
My code is as following. Note that I set the extent (extracted by the input polygon) and additionally the input polygon as in_template_dataset, to make sure that only this extent is kept.
I save the output in the memory workspace as a TIFF file, because it is only an intermediate result. Is the TIFF format appropriate? I need to keep the cell size, exact raster values from the original raster.
extent = arcpy.Describe(inputPoly)
extent_in = str(extent.XMin) + " " + str(extent.YMin) + " " + str(extent.XMax) + " " + str(extent.YMax)
output = r"memory\clip.tif"
arcpy.management.Clip(elevModel, extent_in, output, inputPoly, "", "NONE", "NO_MAINTAIN_EXTENT")


texts: Computing Shortest Path Distance for multiple points in QGIS

texts:
I have been looking for ways on how to compute the shortest distance between one point to multiple points using a specific road network.

I have not been successful in achieving the solutions in this site.

Here's my case:

I have 3 shapefiles:
1. Road Network (lines)
2. Polygon Centroids
3. Source Point (Starting Point)

How do i easily compute the shortest distance between the starting point to the polygon centroids over the road network?


texts: Moving file geodatabases in python without ArcPy

texts:
Please find my script below currently all i'm wanting to do is iterate through a folder and any subfolders extracting anything that is a .gdb directory + its contents to a new folder.

import shutil, errno, re, os

results = []
input_path = r"D:\Profiles\cdixon\Documents\BaseFolder"
output_path = r"D:\Profiles\cdixon\Documents\Extracted_GDBS"
for path, subdirs, files in os.walk(input_path):
    for dir_name in subdirs:
    if os.path.join(path, dir_name).endswith(".gdb") or os.path.join(path, dir_name).endswith(".GDB"):
                    print "Moving " + os.path.join(path, dir_name) + " To " + output_path
                    shutil.copytree(os.path.join(path, dir_name),    output_path)

However when ran I get the error, despite the fact that the Extracted_GDBS folder is empty

  WindowsError: [Error 183] Cannot create a file when that file already exists

If anyone has any ideas how I can fix this that would be great. I've seen some test scripts on here that kind of do what I'm looking for but given I'm not opening up any of the files or doing manipulation, aside from a copy/move, I think this should be ok without using any ESRI Arc libaries.


texts: Making population estimates over large areas for emergency?

texts:
I'm a meteorologist in Florida trying to encourage people to evacuate barrier island areas of south Florida earlier rather than later because of the very real potential of a landfall from a strong hurricane Matthew in under 48 hours, while there are currently no evacuation orders in place in the three big counties of the Miami metropolitan area, and remembering the stories of massive gridlock during Andrew and Floyd evacuations.
I looked for any population estimates of those areas or tools to quickly make one, but found nothing.  I see that Miami Beach itself is 88,000 and looked at a few other cities, but most are split between the mainland and islands.
I have done minor GIS before, including some geocoding with QGIS [though unfortunately don't have the databases installed now or I'd try slogging something together :-(] and have various other limited interaction with GIS overlays, such as utilizing them in GrADS or working with OSM slicers.  But suffice to say my abilities with core GIS programs isn't the best.

I believe the National Weather Service itself must use some similar tool to my needs, as with each warning they can indicate the number of people affected.  But I have no idea of websites or tools available to the public that could produce such products quickly.  And it isn't the first time I've had interest of a quick rough estimate for a certain area.

A county-by-county or city by city estimate for the Florida Atlantic coastline could be very useful information to the public tonight.  I'm certain some have that information, such as emergency managers.  But I am convinced the public is the one needing it in bulk tonight.

I'm not sure it's a perfect SE question, especially on GIS, but I can think of no better place where experts may be able to point me to a site or produce a count in moments when time may be of the essence.


texts: Mapping Disease Incidence in the Northeast in R

texts:
I am trying to map disease incidence by county in states in northeastern states in the US. I have the shapefiles for the entire US but I am having a difficult time deleting the states that I do not need. Here is what I have so far:

soco <- readShapePoly("tl_2014_us_county.shp")
summary(soco)
soco.df <- data.frame(soco)

# Change row name from default to FIPS. I am trying to delete states based on FIPS.
row.names(soco)<-as.character(soco$GEOID)
row.names(soco)

# NAME has a unique 2 number state identifier.
# I am trying to delete the states I do not need based on this number.
ne<-soco[soco$NAME-(31,26 etc), ]


texts: How to keep QGIS plugins always on top of the main window?

texts:
When writing a simple QGIS 3 plugin via Plugin Builder 3, is there a simple way to tell QGIS not to hide the dialog box when the mouse clicks outside it?

I have tried a few options without luck,
How to show QGIS plugin dialog always on top?, https://stackoverflow.com/questions/1925015/pyqt-always-on-top

Or is there an option within QGIS to command plugin displays?

Editing for including my code:

def __init__(self, iface):

    #always on top
    self.dlg = UpdateAttributeDialog()
    self.dlg.setWindowFlags(Qt.WindowStaysOnTopHint)

    """Constructor.

    :param iface: An interface instance that will be passed to this class
        which provides the hook by which you can manipulate the QGIS
        application at run time.
    :type iface: QgsInterface
    """
    # Save reference to the QGIS interface
    self.iface = iface

    # initialize plugin directory
    self.plugin_dir = os.path.dirname(__file__)
    # initialize locale
    locale = QSettings().value('locale/userLocale')[0:2]
    locale_path = os.path.join(
        self.plugin_dir,
        'i18n',
        'UpdateAttribute_{}.qm'.format(locale))

    if os.path.exists(locale_path):
        self.translator = QTranslator()
        self.translator.load(locale_path)

        if qVersion() > '4.3.3':
            QCoreApplication.installTranslator(self.translator)

    # Declare instance attributes
    self.actions = []
    self.menu = self.tr(u'&Update Attribute')

    # Check if plugin was started the first time in current QGIS session
    # Must be set in initGui() to survive plugin reloads
    self.first_start = None

UPDATE

The two suggested lines just needed to be moved to the run(self) section. Replace 'UpdateAttributeDialog()' by the name of your dialog.

 def run(self):
    """Run method that performs all the real work"""

    # Create the dialog with elements (after translation) and keep reference
    # Only create GUI ONCE in callback, so that it will only load when the plugin is started
    #if self.first_start == True:
    #    self.first_start = False
    self.dlg = UpdateAttributeDialog()
    self.dlg.setWindowFlags(Qt.WindowStaysOnTopHint)


texts: Leaflet-ESRI code example that uses Find Task?

texts:
There's an example in the Leaflet-ESRI documentation here:
https://esri.github.io/esri-leaflet/examples/finding-features.html
But the code doesn't seem to actually work, neither on the example page nor when I cut and paste it into jsbin.

I'm trying to create a search/find box for a user to provide text and then find features that include that text in their attributes. I'd love to see examples of this Find task in working action.

EDITED TO ADD: When I load the above page as http:// instead of https://, that does allow the features to be added to the map. However, the Find box still appears not to function properly. I'd like to see an example of the Find tool working properly.


texts: Raster Subtraction with QgsRasterCalculator in PyQGIS

texts:
I'm trying to do a simple subtraction of 2 raster layers using QgsRasterCalculator in PyQGIS. I've searched around and thought I nailed it but the code below didn't yield anything, but there wasn't even a crash.

What could have gone wrong?

Each raster only contains one band.

target_ras_input = 'D:/Users/abc/Desktop/target.tif'
target_ras = QgsRasterLayer(target_ras_input,"target_raster")
QgsMapLayerRegistry.instance().addMapLayer(target_ras)

initial_ras_input = 'D:/Users/abc/Desktop/initial.tif'
initial_ras=QgsRasterLayer(initial_bathy_input,"initial_raster")
QgsMapLayerRegistry.instance().addMapLayer(initial_ras)

target=QgsRasterCalculatorEntry()
target.raster = target_ras
target.bandNumber = 1
target.ref = 'target_raster@1'

initial=QgsRasterCalculatorEntry()
initial.raster = initial_ras
initial.bandNumber = 1
initial.ref = 'initial_raster@1'

entries = [ target , initial ]

final_output = 'D:/Users/abc/Desktop/final_raster.tif'
calc=QgsRasterCalculator ( 'target.ref - initial.ref', final_output , 'GTiff', target_ras.extent(), target_ras.width(), target_ras.height(), entries )
calc.processCalculation()

final_raster=QgsRasterLayer(final_output, "please_work")
QgsMapLayerRegistry.instance().addMapLayer(final_raster)


texts:
I am using WPF to display a map, which is created as MPK by ArcMap.

I am trying to change the DBF table and then update the map, but after I remove and then add new layers, the map won't change, the old layer is still there.

What should I do to refresh the map?

I think it may be because it is using cache or something.


texts:
I'm trying to mosaic the True Marble world imagery 250M geotiffs using gdal_merge.py though one of them doesn't have RGB bands like the rest and gdal_merge fails (same number of bands is listed as a requirement in it's documentation). This is the error I get:

python.exe : ERROR 5: TrueMarble.250m.21600x21600.B4.tif: GDALDataset::GetRasterBand(2) - Illegal band #

This is the tile causing the issue: TrueMarble.250m. 21600x21600.B4

Should I use a different mosaic tool or should I try to coerce the tile into the format gdal_merge requires somehow?


texts: ArcGIS JavaScript Sketch defaultCreateOptions

texts:
I want to change the defaultCreateOptions mode to "click" only. I want to disable adding of point when pointer is dragged. I am trying the following code but it doesn't work.

var sketch = new Sketch({
 view: view,
 layer: graphicsLayer,
 defaultCreateOptions: {
  mode: "click"
 }
});

Adding defaultCreateOptions has no effect. It is still using its default options.


texts: Find intersections between spatial lines and points using R

texts:
This question assumes that I'm using R with the sp and rgdal packages installed. Now that I have a set of spatial segments of class SpatialLinesDataFrame (called 'segmentData') and a set of spatial points of class SpatialPointsDataFrame (called 'poi'). Suppose that the segments and the points intersect each other as illustrated in the example figure below where line segments are in red and POIs are in blue. I'm trying to query using the function 'over' of the sp package in order to find out which lines cross over a given point. I'm having the following 2 problems:

(1) Whenever I try to query using over(poi, segmentData) or over(segmentData, poi), the returned result is always a list of NA's. But apparently, they do intersect each other, as illustrated in the example. How do I get a valid returned list?

(2) The 'over' method only returns one line id that intersects a given POI (if any), but apparently given one POI there could be several segments that cross over it (e.g., POI located at intersection of 2 road segments). How do I query so that for a given POI, I wish to find all the segments that cross over it?


texts:
I have issues trying to write the query for the virtual layer. I am trying to get two layers 'trials' and 'interactions' to merge with the common fields of Longitude and latitude. However I am struggling to get it right without a syntax error occurring. Does anyone know how to correctly write it? The below is what I put

SELECT DISTINCT a.*, d."id"
FROM "Trail" AS a
JOIN "Interactions" AS d ON a."Latitude","Longitude" = d."Latitude","Longitude"


texts: Assigning values to raster obtained form scanned map using QGIS

texts:
I scanned a temperature map for a specific region and I would like that the raster pixels had the original values from the paper map so I can do further analysis (for example, difference of temperature between time periods). I have tried to find an answer in the Google and Qgis documentation, but did not have luck.

This is the map:


texts: Problem Importing CAD features using arcpy.MakeFeatureLayer_management?

texts:
When importing multiple CAD drawings using "arcpy.MakeFeatureLayer_management" and "arcpy.Merge_management" I encountered a problem.  Drawings that have line weights and elevation values assigned to some of the polyline features would not import.  When other polylines that had values of 0 for the line weights and elevations in the same drawing imported just fine.

Why would the "arcpy.MakeFeatureLayer_management" care about line weights or elevation Values?  Has anyone encountered this or a similar problem and found a work around?


texts: How do you manage personal GPS data? (what software solutions)

texts:
I'm interested in creating a personal library of locations which are important to me, sort of like a personal geocache. Maybe spots I've found that I like to swim, or places with good views of the sunset, or places where something significant to me happened. Random things, which I want to be able to keep map reference to in the long term.

I can't figure out the best software approach for this, or understand what's possible and why. Here's some of the requirements:

A phone app (I'm on iPhone) which lets me easily place a marker/record a coordinate when I am somewhere or hear about somewhere. Being able to add it to a specific list/tag it and add descriptions is necessary. If there were some way to link photos that would be cool.

I currently use Maps.me, so maybe that?.

A way to move the GPS coordinates recorded on my phone to my computer easily.
They should be stored in a way which lets me use version control like git to back up and watch the database grow.
A way to view and improving the database on the desktop. Preferably overlaid with a useful map like OpenStreetMap, or better yet a topographical map or even ex OS Maps in the UK.
An easy pathway to return the improved database to a phone-based mapping software, updating the previous versions of that database I had there. Extra points if it looks good and is on a mapping app I might want to use anyway, like Maps.me or Google Maps.
Robust in the very long term, so not a proprietary software data type. Open sources is best, and critical data being saved in a text-like file.

It seems like KML files being imported to and from Maps.me might be a basic solution, but I haven't figured out how to edit/understand these files on a computer, and then put them back on the phone easily, especially if there are many separate lists. Is there an elegant solution, maybe that geocachers like to use? Or should I be going about this another way


texts: How to make PROJ4 on windows including NZ datum shift grids

texts:
Can anyone point me to an example of how to 'make' proj4 including the New Zealand datum shift grids?

I understand I need to include the shift grids in the build, just not sure how to make it on windows.


texts: How to calculate a new raster based on nodata values from original raster, in QGIS?

texts:
I'm using the raster calculator from core in QGIS and I need to generate a new tif raster file (mask) from the nodata values of another tif raster file. I tried isnull(), but no luck. Any one can help?
Thanks


texts:
I have 2 GEOJSON files - 1st contains map of indian states ('India-simple.GEOJSON'), and the second is cartogram based on Gross State Domestic Product of that state ('India-cartogram.GEOJSON'). I tried to transition the normal map to cartogram based on http://metrocosm.com/how-to-make-cartograms-with-animation/ blog by Max Galka -

<!DOCTYPE html>
<meta charset="utf-8">
<head>
<script type="text/javascript">
//<![CDATA[
try{if (!window.CloudFlare) {var CloudFlare=[{verbose:0,p:0,byc:0,owlid:"cf",bag2:1,mirage2:0,oracle:0,paths:{cloudflare:"/cdn-cgi/nexp/dok3v=1613a3a185/"},atok:"a098c04a5388e57588bbd1e6d0499d0f",petok:"225f354af101947ea9756023771c4624bffcb637-1468222727-1800",zone:"metrocosm.com",rocket:"a",apps:{},sha2test:0}];document.write('<script type="text/javascript" src="//ajax.cloudflare.com/cdn-cgi/nexp/dok3v=e982913d31/cloudflare.min.js"><'+'\/script>');}}catch(e){};
//]]>
</script>
<style></style>
</head>
<body>
<script data-rocketsrc="http://metrocosm.com/js/d3.v3.min.js" type="text/rocketscript"></script>
<script data-rocketsrc="http://metrocosm.com/js/d3.geo.projection.min.js" type="text/rocketscript"></script>
<script type="text/rocketscript">

var svg = d3.select("body").append("svg")    //create the SVG panel
      .attr("width", 300)
      .attr("height", 200);

d3.json("http://metrocosm.com/india-simple.geojson", function(usa) {
d3.json("http://metrocosm.com/india-cartogram.geojson", function(usa2) {

var projection = d3.geo.albersUsa()
    .scale(300);

var path = d3.geo.path().projection(projection);

var map = svg.selectAll("path")
   .data(usa.features)
   .enter()
   .append("path")
   .attr("d", path)
   .style("stroke", "rgb(50,50,50)")
   .style("stroke-width", 0.6)
   .style("fill", "rgb(200,200,200)")
   .attr("transform", "translate(-320,-140)");

setInterval(function(){
    map
        .data(usa2.features)
        .transition()
        .duration(3000)
        .attr("d", path);

        setTimeout(function(){
            map
                .data(usa.features)
                .transition()
                .duration(3000)
                .attr("d", path);

        }, 3000)
}, 6000);

});
});

</script>
</body>

Unfortunately, the above html script doesnt do anything, and there are no errors either. How can i achieve the required transition?


texts: ArcMap 10.3 fails to export a clipped image from a mosaic dataset

texts:
ArcMap 10.3 fails to export a clipped image form a mosaic dataset.

I have the following mosaic dataset, I have used Image analysis to clip apart of it based on a polygon, but the issue here is that when I try to export it every time my workstation freezes and stops working

The area of the polygon is about 996,648 Km2 and my workstation is very powerful, I don’t know if it’s a bug or what in this tool.

Now the ArcMAP starts to count but the machine freeze at 1% or more but never reaches 100% , what might be the issue here?


texts: Setting raster values to NoData when raster does not have valid statistics?

texts:
I have a floating point raster with some values of 99999 that I want to set to NoData in Arc 10.2.2. I've tried using the SetNull dialogue box, the SetNull arcpy function, and the Raster Calculator dialogue box as described in Reclassify raster values of 9999 to NoData, leaving all other values alone, in ArcGIS Desktop?, but all of these fail with same error message:

  ERROR 010423: myraster.FLOAT.TIFF does not have valid statistics as required by the operation.

I've tried using the Calculate Statistics tool as part of the Raster Properties toolbox, as suggested by the Help when I click on the error message, but this doesn't seem to change anything.

Any suggestions on how to make it work or alternatives for fixing the error?


texts: Looping over rasters in list with QGIS script?

texts:
Using QGIS 2.14, I want to loop over rasters in a list with a QGIS script.

##Input_raster1=raster
##Input_raster=raster
##resul_rasters= output raster
list=[]
list.append(Input_raster1)
list.append(Input_raster2)
for file in list:
  processing.runalg('saga:rastercalculator', file,[],'ifelse(gt(a,-1),1,0)',False,7,result_raster)

This script just return one result


texts:
Im trying to iterate over features in a vector layer following Iterating over Vector Layer:

from qgis.core import *
import qgis.utils
from qgis.core import QgsVectorLayer

shapefile = '/media/bera/Lagring/GIS/data/ok_riks_Sweref_99_TM_shape/oversikt/riks/ak_riks.shp'

layer = QgsVectorLayer(shapefile, 'borders', 'ogr')
it = layer.getFeatures()

for feature in it:
    geom = feature.geometry()
    print geom.type()
    print 'test'

print 'stop'

It will not go into the for block so all that is printed out is stop

The objects are created:

layer
<qgis._core.QgsVectorLayer object at 0x7fa08d2bcb98>
it
<qgis._core.QgsFeatureIterator object at 0x7fa08d32e218>

What am i doing wrong?


texts: z-Coordinate and interpolation plugin in QGis

texts:
I've a csv file that includes weather stations locations (x,y,z) and precipitation values. I loaded the file by using add delimited text and then converted it into a point shapefile WS84. By using the interpolation plugin I managed to create the different rainfall zones using pseudocolor and later I extracted the contour. However, as I have z values for elevation I'd like to take this into consideration. The interpolation plugin includes the option to select z-Coordinate to load elevation points but I cannot select it. When I open the plugin without any layer loaded the check box works fine, but when the layer is present is not possible to select the check box. The values of elevation column (z) are in mm. Am I  missing anything ?

Qgis 1.8.0, interpolation plugin loaded from the repository, running Windows Vista.

Addition: I believe the answer is below. z data must be within the shapefile not as a feature. Thanks Bernd V!

"Also note that you can use the z-coordinate for the interpolation. This is the elevation of your data that is stored within the geometry itself and not in an attribute column. This is only possible if your GIS data is in 3 dimensions." in gistutor.com/quantum-gis/20-intermediate-quantum-gis-tutorials/… . So it seems you need to generate "real" 3D-data first, but this is a blank spot for me, sorry. – Bernd V


texts: How to fill the gaps and how to export the image

texts:
I would like to work with Landsat 7 images, but most of them have gaps. I would like to fill the gap., Secondly when I import the image in ArcMap, it is shown in black color.
Where am I getting wrong? I am sharing the code link here
https://code.earthengine.google.com/?scriptPath=users%2Fomerreshi90%2FAir-Quality-Parameters%3AKadalundi-2005
Also my code is:
var cloudMaskL457 = function(image) {
  var qa = image.select('pixel_qa');
  // If the cloud bit (5) is set and the cloud confidence (7) is high
  // or the cloud shadow bit is set (3), then it's a bad pixel.
  var cloud = qa.bitwiseAnd(1 << 5)
                  .and(qa.bitwiseAnd(1 << 7))
                  .or(qa.bitwiseAnd(1 << 3));
  // Remove edge pixels that don't occur in all bands
  var mask2 = image. Mask().reduce(ee.Reducer.min());
  return image.updateMask(cloud.not()).update Mask(mask2);
};

var dataset = ee.ImageCollection('LANDSAT/LE07/C01/T1_SR')
                  .filterDate('2005-03-01', '2005-03-31')
                  .filterBounds(Aoi)
                  .map(cloudMaskL457);

var visParams = {
  bands: ['B3', 'B2', 'B1'],
  min: 0,
  max: 3000,
  gamma: 1.4,
};
Map.setCenter(-83.0539, 42.3159, 10);
//Map.addLayer(dataset, visParams);
var test=dataset.median()
Map.addLayer(test, visParams);

//.......exort ......

// Create a 3-band, 8-bit, color-IR composite to export.
var visualization =test.visualize({
  bands: ['B4', 'B3', 'B2'],
  max: 0.4
});

// Create a task that you can launch from the Tasks tab.
Export.image.toDrive({
  image: visualization,
  description: 'Greenest_pixel_composite',
  scale: 30
});


texts: ArcGIS Composite Geocoder Search by Address or Streets Intersection

texts:
I'm struggling with the configuration of a composite geocoder. I have a table of address points consisting of

House number
Street name

and a table of street segments consisting of

Street name
Left cross street name
Right cross street name

Is it possible to create a composite locator that uses the address point as its primary geocoder, then falls back to the street name geocoder to locate the intersection of two streets?


texts: How to enable tracing in QGIS 3 beyond applications for adding polygon purposes?

texts:
I have used the enable tracing feature in QGIS 3 to create new polygons (Add polygon feature), and it works very well. However, i have not been able to use it for other purposes beyond that, for example, to create a ring inside a polygon and similar tasks.

Is it possible to enable tracing for other purposes beyond adding a polygon?

In this moment my only example is creating a ring inside a polygon. I had to create this point by point, and due to the large number of vertex it would be very handy to use the tracing.


texts: Using ArcGIS ModelBuilder to perform In-line variable substitution for input data path?

texts:
Does anyone know how/if you can use in-line variable substitution for an input path name in ModelBuilder for ArcGIS 10?

I used "get field value" to retrieve a code from one feature class, and I want to substitute that text string in two places in another input variable pathname.

For instance, if I used "get field value" to retrieve the code "ABC123", I would use that to call up a raster (make raster layer) with this path name: "D:\Data\Rasters\%ABC123%\el_%ABC123%

Does anyone know if this will work?

ModelBuilder won't recognize this as a valid path name right now.


texts: Update location from text file in Google Earth Pro in real time?

texts:
I'm trying to update the ship's location from a text file in Google Earth Pro in real-time. The Lon/Lat locations of the ship are being written to a text file from the GPS.

I am new to using Google Earth Pro.


texts: Grouping Polygons based on Pipes Upstream from Main Pipeline?

texts:
I am looking to group polygons together that are upstream of a trunk sewer line (Main Pipeline). Basically, the trunk sewers all have other (regular) sewers branching off of them. All of these sewers are currently contained in the same feature class (Sanitary_Pipes). The polygons that I want to group are in their own feature class (Sanitary_Catchments) and overlap with Sanitary_Pipes. The goal is to group together the polygons that overlap with the string any given string of connected regular pipes that branch off of a trunk pipe.

In the screenshot above, the thick green lines are trunk sewers, red lines are regular sewers, the points are manholes along the sewers and the green polygons are catchment areas that are a part of Sanitary_Catchments. The highlighted area is an example what I'm trying to do; all of the polygons that overlap with the set of connected pipes that are upstream of the trunk become grouped.

Is there any sort of tool/collection of tools that can do this?


texts: Dynamic indoor mapping in Native android

texts:
I'm about to work on a project, in which we have to create an indoor map to locate articles inside our store. I have created different layers in QGIS (store, showrooms, articles) shapefiles and published them in geoserver. What I need is an opensource way to get the services from geoserver to show data on an android native application.

The problem is that the store can be changed each 15 days, so I have the DWG file if you could give me any alternative solutions so I don't have to digitize the store again.


texts: How to set identify tool to identify from visible layers as a default when starting ArcGIS

texts: How to set identify tool to identify from visible layers as default when starting ArcGIS without the need to change it from top-most layer to visible layer each time stating arcmap


texts:
We are clipping one raster to another and our teacher gave us this hint for the 2nd parameter "Remember that the Describe function returns an object. The Clip tool will need the string representation of the extent property" - I know I need a string of the extent. Everything that I have tried has given me:

  ERROR 000628: Cannot set input into parameter rectangle.

Here is what I have so far

import arcpy

arcpy.env.workspace = "C:/Student/Newark.gdb"
desc = arcpy.Describe("C:/Student/Newark.gdb/studyarea")
clip = str(desc.extent)
Coord = clip.strip(" NaN")
inRaster = "C:/Student/Newark.gdb/Newark"
outRaster = "C:/Student/Newark.gdb/NewarkClip.jpg"

arcpy.Clip_management("inRaster","Coord,"outRaster","studyarea","","NONE","")


texts:
I am trying to install the package "geographiclib-1.20.tar.gz" Referring to a previous post to install 3rd Party libraries here. I ran OSGeo4W Shell and tried to use pip. This is the error message I am getting. What am I doing wrong ?

How do I install the setuptools ?


texts: Performing simple line-of-sight analysis in QGIS?

texts:
I've just installed QGIS 2.14.3-1 on macOS in an effort to determine if two points in my neighborhood have line-of-sight. I've got a 7.5 arcsecond DEM loaded, and a point layer with my two points (no elevation ascribed yet).

All I want to do is determine if the terrain obstructs line-of-sight from one point to another. What's the easiest way to do this? I don't need a complete viewshed analysis.


texts: Selecting from selection based on attribute table in ArcGIS for Desktop?

texts:
I selected polygons which are less than 313 m2 and from that selection I need to select the ones that are completely within some other polygon.
So, both selections are on the same shapefile.

I tried selection by attribute and by location which seems to respond.


texts: Conversion of Text layers to vector 'Line' Feature Class

texts:
Is there a way to convert Text or Labels to Vectors within Geomedia or utilize the Publish to KML and convert using Google Earth?

I do not have Adobe Illustrator installed here to 'Create Outlines', and obviously not MAPublisher to import them to Illustrator.

I need to create a shapefile with text and other features displayed as lines within the same shapefile.

Note: I do not wish to export to .dgn and convert using Microstation


texts:
I make map with a overview in QGIS. I try to adjust my overview scale in relation with the first map scale. I explain : I want adjust overview scale with a factor 5 or 50, not fixed scale. When my scale is 1:1000 in first map, I want scale in overview is 1:5000, if first map 1:1200, overview 1:6000……

I look in expression but I don't find the good way.


texts: Use a second canvas for drawing an animation in Openlayers

texts:
I want to implement an animation of flowing rivers on top of an Openlayers map. My doubts are if I simply should make a new canvas and position it on top of the map or if there is a more elegant way using a OL layer?

What I try to do is similar to this example using Leaflet. This animation is done using canvas lineTo() function with an interval and works with pixel coordinates drawn on a second canvas. Leaflet offers map.createPane() function which then handles the positioning for you and allows you to link a second canvas to the map.

I have the animation working in OL (jsfiddle) drawing it directly on the map canvas. The problem now is that I have to set the canvas globalAlpha property to simulate the trail effect which makes all other layers disappear. Another possibility would be to draw a semitransparent rectangle over the whole canvas, with the same negative effect of making disappear everything else on this canvas:

context.fillStyle = 'rgba(255, 255, 255, 0.1)';
context.fillRect(0, 0, canvas.width, canvas.height);

I could get it working positioning a new canvas on top of the map, but then I would have to handle all the interactions by myself.

So does there exist something similar to Leaflets map.createPane() function or any other possibility to use a second canvas linked to the same map?

PS: just to make it clear, I want something like that.

PPS: There is a new OL plugin based on windjs which looks really promising.


texts: Preparing dataset to perform co-kriging in R gstat?

texts:
This is the first time I'm using co-kriging in gstat. My problem is that I'm not sure how to prepare the data frame to supply to co-kriging when the variable of interest and auxiliary variables are not measured at the same locations.

Data frame with the variable of interest (z_ar) is voi_gs

head(voi_gs)
             x        y     z_ar
8974  312216.6 530439.8 49.03470
8283  312084.6 530559.8 57.15355
5057  311772.6 530883.8 49.13453
1551  311976.6 531207.8 44.73679
10116 312168.6 530163.8 50.45549
6155  312528.6 530787.8 62.70750

Data frame with the auxiliary variable (z_pts) is aux_gs

head(aux_gs)
            x        y    z_pts
9564 312198.6 530313.8 75.91368
6584 311766.6 530745.8 38.87462
2138 311562.6 531153.8 58.62555
9110 312534.6 530421.8 68.35654
9525 312366.6 530325.8 54.26653
7497 311442.6 530649.8 38.95024

I combined them into one dataframe to supply to variogram() and krige() functions. Since none of the locations are same between voi_gs and aux_gs, I introduced NA values and combined them in the following way:

aux_gs$z_ar=NA
voi_gs$z_pts=NA
comb_gs=rbind(aux_gs,voi_gs)
head(comb_gs)

            x        y    z_pts z_ar
9564 312198.6 530313.8 75.91368   NA
6584 311766.6 530745.8 38.87462   NA
2138 311562.6 531153.8 58.62555   NA
9110 312534.6 530421.8 68.35654   NA
9525 312366.6 530325.8 54.26653   NA
7497 311442.6 530649.8 38.95024   NA

tail(comb_gs)
            x        y z_pts     z_ar
180  312468.6 531363.8    NA 70.54528
8633 312264.6 530511.8    NA 44.34631
7694 312492.6 530631.8    NA 57.30173
1079 312108.6 531255.8    NA 46.96482
2230 311124.6 531135.8    NA 40.36449
2201 312312.6 531147.8    NA 44.85896

and then, constructed cross variogram:

coordinates(comb_gs) = ~x+y
g = gstat(formula=z_pts~z_ar, data=comb_gs)
vg = variogram(g)

but variogram function does not accept NA values and gives me an error. I know I can't have NA values in the data, but I don't know how else I can create the data frame to construct the cross-variogram.


texts: How to Transform Many Projections Into One With SQL Server

texts:
I know very little SQL Server, but am quite familiar with ArcGIS and Python.  The data I want to transform is point data which is in a single table [ID], [X] [Y], [COORDSYS].  There are 5 possible coordinate systems.  NAD83 UTM Zone 13 and 14, NAD27 UTM Zone 13 and 14, WGS84.

Right now I have a Python tool that accesses the table, transforms the points based on the coordinate system and outputs a SHP file that is in WGS84.  It takes about 5 minutes when the user starts it which is less than ideal.

It seems it is possible to do this with SQL Server spatial functions so that the same thing can be done in the background and be ready for the user to just display in ArcGIS.  But I'm not sure how to best change which transformation to use depending on a value in the table using SQL.  I need to know what it's called that I'm trying to do so I can look it up.

Following is the code I'm using to implement this right now.  I suppose that the part which might be taking the longest is creating Polygons out of each distinct set of points, rather than the transformation of each point into WGS84

def addCorner(self, x, y, z):
        pointgeo = arcpy.PointGeometry(arcpy.Point(x,y,z), self.origSR, True)

        try:
            # transforms the points to WGS84 as they are added to the Loop object
            projectedPointGeo = pointgeo.projectAs(self.wgs84SR,self.transformation)
        except Exception as e:
            print("Please ensure that you have all the necessary transformation tables."
                           "You are likely missing the 'sk83-98.gsb' in your ArcGIS installation."
                           "Anyone in the GIS group can help.")

        # Adds each point as an arcpy.Point with Z values
        self.__corners.add(projectedPointGeo.firstPoint)

def CreateLoopsFromAcquire(loopsources, curdir, statusbar):
    if not curdir:
        if not os.path.exists(r"C:\temp"):
            os.mkdir(r"C:\temp")
        curdir = r"C:\temp"

    shapefiles = {}
    for loopsource in loopsources.values():
        print("-- From {0}".format(loopsource.tablename))
        cursor = arcpy.da.SearchCursor(loopsource.tablename,
                                        [loopsource.surveyfield, #0
                                        loopsource.corner_idfield, #1
                                        loopsource.coordtypefield, #2
                                        loopsource.xfield, #3
                                        loopsource.yfield, #4
                                        loopsource.zfield]) #5
        loops = {}

        # Gather Loops (group of corners) into Dict
        i = 0
        for row in cursor: # each row is a loop corner
            loopid = "{0}_{1}".format(row[0],row[2])
            if not loops.has_key(loopid):
                loops[loopid] = Sources.Loop(row[0], row[2])
                i+=1
            loops[loopid].addCorner(row[3],row[4],row[5])
        del cursor

        print("Found {0} loops".format(len(loops.items())))

        newfilepath = os.path.join(curdir, loopsource.outputname)
        if os.path.exists(loopsource.outputname):
            arcpy.Delete_management(newfilepath)

        newfc = arcpy.CreateFeatureclass_management(curdir, loopsource.outputname, "POLYGON", loopsource.templatepath, "SAME_AS_TEMPLATE", "SAME_AS_TEMPLATE", loopsource.templatepath)
        cursor = arcpy.da.InsertCursor(newfc,["SHAPE@","Lp_Survey"])

        for loop in loops.values():
            cursor.insertRow([loop.getPolygonOfLoop(),loop.uniquename])
            # statusbar.StepProgressBar()
        del cursor

        # return a shapefile object so that symbology can be made
        shpobjectname = os.path.splitext(loopsource.outputname)[0]
        shapefileobject = Sources.ShapeFile(
                                    newfilepath,
                                    shpobjectname,
                                    display_expression="Lp_Survey",
                                    symbology=loopsource.symbology)

        shapefiles[shpobjectname] = shapefileobject

    return shapefiles


texts: What is the easiest way to write C extensions to Python for use in QGIS?

texts: I would like to write some functions in C to call from the Python console in QGIS. What is the easiest way to do this? Does there exist anything comparable to Matlab's MEX?


texts: ArcGIS Polygon Volume - Variable height field

texts:
I'm using ArcMap 10.6

I have a TIN of a landscape which has a series of mounds on it. I have isolated polygons for each mound, which are not flat but variable in height as they go around each mound. I would like to calculate the volume and surface area for each mound, with the polygon for the mound treated as the reference plane.

I have tried the 'Polygon Volume' tool but this asks for an absolute height for the reference plane, not one which accounts for the differential height of each vertex in the polygon.


texts: Qgis 2.18 Export TIFF and PDF Fails on Windows 7 not with Ubuntu 16.04

texts:
I have two 64 bit computers.  My work computer is running Windows 7, my home computer is running Ubuntu 16.04.  Both computers have QGIS 2.18 installed on them.  Both computers have similar RAM and graphics cards. I have some QGIS projects stored on an external hard-drive that were largely developed on Ubuntu.

When exporting these projects either PDFs or TIFFs from the Windows machine the hillshade is not included in the export.  On the Ubuntu machine everything seems to export just fine.

I have read some posts about unchecked layers in the layers panel causing rendering issues.  I have ensured that there are no unchecked layers in these projects.

I have saved the Windows project as a separate QGS file and that did not make a difference on the output.  I can confirm the problem exists with two different QGS files although both projects were largely drafted on Ubuntu QGIS.

Why can I not export correctly on the Windows version of QGIS?  I am trying to incorporate QGIS into my workflow so relying on my Ubuntu box at home is not an option.


texts:
I converted an ArcMap VBA code to a Visual Studio dll.  In doing so I had to reference a new namespace "stdole" so my graphic label would show in ArcMap.  The label works on my computer, but does not work on other computers that do not have VS 2010.  I'm guessing it is because the namespace stdole.dll is not on the other computers?  Are there any alternative namespaces that I can reference that the older computers might have?  Below is the VBA pre and raw VB.NET post code conversion (<<<< this defines the code statement where the stdole namespace was needed).

VBA version

'Get the point location

Set pPoint = pFeature.Shape
pTextElement.Text = strLabel 'Add the text string to the Text Symbol

Dim myFont As IFontDisp
Set myFont = New StdFont
myFont.Name = "Courier New"
myFont.Size = 1

pTextElement.Symbol.Size = 5#   'Added Text Symbol Size
pTextElement.Symbol.Font = myFont
pElement.Geometry = pPoint
pGraphicsContainer.AddElement pTextElement, 0

VB.NET version

'Get the point location

pPoint = pFeature.Shape
pTextElement.Text = strLabel 'Add the text string to the Text Symbol

Dim VB6 As Object
VB6 = Nothing

Dim myFont As System.Drawing.Font

myFont = System.Windows.Forms.Control.DefaultFont.Clone()
myFont = VB6.FontChangeName(myFont, "Courier New")
myFont = VB6.FontChangeSize(myFont, 1)

pTextElement.Symbol.Size = 5#
pTextElement.Symbol.Font = VB6.FontChangeName(pTextElement.Symbol.Font, myFont.Name) '<<<< stdole namespace needed
pElement.Geometry = pPoint

pGraphicsContainer.AddElement(pTextElement, 0)


texts:
I am new about Arcgis javascript API. While I am doing a sample application, I confused.
I am analyzing this sample application.

This application is using the http://sampleserver1.arcgisonline.com/ArcGIS/rest/services/Louisville/LOJIC_LandRecords_Louisville/MapServer/0 as FeatureLayer.

And using a proxy

    esri.config.defaults.io.proxyUrl = "/proxy";
    esri.config.defaults.io.alwaysUseProxy = false;

I created a new visual studio project and proxy pages. And I changed FeatureLayer and Base Map.
Program is working. But I opened Mozilla Firebug and traced the traffic. My Application using multiple Http POST data, while zoom in. But Arcgis application Only 2 GET

This is my custom application screencast.
(This is screenCast of my application i.stack.imgur.com/F1C17.png)

This is Arcgis Sample Application Here (help.arcgis.com/en/webapi/javascript/arcgis/samples/widget_measurement/)

(This is Arcsgis Applcation screencast i.stack.imgur.com/muZV9.png)


texts: Calculating Distance and Time on a network without the use of ArcGIS Network Analyst

texts: I'm trying to calculate a distance and time for walking commuters (about 30,000 records) on a street network (NYC Lion).  I have a starting and ending address for each record that have been geocoded.  As far as I know I don't have an ArcGIS Network Analyst extension license at the moment (maybe there's a free trial?) so I'm curious if anyone knows how to do this without Network Analyst.  I have QGIS but I'm not sure if you can batch process records in this capacity.


texts: Use gdal installed with QGIS in Python Windows 7

texts: On my Windows 7 PC, I installed QGIS 2.18 and gdal is part of this Installation.   How do I connect Python 3.6 to gdal?


texts: Hiding previously submitted survey points on Crowdsource Reporter?

texts:
Is there a way to hide previously submitted survey points on the maps in Crowdsource reporter?

My project involves having users point locations on the map where infrastructure may impact biodiversity and have them fill in survey questions as attributes. Just to provide more context, my initial idea was to have a survey form with several sections each containing a geopoint function. But since survey123 only supports one geopoint per form, I created several surveys based on each section and incorporated them into CR to allow pointing locations for each survey section in one platform as opposed to having several separate forms to fill.

But the problem arises as we would like to treat this platform as similar to a normal survey as possible where each survey/response is a unique entry while keeping them confidential at the same time (due to demographic related questions). Hence it's best if we could just hide previously submitted survey points. I have tried seeing if crowdsource reporter is configurable but so far there are no options.

Is there a function that can achieve this? Or any ideas which ArcGIS online platform or web apps that can support this?


texts: Can a CartoDB-Torque map be used for polygon instead of a point data?

texts: Can you use Torque for country data - polygons instead of points? Make animated maps with polygon data instead of point data? I would like the color of countries to change over time as a particular variable gets smaller and larger.


texts: Windows' Offline Files for GIS projects?

texts:
I'm about to start yet another project which will consume and produce several to tens of GB of data. For efficiency, performance, and source code control* I want to do everything on the local machine. For backup, sharing with colleagues, and future proofing everything needs to be on the server (in our organization only network resources are backed up).

Windows has a built in feature for scenarios like this called Offline Files which purports to be more or less fire and forget. Has anyone used it for large data sets like we work with? Is it worth attempting/exploring? (And is there any control over where the data is synced on the local machine?)

It's not that difficult to setup up a batch file and a scheduled task to robocopy or whatever from the local machine to the server, but I'd rather not write, debug and maintain that in addition to the actual project (times N projects). If Offline Folders work better, great! If not, well... ;-)

* A portion of the project involves code and version control (e.g. git/mecurial), and it is ill-advised to run VC within remote folders.


texts: Is there a way to improve the time for opening a QGIS project?

texts:
We have a QGIS project using postgresql database "made with" osm2pgsql.
The database is about 68GB.

Our Configuration is :
AMD-FX 8 cores @ 4.7Ghz
32 Gb RAM
500 Gb SSD (read/write 600Mb/s)
Ubuntu 14.04 LTS
QGIS 2.8.2

When we open the project, it takes about 10mn.
After we can navigate, zoom in/out, show hide layers... without any problem.

My question is there a way to improve the time for opening the project ?


texts: Optimizing my script; taking too long; distance between many objects PostGIS

texts:
I have many many points and want to calculate the shortest distance from each point to the closest line of a table with many many lines. I think my underneath code is correct, I however never get an answer because pgadmin just keeps 'waiting for query execution to complete'. How can I optimize this underneath script?

SELECT
    max(br.id), max(to_timestamp(br.time)), max(at.linknummer),
    min(ST_Distance(st_transform(br.geom,3857), at.geometry))
FROM
    schema1.data1 AS br,
    schema2.data2 AS at
GROUP BY
    br.id, br.time


texts:
I know that QGIS can add values from a raster to a set of points on top of the raster through the SAGA "Add raster values to points" tool. However, I would like to generate the points themselves based on some threshold of difference in the raster values.

For example, over a body of water there would be few points generated because the change in elevation is small, but over a rugged mountainous terrain there would be more points to match the changes in elevation. Is there such a tool that QGIS has or can use through its Processing Toolbox?


texts: Any case when one must use older ArcPy cursors (not da module based)?

texts:
I was using the arcpy.da cursors in ArcGIS 10.2 to iterate for quite some time now yet didn't use "classic" cursors that became available in 10.0 a single time. Is there any situation where one would must (due to a limitation or a similar reason) use arcpy.%type%Cursor instead of arcpy.da.%type%Cursor?

Esri help doesn't seem to outline any case when one would need to use classic one:

  ArcGIS 10.1 added a new data access module (arcpy.da)  The previously
  existing cursors (that are still listed under arcpy) are still
  functional and valid; however, the new arcpy.da cursors include
  significantly faster performance.

So, I am just wondering why Esri keeps updating the older arcpy cursors in 10.2 when the Data Access module is available since 10.1 and provides "better" cursors? Is there any other reason apart from keeping 10.0 scripts that use arcpy (not da module based) cursors running?


texts: I know when you zoom into Google Earth and the icons come up for places, car washes are included in places. Is there a way to extract them and have them as a KMZ file IN A SERTERATE machine?


texts: Writing Python script tool for Data Driven Pages?

texts:
I'm very new to writing my own scripts.

What I am trying to accomplish is exporting JPEGs with a world file that follow data driven pages.  I would like to automate this process as much as possible, creating a tool for it.

The following script I was able to get to work if I am in the mxd file and create a folder previously that I add to the script (this test run I used C:\Data\test\PIC), then run it .  Here it is:

    import arcpy
mxd = arcpy.mapping.MapDocument("CURRENT")
df = arcpy.mapping.ListDataFrames(mxd)[0]
for pageNum in range(1, mxd.dataDrivenPages.pageCount + 1):
...     mxd.dataDrivenPages.currentPageID = pageNum
...     print "Exporting page {0} of {1}".format(str(mxd.dataDrivenPages.currentPageID), str(mxd.dataDrivenPages.pageCount))
...     arcpy.mapping.ExportToJPEG(mxd, r"C:\Data\test\PIC" + str(pageNum) + ".jpeg",df,df_export_width=1600,df_export_height=1200,world_file=True)
del mxd
del df

Now I'd like to take it a step further, creating a tool.  Something with parameters that allow me to pick the folder, the file name,, the width, height, etc.  So I did a test with an input for folder and file name.  It didn't work, all it did was open up the script in Notepad when it was complete.

What am I doing wrong?

Am I going about this in the right way at least?

The parameters I added were Folder for 0 and String for 1.  Here's the new script I tried when I made it into a tool:

import arcpy

mxd = arcpy.mapping.MapDocument("CURRENT")
df = arcpy.mapping.ListDataFrames(mxd)[0]

in_folder = arcpy.GetParameterAsText(0)
file_name = arcpy.GetParameterAsText(1)

for pageNum in range(1, mxd.dataDrivenPages.pageCount + 1):
...     mxd.dataDrivenPages.currentPageID = pageNum
...     print "Exporting page {0} of {1}".format(str(mxd.dataDrivenPages.currentPageID), str(mxd.dataDrivenPages.pageCount))
...     arcpy.mapping.ExportToJPEG(mxd, in_folder + file_name + str(pageNum) + ".jpeg",df,df_export_width=1600,df_export_height=1200,world_file=True)
del mxd
del df


texts:
I try to use pgr_Dijkstra through Postgis and Pgadmin to find the shortest path in a network between a lot of nodes.
So I thought I could to "automate" the process with a PLSQL query. Here is a query I simplified as example. It's a loop where I try to use pgr_Dijkstra between start nodes and end nodes having respectively _ids and _ide as id in the table network:

DO $$ DECLARE

_ids integer;
_ide integer;

BEGIN
    _ids:=161794;
    _ide:=161794;
    DROP TABLE IF EXISTS test_dijkstra;
    CREATE TABLE test_dijkstra(seq integer, node integer, edge integer, cost double precision);
    WHILE _ids <=161797 LOOP
        WHILE _ide <=161797 LOOP
            INSERT INTO test_dijkstra(seq, node, edge, cost)
            SELECT seq, id1 AS node, id2 AS edge, cost
            FROM pgr_Dijkstra('SELECT id, source, target, cost, x1, y1, x2, y2 FROM network', _ids, _ide, FALSE, FALSE);
            _ide:=_ide+3;
        END LOOP;
        _ids:=_ids+3;
    END LOOP;
END$$;

My problem is that even with this simple example, the NOTICE tells me :

NOTICE:  Deprecated function
CONTEXT:  instruction SQL « INSERT INTO test_dijkstra(seq, node, edge, cost)
            SELECT seq, id1 AS node, id2 AS edge, cost
            FROM pgr_Dijkstra('SELECT id, source, target, cost, x1, y1, x2, y2 FROM network', _ids, _ide, FALSE, FALSE) »
fonction PL/pgsql inline_code_block, ligne 13 à instruction SQL

I think the function doesn't want to have variables (_ids and _ide) as parameters, but I can't figure out why.

Is there a way to make it work? I have a lot of nodes, so I can't process pgr_Dijsktra manually for every couple!


texts: Opening PDF generated by Python AddIn using Report (*.rlf) file automatcally?

texts:
As a test I have:

Authored a test.mxd with one layer called Localities in a single data frame called Layers
Authored a test.rlf using View | Reports > Create Report of ArcGIS 10.2.2 for Desktop
Created a Python AddIn which creates the expected PDF file from the Report (*.rlf).

This is the code which all works to create the PDF report from the Localities features within the current extent:

import arcpy
import pythonaddins

class ReportForExtent(object):
    def __init__(self):
        self.enabled = True
        self.checked = False

    def onClick(self):
        mxd = arcpy.mapping.MapDocument("CURRENT")
        df = arcpy.mapping.ListDataFrames(mxd,"Layers")[0]
        lyr = arcpy.mapping.ListLayers(mxd,"Localities",df)[0]
        arcpy.mapping.ExportReport(lyr,r"C:\temp\test.rlf",r"C:\temp\test.pdf","EXTENT",extent=df.extent)
        del mxd

Is there a way to have the Python AddIn not just create the PDF file but also open it in Adobe Acrobat Reader?


texts: QGIS - Struggling with r.reclass, North must be larger than south, etc - error log in post

texts:
I'm currently in between ArcMap licenses so I'm working through QGIS and trying to reclassify a land use Raster from the USDA with 256 values. I have a text file I'm trying to use as the input:

Here is the symbology tab of my raster, let me know if you need to view different aspects of it. I should also note I can't find an attribute table for the raster, so for all I know, it has to do with that.

Here are my parameters, pretty basic stuff:

I'm trying to run this in the r.reclass tool in GRASSGIS but I'm getting some errors shown here:

From what I can tell, it accepts my rules when I enter them as text rather than upload them as a file, but I'm still getting an error when I try to run the script.

EDIT: After fixing the north larger than south, I'm still getting other errors - log  here


texts: Changing drawing order of excluded symbols in ArcGIS Pro

texts: I have a point feature class with time enabled. I have a time series where I am showing values per day of this point feature class using proportional symbols and I have excluded 0 values and given them unique symbology. When I show more than one day at a time in the time series, the symbology of the excluded values appears on top and hides the symbology of the point feature class. Is there a way to change the drawing order of excluded values?


texts: What determines Zoom Extent in MapBox Studio

texts:
When uploading data to Mapbox Studio, a zoom extent is attributed to the dataset/tileset after processing. What controls this attribute?

I've uploaded datasets that I want to use at Zoom Extent different than those that are assigned and I don't see any way to edit this.

I've seen this behaviour when adding GeoJSON or shp


texts:
I am very new to R and to raster package and I was hoping that someone could provide some guidance. I am trying to run a function that averages a series of rasters, 300 to be exact, they are in geotif format and each is of them has 44 megabytes. For this I am using a raster stack and some parallelization tool in R to make it faster. The problem I run is that when I try this on a stack that has 70 or more rasters R stops responding, I left it overnight and the still no response. Are there any limitations of the number of rasters/size I can include in a stack? The computer I am using is a Windows server with 40 cores and 384 gigabytes of ram.

Below the code I am using:

library(spatial.tools)
library(snow)
library(parallel)
library(gdalUtils)
library(raster)

# Funciones
ff <- function(x) mean(x, na.rm=T)

# ......

setwd('S:/Miguel_Fernandez/MODIS_GEOTIFF')

lista_1.tif <- list.files(pattern="_1.tif$")

# test with tiffs
beginCluster(n=40)
for (i in 1:12){
  nn <- (i*25)
  test50_stack_1 <- stack(lista_1.tif[1:nn])
  t.ini <- Sys.time()

  res <- clusterR(test50_stack_1, fun = calc, args=list(fun=ff))

  print(paste ('loop', nn, Sys.time()-t.ini))
}
endCluster()


texts: How to reproject world raster from wgs84 to World Miller Cylindrical (EPSG 54003) and keep New Zealand?

texts: I have a world raster in WGS84, downloaded from Natural Earth. When I set projection to World Miller Cylindrical (EPSG 54003) I lose parts of New Zealand as you can see in the attached pictures. Is there any way to solve this problem in QGIS (version 2.8 or 3)?


texts: How to display a PostGIS table on Google maps?

texts:
I am attempting to take a shape file and display it on a map:

I can do something like this:

SELECT ST_AsText((ST_Dump(the_geom)).geom) as the_geom FROM shapes;

But the problem is I still have to do a lot of parsing of the data to get it to display on a google map. Is there a way to do this easily via google's javascript api?


texts: How to get the georeferencing mean error?

texts: I'm georeferencing pictures and want to see the error. When using the affine transformation the mean error is displayed in meters in the bottom right corner of the georeferencing window. For the transformation I use the projective one and the error is displyed in pixels. I want to know the error in meters. Is there a way to change the unit?


texts: Getting a WMS-like bounding box from WMTS services

texts:
I'm new to GIS and confused by all the systems that are used to request data. Every provider seems to use a different system.
What I need is to request an area between coordinates (WGS-84) as a picture. This is very easy by using bbox=...... in the URL on WMS-servers but those seem to be slow and rare.

Now I wanted to use bounding boxes (or similar) with WMTS (MapBox, Google Maps etc) but as far as I can see, the service doesn't support bounding boxes anymore.

Is there a way to simply get this functionality back? I'm using it in a program using C#, so if there is a library for this it would also do the job.


texts: How can I use a raster layer over a vector layer to read in the next UTM area?

texts:
I'm still a beginner with QGIS and now have a problem that the numerous training videos don't show me and I can't find anything in the manual.
I received two layers from Swisstopo. One with raster data (DEM) and one with vector data (waters). These are now in completely different display areas (see picture).
I would now like to lay the waters over the heights. To do this, I have to move everything down (or up?) By about 300,000 units. West-East is probably fine. Is there an easy way? (I suspect that there are, but I just can't find them or don't know.

Ich bin noch Anfänger mit QGIS und habe nun ein Problem, das mir auch die zahlreichen Schulungsvideos nicht zeigen und ich finde nichts im Handbuch dazu.
Ich habe von Swisstopo zwei Layer erhalten. Einen mit Rasterdaten (DHM) und einen mit Vektordaten (Gewässer). Diese liegen nun aber in ganz unterschiedlichen Darstellungsbereichen (siehe Bild).
Ich möchte nun die Gewässer über die Höhen legen. Dazu muss ich aber alles ungefähr 300000 Einheiten nach unten  (oder nach oben?) verschieben. West-Ost ist vermutlich in Ordnung. Gibt es da eine einfach Möglichkeit? (ich vermute, dass es die gibt, ich sie aber einfach nicht finde bzw. nicht weiß. Danke schon im Voraus! :-) Gerhard


texts: Looking for webservice for distances (road) between UK towns

texts:
I am thinking about writing a program to check some mileage claims (about 45'000 of them actually) made by British MPs recently.

The data I have is quite course - the origin and destination of the trips is generally provided only at the Town level. What I am looking to do is basically look at the mileage figures they have submitted, look at the "calculated" mileage obtained from a web service, and then compare / contrast. (I'm thinking of perhaps being interested in cases where the claimed mileage exceeds 40 miles, and the claimed mileage exceeds the website figure by a factor of at least 1.4 for instance).

Code wise I'm fine, but what I don't have is a webservice I can use to get the distances. Google Directions is out (has to be used with Google Maps) - and I am looking at writing this as a console app in C#.

Can anyone suggest any suitable service?

Martin.


texts: Scale-dependent rendering in QGIS PDF output

texts:
If I export my QGIS layout to PDF with "Always export as vectors" checked the resulting PDF appears to render vector data as expected. Line work appears crisp at any zoom level.

Given that the PDF is georeferenced, and contains vector data for the map's features, is there any way to have my map's scale-dependent rendering rules apply in the PDF? I'd like to have some labels only appear as the user zooms in to the PDF. Presumably the resulting PDF would not be compatible with all PDF viewers.


texts: Open-source reader of Leica C10 terrestrial LiDAR proprietary file format

texts: I have a Leica C10 terrestrial LiDAR (laser scanner) system. The system creates and stores point cloud data in a proprietary .bin file. Currently I have to use Leica's proprietary LiDAR software, Cyclone, to convert the files generated by the system to a format (e.g. LAS, e57, or ASCII) that can be read by other open-source software for processing point cloud data. Does anybody know of an open-source software that can read and convert the C10 .bin format? There does not appear to be a public description of the file format such that I can write my own translation tool.


texts: Is it possible to edit a layer in QGIS while there is an active query?

texts: I am new to Quantum and I am digitizing polygon layers using the editing features. I would like to apply simple queries ('date > 2012/08/01'), however while the query is active I cannot edit the layer. Is this possible?


texts: Filtering out a subpart of a polygon shapefile

texts:
I have ArcGIS 10.1 for Desktop software.

I am working with a polygon shapefile (having n rows) in ArcMap. The project demands splitting this file into two polygon shapefiles having n1 and n2 rows such that n1+n2=n.

I can filter out the data using a where statement on the attribute table, but the resultant data does not have the field called "Shape", due to which I cannot plot the polygons.

I am new to ArcMap software.


texts: Unusual behavior of coordinate grid in QGis Composer

texts:
The coordinates in the grid in QGis composer display incorrectly when using the "degrees, minutes, seconds" option (or even decimal degrees). It displays West and East in the longitude and North and South in the latitude (or positive and negative) and sometimes it even mixes UTM and degrees.

It doesn't seem to have any problem with decimal format (UTM) as shown below, but no matter what I do it still has the same trouble with degrees.

I have done multiple maps in QGis 2.14 and 2.12 having the same problem and I remember that in earlier versions this did not happen. Is there any way to fix this?

June, 2018: It seems new versions of QGis don't have this problem.


texts: How to move a point on a touch screen using esri map, dojo and javascript

texts: I am new to map development. For touch screen devices, I need to add a point on an ESRI map and then move this point to other location on the map and get its coordinates.


texts:
Possible Duplicate:
  How to move SHP to gdb with ArcView license?

Let's cut it short, I have numerous mxd with layer coming from .SHP and various servers WMS that I need to move to a gdb.
What I want is to open the mxd and run a script that select all .SHP layers in the TOC create a gdb, export the .SHP, and create a new MXD with the new path pointing to gdb.

Following is the code I wrote, but I can't have it to select only the .SHP, I tried with describe dataType but they all come out as "FeatureLayer"...need help!

# Import system modules
import arcpy, datetime, os, traceback
from arcpy import env

# Load required toolboxes...
arcpy.AddToolbox("C:\\Program Files (x86)\\ArcGIS\\Desktop10.0\\ArcToolbox\\Toolboxes\\Data Management Tools.tbx")

# Get the active map document
import arcpy.mapping
mxd = arcpy.mapping.MapDocument ("CURRENT")

# List broken links
# arcpy.mapping.ListBrokenDataSources(mxd)

# Set environment settings
ws = "C:\\PARK"
arcpy.env.workspace = ws
print("env.workspace completed successfully")
installdir = arcpy.GetInstallInfo("desktop")

#check GDB exist
gdb_nam = "test.gdb"
gdb_full_path = os.path.join(ws,gdb_nam)
if os.path.exists(gdb_full_path):
    arcpy.Delete_management(gdb_full_path)
print("GDB checked")

# Execute CreateFileGDB
arcpy.CreateFileGDB_management(ws, gdb_nam)
print("CreateFileGDB completed successfully")
outWorkspace = gdb_full_path

# Iterate

list = []

for df in arcpy.mapping.ListDataFrames (mxd):
    for fc in arcpy.mapping.ListLayers (mxd, "", df):
        if (fc.datasetName.endswith(".shp")):
                    list.append(fc)
for fc in fclist:
    arcpy.FeatureClassToGeodatabase_conversion(fc, gdb_full_path)


texts:
I created an own map for my p&p rpg group and now I want to transfer it into QGis.
That's easy to make, but I'd like to have an own CRS for that map but don't know
how to use the right proj.4 definitions.
Maybe anyone here can help me with that?
The map will be 30.000km from West to East and 21.000km from north to South.

David


texts: How to auto-add new values to unique values classification values range?

texts:
I am using ArcGIS Desktop 10.1, ArcGIS Server 10.1 & SQL Server 2008 software.

I am adding the query layer and directly publishing to ArcGIS Server and accessing the map in a JavaScript application. The database is auto-updating in every 2 hour.

I am changing the layer properties (Symbol ----category ----unique values) of the layer but it is not auto-updating in ArcMap for new records and in the REST service.

What I tried:

I have tried using geo-processing task but not successful.
Auto-rendering option in ArcGIS JavaScript API

So is there any workaround for this? (Should I need to reset the parameters in Geo-processing task?)


texts: Using SAGA Cluster Analysis for grids projection

texts:
I am having an issue with the SAGA cluster analysis for grids tool within QGIS 2.4.  I am actually taking this as an online class through Del Mar and we are going over supervised and unsupervised classifications of aerial imagery.

I have my CRS set as EPSG 32616 for the shapefile I want to classify, however, when using the tool I keep getting an error that the CRS is undefined and it defaults to 4326.

On top of that, the resulting output is upside-down.  When looking in the layer properties the CRS is 32616 for the Clusters but the Statistics table is defined as 4326.  My Project Properties are set to on the fly and I set the CRS as 32616.  I have attached an image of the issue as well.

Anyone know what may be causing this to happen?

I have gone through this 3 or 4 times with the same result.


texts:
In my web application, I am retrieving geospatial data from a server backend (Java) to the browser frontend (JavaScript). I am considering two options of how to serialize the data in this communication - XML and JSON.

XML seems to be a faster option as it can be deserialized by browser itself (into a documentElement), while JSON must be deserialized by a JavaScript library. On the other hand, there is a standard for storing geospatial data in JSON (GeoJSON), so using JSON would make the application more standard-compliant.

So my question is, is there an XML equivalent to GeoJSON? I found out that the term GeoXML has already been coined, but it is the name of a library, not a format standard.


texts: Apply a mask to multiple images(for same geo location) in a image collection

texts:
I have a image collection which holds multiple images of the same area. I can apply mask individually to those images but how do I apply a mask to all the images in the collection. Say the mask i have to apply is to filter out all the pixels with value more than a given value. Since i have multiple images for the same location, I want to 'and' the logic so that i get a single image with pixels displayed which satisfy the condition for all the images. Here is some code i wrote which works for fine when i have single image:

    Map.setOptions('SATELLITE')
    var image = ee.ImageCollection("imagesForIrvineCalifornia");

    var min = image.aggregate_min('min');
    var max = image.aggregate_max('max');

    var VIS_OPTIONS = {
      'b1': {
        description: 'b1' +
                     'b1 descripion',
        visParams1: {min:min.getInfo(), max:max.getInfo(), palette: ['DCF5E9', 'BDDEC9', '9DC7AA', '82B38E','689E75', '4E8A5B', '387847'], bands: ['b1']}
      }
    };

    var maskingFunc = function(val){
      return function(i){
        var band = i.select('b1');
        var mask = band.lte(val);
        return i.updateMask(mask);
      }
    }

    var maskedImage = image;
    var cutoffVal=3000 ;
    maskedImage = image.map(maskingFunc(cutoffVal));

    Map.addLayer(maskedImage.select('b1'), VIS_OPTIONS['b1'].visParams1, '(b1)', true);
    Map.centerObject(Map.layers().get(0).get('eeObject'));