texts: ERROR 000732 in code to Save shape as .lyr?

texts:
I am using a code (modified from Midavalo) to save .shp files as .lry files. The code runs well and displays the .lr files in the table of Contents. However, the Python displays Runtime error. On the other hand, I don’t know how to save .lyr
I need the .lyr files to export as .kml files

import arcpy
from arcpy import env
env.workspace = r"D:\ss3\polygons"
buf= r"D:\ss3\polygons\001.shp"
ou2= r"D:\ss3\layers\001.lyr"

arcpy.MakeFeatureLayer_management(buf, 'buf_layer')
arcpy.SaveToLayerFile_management('buf_layer', ou2, "ABSOLUTE")

This the error:

  Runtime error  Traceback (most recent call last):   File "",
  line 7, in    File "c:\program files
  (x86)\arcgis\desktop10.5\arcpy\arcpy\management.py", line 6965, in
  MakeFeatureLayer     raise e ExecuteError: ERROR 000732: Input
  Features: Dataset D:\ss3\polygons\001.shp does not exist or is not
  supported


texts: I have a raster file which has a range (min: 80 max: 120). I would like to reclass this raster by using quantiles. Is that possible?


texts: Changing pixel size in 'asc' file using QGIS?

texts:
I have soil data and BioClim data from Australia

soil data:

rows 2200
columns 1400
pixel size 0.0025, -0.0025 255

BioClim data:

rows 660
columsn 420
pixel size 0.00833333, -0.00833333

I am using QGIS but happy to use commandline or Python. How can I resample the soil data to match the pixel size of the BioClim data? I would like to add this to my Maxent model and I am sure it needs to be the same resolution.

I have looked for duplicate questions but could not see anything that specifically did this.


texts: Openlayers 3: Draw Polygon (convex hull) around a cluster of points

texts:
I'm grouping a set of points as a cluster. Then I'm trying to draw a polygon around the points.

Please see this link: Draw Poly

I tried by sorting the points based on co-ordinates. But still I'm not able to draw the polygon correctly.

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

var map = new ol.Map({
  layers: [
      new ol.layer.Tile({
      source: new ol.source.OSM()
      }),
      new ol.layer.Vector({
          source: vectorSource
      })
  ],
  target: 'map',
  view: new ol.View({
    center: [-16, 75],
    zoom:2
  })
});

//clusters
var min=1;
var max=20
for (j = 1; j < 3; j++) {
var polyCoords = [];
for (i = 0; i < 10; i++) {
    var x = Math.random() * (max - min) + min;
    var y = Math.random() * (max*2 - min*2) + min;
  polyCoords.push( ol.proj.transform([-x,y],'EPSG:4326', 'EPSG:3857'));

  var point = new ol.geom.Point([-x, y]).transform('EPSG:4326', 'EPSG:3857');
  var fea = new ol.Feature({geometry:point});
  vectorSource.addFeature(fea);
}

var thing = new ol.geom.Polygon( [polyCoords]);
var featurething = new ol.Feature({
    name: "Thing",
    geometry: thing
});
vectorSource.addFeature( featurething );
min=max+1;
max=max*2;
}

How can I draw a polygon around all the points?


texts: Copying only certain Fields/Columns from shapefile into new shapefile using Modelbuilder in ArcCatalog?

texts:
I am trying to take a shapefile, with way too many fields/columns, and copy over only certain fields into a new shapefile.

Normally, I would just delete the fields I don't need from the original shapefile, but I don't want to come across an instance where the GIS operator decides to remove a field that my automated .py scripts is trying to delete, causing the program to fail.

At this point, I have not physically tried anything because every tool I have come across I have been able to find documentation proving it's lack of usefulness (Fields -> Transpose Field and General -> Copy).

Fields -> Transpose Field does copy individual fields but it transposes them from a column into a row appendage, which does not help.

General -> Copy copies the entire shapefile and will not let me choose only certain columns. What I want is something similar to this tool but one that will allow me to choose only certain fields to include in the copy.

I use "Feature Class To Shapefile (multiple)" to convert my features into shapefiles, and I tried looking into a "Feature to Shapefile" solo Tool that would let me pick the columns/attributes of the feature to convert. That has so far turned up empty as well.

I am using the ArcCatalog ModelBuilder to perform this task, not raw Python scripts. I know enough Python to go that route if needed, but a Toolbox Tool is the preferred goal. This task will be part of a larger whole that will be run via a batch file calling the exported Python code from the Models. It has to be something that can be run automatically, so manually copying/pasting is not applicable.

As for a screenshot of the Model, my current setup has each task as individual models: one model copies rows, one creates shapefile from feature class, one adds a field, one populates that field, etc. I'm sure this is probably inefficient but this task has been passed on to me from someone else and this was how it was done. I can add screenshots of the other models I have but I'm not sure what use that will provide, as it doesn't give much of a context.

Basically, I'm looking to do something like Copying data from one attribute field to another using ArcGIS Field Calculator? but with about 10 different fields. I don't want to do them individually, but using one or maybe two tools.


texts:
I downloaded forest cover rasters from Google Earth Engine with code from MapBiomas in 30 x 30m resolution. However, when operating with the rasters, I saw that the resolution came as degrees: 0.0002694946, 0.0002694946  (x, y). From this question: How to convert resolution in degrees to resolution in meters? I understood it's not a good idea to convert degrees to meters, and I don't know how what else to do to get that data in meters since I already selected the resolution as 30 x 30m. I need to calculate the area for a cropped raster.
Is that still possible without having to redownload the data and, if I do need to redownload the data, how do I get GEE to actually give me the resolution in 30 x 30m?
This is the info of the first raster (I have a stack 2010-2019):
class      : RasterLayer
dimensions : 20700, 33207, 687384900  (nrow, ncol, ncell)
resolution : 0.0002694946, 0.0002694946  (x, y)
extent     : -53.11038, -44.16127, -25.35809, -19.77956  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
names      : mapbiomas.brazil.collection.50.saopaulo.2019
values     : 0, 255  (min, max)

And this is the GEE source code I used: https://code.earthengine.google.com/c94f49b3d7037af410daad3ec6b85b66?noload=1


texts: Error executing grass - v.net.steiner within R

texts:
I'm trying to use the v.net.steiner function within R, using the grass tools. But I get this error:

Error in execGRASS("v.net.steiner", input = "lx", output = "st", node_layer = "px",  :
  The command:
v.net.steiner --overwrite --verbose input=lx output=st node_layer=px terminal_cats=0-1000000
produced an error (1) during execution:
Number of terminals: 0
ERROR: Terminales insuficientes (<2)

Using this code (where all seems to work just fine until the steiner process):

library(rgrass7)
library(rgdal)
library(rgeos)
library(sf)
library(raster)
library(mapview)

# PARAMETERS
lpath <- "./ls.gpkg"
ppath <- "./ps.gpkg"
output <- "./Output/"
crs <- "+init=epsg:25829"

#-------------------------------------------------------------------------------
# GRASS config
#-------------------------------------------------------------------------------
G <- initGRASS(gisBase = "/usr/lib/grass76",
               home = getwd(),
               gisDbase = getwd(),
               location = 'grassdata',
               mapset = "PERMANENT",
               override = TRUE)
G
use_sf()
execGRASS("g.proj", flags = "c", epsg = 25829)

#-------------------------------------------------------------------------------
# PROCESS
#-------------------------------------------------------------------------------
l <- readOGR(lpath)
p <- readOGR(ppath)

# MAP OF LINES, NODES AND POINTS OF INTEREST
mapview(l, cex = 1, map.types = c("OpenStreetMap")) +
    mapview(p, color="red")+
    mapview(st_combine(st_as_sf(l)), color= "blue")

#-------------------------------------------------------------------------------
# DISSOLVE
#-------------------------------------------------------------------------------
lx <- aggregate(l, dissolve = TRUE)
# or
lx <- st_cast(st_as_sf(lx),'MULTILINESTRING')
#or
lx <- gLineMerge(lx)

# save layer to test it within qgis
st_write(st_as_sf(lx),"lx.gpkg")

#-------------------------------------------------------------------------------
# GRASS
#-------------------------------------------------------------------------------
# import clipped layers into grass
writeVECT(st_as_sf(lx), 'lx', v.in.ogr_flags = c('o','overwrite'))
writeVECT(st_as_sf(p), 'px', v.in.ogr_flags = c('o','overwrite'))
execGRASS("v.info", map = 'l')
execGRASS("v.info", map = 'p')

# create lines map connecting points to network (on layer 2)
execGRASS("v.net", input = 'lx', points = 'p', output = 'lxc',
          operation = 'connect', threshold = 500,flags = 'overwrite' )

# clean lines
execGRASS("v.clean",
          input = 'lx',
          output = 'l_clean',
          tool = c('break','snap'),
          threshold = c(10,10),
          flags = c('c','overwrite'))

# steiner tree
execGRASS("v.net.steiner",
          input = 'l_clean',
          output = 'st',
          node_layer = 'px',
          terminal_cats = '0-1000000',
          flags = c('overwrite','verbose'))

This reads two layers (lines and points), dissolves and cleans the lines, import into grass and perform an steiner tree.

¿Why do I do that dissolve part? Because, executing this v.net.steiner tool within QGIS it also fails, unless I do just a simple dissolve over lines.

Here is the error I get if I use just those lines and v.net.steiner:

...And here is the successful tool if I just do a simple dissolve over lines before:

The final goal is to execute this succesful v.net.steiner tree within R, by making the correct dissolve over lines or including another tool that avoids the error.

DOWNLOAD THE LAYERS

PD: I use QGIS 3.10 in Ubuntu 19.10.


texts:
Is there a way to do a screenshot of the composer in QGIS?

I am regularly doing maps in a rush to insert in reports or presentation and it could save time if I would not have to save them and then import them in the report/presentation.

In addition it would allow not to store plenty of draft maps.


texts: gpspipe returns correct timestamp, but wrong date (21 dec 1995)


texts: Avoiding intersecting lines in transect generation

texts:
I'm trying to generate transects using the QGIS transect tool under vector geometry as a step in generating an interpolation layer.  Unfortunately, the transects are intersecting, which will lead to problems in subsequent steps.  How can I trim the resulting transects so they don't intersect, or otherwise end up with the same result?

I can generate transects perpendicular to my source layer or a trend. Here I've used a smoothing and generalization routine to get transects following the trend:

But many of the transects cross if they are near any bends.  How can I avoid these intersections in the transects, to get something similar to this:

I could go through and manually trim each of the intersecting lines, but I'm looking for a more elegant/ less labor intensive solution.


texts: Change color for different string in geoserver

texts:
I have a DB in pgAdmin and in one table i have an field named "description" that only can have 2 descriptions: "car" or "bike"
in geoserver i set this label and work well:

<Label>
 <ogc:PropertyName>description</ogc:PropertyName>
</Label>
<Font>
 <CssParameter name="font-family">Arial</CssParameter>
 <CssParameter name="font-style">Normal</CssParameter>
 <CssParameter name="font-size">12</CssParameter>
</Font>
<Fill>
 <CssParameter name="fill">#b7d5ff</CssParameter>
</Fill>

but i want different colors for the 2 possibles descriptions. who is possible in geoserver?


texts: Extract highest point in raster and convert to point vector

texts:
I have a raster DEM and I want to extract from it the highest point found in every polygon of a polygon shapefile.

The result that I would like to obtain would be a point shapefile with 2 collumns: polygon ID and height.

I'm working in QGIS.

I tried SAGA raster statistics to polygon, but the result is just a polygon with the height value. I am trying to extract the point to a point vector layer.


texts:
I want a code or an option in QGIS or in a web map using Qgis2web (in leaflet or openlayers), when I click on a feature on the map, automatically it should register the ID (from attribute table) of the feature in an text document  file or in a csv file.

Is it possible?


texts: How to use Draw.js in a Leaflet offline map

texts:
I'm trying to build a data-collecting map using Leaflet. To do this I have to download the Leaflet JavaScript files but I can't find out how to download a working draw.js.

When I download draw.js I get this error

index.html:138 Uncaught TypeError: L.Control.Draw is not a constructor

I'm using this code for the map variable

var map = L.map('map', {drawControl: true, center: [42.381899, -71.122499], zoom: 13});

and in my code block

var drawControl = new L.Control.Draw({
  draw : {
    polygon : false,
    polyline : false,
    rectangle : false,
    circle : false
  },
  edit : false,
  remove: false
});


texts:
I am presently working with some open data in an original .csv file format with some lat, long columns.

I want to do a "Select by expression for each mouth" :
if I do for example : "DATE" > 2016-01-01 it select all the rows in the table.

I believe it could be a problem of data format in the column, but I am not sure.

Does someone know how I could confirm my query is correct and it is a way to confirm the data type in the column "DATE" ?


texts: Counting number of points differing in specific attribute within polygon in ArcGIS Desktop?

texts:
I have a grid of polygons and a shapefile with several types of points (POIs). Now I would like to count the number of POI classes lying within each polygon?

With spatial join I only managed to get the total number of POIs within one polygon.

Does anyone have an idea? I use ArcGIS 10.5


texts: How to calculate the grid point precipitation from the station data?

texts:
I have a problem to convert my hourly station precipitation data to the surrounding grid points.

I need to convert 4 station hourly rainfall data to 18, 60, and 240 grid points for 0.0625,0.03125, and 0.015625 respectively.

I know how to do it in ArcGIS.

I have 20 years hourly precipitation data that need to be interpolated to those grid points.

Making 20 years hourly raster data seemed impossible to me to start with.


texts: How to get a list of all objects contained in my different QGis layers

texts:
I would like to know how we could do in QGis to get a list of object contained in my different opened layers:
We have n .SHP files or PostgreSQL Views (they all contain at least an "object type" column and an "object id" column) opened in a QGis Project to make a map.
When the map is done, I want to save my .QGS file on a server and then I want (by a Python program ?) to get a .CSV file with 3 columns (object type, object id, path of my .QGS file) - or even better to write the data directly into PostgreSQL.

Does anybody ses how to make that easily on QGis ?


texts: Coordinates heat map results are not the same as the shapefile

texts: I have done with a heat map shapefiles as shown below. But the results of heat map coordinates are not the same as the coordinates of shapefiles. Is there a procedure that is wrong?


texts: Creating integer field from string with characters in QGIS

texts:
I am trying to create a new attribute field in a layer that converts a string with characters into an integer.
EXAMPLE:
234' -> 234
or
234'* -> 234
Any ideas on how to accomplish this?


texts: Is Orbit EOS Lite Free GIS Map Server an OpenSource software?

texts: I'm writing an essay about web mapping software and I can't find out is Orbit EOS Lite Free GIS Map Server an OpenSource software or not? Also, I would be thankful if someone gives me some links with some list of web mapping programs. Thank you.


texts: Creating Vector Attribute Field with Raster Value in QGIS

texts: I want to create a new attribute field contains data coming from a raster in a vector file. I tried the GRASS v.drape tool but the analysis did not create any new field in shapefile. Also, I tried the drape tool in QGIS and the result was the same.


texts: Is there a tool to track edits that were made?

texts:
Is there a tool/report that keeps track of what edits were made that I could run prior to or after reconciling my versions?

This would help me to know if there are any changes I need to make to any maps or databases that are disconnected from our SDE.

I'm using ArcGIS 10.1


texts: Is it possible to sync the extents of two open mxds?

texts: I've never seen it, and it's not vital, but I'm wondering if anyone knows of a way to sync the extents of two open mxds?  I'm pulling an old ims service into one mxd, and have another mxd open that i will use as the basis of an ArcServer service.  I'me trying to match up the layers and symbology, and it would be cool if I could change the extent in one mxd and have it sync with the other.  Is it possible?  Thanks in advance for any help.


texts: Open MapGuide MrSID rasters with FDO GDAL. invalid extent. Not drawing

texts:
When I try and save the resource I get: Spatial context extent appears to be invalid (or Default).

The XML:

<?xml version="1.0"?>

  OSGeo.Gdal.3.4

    Default
    This coordinate system has been overridden.
    LOCAL_CS["XY-MT",LOCAL_DATUM["X-Y",10000],UNIT["Meter", 1],AXIS["X",EAST],AXIS["Y",NORTH]]
    LOCAL_CS["XY-MT",LOCAL_DATUM["X-Y",10000],UNIT["Meter", 1],AXIS["X",EAST],AXIS["Y",NORTH]]
    Dynamic

        -10000000
        -10000000

        10000000
        10000000

    0.001
    0.001

#

These sids are 205mb a piece and theirs 53 of them (sat images). I decoded them to geotiff and they came out 2+gig ea.

What do I need to do to get these to display in MG? Whats the smartest method to get them web ready?


texts:
I am working on a volunteer project and we are in desperate need of a baseline layer. I have been unable to find it anywhere (there are plenty of shoreline and coastline layers out there) does anyone have a source for this data?

Thanks,
Dave
Added:
The baseline from which the territorial sea is measured is mainly the low-water line along the coast as marked on large-scale charts officially recognized by the coastal state.


texts: How to calculate false easting/northing for Lambert Conformal Conic projection

texts:
I'm defining a new GRASS location using the Lambert Conformal Conic projection. My understanding is that the values for the standard parallels and central meridian are dependent on my area of interest; I'm mapping Canada, for which common values are:

1st parallel: 49
2nd parallel: 77
Central meridian: -91.86667 (ie, 91 52)

I'm not sure how to determine where to put the central parallel (or if it matters) and what values to use for false easting/northing.

I have found lots of definitions of these parameters online, but no explanation of how to calculate them. Usually this means it would be self-evident if I knew what I was doing!

How do I determine the false easting/northing values for a LCC projection?

UPDATE

After some experimentation, I realize that the false easting/northing values, and the central meridian/parallel, don't change the projection at all; rather, they just shift the location of the origin (i.e., the point considered (0, 0)) for the map. That being the case, is there any convention for setting the origin and/or false easting/northing, or do we just set them to any convenient values for our particular project?


texts: Changing default ArcMap polygon selection symbology

texts:
The default ArcMap polygon selection symbol is a cyan outline.  This symbol can be manually changed for each layer within each map document using Layer Properties > Selection.
However, I would like to change the default selection symbol for all polygon layers in all map documents (including currently existing .mxds as well as future) from the cyan outline to some other symbol of my choosing, without having to cumbersomely modify each layer in each .mxd one-at-a-time.
For example, I may wish to change the selection symbology to a marker symbol, as shown in the screenshot below.  While I can do this per layer, per .mxd, I would like make such a symbol the application-wide default:

Customize > ArcMap Options does not offer this ability.  A search through the ESRI help documents, as well as an internet search turned up nothing.
How can I set the default polygon selection symbol from a cyan outline to some other symbol?


texts: UpdateCursor on two tables with value from SearchCursor

texts:
I try to keep myself busy and learn some Python.
I have two tables:

test2.shp (for searchcursor)

test.shp (needs to be updated)

I want to update test.shp with a UpdateCursor. The column which needs to be updated is column2 if column1 of test.shp satisfies a certain criteria (column1 = aa).
So, if column1 = aa, insert value in column2.

I try to take the values for the update from test2.shp with a SearchCursor.

That's my code:

import arcpy

arcpy.env.workspace = "c:\\Users\\python_ex"

scursor = arcpy.SearchCursor("test2.shp", "column1" = 'aa', "", ["column1", "column2"], "")
ucursor = arcpy.UpdateCursor("test.shp", "", "", ["column1", "column2"], "")

for srow in scursor:
    if srow == "aa":
        a = srow.getValue("column2")
    else:
        continue
        for urow in ucursor:
            if urow == "aa":
                urow.setValue("column2", a)
                cursor.updateRow(urow)
            else:
                continue
        del ucursor, urow
del urow, srow

I get the error: Keyword can't be an expression where I define the scursor.
But I guess that is not my only mistake...

I am working on ArcGIS for Desktop 10.2.2
Any ideas?

EDIT 8/7/14 12:20pm

Alright, this is what I came up with:

import arcpy

#set to folder where features are located
arcpy.env.workspace = "c:\\Users\\python_ex" #on windows use \\ instead of /

#---------------------------
#define variables for cursor
#---------------------------

"""SEARCHCURSOR"""
"""SearchCursor (in_table, field_names, {where_clause}, {spatial_reference}, {explode_to_points}, {sql_clause})"""
#field taken as "identifier". i.e. if column1 == a, take the value from the same row of column2
stable = "test2.shp"
sfield = ["column1", "column2"]
#where clause will be written in SQL and not assigned to variable

"""UPDATECURSOR"""
"""UpdateCursor (in_table, field_names, {where_clause}, {spatial_reference}, {explode_to_points}, {sql_clause})"""
#field taken as "identifier". i.e. if column 1 == a, take the value of the SearchCursor and insert into same row, but column2
utable = "test.shp"
ufield = ["column1", "column2"]
#where clause will be written in SQL and not assigned to variable

#--------------------------
#start the loop
#--------------------------

with arcpy.da.SearchCursor(stable, sfield, """"column1" = 'aa'""") as scursor: #SQL queries in python are build that way http://resources.arcgis.com/en/help/main/10.2/index.html#//002z0000001r000000
    for srow in scursor:
        print srow
        with arcpy.da.UpdateCursor(utable, ufield, """"column1" = 'aa'""") as ucursor: #if used with statements, no need to close all the cursors
            for urow in ucursor:
                ucursor.updateRow(srow)

BUT, that only takes the last the last found value of aa in test2.shp and updates all the other aa in test.shp with this value!


texts: How to address the new "Task-Completed" QgsMessageBar in Python?

texts:
I am really amazed what huge steps the QGIS development has taken in the last months. The new iconset, the improved Print composer and many more changes to be seen in QGIS 2.0 are really awesome!

Among many other cool improvements there is also this nifty blue task-bar, which I find very appealing.

Is it possible to address this bar in python for scripts and plugins (for instance to show a small message after finished computations)?
I am asking this question here, because i believe that many QGIS developers also visit this Q&A site.

PS: The bar is only visible in QGIS dev and only for some operations.


texts: Finding N nearest points to point using JTS?

texts:
I'm trying to build a mobile app, and find the x number of points closest to a users location that satisfy certain conditions that the user chooses on the front end.

My current solution is to incrementally search a small area around the user and expand the search area until the x items are found or the maximum bounds are reached, using a JTS QuadTree.

public static Set<String> searchQuadTree(Quadtree quadTree,
                                         final MyBoundary, data_bounds,
                                         final ArrayList<String> types) {

    //used to search QuadTree
    double shell = SHELL;
    final double EXPAND_SHELL = 5000;
    final double MAX_SHELL_SIZE = 70000;
    final int MAX_SEARCH_ITEMS = 20;

    Set<String> resultSet = new HashSet<String>();

    //create polygon to search for points that intersect with it
    Polygon polygon = data_bounds.searchAreaPoly(shell, 0);

    List<MyQuadNode> items = quadTree.query(polygon.getEnvelopeInternal());

    if(items != null) {

        GeometryFactory factory = new GeometryFactory();
        //search until the min search items are found using search criteria
        while ((resultSet.size() < MAX_SEARCH_ITEMS) && (shell < MAX_SHELL_SIZE)) {

            Iterator<MyQuadNode> iterator = items.iterator();
            while (iterator.hasNext()) {
                MyQuadNode item = iterator.next();

                Point point = factory.createPoint(getWGSCoord(item.getLongitude(), item.getLatitude()));

                    if (types.contains(item.getType()) && polygon.contains(point) && resultSet.size() < MAX_SEARCH_ITEMS) {
                        resultSet.add(item.getImageUrl());
                    }
                }
            }

             if(resultSet.size() < MAX_SEARCH_ITEMS) {
                //expand shell slowly at first;
                if (shell < 1000) {
                    shell *= 5;
                } else if (shell < 5000) {
                    shell *= 2;
                } else {
                    shell += EXPAND_SHELL;
                }

                //search a larger area
                polygon = data_bounds.searchAreaPoly(shell, 0);
                List<MyQuadNode> allItems = quadTree.query(polygon.getEnvelopeInternal());

                //avoid filtering the same items twice
                allItems.removeAll(items);

                items = allItems;
            }
        }
    }
    return resultSet;
}

This will work fine (although it doesn't really find the nearest neighbors) for a small data set, if the data set grows large however this seems inefficient, i'm estimating efficiency is around x(nlogn + n + n), where x = MAX_SEARCH_ITEMS.

Has this problem been solved previously?


texts: Installing GDAL on Ubuntu 20.04 - unmet dependencies

texts:
I'm trying to install GDAL on my Ubuntu 20.04 with the following commands:
sudo add-apt-repository ppa:ubuntugis/ppa && sudo apt-get update
sudo apt-get update
sudo apt-get install gdal-bin
sudo apt-get install libgdal-dev

After running the last line, I get the following error:
Reading package lists... Done
Building dependency tree
Reading state information... Done
Some packages could not be installed. This may mean that you have
requested an impossible situation or if you are using the unstable
distribution that some required packages have not yet been created
or been moved out of Incoming.
The following information may help to resolve the situation:

The following packages have unmet dependencies:
 libgdal-dev : Depends: libpq-dev but it is not going to be installed
E: Unable to correct problems, you have held broken packages.

What is the problem here?


texts: Why won't lines from GPX simplify/generalize/smooth?

texts:
I am working with GPX data. I imported them to ArcGIS by 'GPX to feature' and then I created lines from these points. I am trying to smooth (or simplify) the lines now, but the tools 'smooth line' and 'simplify line' are not working. (Maybe it is because of the big amount of points..)

Do you know any tools or hints how to smooth these lines created from GPX?

Thank you for your answers


texts: Why won't Mapsurfer OSM layers open anymore from the QuickMapServices plugin in QGIS version 3


texts: Is there a way to describe what side of the polygon is "inside"

texts: Let's say I want to create a polygon of northern half of earth, the line string of the polygon will be on the equator. In JTS (Java) with SRID 4326 is there a way to say if the north side of Earth or south side of Earth is the "inside" part of the polygon?


texts: How to convert GeoJSON with mixed and NULL geometries to SQLite?

texts:
I am trying to convert a GeoJSON file into SpatiaLite, so I can load it as a SpatiaLite layer into QGIS, as described here.

Converting the file in QGIS with "save as.." did not work (see also here), which is why I tried ogr2ogr in the command line to do the trick. After a few attempts I managed to convert it as follows:

/Library/Frameworks/GDAL.framework/Programs/ogr2ogr -f sqlite -dsco spatialite=yes mydata.sqlite mydata.geojson -nlt PROMOTE_TO_MULTI

This works fine for an input file with mixed geometries but does not seem to work with my file, which also contains NULL geometries. In this case, I end up with a .sqlite file that I can open in QGIS as a vector layer, but not as a SpatiaLite layer. Here it only loads the attribute table and no geometries, which are listed as "Unknown Geometries".

Does anyone have an idea how I could resolve this?


texts: Filtering dates by validity in QGIS 3.14

texts: I have a WFS data, which I'm trying to filter by DATE. I'm using query builder. I need to filter valid dates which are max three years older than current date. I can do it if I calculate three years my self, like "DATE" >= '2018-05-10T00:00:00.000'. Is there any expression that calculates that automatically every day, and shows only features that are valid?


texts: Change attribute field from whole number to decimal in QGIS

texts: Is there a way to change the number type once a field has been created in the attribute table in QGIS? I selected 'Whole number (integer)' initially but now need to switch to 'Decimal number (real)' in order to limit the number of decimal places to 1. I've already filled in all the rows of the column and so I'm wondering can I change the type and covert the existing numbers to a single decimal place or do I have to create a new field column selecting decimal from the start?


texts: OSGB Coordinates not showing correctly in QGIS 2.14

texts: The coordinate box at the foot of my map project is not showing correctly for OSGB1936 British National Grid EPSG 27700 in QGIS 2.14. This was working perfectly before I upgraded. Screen shot attached showing how they are represented. Does anyone know how to fix this?


texts:
I want to be able to extract two sets of coordinates which are the midpoints of two shorter sides of a rectangular polygon (purple points).

I know there is an intersection tool in QGIS which can get me those, however I can't figure out how to automatically generate the midlines which are used in the intersection tool. Is there a way of doing that in QGIS/GRASS? I would like to automate the process as much as possible so manual line creation is not an option.


texts: Extent of shapefile in meters instead of decimal degrees QGIS

texts:
i'd be thankful for your help.

My QGIS-project has the CRS WGS 84 EPSG 4326 and so does the shapefile I created. The coordinates on the edge of the window however shows, I think, meters, although I defined the displayed coordinates should be decimal degrees.

The other problem is:
I want to upload the shapefile on a google maps plug-in on a specific website. It doesn't work and I guess it has to do with the wrong extent of the shapefile are 7-digit numbers, but the units are defined as degrees.

Can somebody help here?


texts: Connection property set error on registering PostgreSQL geodatabase with ArcGIS for Server?

texts:
We're in the middle of a migrating process.

We've imported a geodatabase on a new Linux machine and we're trying to publish some map services to this new ArcGIS Server (that resides on the same machine).

We've upgraded the 10.1 geodatabase to a 10.2 one.

The geodatabase was imported via pgAdmin backup- restore procedure then we used the Upgrade Geodatabase tool from ArcToolbox to go to 10.2.

We've encountered problems when trying to register the geodatabase with the server.
The error we get is:

  The connection property set was missing a required property or the
  property value was unrecognized. This release of the GeoDatabase is
  either invalid or out of date. DBMS table not found

Same error occurs with a new created geodatabase as well.

Details:

ArcGIS Server 10.2 on a Linux Machine
PostgreSQL Geodatabase 9.2 on the same Linux Machine


texts: Adding VBA into a Python script generated by ModelBuilder

texts:
I have created a model that does a calculate field and have used VBA in the code block. But when I export it to Python, it doesn't work.
Here is the VBA code:

new_value = [HDATA_N_LAREF_ACTIVITY_TEXT]
if len(new_value) > 100 then
    desc = Left(new_value,100)
    desc = Left( desc, InStrRev( desc, " ") - 1) & " " & "..."
end if

When i exported it it came out as this....

# Process: Calculate Field...
gp.CalculateField_management(DWWEB_FINAL, "DESCRIPTION", "desc", "VB", "new_value = [HDATA_N_LAREF_ACTIVITY_TEXT]
if len(new_value) > 100 then
    desc = Left(new_value,100)
    desc = Left( desc, InStrRev( desc, \" \") - 1) & \" \" & \"...\"
end if
")

EXTRA INFO:

I have had a look at a previous script that was created (before I started this role), where within the Python script the VBA code is all on one line and contains "\n", where there is a new line within the code, so it should end up like this:

gp.CalculateField_management(myTempWorkspace + os.sep + "DWWEB_Joined", "DESCRIPTION", "desc","VB","new_value = [HDATA_N_LAREF_ACTIVITY_TEXT]\\nIf len(new_value) > 100 then \\ndesc = Left(new_value,100)\\ndesc = Left( desc, InStrRev( desc, " ") - 1) & " " & ""...""\\nend if\\n\\n\\n")

but in PythonWin it gives me the following error:

>>> Unhandled exception while debugging...
Traceback (most recent call last):
  File "D:\ScheduledTasks\Scripts\DES_HYMS_Streetworks_v1.py", line 255, in ?
    gp.CalculateField_management(myTempWorkspace + os.sep + "DWWEB_Joined", "DESCRIPTION", "desc","VB","new_value = [HDATA_N_LAREF_ACTIVITY_TEXT]\\nIf len(new_value) > 100 then \\ndesc = Left(new_value,100)\\ndesc = Left( desc, InStrRev( desc, " ") - 1) & " " & ""...""\\nend if\\n\\n\\n")
RuntimeError:
(Microsoft VBScript compilation error) Syntax error
Failed to execute (CalculateField).

But that doesn't run, so what do i have to do to get the code block to work?

Thanks


texts: Applying machine learning methods in orfeo toolbox into DEM raster

texts:
I would like to apply machine learning methods in orfeo toolbox to classify some landscape features such as shrub and coarse woody debris in woodland. I created DEM raster image from very dense ground LiDAR point cloud (1 cm resolution).

First, I vectorized three types of 71 training samples.

Then I calculated image statistics with CalculateImageStatistics and trained the image with TrainImagesClassifier libsvm classifier successfully.

Finally, I run ImageClassifier to classify the image. However, the result image has pixels with only 0 value. I tried with random forest classifier as well, but got the same result.

Is it not working because orfeo toolbox cannot classify a DEM raster dataset?


texts: PyQGIS custom script writes file as optional

texts:
I had this problem in my custom script in PyQGIS.

I create the variable as input in this way:

##lista_formaciones_forestales_con_estrategia = file

In theory this is the correct way to write a non-optional parameter as input, but when I open my script it is writing as optional:

Is it a problem with my QGIS version? My QGIS version is 2.16.1.


texts: ArcGIS Server cache from outside wms service

texts:
I have an arcgis server web application (version 10.0) that utilizes an external wms service for displaying some information layers on top of a local basemap. This external service has a tendency to be rather slow and buggy so I am wondering if it would be posible to use a remote wms service as a source for creating a local cache from within arcgis server, following its tile scheme.

The method would be something like: a) connect to thewms service in ArcGIS Desktop b) save the mapfile and publish it to our local server c) set the appropriate scales via ArcCatalog d) creating the cache from ArcToolbox

Has anybody tried this approach and are willing to share their experience?


texts: Overlay analysis of two raster by attribute?

texts:
I have two rasters (one showing land use and the other showing soil texture). I'd like to find the areas where there is a specific land use on a specific soil texture.

Is there a way to overlay two rasters but with a boolean statement or some kind of qualifier?

I have tried zonal statistics but each attribute table for my rasters have multiple columns and I only want to compare values between two of the columns (soil texture in one table and land use type in the other table).


texts: PostGIS uuid type not visible in QGIS or GeoServer

texts:
This is probably a simple question, but nothing comes to my mind right now.

I have created a view in PostGIS containing multiple columns, including one column of data type uuid.

When I try to visualize the view as a layer in either QGIS or GeoServer the column of data type uuid is simply omitted while the rest of the columns looks fine.

Is the data type uuid not supported by neither QGIS or GeoServer or is it a simple mistake on my part?

Versions:

QGIS: 1.7
PostGIS: 1.5.2
PostgreSQL: 9.0.1
GeoServer: 2.0.2


texts: Adding HERE map to ArcGIS Desktop using WMTS?

texts:
I need here maps on ArcGIS Desktop so can ı do it using WMTS or API?

I tried add wmts services with Portable Basemap but did not work

I searched wms or wmts services for here map but ı did not find it

I am using Arcgis version 10.4.1


texts:
In the 10.3 REST API they introduced:

Supports pagination in a query layer. Use resultOffset and resultRecordCountparameters to page through a query result.

I don't see any options to take advantage of these features in the ArcGIS .NET Runtime SDK.

Am I missing it? Is it a planned feature? Is there way I could add this functionality?


texts: Create Table with postGIS geometry in Python

texts:
I'm creating the following table with a postGIS geometry column in python but get  an error psycopg2.ProgrammingError: type "geometry" does not exist
 What should I use to be able to use postGIS in python?

import psycopg2

def main():
        conn_string = "host='localhost' dbname='postgres' user='postgres' password='postgres'"
        conn = psycopg2.connect(conn_string)
        conn.set_isolation_level(ISOLATION_LEVEL_AUTOCOMMIT)
        cursor = conn.cursor()

        createTable = """  CREATE TABLE map
                        (
                            id INT PRIMARY KEY DEFAULT NULL,
                            geom GEOMETRY DEFAULT NULL,
                        );
                        """
        cursor.execute(createTable)
        conn.commit()
        cursor.close()
        conn.close()
if __name__ == "__main__":
    main()


texts: Layers created by GeoServer do not display in Mapbender3

texts:
I have installed GeoServer 2.7.2 on my local windows 7 machine, which retrieves data from a remote PostgreSQL server. The data is represented in layers just fine!

I also have installed Mapbender3 locally and I have introduced the local GeoServer as the layer source. The problem is that the layers are not displayed in Mapbender application while their names are listed there in Mapbender's layer tree. I have set the map's SRS to ESPG:3395 as it is declared originally in GeoServer.

Please help me with why the layers are not loaded in the Mapbender application!


texts: how to get IVersionedTable from IVersionedWorkspace

texts:
I'm creating a custom Audit History component that will show old values (base values) and the new values (delta versioned tables).

I have been able to get the ObjectIDs of different change types between the base and delta tables changes such as:

esriDifferenceType.esriDifferenceTypeInsert

esriDifferenceType.esriDifferenceTypeDeleteNoChange

esriDifferenceType.esriDifferenceTypeUpdateNoChange

etc.

but I can't figure out how to access the actual versioned table (IVersionedTable - I'm assuming) so that I can compare each field's value to the base table (If there was an update, I need to know which field was updated, and the difference between the 2 values). Or, is there a better way of doing this?
Am I on the right track?

private static void getVersionedFeature(string parentVersion, string specificVersion)
    {
        IWorkspace pWorkspace = connectToAtlasEditSDE();
        IVersionedWorkspace pVersionedWorkspace = (IVersionedWorkspace)pWorkspace;
        //get a enumeration of all the versions on the versioned workspace
        IEnumVersionInfo enumVersionInfo = pVersionedWorkspace.Versions;
        enumVersionInfo.Reset();
        IVersionInfo versionInfo = enumVersionInfo.Next();
        while (versionInfo != null)
        {

            string versionName = versionInfo.VersionName;
            if (versionName == parentVersion)
            {
                IEnumVersionInfo childrenVersionInfo = versionInfo.Children;
                IVersionInfo childrenInfo = childrenVersionInfo.Next();
                while (childrenInfo != null)
                {
                    string childName = childrenInfo.VersionName;
                    if (childrenInfo.VersionName.Equals(specificVersion))
                    {
                        // get table here?!
                    }
                    childrenInfo = childrenVersionInfo.Next();
                }
            }
            versionInfo = enumVersionInfo.Next();
        }
    }


texts: QGIS/GeoJSON: Dissolve featurecollection consisting of polygons and multipolygons into a large multipolygon

texts:
I've read this question. I would like to know how to do this in QGIS.

I have a Geojson file with a featurecollection, starting with:

{"type":"FeatureCollection","features":[{"type":"Feature","geometry":{"type":"Polygon","coordinates":[[[5.21846588339...

But there are also quite a few MultiPolygons in the collection.

When I do the following:

Import the GeoJSON file in QGIS (add layer --> vector layer)
Select all of it
Do a dissolve (vector --> geoprocessing tools --> dissolve)

Then all the normal polygons get merged into a giant multipolygon, but not the multipolygons in the original collection.

What is the best way to include them? I've thought about manually "splitting" them into polygons in the geojson file with a text editor (because individual feature metadata is not important/removed after the dissolve anyway). But I have to do this more often in the future, and it seems like this should be doable automatically :)

Before Dissolve:

After Dissolve:


texts: Expression-based labeling in QFIELD 0.8.2 not working

texts:
It seems that labeling is not possible in the actual versions of QFIELD any longer (I'am working with 0.8.2 now) when there is a scale-based expression There is no labeling in this case although the new "busy-indicator" shows after disappearing - that rendering /labeling came to an end.
Working for example with a simple ASE WHEN $scale <2000 THEN "fieldname" expression does work in QGIS 2.14/2.16.2 but not in QFIELD (its the same with using variables instead of numbers for the scale). Simple labeling without expression works without problems.

(I am using a sqlite-database with a multipolygon-geometry for the test with Android 5.0.1.)
(In my regular project I use more complicated conditions to build different contents by a scale-based expression - this very simple example might have worked by other settings in the labeling-gui - but I tried to find out where the problem lies)


texts: How to fix incorrect land cover pixel data in GeoTIFF

texts:
I have downloaded a large dataset of 7,046 global land cover TIFF files from here:
http://data.ess.tsinghua.edu.cn/fromglc10_2017v01.html
The dataset is 10m per pixel so is very high resolution, but the problem is that it is not accurate in some areas. The biggest problem is in cities with tall buildings, where the building's shadow is interpreted as "water" instead of "artificial surface". Note that this issue also happens in some other places as well (e.g. coal fields) but this is less important for me than correctly showing the land cover in major cities.
This image shows an example of Sydney, Australia which highlights the problem:

The blue "water" in the middle of the image is tall building shadows which should be red (artificial surface).
I would like to correct this situation but the only way I can think of doing it is manually using a process as follows:

Start working in descending order from the most populated cities in the world with tall buildings, e.g. from this link: https://en.wikipedia.org/wiki/List_of_largest_cities. This is simply a prioritisation exercise so that I'm focusing on the largest cities first and leaving smaller cities until later (if ever)
Using GIS software, overlay satellite imagery over the GeoTIFF at 50% transparency and visually compare the GeoTIFF data with actual land information (i.e. a visual check) to make sure it really isn’t water. For example, many cities have rivers, lakes or ponds in the middle of them that could be correctly classified
Re-colour any inaccurate pixels as the correct land type (e.g. artificial surface in most cities), taking care if the city borders a genuine water body since the land border should be drawn as accurately as the rest of the GeoTIFF file. This is easier said than done since the software I use (Global Mapper) doesn't seem to allow you to change raster data, so this could be a laborious exercise of exporting to image software (e.g. PhotoShop) to change the colour values of specific pixels

Has anyone ever fixed a similar issue to this before? Is there a more efficient way of solving this problem than the manual process I've outlined above? Even better, does anyone know of an automated solution?


texts: Choosing polar projection and tiling scheme for web mapping of polar regions?

texts:
I am looking for examples of web maps focusing on polar regions (Arctic or Antarctic) preferably using a polar projection.

What polar projection and tiling scheme could be used for such web mapping applications?


texts: Memory leak while using in_memory ArcGIS 10.1

texts:
I am running a Python script to make point feature layer from a series of netCDF files and then perform a spatial join with a number of polygons and save the results as a table. I am using arcpy.MakeNetCDFFeatureLayer_md and address its output to a variable in 'in_memory' workspace. My problem is that even though I am deleting the 'in_memory' after each loop using arcpy.Delete_management('in_memory'), it does not free up my RAM and after a couple of iterations I run into Memory Error. The function that is being called for each netCDF file looks like this:

def spatialJoin(nc):
    arcpy.env.overwriteOutput = True
    arcpy.env.workspace = 'in_memory'
    arcpy.MakeNetCDFFeatureLayer_md(nc, "SM", "lon", "lat", "SMPoints", "time;lat;lon", "", "", "", "BY_VALUE")

    <does the spatial join>

    print arcpy.Exists('SMPoints')
    arcpy.Delete_management('SMPoints')
    print arcpy.Exists('SMPoints')

And the console prints:

True
False

for each iteration, which means Delete_management is working. But checking the python process in windows task manager shows that the data is piling up in RAM! I have tried deleting individual features inside 'in_memory' as well as deleting the whole 'in_memory' but nothing seems to work. Am I missing something here?

P.S. Spatial join does not affect the problem here, problem exists even when I completely comment out the spatial join section.


texts: Creating Info Window in CartoDB to scroll through features for different points with same geometry?

texts: How do I create an info window where you can scroll through and see features for different points that have the same geometry? I am making a map of newspapers. For example, there's a cluster of 6 at New Orleans and 272 at New York. My info window only shows info for one newspaper even though there are obviously many more. How do I fix this?


texts: My script in pyqgis run (grass7:I.pca) runs fine but output folder empties

texts:
this is my PyQGIS scrpit (QGIS 3.10).  this code runs fine but the output folder doesn't show results.  the MyFeedBack function shows me that it runs correctly and matches the Registry that is generated when the action is performed from the Grass7 i.pca module.  please i need help i don't know what else to try

import os
import sys
import gdal

path=r'D:\ESCUELA\clip_L8'
output= r'D:/ESCUELA/PCA'

bandList = [band for band in os.listdir(path)if band[-4:]=='.tif']
input = [os.path.join(path,band) for band in bandList]
print(input)

class MyFeedBack(QgsProcessingFeedback):
   def setProgressText(self,text):
       print(text)
   def pushInfo(self,info):
       print(info)
   def pushCommandInfo(self,info):
       print(info)
   def pushDebugInfo(self,info):
      print(info)
   def pushConsoleInfo(self,info):
      print(info)
   def reportError(self,error, fatalError=False):
      print(error)

def pca(input,output):
params ={'input':input,\
    'rescale':[0.300],\
    'percent':99,\
    '-n':True,\
    '-f':False,\
    'output':output,\
    'GRASS_REGION_PARAMETER':None,\
    'GRASS_REGION_CELLSIZE_PARAMETER':0}
feedback = QgsProcessingFeedback()
res = processing.run("grass7:i.pca",params,feedback=MyFeedBack())
print(res)
return

pca(input,output)


texts: Geographic Transformation in Alaska Between NAD83 (NSRS 2007) and NAD83

texts:
What is the best geographic transformation for Alaska to go from NAD83 to NAD83 (NSRS2007)?

Esri available transformation for GCS_NAD_1983_NSRS2007 to NAD83:

WGS_1984_(ITRF00)_To_NAD_1983 + NAD_1983_NSRS2007_To_WGS_1984_1
NAD_1983_To_WGS_1984_5 + NAD_1983_NSRS2007_To_WGS_1984_1
NAD_1983_To_WGS_1984_1 + NAD_1983_NSRS2007_To_WGS_1984_1

The suggestion from the below post is to use WGS_1984_(ITRF00)_To_NAD_1983.
 Geographic Transformations in Alaska Between NAD83 and WGS84

Looks like this is still the best option for Alaska for NAD83 to WGS84. I am wondering about the second half of the transformation "chain." Because NAD_1983_NSRS2007_To_WGS_1984_1 is a null transformation it is putting the NAD83 into the WGS84 realization that is tied to ITRF 2000 and then saying that this is equal to NAD83 (NSRS2007). I would think that we would need a transformation "chian" like WGS_1984_(ITRF00)_To_NAD_1983 + WGS_1984_(ITRF00)_To_NAD_1983_CORS96 to account for the change from NAD83 to NAD83 (NSRS2007).

(Transformation parameters can be found in this Esri document:
http://downloads2.esri.com/support/TechArticles/geographic_transformations_1040.pdf)


texts: Polygons generated using python-script appear stretched?

texts:
We are trying to generate a grid of hexagons that is supposed to span the whole of Germany. The code we are using to generate this grid is as follows (excluding some sqlalchemy related parts):

from shapely.geometry import Polygon, Point
from pyproj import Proj, transform

class HexGenImpl(object):
def __init__(self, size):
    self.size = size

@property
def width(self):
    return self.size * 2.0

@property
def height(self):
    return self.width * ( math.sqrt(3) / 2.0 )

def colFactor(self, col):
    return (col * self.width * (3.0 / 4.0))

def rowFactor(self, row):
    return (row * self.height)

def __call__(self, col, row):
    for angle in range(0, 360+60, 60):
        theta = math.radians(angle)
        x = self.colFactor(col)
        y = self.rowFactor(row)
        # Offsetting every other column to fit the grid
        y += (col & 1) * 0.5 * self.height
        x += math.cos(theta) * self.size
        y += math.sin(theta) * self.size
        x = round(x, 10)
        y = round(y, 10)
        if x == -0: x = 0
        if y == -0: y = 0
        yield (x, y)

def generateHexGrid():
    session = Session()
    hexagon_generator = HexGenImpl(.1)

    ap = Proj("+proj=utm +zone=32 +ellps=WGS84 +datum=WGS84 +units=m +no_defs")
    # Formerly just Proj("+proj=EPSG:4326")

    for a in range(100):
        for b in range(400):
            hexagon = hexagon_generator(a, b)
            hList = list(hexagon)
            glist = [ap(p[0], p[1]) for p in hList]
            pol = Polygon([[p[0], p[1]] for p in glist])
            tHex = Hexagon(geom=pol.wkt)
            session.add(tHex)
        print(a)
        session.commit()

This code works in that it generates a grid that can be rendered using our Geoserver without any problems. Sadly, whenever we try to use the layer in our front-end-application (which uses Leaflet) the hexagons appear stretched. This seems to be a problem with the projection, yet I haven't been able to figure out how one would go about projecting our coordinates (which are basically points on a cartesian plane) to something that gets rendered properly.
The closer the hexagons are to the equator, the less they appear to be stretched out.

I tried multiple projections, offsetting the coordinates, changing settings on the GeoServer etc. yet the problem persists.

How do I go about this?


texts: Raster Calculator or Local-Combine tool? (in ArcPy)

texts:
i have 3 classified rasters and i wanted to use the following query from r. calculator:

(rclslope>=45) & (rclmosaic>=2000) & (rclfocal>15)

to arcpy. How can i write it?
Is it better to use local - combine tool with arcpy to do it?

I am using ArcGis 10.2


texts: Warping radar maps for overlay into Google Maps

texts:
For a small project I'm looking to add various weather maps to overlay on Google Maps.  I've can pull several maps from various sources, however I also know that they are not in the right projection to just dump into Google Maps.
Unfortunately I'm not an expert in GIS.
From what I've been able to read, I can assume that the images I get from the radars are Gnomonic projection, and I can get a reasonably accurate Latitude and Longitude of they centre pixel of the map.  The scale of the map (in terms of distance for each concentric circle) is also published.
So far I've pieced together that I need to use goal-translate and gdalwarp to first locate the image, and then warp it into the EPSG:3857 projection that Google Maps uses.
I read someone else's question and the answers in:
Convert gnomonic projection to Mercator
Form there I see the steps as:
gdal_translate -a_srs "+proj=gnom +lat_0=[Centre Lat] +lon_0=[Centre Lon]" -a_ullr [EXTENTS] infile.png translated.png

to actually locate the image with a projection, and from there I can use
gdalwarp -overwrite -s_srs "+proj=gnom +lat_0=[Centre Lat] +lon_0=[Centre Lon]" -t_srs EPSG:3857 translated.png warped.png

My issue comes from that while I understand the basics of what I'm doing, I'm clueless as to what to use in the [EXTENTS].  I understand that they represent the coordinate edges of the image, but I don't know what units they are or how to get them.
I see that the person asking the question referenced cs2cs which looks useful and I can get information out of it, however I'm unsure of what values to put into it to get what I want.
To provide a concrete example:
I can fetch a radar image from Environment Canada.

From there I crop out the border and scale to the side, and get an image that is 478x478 px.  Form the published data I know that the centre pixel is at -113.3881 51.0264, and that the scale is 1 px/km.  Meaning that the mid points of the edges of the image are 239km from the centre.
From here I am reasonably sure that I need to go through gdal_translate to first locate the image, but I'm not sure what to use for extents.  I think cs2cs could help me, but I'm not sure what to provide it for the coordinates.
Can I generate these extents from the knowledge that the edges are 283km from the centre? (at least in the mid points, I know that doesn't apply to the corners)
Am I going down the correct rabbit hole at all?
Ultimately I would love to generate a pipeline that takes in:

The image
The lat/lon of the centre of the image
The scale in km/px of the image

Warp the image to be suitable to display on Google Maps, preferably with the pipeline telling me what the lat/lon extents I need to use for Google Maps are.
Could anyone nudge me in the right direction?


texts: Layer group with large number of layers vs large vector layer in GeoServer?

texts:
In my Geoserver project, I need to put each feature in a separate layer, and combine all these layers in one layer group.

The problem is when the number of features increases ( > 1,000 feature) the performance of displaying and querying the layer group degrades significantly compared to putting all the features in a single layer. I'm using Openlayers for the display.

The expected number of features is > 10,000 feature, so it will be a big performance issue.

What is the reason, and what is the solution?


texts: (ms4w) oracle spatilal problem of drawing compound polygon for Oracle 10g

texts:
When i draw a compound polygon

SDO_GEOMETRY(
    2003,  -- two-dimensional polygon
    NULL,
    NULL,
    SDO_ELEM_INFO_ARRAY(1,1005,2, 1,2,1, 5,2,2), -- compound polygon
    SDO_ORDINATE_ARRAY(6,10, 10,1, 14,10, 10,14, 6,10)
  )

by using msw4 accompanied with mapserver 5.4 or the previous version, they are successful to create an image.

But when i change to msw4 3.0.6 with mapserver 6.0.3 (libmap.dll of Oracle 11g), i always get the problem :

msDrawShape(): General error message. Only polygon shapes can be drawn using a POLYGON layer definition.

I think it comes from the file libmap.dll, maybe libmap.dll for Orl 11g doesn't work welll with Orc 10g

Can you suggest me what to do ?

Iam using Oracle client 10.2.0.3, and the underlying database is Oracle db 10g (10.2.0.4)


texts: reload ServerVector with different zoom levels in openlaers3

texts: how to reload ServerVector with different zoom levels ? I don't want to use ol.loadingstrategy.createTile strategy because it sends many requests with each zoom -one request with each tile - which takes some time. I want to force it to send one request in each zoom level because WFS parameters are changed in every zoom level


texts: How can I calculate different polygons by different attributes?

texts:
I'm doing my master thesis and want to calculate different polygons (in QGIS 1.8 for Mac) in four different columns which I defined by different attributes. i.e. some of the attributes in the first column are named "1112" and some "1113",.... So how can I calculate the area of the polygons in column "1" with the attributes "1112"?

The numbers here are only examples, my thesis is about greenspaces, so it would be to complicate to explain here...all of the numbers have a different meaning and have to be joined.
I need to sum these attributes to make a statement for the thesis.

Can anyone please help me?
Thanks a lot
 Markus

@M_Gebhart


texts: ArcMap crashes when opening extension folder?

texts:
I need to activate the Spatial Analyst extension, however each time I select "extension" ArcMap 10.3.1 crashes.

This is what I'm doing: while in ArcMap, I navigate to Customize/extension then the program crashes, everything single time. I can access the extension folder through ArcCatalog, but that is obviously not what I need.

Has this happen to anyone?


texts: Unable to connect Postgis from QGIS 1.8.0

texts:
I have  installed PostgreSQL and PostGIS  on my laptop (Win 7 64 bit OS). I am running QGIS 1.8.0 on my machine .

When I am trying to connect Postgis I am getting connection failed error then I tried following things :

Details:

Port :5433 User name : postgres  Password :postgres

I checked the pg_hba.conf file (All default settings)
Restart the PostgreSQL service and tried
I have reinstalled PostsreSQL and getting same error.
Tried using empty fields while connecting still getting an error. (Not sure what to put in service field..tried using postgres and postgresql)
My firewall is disabled.
The database PostgreSQL (version  9.2) is installed on my local machine also postgis (version 2.0.3-2) installed on my machine and both are for 64 bit setups.

7.Should I try lower version of Postgresql with qGIS ??

Researching on Google also and found similar  bugs related to connection issue(not sure)


texts: Accessing single scene image file in Google Earth Engine

texts:
I am reading the example in Google Earth Engine was many using the Image Collection but I can't find a single scene reading the image so I am trying to calculate for single scene using Landsat 8 OLI for 09 January 2019 but have a problem reading the scene then calculate NDVI from single scene, this is the script :
var tuban = ee.FeatureCollection("users/iqbal/Tuban");

var image = ee.Image('LANDSAT/LC08/C01/T1/LC08_L1TP_119065_20190109_20190130_01_T1')
// Get information about the bands as a list.
var bandNames = image.bandNames();
print('Band names: ', bandNames); // ee.List of band names

// Get projection information from band 1.
var b1proj = image.select('B1').projection();
print('Band 1 projection: ', b1proj); // ee.Projection object

// Get scale (in meters) information from band 1.
//var b1scale = image.select('B1').projection().nominalScale();
//print('Band 1 scale: ', b1scale); // ee.Number

// create function to add NDVI using NIR (B5) and the red band (B4)
var getNDVI = function(img){
  return img.addBands(img.normalizedDifference(['B5','B4']).rename('NDVI'));
};

// map over image collection
var l8ndvi = image.map(getNDVI);

var composite = l8ndvi.qualityMosaic('NDVI').clip(tuban);
print(composite);

// Visualize NDVI
var ndviPalette = ['FFFFFF', 'CE7E45', 'DF923D', 'F1B555', 'FCD163', '99B718',
               '74A901', '66A000', '529400', '3E8601', '207401', '056201',
               '004C00', '023B01', '012E01', '011D01', '011301'];
Map.addLayer(composite.select('NDVI'),
            {min:0, max: 1, palette: ndviPalette}, 'ndvi');

There is some error for accessing the single scene like this:
Band names:
List (Error)
Image.load: Image asset 'LANDSAT/LC08/C01/T1/LC08_L1TP_119065_20190109_20190130_01_T1' not found.
Band 1 projection:
Projection (Error)
Image.load: Image asset 'LANDSAT/LC08/C01/T1/LC08_L1TP_119065_20190109_20190130_01_T1' not found.
image.map is not a function
``
Any suggestions as a solution for this?


texts: Conditional operation with two raster files

texts:
I need to perform an operation on a raster file. My first attempt was with mapcalc in GRASS, but I got lost with the way it handles the 0 values, I am around 0, there migth be many zero values. I am working with double numbers. Maybe, it can be done in QGIS, which will be nice.

In simple words, this is like adding a building on a surface, but limiting it to 2 m and modifying anything if there no building or it is below surface.

The input:

Raster1: a regular dtm with some null values in the bounduaries (I could have it in many format)
Raster2: a same size and resolution but with a small surface within non null values of the previous raster

The logical operation, raster3 is the result:

if raster1 = null then raster3:=null
else
  if raster2 = null then raster3:= raster1
  else
     if raster2 <> null then
       if raster2-raster1>0 and raster2-raster1<2 then
          raster3:=raster2
       else
          if raster2-raster1>=2 then raster3:=raster1+2
          else raster3:= raster1

It looks that I have a code, but I have never coded anything in GRASS or QGIS


texts: Interrogating and displaying polygon created from points in PyQGIS

texts:
I've been working on this problem for several questions now - the most recent one here.  I believe the code below had created the polygon from the points, and added this polygon to a layer as a polygon in a feature object.  The layer is also added to the table of contents. However, for some reason I don't understand the polygon does not display on the canvas?
# Attach libraries
from qgis.core import * # attach main QGIS library
from qgis.utils import * # attach main python library
import os # attach operating system library

# Set a object for the current project
Prj = QgsProject.instance() # Object for current project

# Create an array [] object with the polygon vertices
vrtcs = []
vrtcs.append(QgsPointXY(396100,8969000))
vrtcs.append(QgsPointXY(396100,8973900))
vrtcs.append(QgsPointXY(397900,8973900))
vrtcs.append(QgsPointXY(397900,8969000))

# Create a polygon from the coordinates
ply_01 = QgsGeometry.fromPolygonXY([vrtcs])

# Create a feature object then put the polygon into the feature
ftr = QgsFeature()
ftr.setGeometry(ply_01)
print(ftr.geometry())

# Create a scratch layer for the feature, in the desired CRS
lyr = QgsVectorLayer('Polygon?crs=epsg:29194', '200909_Bdy',"memory")
Prj.addMapLayers([lyr])

# Set an object for the data provider for the layer (not sure why I need this?)
prv = lyr.dataProvider()

# Make the layer editable
lyr.startEditing()

# Add the feature with the polygon to the layer (this fails)
lyr.addFeatures(ftr)

# Save the layer changes
lyr.commitChanges()

# refresh the canvas to show the polygon ?
iface.mapCanvas().refreshAllLayers()

I also expect there is method to interrogate each object to determine what objects have been appended as I go like the line print(ftr.geometry()) which reports the that the object contains a polygon with listed coordinates.
Is there a similar function, which reports what is in the layer?


texts: Splitting one muliple point layer to many layers of one point in ArcMap?

texts:
I need to do a cost distance calculation on multiple points in a layer. If I am not mistaken, in order for the calculation to be done on every point separately, I need to run the calculation separately.

How can I divide the original point layer into many one point layers which will be placed in a geodatabase and used for calculation in a module?


texts: Excel to Table Conversion Tool not working

texts:
When I input the table in the dialog box, I get this rather long error:

ERROR  updateParameters Execution Error: Runtime error Traceback (most recent call last):
File "c:\program files (x86)\arcgis\desktop10.3\ArcToolbox\Toolboxes\Conversion Tools.tbx#ExcelToTable_conversion.UpdateParameters.py",
    line 4, in
File "c:\program files (x86)\arcgis\desktop10.3\ArcToolbox\Toolboxes\Conversion Tools.tbx#ExcelToTable_conversion.UpdateParameters.py",
    line 40, in updateParameters
File "C:\Program Files (x86)\ArcGIS\Desktop10.3\ArcToolBox\Scripts\ExcelToTable.py",
    line 125, in get_sheet_names workbook = xlrd.open_workbook(in_excel)
File "C:\Python27\ArcGIS10.3\lib\site-packages\xlrd\__init__.py",
    line 416, in open_workbook ragged_rows=ragged_rows,
File "C:\Python27\ArcGIS10.3\lib\site-packages\xlrd\xlsx.py",
    line 705, in open_workbook_2007_xml ensure_elementtree_imported(verbosity, logfile)
File "C:\Python27\ArcGIS10.3\lib\site-packages\xlrd\xlsx.py",
    line 46, in ensure_elementtree_imported ET.iterparse(_dummy_stream)
File "C:\Python27\ArcGIS10.3\Lib\xml\etree\ElementTree.py",
    line 1194, in iterparse parser = XMLParser(target=TreeBuilder())
File "C:\Python27\ArcGIS10.3\Lib\xml\etree\ElementTree.py",
    line 1446, in __init__
"No module named expat; use SimpleXMLTreeBuilder instead" ImportError: No module named expat; use SimpleXMLTreeBuilder instead

Doing some research, I found here that I might need to upgrade the python.  I am new to Python, and I installed 2.7.8 and the tool worked previously.  I use Windows 8.1 and it's a 64-bit system.  I tried installing the newer python and I got "missing VCRUNTIME140.dll" and I have not found where to fix that yet.


texts: Create Fishnet issues with measurement distance

texts:
I am trying to create a regular soil sampling grid based upon digitised field extents.

I had posted a previous question today relating to issues with coordinate systems and displaying the sample points accurately in TM75 Irish Grid, I've now resolved this issue.

I want my soil sample points to be on a 35m grid basis.

However, when I now perform create fishnets and specify 35 for the number of rows and columns, when the points are projected into TM75 these points are now 20m apart (when measured in WGS84 these points are 35m apart).  I have been digitising and producing the gridded points (through fishnet and clip) using WGS 84 before projecting these points to TM75 in a new map document. Does the issue with distance spacing lie with WGS 84? What should I be specifying in the number of rows and columns, I can't seem to find any information on WGS 84 and fishnets.


texts: OpenLayers: Smoothly change tile source on zoom

texts:
I'm trying to load tiles from different subdomains depending on my current viewZoom level (i. e. all possible zoom steps). The resolutions for each step are calculated from the user's screen size. They might have identical gridZoom but different URLs.

viewZoom zoom   subdomain
0        1.0    subdomain1
1        1.2    subdomain1
2        1.4    subdomain2
3        1.6    subdomain3

For simplicity, I'm giving the actual zoom values for each viewZoom step instead of the resolutions. Therefore, gridZoom = Math.round(zoom). This is the {z} value that is loaded from my tile source.

I've tried to calculate the URL from which the tiles are loaded depending on the current viewZoom level using tileUrlFunction:

var tilesource = new XYZ({
// ...
    tileUrlFunction: function(tileCoord) {
      var z = tileCoord[0];
      var x = tileCoord[1];
      var y = tileCoord[2];
      var viewZoom = map.getView().getZoom();

      if(viewZoom <= 1) {
       path = 'subdomain1';
      } else if (viewZoom <= 2) {
       path = 'subdomain2';
      } else {
       path = 'subdomain3';
      }

      return 'http://mytileserver.net/' + path + '/' + z + '/' + x + '/' + y + '.png';
  }
}

With this approach, sometimes tiles from previous subdomains will be shown, probably due to caching. I've tried to force a reload whenever the subdomain changes:

map.getView().on('change:resolution', refreshMap);

function refreshMap() {
  var newZoom = map.getView().getZoom();

  var old_subdomain = getSubdomain(currZoom);
  var new_subdomain = getSubdomain(newZoom);

  if(old_subdomain != new_subdomain) {
    tilesource.refresh();
    tilesource.setTileUrlFunction(tilesource.getTileUrlFunction());
  }
  currZoom = newZoom;
}

In this case, no tiles from wrong subdomains will be loaded. However, I'm experiencing several issues with this approach:

If I zoom out, in some cases, there won't be any tiles at all for some areas. I guess, it's the same problem with caching again.
The switch between subdomains is very abrupt. Unlike "normal" zoom operations, old tiles are deleted before new ones are loaded. I would rather like to have a smooth transition with  old tiles to be removed only when the new tiles are ready.

If I'm using refresh on every viewZoom change, the first problem will be gone. However, in that case I'm losing all the beneficial effects of caching and the second problem becomes even more severe, occuring on any zoom operation.

I hope there is a solution for smooth loading of tiles from different sources depending on the viewZoom.


texts: raster2pgsql not in C:\Program Files\PostgreSQL\12\bin

texts:
After adding both

C:\Program Files\PostgreSQL\12\bin
C:\Program Files\PostgreSQL\12\lib

to my windows environmental variables PATH thing, I can't access raster2pgsql:

'raster2pgsql' is not recognized as an internal or external command, operable program or batch file.

I had a look in C:\Program Files\PostgreSQL\12\bin, but I can't actually find raster2pgsql there? Isn't this strange? I would expect it there?

I have postgresql (12.0) and postgis (3.0.0) installed


texts:
I'm working with TextSymbolizer in Mapnik. I have a shapefile with a attribute which is a real. I wish to print the value of this attribute as integer and I don't find how I can do it in Mapnik.

I must not change the source shapefile. Can someone advise of the recommended way to do this?


texts:
I want to get real line segment dataset to do an experiment. So I searched it, and found the site that includes the real data: http://www.wsdot.wa.gov/mapsdata/geodatacatalog/
From this, I downloaded railroad data and opened it with QGIS. I want to get the end points of each line segment (or feature) in the data. When I click the one feature in vector layer, I can see the firstX, firstY, lastX and lastY in 'Identify Results' window. I think this information is what I needed.
However, I cannot find the way to export this information to .csv-file or text file. I searched, but all of these solutions are not working for me:

Getting list of coordinates for points in layer using QGIS?
How to get X and Y coordinate values of the map?

How can I export end points of each line segment in such dataset?


texts:
I have tried to put cluster plugin in geoserver\WEB-INF\lib.whenever i tried to put those plugins,geoserver failed to load, but i remove those plugins from lib geoserver works.
Can anyone please help me..
i have attached the screenshots of plugin which i used.


texts: ArcGIS Server - Feature Access from view

texts:
I want to create a Feature Access Service using an Oracle View as DataSource.

In ArcMAP I can add the layer without problems, but when I try to publish the service, I get a "00134: Layer's data source is not supported" error.

Any ideas why?

I am using ArcGIS Desktop and Server, both version 10.3, and Oracle 12c


texts: Persistent datum conflict error. Datums are the same

texts:
I'm clipping a number of features using other feature layers, and was careful to get everything into WGS84 beforehand, using what I thought were the proper datum transformations.  Now I keep getting Error 001003, Datum conflict between input and output.  I'm not sure if this is referring to my input shapefiles, or if it's a problem with output coordinate environment settings, or none of the above.

Contents of the .prj file of the input feature:

  GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],VERTCS["Unknown
  VCS",VDATUM["Unknown"],PARAMETER["Vertical_Shift",0.0],PARAMETER["Direction",1.0],UNIT["Meter",1.0]]

Contents of the .prj file of the clip feature:

  GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]]

Any ideas?


texts: Default values in shapefile - point sampling tool

texts:
Is there a way in QGIS to get values in my point shapefile from different layers via point sampling tool as "default values" ?

E.g.

When I add a point in the shapefile the altitude (raster) and county (polygons) are auto-completed and I do not have to run the point sampling tool after every session to complete these values.


texts: Is it possible to simplify geometry only for polygons with certain attributes?

texts: I have a map of electoral districts. Urban areas have more densely packed districts whose borders are more important to keep at high resolution. More rural districts have simpler geometries and have border details which are less important. I'd like to reduce the number of vertices where I can but keep a higher resolution in urban areas. Is this possible?


texts: IMERG precipitation from mm/hr to mm/month using NCO

texts:
I have IMERG monthly rainfall data with unit mm/hr, and in nc4 format. Downloaded from https://disc.gsfc.nasa.gov/datasets/GPM_3IMERGM_06/summary?keywords=IMERG
And I want to use the data to calculate SPI and required precipitation in mm/month. After reading various thread, NCO is the best software to do the conversion. The problem, I never used NCO before and not familiar with it.
After read the documentation, I found it could be done using ncap2 command but I am still confuse with it.
Case: for month with 31days, I need to multiply the data with 744 (24h*31day) so I can write the code
ncap2 -s 'precipitation=744*precipitation' in.nc4 out.nc4

And its worked.
The problem, I have merged all the nc4 into single nc4. How to apply above code with that condition?


texts: Small sample size for geographically weighted regression; limitations

texts:
I am studying West Nile Virus risk factors in a specific county using mosquito surveillance data over a ten year study period. The highest resolution of the mosquito data that could be obtained was at the city level. There were 2,040 traps set over the study period, but since the data standardization was at the city level, this resulted in my sample size being N = 57. I have performed OLS regression in R and in ArcMap 10.2. I selected my explanatory variables based on known predictors from the literature and then used the leaps function in R to select the best fit model. Due to the low number of observations, I did not want to run the risk of overfitting the model so I decided to keep the number of predictors under 6. The best subset of variables from the leaps function was y = β0 + β1x1…+ β4x4 + e; where y is the percentage of WNV positive traps in a city, and the explanatory variables are mosquito abundance, population density, percent under federal poverty level, and percent Hispanic. The OLS regression in ArcMap produced a significant Koenker statistic, indicating spatial non-stationarity. I want to run a geographically weighted regression but I suspect my sample size is too small.

I am self-taught ArcMap and this is the first time I’ve done regression using this software. I’ve done some tutorials and it is suggested that a small sample size should NOT be used for GWR; however, I wonder if it is dependent on the application. My analysis is more exploratory than anything, and I would still like to present a GWR regression but also present the limitations associated with the small sample size and encourage mosquito control programs to geocode their trap locations.

What limitations are associated with a small sample size in GWR and is there any way to quantify these limitations using the GWR output?


texts: Calculating Area with GSSURGO vector data

texts:
I'm attempting to calculate acreage of the various "prime farmland" classifications in the GSSURGO dataset for the Omaha Reservation in Nebraska using QGIS but when I sum the total it's about 86 times more acreage than the total area of the reservation. The steps I followed to get to this calculation was:

Load the GSSURGO MUPOLYGON layer into QGIS and reproject to EPSG 5070.
Join with MAPUNIT Table (that contains prime farmland classifications for each soil unit in MUPOLYGON) using shared MUSYM field and export as joined shapefile.
Clip joined GSSURGO data to include data within exterior boundaries of Omaha Reservation.
Add new Decimal String field and Calculate acres for each polygon within the clipped layer using the Calculate Geometry polygon.
Generate Summary Statistics using Vector>Analysis Tools>Basic Statistic for Fields.

The result it gives me is 17,0120,077 acres. However, the area of Omaha Reservation (calculated using boundary file used in clip operation in step 3) is only 198,395 acres. I've double checked I was calculating acres instead of square meters or square feet).

Why are these totals so wildly off?


texts: ArcPy function for changing data source server properties inside mxds?

texts:
Is it possible to change the parameters highlighted in the screenshot below? If so which ArcPy function would do this? I know about and do not want to use the function that references .sde files on C:\
In the example below AUWEISQL1 will change I want to AUWEISQL3.


texts: Mouse wheel events not fired when focus on map element INSIDE winform web browser

texts:
i'm displaying a openlayers map inside a winform web browser,
since i'm getting focused on map element(with a click on the map), i have not getting mouse wheel events any more,even not document mousewheel event, to start getting the events again i need to click on some other element rather than a map element.
when the app is running outside of the web browser(in IE,chrome,Firefox) everything working fine.
Has anyone seen a similar problem?


texts: Instalación masiva de QGIS - Automated installation of QGIS

texts:
Solicito lo siguiente requiere la línea de comandos que permitan realizar la instalación automatizada o desatendida para poder enviarla desde SCCM para la Instalación de QGIS

Google Translate - Spanish to English

I need command line syntax that allows for automated/unattended installation of QGIS via System Center Configuration Manager SCCM.


texts: Invalid field name by using GetFieldAsString from ogr.py

texts:
I use a little algorithm to transform values from an openstreetmap-file into a table (.dbf).

I use the function: GetFieldAsString("string") from ogr.py and everything is fine, but when the "string" contains a char like "_" or ":" the error:

ERROR1: Invalid field name: 'string_string'

appear.

(The function I use comes from there)


texts: QGIS how to make options in valuemap invisible when outside constraint

texts:
I am making a geopackage to plan aids to navigation in a maritime environment.

I have made the design of my attributes form with the drag/drop designer.
The constraints I put in work i.e. If i choose a lateral portside-hand mark i can not press o.k. unless the color red is chosen, although all the other colors are still visible in the valuemap i've set.

Is it possible to blank out the values in the valuemap when they fall outside of the limits of the constraints i've set?
I would like to make the form idiot proof so the end-user can not choose an option that is outside the rules en regulations in the field.


texts: QGIS 2.18 crashes every time I zoom to a layer

texts:
When I zoom to any layer, QGIS 2.18 crashes without an error code. I also tried completely changing the CRS, removing all map services. What else can I try?

The CRS was WGS 84 / Pseudo Mercator. I also tried changing it all to NAD83HARN Mercator and the same problem persists. The "zoom to layer" seems to be the only trigger to this crash. The zoom starts to happen, then the moving blue "thinking" indicator at the bottom of the screen moves, then it stops, and then it crashes. The only layers are shapefiles.


texts: Join or intersect fields without creating new fields

texts:
GIS is not my full time gig and I’m a bit stuck.
I have 2 layers, both with the same fields with the exception that my primary layer has extra fields and is based on a different set of polygons.
How do I intersect or join (specific fields from) my primary layer and my secondary layer without having new fields (tacked on the end) and a new output layer?
I have tried intersect and join/spatially join and they both work in that the data from from secondary layer is joined or intersected but in a new output.
What I want to do is be able to mash the identical fields from both layers and suck the attribute data from my secondary layer into my primary layer (as it does in intersect) without creating a new output and duplicating fields.
I’d also like to be able to do this after I add a new polygon to my primary layer.
I use both ArcMap and ArcGIS Pro and know nothing about scripts if that helps.


texts: Ordinary Least Squares Point Aggregation

texts:
I want to run an ordinary least squares regression in ArcGIS Pro on my dataset which is composed of about 400 fire occurrences (point data) and demographic variables (100x100m polygons of income, education, etc). I plan to use the fire occurrences for my dependent variable and the demographics for my explanatory variables.

The problem I have is with aggregating my fire incident point data. Essentially, I don't want to aggregate my incidents into larger polygons (larger fishnet, census districts, postal districts, etc) and lose the accuracy of my 100x100m demographic variables (which when spatially joined with the point layer only yield maximum 3 join count).

I have tried the Integrate/Collect Events method of aggregation suggested in the help documentation (and ESRI seminar videos), however because I do not have so many coincident points I am finding myself with a low ICOUNT as a result (again maximum of 3 events aggregated). I know that the dependent variable input for OLS requires more variation in the weighted value.

Does anyone know of a method that I can both aggregate my points properly and preserve the high detail of my 100x100m demographics for my regression?


texts: Preventing QGIS from losing information in attribute table during KML to SHP conversion?

texts: I have a Shapefile with certain information in the attribute table. I convert it with QGIS 2.14 to KML. I convert the KML back to SHP. When I open the attribute table, all the infromation is lost. How can I prevent QGIS doing this?


texts:
What is the best way to use routing table generated by osm2po? When I apply osm2po on my *.osm.pbf file, the generated table looks like this:

Here is my routing script:

SELECT d.seq, d.path_seq, d.edge, d.cost, h.osm_name, h.kmh, h.x1, h.y1, h.x2, h.y2, h.geom_way  FROM pgr_dijkstra('
    SELECT id,
         source,
         target,
     cost,
         reverse_cost
        FROM hh_2po_4pgr',
    (select source from hh_2po_4pgr h order by st_makepoint(h.x1, h.y1) <-> (select st_makepoint(lon, lat) from fnc_geocode('my source address')) limit 1),
    (select target from hh_2po_4pgr h order by st_makepoint(h.x2, h.y2) <-> (select st_makepoint(lon, lat) from fnc_geocode('my target address')) limit 1),
    directed := false) d
left join hh_2po_4pgr h on d.edge = h.id

When I display route in Geoserver target is not correct, as you can see on the picture:

Can you tell me what I am doing wrong? I need my route to end on the blue point. I am using osm2po 5.1.0. Function fnc_geocode('xxx'). It has a function that returns coordinates of the entered address. With osm2pgrouting it works fine. Now we have more data, so this is the reason why I decided to use the osm2po tool instead.


texts: Extracting the pixelwise proportion of classes of a finer raster in a coarser one

texts:
I have two rasters, one has a 1 km resolution. The other has a spatial resolution of 30 m. This latter pixelwise is a classification. I want to overlay the 1 km raster on the classification raster and extract for each 1 km pixel the distribution of classes corresponding to the 30 m classification.

The 30 m classification includes 7 classes so I would want this process to produce 7 new layers at 1 km. One for each class. Each pixel of 1 km will now have associated a corresponding % of the presence of each class present in the 30 m classification. If one made a table where each row is a pixel then one would have a table like the following:

One additional difficulty is that the rasters are quite big. I would really like to achieve this using python.

EDIT:

I have added subsets of my two datasets a this link (at 30m and 1000m resolutions respectively), in case anyone wants to propose a proof of concept; would be amazing.


texts:
I have a problem with the "export map options" of ArcMap 10.2
There are only two parameters in the option window that are: general and format.
I set high dpi at 300 but the map exported has not a good resolution as I parametrized and as I wish.

Why aren't there the other parameters like pages, security, advanced?

I tried also using ArcPy but I didn't find the correct Python code.


texts:
I want to convert TIGER data

GEOGCS["GCS_North_American_1983",
    DATUM["North_American_Datum_1983",
        SPHEROID["GRS_1980",6378137,298.257222101]],
    PRIMEM["Greenwich",0],
    UNIT["Degree",0.017453292519943295]]

to WGS84. How does ogr2ogr determine which transformation to use and how big would the introduced error be?

Study area is the county of Kauai in Hawaii, USA.


texts: Feature Class to Feature Class error: 999999 (ArcGIS Pro 2.5)

texts:
I've tried this with the new ArcGIS Pro 2.5 Python Notebook and the Python compiler from within Pro.

I have about 15 layers that I'm reading in from a couple different external servers. I can read in test layers from all unique servers just fine, but there are three layers that after I read them into my table of contents, they fail when using the Feature Class to Feature Class GP tool.

Error:

ExecuteError: ERROR 999999: Something unexpected caused the tool to fail. Contact Esri Technical Support (http://esriurl.com/support) to Report a Bug, and refer to the error help for potential solutions or workarounds.
ERROR: code:500, Error performing query operation, Internal server error.
Failed to execute (FeatureClassToFeatureClass).

What I have in Notebook and what I've tried in the compiler:

aprx = arcpy.mp.ArcGISProject(r"C:\Project\Pro_Points2018.aprx")
mapx = aprx.listMaps("Map")[0]

env.overwriteOutput = True

arcpy.env.workspace = r"C:\Project\Pro_Points.gdb"
out_loc = r"C:\Project\Pro_Points.gdb"

mapx.addDataFromPath(r"https://gisportal.hrsa.gov/server/rest/services/Shortage/MedicallyUnderservedAreas_FS/FeatureServer/0")

arcpy.conversion.FeatureClassToFeatureClass("Medically Underserved Areas /Populations – perimeter polygons", out_loc, "MUA_SHP_DET_CUR_VX", '')

Other posts point to the error being caused by import/export name being the same, but I've already verified that.

All I'm trying to do is read a FC in from someone's server (that I have access to), and save it off locally in my GDB.


texts: Use existing property pages with new IFeatureRenderer

texts:
I am trying to override the drawing of a feature within a custom plugin data source.  I believe this entails making a custom IFeatureRenderer.

However, the "renderer" I'd like to design would be 99% the same as a ClassBreaksRenderer, but with a custom Draw override.

What would be the simplest manner to be able to reuse the existing functionality of the ClassBreaksRenderer (including and especially the UI property pages) while still being able to provide a custom Draw implementation?


texts: QGIS Merge rasters without projection edges

texts:
I have been trying with no luck to merge forest canopy rasters available through the NLCD (e.g. Grids 102 and 105) using QGIS in Mac (as well as R, and SAGA in Windows).

The GeoTIFF grids import and display fine.  However, when I try to merge the two, triangular black spaces appear around the edges of each, and where they overlap the zero/null values of one are superimposed on the actual values of the other.  It looks like this:

I suspect this is a projection problem, as the grids do not align with their mutual boundary.  The files import with a "user defined" projection that is similar to NAD83/Conus Albers (EPSG:5070).

I have tried straight merging, and re-projecting before merging to various projections, but nothing seems to work.  Strangely enough, however, the "land cover" (as opposed to "forest canopy") GeoTIFF rasters provided by the NLCD corresponding to the same areas - and with the same projection - merge without issue.

How can I perform a proper merge with these rasters?


texts: How do I use a filter on an uncached WFS layer in QGIS?

texts:
I am trying to create a filter to show only specific features based on an attribute query within a WFS using QGIS. This works fine when using the standard WFS loader and creating a filter (based on the type "attribute"='value'). However this only works when the features are cached.

I want to be able to use a filter like this for an uncached WFS layer. Does anyone have any idea how to do this?

Example using a BGS public WFS.

A standard cached query would look like this:

http://ogc.bgs.ac.uk/cgi-bin/BGS_Bedrock_and_Superficial_Geology/wms?SERVICE=WFS&VERSION=1.0.0&REQUEST=GetFeature&TYPENAME=GBR_BGS_625k_SLS&SRSNAME=CRS:84&FILTER=<ogc:Filter xmlns:ogc="http://www.opengis.net/ogc">
 <ogc:PropertyIsEqualTo>
  <ogc:PropertyName>LEX_ROCK_D</ogc:PropertyName>
  <ogc:Literal>PEAT</ogc:Literal>
 </ogc:PropertyIsEqualTo>
</ogc:Filter>

When I untick the cached option, the layer source is as follows:

http://ogc.bgs.ac.uk/cgi-bin/BGS_Bedrock_and_Superficial_Geology/wms?SERVICE=WFS&VERSION=1.0.0&REQUEST=GetFeature&TYPENAME=GBR_BGS_625k_SLS&SRSNAME=CRS:84&BBOX=-3.4085901736001083,57.05755199511002473,-2.94851240830698069,57.42328922738386154&FILTER=<ogc:Filter xmlns:ogc="http://www.opengis.net/ogc">
 <ogc:PropertyIsEqualTo>
  <ogc:PropertyName>LEX_ROCK_D</ogc:PropertyName>
  <ogc:Literal>PEAT</ogc:Literal>
 </ogc:PropertyIsEqualTo>
</ogc:Filter>

The main difference is the bounding box. What happens with the second layer definition is that the layer loads but has no features within it.

I am able to create a GeoServer layer within my WFS that meets this query, but I am trying to find a way that users will be able to replicate without intervention at the GeoServer end.


texts: How can I convert geoJSON data to GPX format

texts: I want to know how I can convert geoJSON data to GPX format.


texts:
I'm trying to hide some columns in my vector data attribute table: Layer properties - Fields - Widget Properties - Hidden. This should result in the fields being invisible in the attribute table.
However, it doesn't work in my case - the hidden fields still show up (with their values visible), even after closing and re-opening the attribute table. Has anyone had this problem and knows how to solve it?


texts: Joining attributes from start and end point to line feature in QGIS?

texts:
I have a line shapefile and a point shapefile that are start and end point of these line features. Now I want to join ALL attributes from both start and end points to the line shapefile.

How do I do this using QGIS?


texts: I am trying, in ArcGIS Online specifically, to reference a map service layer by its name as opposed to by its reference ID.  Does anyone know how to do this in the URL?


texts: Why are some OSM GPX line string endpoints connected?

texts:
I am working with GPX files from OpenStreetMap to render the public transportation lines of Berlin. I converted them to GeoJSON using the togeojson script. For some reason some of the lines are rendered as a "circle". The first and the last point are connected for some reason. The effect can be seen at the GitHub project I linked.

What am I doing wrong?


texts:
I have PostGIS table with 2 geometries ( "GEOMETRY(MULTILINESTRINGZ)" and "GEOMETRY(LINESTRINGZ)" ). There are lines made just from 2 points - XY are WGS and both of them have random elevation Z (I hope default units are meters).

Example from ST_AsTExt(geom):
1) LINESTRINGZ (2.17766 41.3885803 0,2.18333 41.3833 754)
2) MULTILINESTRING Z ((2.17766 41.3885803 0,2.18333 41.3833 754))

I tried to visualize them in ArcScene, ArcGlobe, qgis2three and I can't see them as 3D nowhere. They act as 2D...

I also tried export to SHP with ogr:

ogr2ogr -f "ESRI Shapefile" "~/myshape" PG:"user={user} password={pass} dbname={dbname} host={host}" {table} -geomfield GeomColumn

... and checked with "ogrinfo file.SHP" which just says it's Linestring for both files (nothing about multilinestringz and Z coordinate).

None of two SHP's were displayed as 3D.

So, is there any error in my geometry types?

Should I specify CRS in data type?

Should I add Z coordinate CRS/unit somehow?

Are any other options required to set them as 3D in Arc* products?

Is there any other easy option to export and check 3D of at least few lines somehow (I am now going to explore possibilities of 3D in ST_AsGeoJson, later maybe Grass GIS)?


texts: How to create a 9.3 file geodatabase with ArcObjects SDK for ArcGIS 10?

texts:
In ArcGIS 10 I could use the Create File GDB geoprocessing tool to create a 9.3 file geodatabase or the Create Personal GDB geoprocessing tool to create a 9.3 personal geodatabase.
But I cant't find a similar method in the ArcObjects SDK for ArcGIS 10. The IWorkspaceFactory2.Create method seems only to create geodatabases compatible with ArcGIS 10.

Possible workarounds are:

using the geoprocessing interface
centrally storing one empty 9.3 geodatabase and using that as
template geodatabase

But is there really no possibility to create a 9.3 Geodatabase with the ArcObjects SDK for ArcGIS 10?


texts: How to find cities without nearby bigger city in SpatiaLite?

texts:
I have a table cities, each with a population INTEGER and a geom POINT location (WGS84). I'm trying to select all cities where there is not a more populous city within 50 km using the latest SpatiaLite (4.3.0).
This has to use a spatial index because there are almost 200000 cities, so doing a full table scan for each of them would be too slow. I understand that (unlike PostGIS) SpatiaLite won't automatically use my spatial index, but I found this question that I thought should help.
This is the setup (slightly simplified):
sqlite> .load mod_spatialite.so
sqlite> SELECT InitSpatialMetaData();
1
sqlite> CREATE TABLE cities (
          id INTEGER NOT NULL PRIMARY KEY,
          name STRING NOT NULL,
          population INTEGER NOT NULL
        );
sqlite> SELECT AddGeometryColumn('cities', 'geom', 4326, 'POINT', 'XY', 1);
1
sqlite> INSERT INTO cities VALUES (1, 'One', 500, MakePoint(0.0, 0.0, 4326));
sqlite> INSERT INTO cities VALUES (2, 'Two', 300, MakePoint(1.0, 0.0, 4326));

Here's what I thought should work:
sqlite> SELECT * FROM cities AS a
        WHERE NOT EXISTS (
          SELECT i.ROWID
          FROM SpatialIndex AS i, cities AS b
          WHERE
            i.ROWID = b.ROWID AND
            f_table_name = 'cities' AND
            search_frame = ST_Expand(a.geom, 50000.0) AND
            b.ROWID <> a.ROWID AND
            ST_Distance(a.geom, b.geom, FALSE) < 50000.0 AND
            b.population > a.population
        );
1|One|500|
2|Two|300|

But no matter how large I set the distance, this query always returns both cities. Maybe this is because the spatial index seems to be empty:
sqlite> SELECT * FROM SpatialIndex;
sqlite>

So then obviously, no matter what I put in the innermost WHERE clause, it's going to come back empty.
What am I doing wrong?
Also, is there any up to date documentation? The SpatiaLite cookbook seems to be outdated, because it doesn't use this SpatialIndex virtual table.