texts: Parts of the lines aren't split when using Split Line At Point

texts:
The line shapefile is obtained from the Dissolve toolbox of all routes, so there is only one polyline in this shapefile. The point shapefile is already there. The Split Line At Point should allow me to get another shapefile which will be split according to the points. I notice that some of the lines are split, but some of them are not. The common feature for those lines don't get split is there is no turning line at that point, as shown below. This means the line is split perfectly when the point is the intersection, but for those points are just in the middle of a line, the split function will not work. Also, the points used to split the line are obtained from this line shapefile, so they are perfectly located along the line. This is a new question, which is not the same as the similar questions asked before since the problematic points are all the middle points along the line without intersecting with other lines. Any comments on this issue?

The following picture shows the problematic points where the line can't be split.

I am using ArcGIS Desktop Advanced License 10.3.1.


texts:
I have a line layer in which the lines have several gaps.

How can I quickly close the gaps in QGIS 3.4?

Is there a tool or plugin available I haven't found yet?


texts: SRTM data represented in 3D using Unity. Problem texturing with sat images

texts:
I'm currently working on a project that represents geographic data in 3D using Unity. It basically loads, when given a GPS position, the area around it. For example, given the coordinate 50.5N 6.5E it will load the terrain from 50N to 51N and 6E to 7E. Once the terrain is loaded a texture is applied to the 3D mesh, this texture can be a Landsat image or a Openstreetmap map for example.

I downloaded elevation data (*.hgt files) from NASA’s Shuttle Radar Topography Mission (SRTM). These files have a resolution of 1 arc-second (30m) and contain a 3601 x 3601 matrix of points that represent the height for that coordinate, being the first value of the file the bottom-left coordinate (50N, 6E in my previous example).

To calculate the distance along the longitude between the height points I represent, I apply Haversine's formula (in the file distances between points are not specified). With it I obtain the distance between longitudes for the given latitude, which I can then divide by 3600, obtaining the distance between points.

For Haversine's formula i don't use a fixed value for earth's radius. Depending on the latitude I adjust it between the value at equator and the value at the poles.

For the distance along the latitude its easier since I only have to divide the constant distance (111.2 Km) by the number of points (3600).

Up to here I think that everything is fine and that I have done everything properly. The problem comes now, when applying a texture to the 3D mesh.

To obtain an aerial image of the area I am representing I decided to use Universal Maps Downloader. It allows you to introduce a set of latitudes and longitudes, choose a source (Google Street Maps, ArcGis Online Topography, Open Street Maps...), a level of zoom, and downloads the images of the area as a set of tiles that can be merged into a single image file using the included Map Combiner.

I import this image into Unity and apply it as a texture to the mesh (I don't think that talking about how I adjust the UVs is relevant), and it works, but not perfectly.

The problem is that the texture (call it aerial image or map, as you prefer) doesn't match the topographic data well. Depending on the location I load the image has an offset to the 3D mesh, the direction in which the image is displaced follows no pattern apparently (sometimes a bit to the left sometimes to the right, sometimes a lot, sometimes not much, sometimes up, sometimes down....).

I can adjust it manually adding an offset to the UV mapping, but that is only possible if I have visual references like rivers or coastlines. And since the offset doesn't follow a pattern I can't program my own offset adjustment.

The question is basically: does anyone know where is my mistake or how could I improve my results?

Maybe I'm downloading aerial imagery from the wrong place (and it is not precise) and I should use a different source. So far this option is the only one that lets me obtain images of quality for a given area defined by a set of coordinates in a quite simple and fast way. I couldn't find anything similar.

Maybe I'm doing something wrong when representing the SRTM files.

Perhaps I took a wrong approach from the beginning...

28N, 17W Tenerife and 50N, 6E NRW-Cologne

0N, 6E Sao Tomé and Principe


texts:
I'm new in QGIS and I'm trying to do very simple thing:

Filter OSM data: natural=water
Add Vector Layer, select multipolygons

Here is what I got with QGIS:

Some islands are filled with blue (why?), but others are not (as expected).
I've checked with JOSM that all islands and the lake mapped correctly (inner/outer roles), and the only difference is that the 1st ones are "closed ways" and the 2nd ones are "relations" with type=multipolygon. But, as far as I understand, there is nothing wrong with it.

So, what I'm doing wrong and how can I get the proper result (all islands are not filled) with QGIS?


texts: Sort User Defined Features Classes Into Premade Feature Dataset Using ModelBuilder?

texts:
I have a large Model where the user selects the text file, it creates the xy file, then decimates the file, thus creating a feature class.

I would like this decimated feature class to get sorted into my already predefined datasets.

I have the user pick the type animal, in my example if someone picks cat from the drop down I want it to go to the "Cat_101" feature dataset based on the fact they picked keyword "Cat".

I am pretty new to python and model building.

I do understand that this will likely take some python if / then / else statements, but I honestly have no idea where to start with this.

My goal is make this happen behind the scenes.  I do understand the user can always select where to save the file, but I want to remove that from them and save it based off their choice of "Cat", "Dog", "Bird".


texts: When using arcpy can I only calculate fields of a subtype and not the entire feature class?

texts: I'm appending two sepearte shapefiles into one feature class with different subtypes.  Is there a way to only calculate fields from a certain subtype?


texts: Automatically create polygons from top right and bottom left coordinates?

texts:
I have a CSV file with hundreds of records. Each record contains 2 matching coordinates, 1) top right corner of a rectangle, 2) bottom left corner of a rectangle.

Polygon ID, Lat 1,       Lat 2,       Long 1,       Long 2
1,          27o 28' 28", 27o 33' 30", 86o 06' 04",  86o 15' 00"

Is there a method which I can use to automate the creation of these polygons? I anticipate there being a few steps.


texts: ogr2ogr srs options: where is the well known definition determined?

texts:
In the ogr2ogr documentation it says that when using -a_srs srs_def:

  Srs_def can be a full WKT definition (hard to escape properly), or a well known definition (ie. EPSG:4326) or a file with a WKT definition.

Where are the default well-known definitions listed?
Do they include ESRI as well as EPSG?
If not, why not, and can I add all the ESRI definitions listed at spatialreference.org?
Is there an easy way to load an srs definition using an SRID from a PostGIS spatial_ref_sys table, aside from just running some sql to get the text?

Thanks!


texts: Spatialite as a DBMS-server with GeoDjango?

texts:
I would like to add a spatial component to my Django application. I turned out to GeoDjango immediatly. I am wondering now if spatialite could be used as DBMS-server for an application in Prod. Grosso modo, the app will allow to store shapefiles, kml and others spatial data on server database an allow spatial queries on these files. Can spatialite be used for this?

Can spatialite achieve multi-task simultaneously, Many queries in the same time? I have read this http://www.gaia-gis.it/spatialite-2.3.1/spatialite-arch-2.3.1.html and they say spatialite is for personal use only. Is it a limitation for a multi-users inline GeoDjango application?

I would prefer spatialite better than PostGIS for its simplicity and lightweight.

Do you have advices about feasability and performances issues, or else?


texts: Problem on exporting .xyz from raster layer

texts:
I am new in QGIS, and I  like to export xyz file from Raster(.tif). I saw last posts related to this object, but my xyz file only have the latitude and longitude.

How can I add the Elevation to my xyz file?

Perhaps I need to make a DEM and then join the elevation tab to my current xyz file.


texts: ModelBuilder - Automatize names output files

texts:
I'm making a ModelBuilder to facilitate some geographical treatments which concern many CSV files in a same workspace.

For that, I made an iteration on that workspace, it works but I have a problem with the output files' names.

The csv files which are input looks like that : "456789_V1_58412_FTMP_68513_PMSKD.CSV" or "123456_V2_58412_FTMP_68513_PMSKD.CSV" .

For the output names, I need the 9th first characters of the originals CSV files with "Points_" at the beginning. To resume, I need the output names looks like that : "\Défault.gdb\Points_456789_V1" or "\Défault.gdb\Points_123456_V2".

I have made some tests by writing "Points_%Name%" in the output names in my Model builder but the results are not good. I have that :"\Défault.gdb\Points_456789_V1_58412_FTMP_68513_PMSKD" .

I have 2 problems with this :
 1 - It doesn't match with what I want to have as a result name
 2 - The name is too long for a geodatabase so my names results are heterogeneous

Please, do you have some ideas to help me?


texts: I want to transform Gauss-Kruger coordinates into WGS84 coordinates by applying the Bursa-Wolf transformation. I know the Bursa-Wolf parameters for the desired transformation. What confuses me is that the Bursa-Wolf formula expects 3D coordinates (X, Y, Z), whereas my input Gauss-Kruger coordinates are only 2D (northing, easting), and my output WGS84 coordinates are also only 2D (latitude, longitude). How do I map my 2D input/output coordinates to the 3D coordinates in the Bursa-Wolf formula?


texts: Clean local OSM website database but only the geographical data

texts:
I also posted this question at help.openstreetmap 2 weeks ago but got no answer, so I hope someone here have the answer.

I ran ./osmosis --truncate-apidb in my local OSM website to clean the data, so I could import new geofabrik data, and it worked perfectly.

But I noticed that, beyond the geospatial data, I wiped also the users accounts, iD configuration, and I think even the pre-rendered tiles were gone, because I feel the rendering is slower now.

So, there is a way to clean only the spatial data?


texts: Can I link to outside website within ArcGIS Online dashboard?

texts:
I have three mini-dashboards embedded within a single dashboard, and I would like a 4th (default landing) tab (circled in red in pic) to embed this external URL.
I dont want the user to have to click on the embedded URL; I want them to see the website the URL links to by default (without having to click on a link).
I've tried many ways to add as embedded content but it shows the sad face and says that URL refused to connect.  How to fix this?


texts: Projecting Shapefile from WGS to Web Mercator changes output Field Type

texts:
Using ArcGIS 10.6 I am trying to project a shapefile data from Here street data-set in GCS WGS 1984 to WGS 1984 Web Mercator Auxiliary Sphere. While the original data has an attribute like

it changes the POINT_NUM and Z_LEVEL data type from Short Integer to Long Integer

Can you please let me know why this is happening and how I can prevent it?


texts:
I've recently tried using Shapely to dissolve and intersect large sets of 100'000+ geometries, some of them quite complex. This would take me at least one hour, at which point I would usually shut the process down.

I've been recommended PostGIS a few times as a super fast alternative, and I am wondering if it would indeed be faster? From my understanding both PostGIS and Shapely use GEOS, so would there really be such a speed difference between the two?


texts: Inserting %20 instead of a space as a string

texts: I have 600 records that I need to be linked on a website.  As clickable urls. But when I have the URL's, they only link to the 1st set of text before a space.  So is there a code in QGIS that I can put in to alleviate this.  I have tried:   replace( "url"," ") with "%20" But this requires some more work.


texts: How to add a vector layers to my map using OpenLayers 3?

texts:
I'm trying to add a vector layer manually to my map but i receive this error: Uncaught TypeError: b.attachEvent is not a function
the code is :

// adding data
// Define Geometries

var point = new ol.geom.Point(
    ol.proj.transform([10.8,36.7],'EPSG:4326','EPSG:3857')
);
var circle = new ol.geom.Circle(
    ol.proj.transform([10.8,36.7],'EPSG:4326','EPSG:3857'),
    600000
);
// features
var pointFeature = new ol.Feature(point);
var circleFeature = new ol.Feature(circle);
// source
var vectorSource = new ol.source.Vector({
    projection: 'EPSG: 4326'
});
vectorSource.addFeature([pointFeature,circleFeature]);
//vector Layer
var vectorLayer = new ol.layer.Vector({
    source: vectorSource
});
// add layer to the map
map.addLayer(vectorLayer);

I can't resolve it.


texts: my module dialogs are not appearing in GRASS 78 (for 76 it is working) e.g. I want to change the colour of a vectorlayer and click on the properties button, but no module dialog will open.


texts: Creating point feature at intersection of two lines using ArcGIS Desktop?

texts:
How do I create a point feature at an intersection of highways and trails?

I am using ArcGIS Desktop 10.5.1.


texts:
I'm trying to build an interactive map, without any knowledge in javascript, and i found myself a bit lost in this large world.
I've already designed my style in mapbox studio. And now, trying to make it work on the web.
I'd like to display/hide multiple layers, but with one button and the default visibility set as none (tried to set "visibility" as none, but doesn't seem to work). (user must click on 1 button to load all layers in a group)

Followed this tutorial : https://www.mapbox.com/mapbox-gl-js/example/toggle-layers/
I could accomplish it, but i got 1 button for each layer.
So i tried to group them, following this tutorial: http://leafletjs.com/examples/layers-control.html

But things aren't that easy, I'm trying to mix things every way I can, but can't figure how to make it work...

Here is my basic code:

mapboxgl.accessToken = 'pk.eyJ1IjoiY2RyaWNlIiwiYSI6ImNpcDZvMTU2cDAwMGR2MW03aXZ1Y2xpdXAifQ.DIB1RcA8uDyXSVvyz7Skcg';

// Set bounds
var bounds = [
    [5.336040, 47.678506], // Southwest coordinates
    [7.202343, 48.697659]  // Northeast coordinates
];

var map = new mapboxgl.Map({
    container: 'map',
    style: 'mapbox://styles/cdrice/cip6vwuus000wdkm3q2adkr24',
    center: [6.95, 48.28],
    zoom: 8.5,
    maxBounds: bounds // Sets bounds as max
});

map.on('load', function () {

map.addLayer({
        'id': 'reciprocite',
        'type': 'fill',
        'source': 'reciprocite',
        'layout': {
            'visibility': 'none'
        },
        'source-layer': 'reciprocite-gap'
    });

map.addLayer({
        'id': 'reciprocite',
        'type': 'fill',
        'source': 'reciprocite',
        'layout': {
            'visibility': 'none'
        },
        'source-layer': 'reciprocite-gpv'
    });
 });
addLayer('Reciprocite', 'reciprocite-gpv');
addLayer('Reciprocite', 'reciprocite-gap');
function addLayer(name, id) {
    var link = document.createElement('a');
    link.href = '#';
    link.className = 'active';
    link.textContent = name;

    link.onclick = function (e) {
        e.preventDefault();
        e.stopPropagation();

        var visibility = map.getLayoutProperty(id, 'visibility');

        if (visibility === 'visible') {
            map.setLayoutProperty(id, 'visibility', 'none');
            this.className = '';
        } else {
            this.className = 'active';
            map.setLayoutProperty(id, 'visibility', 'visible');
        }
    };

    var layers = document.getElementById('menu');
    layers.appendChild(link);
}


texts: Check if polygon intersects with coordinates in a .csv file with PyQGIS standalone script

texts:
I have a bunch of .csv files with coordinates and a polygon. In my standalone script I would like to check if the polygon intersects with the coordinates.
So, when I write:
            for feature in polygon.getFeatures():
                if feature.geometry().intersects(file.geometry()):

It doesn't work because the coordinates in my 'file' are of type string and not a QGIS Rectangle.
I can use numpy to sort of get the extent of the coordinates as string by writing:
        data = np.genfromtxt(file, delimiter=" ")
        x = data[:, 0]
        y = data[:, 1]

        xmin = min(x)
        xmax = max(x)
        ymin = min(y)
        ymax = max(y)

I feel like I have all the extents but now I need to find a way to check if they intersect. I'm stuck at this point. How can I do this?
I'm using QGIS 3.10.


texts: Loading one raster into QGIS from a table of rasters in PostGIS

texts:
I have a question about loading raster layers into QGIS from PostGIS, it's a an extension from the question asked [here] : Loading PostGIS raster layer into QGIS
'DB Manager', 'add to canvas' works fine, when my PostGIS table contains only 1 raster.

However, I have set up a table where each row has it's own raster, and I am trying to load just the raster that results from a query. So for example, I run the query:

SELECT * FROM raster_table WHERE rid = 8;

This returns the row I am after - which contains a value in the rid column and a value in the raster column relating to the raster, however when I try to load the raster by checking "Load as new layer" and completing the dialog box, I can only get a table to load containing just the rid. I can't get a raster to load from the query.

Could anyone explain if it is possible to do it this way, or does a raster need to be saved in PostGIS individually?


texts:
I'm using the opengeo stack, and want to double check my understanding of each element.

PostGIS - The database for storing my data
GeoServer - The link between the online part / the database

Here's where I am not fully sure:

Geoexplorer, viewing, editing data visually in the POST database
Client SDK - The same as geoexploer but I have to add the functionality?

If I want to add Graphs to the map, i.e. the user clicks on a point, how would I go about that? using SDK and the GeoEXT?


texts: Parse attribute table value to export filename using data driven pages and arcpy.mapping

texts:
I'm trying to customize a data driven page looping statement workflow and want to set a better filename to the final export.  Currently the file name uses the page index (pgIndex), but I would like it to use a specific field in the attribute table.

Here is my file naming looping statement:

for pgIndex in range(1, tempDDP.pageCount + 1, 1):

  # Create a name for pdf file
  temp_filename = r"C:\Desktop\Test\temp_pdfs\MB_" + \
                            str(pgIndex) + ".pdf"
  if os.path.exists(temp_filename):
    os.remove(temp_filename)

Any idea how I can get the result I want?  I know this is possible if you use the built in export tool for arc 10, but I need more functionality, and therefore a customized script.


texts: Getting “Segmentation fault: 11” with Shapely in Python on OSX 10.12.6

texts:
I have a number of data processing scripts that worked a few months ago but now getting seg faults.  In one script I had to pull out code from a method and embed in the calling method.  In this case, I am performing the simplest of tasks.  I am trying to use OGR and Shapely to find the nearest street to a point.

OSX 10.12.6
Python 3.5.2 :: Anaconda custom (x86_64)
Shapely 1.6.1

from osgeo import ogr, osr
from shapely.geometry import Point
...
def FilterNearbyStreets (self, logLevel, pntSource):
    shpPoint = Point(pntSource.GetX(), pntSource.GetY())
    enoughSegments = False
    bufferSize = 0.001
    while (enoughSegments != True):
        print ("About to buffer {}".format(bufferSize))
        buffer = shpPoint.buffer(bufferSize)
        print ("Buffered!!!")

We never get to Buffered!!!

Segmentation fault: 11


texts: Installing Python for ArcPy with ArcGIS Desktop and ArcGIS Pro?

texts:
I am trying to learn Python here. I have ArcGIS Desktop 10.4 and ArcGIS Pro 1.3 installed in my machine. I understand that products install Python 2.7 and 3.2. Do I still have to download Python in order to run Python? I tried to use the command prompt and type 'python' and it is not recognized.

I have tried to install Python in my old machine, and I ended up messing up a lot of things which I did not know why. In the end, I was not able to run Python. I am trying to learn programming here.


texts: Creating and manipulating an HTML frame in QGIS 3.2.0 Print Composer using PyQGIS

texts:
I'm trying to create, reposition and resize an HTML frame in a QGIS 3.2.0 Print Composer layout using PyQGIS:

proj = QgsProject.instance()
comp = QgsPrintLayout(proj)
blockStats = QgsLayoutItemHtml.create(comp)
blockStats.setFixedSize(QgsLayoutSize(34.3,25, QgsUnitTypes.LayoutMillimeters))
blockStats.attemptMove(QgsLayoutPoint(200,89, QgsUnitTypes.LayoutMillimeters))
blockStats.setFrameEnabled(True)
blockStats.setFrameStrokeWidth(QgsLayoutMeasurement(0.5,QgsUnitTypes.LayoutMillimeters))
comp.addLayoutItem(blockStats)

Lines 4 and 5 above give me "object has no attribute" errors for 'setFixedSize' and 'attemptMove' - how can I reposition and resize this object?

Using the answer by @ndawson I've revised my code as follows:

proj = QgsProject.instance()
comp = QgsPrintLayout(proj)
blockStatsHTML = QgsLayoutItemHtml.create(comp)
blockStatsHTML.setHtml("test
test")
blockStatsHTML.loadHtml()
blockStatsFrame = QgsLayoutFrame(comp, blockStatsHTML)
blockStatsHTML.addFrame(blockStatsFrame)
blockStatsFrame.setFixedSize(QgsLayoutSize(34.3,25, QgsUnitTypes.LayoutMillimeters))
blockStatsFrame.attemptMove(QgsLayoutPoint(100,100, QgsUnitTypes.LayoutMillimeters))
blockStatsFrame.setFrameEnabled(True)
blockStatsFrame.setFrameStrokeWidth(QgsLayoutMeasurement(0.5, QgsUnitTypes.LayoutMillimeters))
comp.addLayoutItem(blockStatsFrame)
comp.moveItemToTop(blockStatsFrame)

I'm no longer getting Python errors, however the frame still isn't appearing in the print layout. I've tried changing the position, using positive and negative coordinates, but I can't make it show up.

Is there another step required to make it visible?

Just an update and clarification - I've added the following line at the end of the code:

lmgr.addLayout(comp)

When I open the layout in the print composer, the Items Panel contains only "Map 1" - there is no listing for the blockStatsFrame or blockStatsHTML item. So this tells me any other attempts at removing or resizing the item would be futile, since it doesn't exist in the layout as far as the print composer is concerned.

Apart from comp.addLayoutItem(blockStatsFrame), should anything else be necessary in order to add the item to the layout?

In case anyone is facing the same issue, the solution for me was in using attemptSetSceneRect to set the dimensions of the frame. Without this, the frame didn't appear in the layout at all. Note however that after using this method, setFixedSize doesn't appear to work. Relevant parts of my working code are below:

proj = QgsProject.instance()
lmgr = proj.layoutManager()
comp = QgsPrintLayout(proj)

# create block stats HTML multiframe
blockStatsHTML = QgsLayoutItemHtml.create(comp)

# create frame to show content from blockStatsHTML
blockStatsFrame = QgsLayoutFrame(comp, blockStatsHTML)
blockStatsFrame.attemptSetSceneRect(QRectF(10, 10, 30, 20))
blockStatsFrame.setFrameEnabled(True)
blockStatsHTML.addFrame(blockStatsFrame)

# set HTML contents
blockStatsHTML.setContentMode(QgsLayoutItemHtml.ManualHtml)
blockStatsHTML.setHtml("test
test")
blockStatsHTML.loadHtml()

# reposition blockStatsFrame
blockStatsFrame.attemptMove(QgsLayoutPoint(10,10, QgsUnitTypes.LayoutMillimeters))

# set frame outline
blockStatsFrame.setFrameStrokeColor(QColor(255, 0, 255))
blockStatsFrame.setFrameStrokeWidth(QgsLayoutMeasurement(1, QgsUnitTypes.LayoutMillimeters))

# set frame background
blockStatsFrame.setBackgroundEnabled(True)
blockStatsFrame.setBackgroundColor(QColor(255, 255, 255))

# add frame to layout
comp.addLayoutItem(blockStatsFrame)

# add layout to layout manager
lmgr.addLayout(comp)


texts: How to create a distance from road raster in QGIS? (Maxent Bias)

texts:
I am trying to create a raster to show the distance from a road for an area in Australia.

I have the Open Street Map data downloaded and converted to a *.SHP file.

Is there a way to produce a raster showing a value relating to the distance to the nearest road?


texts:
Within a QGIS-course I would like to tell a bit about the proj.4 parameters. However, I'm not sure about the "+no_defs" parameter.
The proj.4 website (http://proj4.org/parameters.html#parameter-list) says:

  +no_defs   Don't use the /usr/share/proj/proj_def.dat defaults file

Can anyone explain me what this means? (What would it mean to use the /usr/share/proj/proj_def.dat defaults file?)


texts: I used the ModelBuilder to create a module that can automatically create a database to store the intermediate outputs and save the final output (one shapefile) in a folder. User can choose which folder to store the final output and the database will be created within that folder too. To avoid any confusions or overwrite issues, I'd like this module to automatically delete the database after the final output being generated. So I added the "Delete" function and used the final output as the precondition. However, it leads to an ERROR 000601 that saying the database may be locked by another application so it can't be deleted. I'm wondering if there is anyway that I could run the "delete" database function in the same module that creates the database.


texts:
Code below is working, but after opening  - there are No features   in  shp file.
What can be a problem?

import osgeo, os
from osgeo import ogr, gdal

shpFile = r'C:\Users\ReBurns.shp'

driver = ogr.GetDriverByName('ESRI Shapefile')

ptsDS = driver.Open(shpFile)
ptsLayer = ptsDS.GetLayer()

#Passing criteria  to select polygons:
myQuery = 'Area_ha > 39000.00'
ptsLayer.SetAttributeFilter(myQuery)
print ptsLayer
ptsLayer.GetFeatureCount()

#Save selected features to disk
outds = driver.CreateDataSource(r'C:\Users\ReBurns111.shp')
outlyr = outds.CopyLayer(ptsLayer,'ReBurns111.shp')

Got error:

Traceback (most recent call last):

  File "<ipython-input-167-1e169b0c7c9c>", line 44, in <module>
    outlyr = outds.CopyLayer(ptsLayer,'ReBurns111.shp')

AttributeError: 'NoneType' object has no attribute 'CopyLayer'

No selected files, but how to pass  selected by criterias?


texts: Layers with same spatial reference do not overlay?

texts:
I have two shapefiles (Parks and VMCC), and both the shapefiles have the same spatial reference (NAD 1983 HARN Adj MN Dakota (US Feet)).

The Parks.shp falls at the exact location where it should, the issue occurs with VMCC and it is falling at an inaccurate location.

Extent information VMCC

What all I have tried:-

Used the project tool to set the spatial reference for both the layers to Web Mercator
Deleted the .prj file and used the define projection tool to give the required projection

The behavior is reproducible in ArcMap and ArcGIS Pro.
What is a possible reason for this?


texts: Identification of cloudy days based on Net Radiation data

texts: How can I identify that a particular day is cloudy or clear sky based on the daily time series of Net Radiation (calculated at FluxNet/AMERIFlux sites)? I read somewhere that clear sky radiation has sinusoidal behavior. If we get different behavior that means it is induced by clouds. I don't know how exactly approach this problem? Here is the sample data: https://drive.google.com/file/d/1f1YfWgPbg3fxhFG1DTNH9Ex3lTDqLBZp/view?usp=sharing


texts: Intersecting geotiff with another one using alpha-channel (or nodata value) in GDAL?

texts:
I have two geotiffs related to one location but with different areas covered:

1.

2.

How can I intersect their opaque parts into single geotiff like that?


texts:
How can I create a point that identifies the outer most point/corner of a polygon (one point per polygon)?.

It is easy to identify the centre point, but I cannot find options for the outer most corner.


texts:
I'm trying to see what color options I have for formatted text for arcgis desktop 10.5. Mostly I'm trying to get these colors to work in a legend text string. I know the basic color options:

<CLR BLUE ='255'> TEXT </CLR>
<CLR GREEN ='255'> TEXT </CLR>
<CLR RED ='255'> TEXT </CLR>
<CLR YELLOW ='255'> TEXT </CLR>
<CLR MAGENTA ='255'> TEXT </CLR>

I know that you can lighten or darken the text color by adjusting the numerical value.

Are you somehow able get colors between the primary RGB and CMYK bands? For example, orange?


texts: unable to draw polygon without attribute window popping up?

texts: I'm fairly new to QGIS and an trying to draw a polygon via "Adding Feature". When I make the first left click to start, the attribute window pops up instead of letting me continue on with additional points. These instructions don't seem to offer a solution?


texts:
I am currently using gstat to perform ordinary kriging on daily fields of some variable - that is a new variogram is fit for every daily set of observations.

I would like to compare this with pooling the data over a month (say) and computing a single variogram model from the aggregated data. Note that the locations are fixed in time.

It seems to me that the variogram should be binned only from correlations within a single day, but that aggregating these terms over all days will improve the smoothness (maybe) of the resulting variogram.

My question is how do I achieve this in gstat? I have found this (unanswered) question: https://stats.stackexchange.com/questions/202140/kriging-fit-method-and-dx-argument-in-variogram

This mentions using dX=0 but this seems to only ignore contributions between the same spatial locations.


texts: "There is not enough memory" Error when accessing IArea.LabelPoint

texts:
I have an intermittent error that I'm having a terrible time tracking down so I thought I'd see if anyone has encountered something similar.

I am trying to access the label point of a polygon feature to put a text element there.  Using ArcGIS 10.0 SP5, VB.net.

When I attempt to set the element's geometry to IArea.LabelPoint, I alternatively get "There is not enough memory," "Operation attempted on empty geometry," or "The operation could not be performed because the object was geometrically degenerate."  So far this is only happening for one feature - other features in the same feature class work fine.  The feature class is in an SDE so I can't do a Check Geometry.

The feature geometry is NOT empty (checked with IFeature.IsEmpty).  I've checked the IPolygon4.IsEmpty property and it is False...until I try to do something with IArea, at which point sometimes changes to True (but not always)...regardless, when I check for empty geometry, the check passes fine, so doesn't stop the error.  I've included the relevant bits of my code below - any ideas?  At this point I don't even care about making the element for this particular feature, I'd be happy just to know what condition to check so that I don't get a fatal error when I try to access its LabelPoint.

(pFeature is an IFeature determined ahead of time)

Dim pTextElement As ESRI.ArcGIS.Carto.ITextElement
Dim pElement As ESRI.ArcGIS.Carto.IElement
Dim pPolygon As ESRI.ArcGIS.Geometry.IPolygon4
Dim pArea As ESRI.ArcGIS.Geometry.IArea
 If pFeature IsNot Nothing Then
     If pFeature.Shape IsNot Nothing AndAlso Not pFeature.Shape.IsEmpty Then
          'get the geometry
          pPolygon = pFeature.Shape
          'get the pointer to the IArea interface so we can get the labelpoint
          pArea = CType(pPolygon, ESRI.ArcGIS.Geometry.IArea)
          'new TextElement
          pTextElement = New ESRI.ArcGIS.Carto.TextElement
          'Query Interface (QI) to an IElement pointer and set the geometry
          pElement = pTextElement
          pElement.Geometry = pArea.LabelPoint
      End If
  End If

The error occurs on the pArea.LabelPoint, and debugging shows that the error is with the LabelPoint; pArea is not Nothing and has 0 area and one of the error messages listed earlier for Centroid and LabelPoint.

Edit:
In case it helps, the various properties of the 3 main objects during debugging:

Edit #2:
Because it didn't paste well in the comments, here was my workaround:

        Dim pSafeFeature As ESRI.ArcGIS.Geodatabase.IFeature
        Dim pParcelFC As ESRI.ArcGIS.Geodatabase.IFeatureClass
        If pFeature IsNot Nothing Then
            pParcelFC = CType(pFeature.Class, ESRI.ArcGIS.Geodatabase.IFeatureClass)
            pSafeFeature = pParcelFC.GetFeature(pFeature.OID)
...

Proceeding with pSafeFeature instead of pFeature, everything works fine.


texts: How to calculate color distances in raster data in QGIS?

texts:
I have to compare aerial photographies with geomagnetic data. The photos are georeferenced scanned slides showing structures on the surface. The magnetic data have a certain unit, but the photos do not. For a comparision I need a unit for the structures in the photos (like Delta E colour difference). I vectorized the structures of interest. Now I would like to use some kind of function to calculate the colour difference between a structure within the polygon and the surrounding soil.

Does anyone know such a function in QGIS?

I thought about converting the photos into greyscale images first, to compare the greyscale cell-values. But it would be much better to use the color images.


texts: Overlaying of raster in GRASS for visualization?

texts:
I have generated  insolation time rasters map in grass using r.sun module.

How to overlay insolation time rasters over DSM depicting the color table for visualization?


texts: error downloading from OpenStreetMap ("..server replied: Gateway Time-Out")

texts:
I've been trying to download OpenStreetMap data (vectorial > OpenStreetMap > Download data), but I get the error shown in the image below. Anyone know how to fix it? I have tried in different computers (laptops and pc), this error started to appearing today I used this plugin weeks ago.

I'm using QGIS 2.8.3, Windows 7, 64 bits.

Sorry for the words in Spanish


texts: Georeferencing uneven/irregularly gridded rasters using GDAL?

texts:
I'm trying to take a set of NetCDFs and turn them into georeferenced datasets using GDAL.

However, the only Geospatial data is a set of 2-D lat/lon arrays. The projection type is not lat/lon. THe arrays give the lat/lon of every individual cell.

I can't use a geotransform to reference the data because of the irregular gridding.

How can I create a dataset (such as a GeoTiff) that properly references the data?

Ideally I would then warp to some even grid.

I tried using an even spread of Ground Control Points, but the resulting GeoTiff did not display correctly... Despite giving GCPs up to 90 deg lat, ArcMap claims the GTiff has an extent ending at 70 deg latitude...

I was able to solve the issue using gdalwarp -geoloc as  suggested. I created 3 VRTs (1 with cell values, 1 with latitude, 1 with longitude), and wrote the lat/lon vrts in a geoarrays in the data vrt metadata. Then using gdalwarp -geoloc worked like a charm.

The issue with extent was unrelated and due simply to the large difference between the original coordinate system and lat/lon. It was solved by explicitely stating the extent instead of having GDAL try to guess what it should be. (i.e. I added -te -180 0 180 90 to the gdalwarp line)


texts: Creating weight raster for flow accumulation in ArcGIS for Desktop?

texts:
I am studying flood flow using ArcGIS for Desktop.

When creating flow accumulation in hydrology, there is a weight raster.

However, İf you do not have it, the programme sets its value to 1 which is not a suitable approximation for the simulations.

How can I create a weight raster in ArcGIS for Desktop?


texts: QGIS show gdal_polygonized file in EPSG:54004 instead of EPSG:3857

texts:
I am polygonizing a TIFF with the gdal_polygonize.py and the resulting Shape is shown in QGIS in the wrong projection - EPSG:54004 instead of EPSG:3857.

I have a TIFF in EPSG:3857 and polygonize it with

gdal_polygonize.py input.tif -f "ESRI Shapefile" output.shp

After doing this, the corresponding .prj file looks ok to me and I also checked it with various other sources like this

PROJCS["WGS_84_Pseudo_Mercator",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Mercator"],PARAMETER["central_meridian",0],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1],PARAMETER["standard_parallel_1",0.0]]

When I load the output.shp into QGIS it automatically detects the file in EPSG:54004 which results in an offset compared with the source. After changing the Shapes projection in QGIS manually to EPSG:3857 both files, the TIFF and the Shape, do align.

Could this be a problem with QGIS or a problem with GDAL?

Offset between TIFF and Shape in the red rectangle

Update

gdalsrsinfo on the original TIFF shows

PROJ.4 : '+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs'

OGC WKT :
PROJCS["WGS 84 / Pseudo-Mercator",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Mercator_1SP"],
    PARAMETER["central_meridian",0],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs"],
    AUTHORITY["EPSG","3857"]]


texts: A reproducible example for the use of separability measures in improving selection of spectral variables

texts:
I have been working with the tutorials in GEE on supervised and unsupervised classification, but I haven't come across any examples on the use of separability measures such as Jeffreys-Matusita, transformed divergence or B-distance for instance.

Here's my code to run a classification using CART()

// Add a clipping function to clip the image collection

var clipper = function(image){
  return image.clip(WG);
};

// Import raw Landsat 5 scenes and clip them over your study area
var l5filtered = ee.ImageCollection('LANDSAT/LT05/C01/T1').filterDate('2008-01-01','2012-12-31').map(clipper);

// Obtain a cloud-free composite of the scenes filtered
var cloudfree = ee.Algorithms.Landsat.simpleComposite({
  collection: l5filtered,
  asFloat: true
  });
Map.addLayer(cloudfree, {bands: ['B4','B3','B2'], max: 0.3}, 'composite');

// Specify the bands of the Landsat composite to be used as predictors (p)
var predictionBands = ['B1','B2','B3','B4','B5','B6','B7'];

// Step 2: Creating training data manually
// Draw a polygon around a bare ground area and import it as a featureCollection
// Click New property and import the new 'class' and give it a value of 0
// Similarly vegetation and water features were created

// Merging the three featureCollections to obtain a single featureCollection
var trainingFeatures = bare.merge(vegetation).merge(water);

// Note - we use .sampleRegions, if we have a featureCollection and we want to sample an image
// given a featureCollection

var classifierTraining = cloudfree.select(predictionBands)
                                  .sampleRegions({
                                    collection: trainingFeatures,
                                    properties : ['class'],
                                    scale: 30
                                  });

// Let's instantiate a CART and train it
var classifier = ee.Classifier.cart().train({
                   features: classifierTraining,
                   classProperty: 'class',
                   inputProperties: predictionBands
});

// Classify the image
var classified = cloudfree.select(predictionBands).classify(classifier);
Map.addLayer(classified, {min:0, max:2, palette: ['red','green','blue']}, 'classified');

// ACCURACY ASSESSMENT
// Step 1: Partition the set of known values into training and testing sets

// Adding a randomColumn of values ranging from 0 to 1
var trainingTesting = classifierTraining.randomColumn();

var trainingSet = trainingTesting.filter(ee.Filter.lt('random',0.6));
var testingSet = trainingTesting.filter(ee.Filter.gte('random',0.6));

// Now run the classifier only with the trainingSet
var trained = ee.Classifier.cart().train({
  features: trainingSet,
  classProperty: 'class',
  inputProperties: predictionBands
});

// Now classify the testData and obtain a Confusion matrix
var confusionMatrix = ee.ConfusionMatrix(testingSet.classify(trained)
                                                   .errorMatrix({
                                                     actual: 'class',
                                                     predicted: 'classification'
                                                   }));

// Now print the ConfusionMatrix and expand the object to inspect the matrix()
// One can also obtain basic descriptive statistics from the confusionMatrix
print('Confusion matrix:', confusionMatrix);
print('Overall Accuracy:', confusionMatrix.accuracy());
print('Producers Accuracy:', confusionMatrix.producersAccuracy());
print('Consumers Accuracy:', confusionMatrix.consumersAccuracy());

Unsure of how to calculate separability measures. I have shared the assets necessary. Here's the code link as well - https://code.earthengine.google.com/cb156ba6f71b6448f08bf106162c6373


texts: Address locator alternative for ArcGIS Online

texts:
I do not have a license for ArcGIS Server.

I created an address locator using my county's centerline data. I edited my school boundaries so they "snapped" to those centerlines. Now when I do population analysis I can know that my selections of geolocated students are in the right boundary. Great! Now lets make a webmap for families to search and find their home schools! Found Here

ISSUE:  The Esri geocoder is different "off" from my centerlines in a few places.  On a few streets, where the attendance boundaries should go down the middle of the streets, it is telling families the wrong school because the (polygon) boundary is a few feet over the centerline.

Example: 10741 FAIR OAKS BLVD, FAIR OAKS 95628 was an address of concern. The app pop up says the wrong school because of variation from county centerlines used to create the boundaries and the ArcGIS map.

WHAT I WANT TO DO:  Use a "free" or very cheap locator that I can create using my already made address locator or something open I can edit when issues on a few street arise.  I want to add this as the address search bar in my app.

Is this possible? If so, how?  Details needed, I have very basic Python experience and am the only GIS professional at my job.


texts: ArcGIS Pro (2.3.0) Forest-based Classification and Regression Issue

texts:
I keep getting this error when I'm using the 'Predict to raster' Prediction Type in the Forest-based Classification and Regression tool:

  'Problem extracting values from a.img using origin location [x,y]'

The given origin location is not on a NoData area, it's value is: Stretch.Pixel Value 0.016262

Is the issue that it's in the ERDAS Imagine format? Should I make it a GRID or .tif? Every raster is defined in the same GCS and the same PCS. Additionally, the data frame is in the same GCS & PCS.

The Forest-based Classification and Regression tool works fine in the 'Train only' Prediction Type mode.

Any help is greatly appreciated.


texts:
I've always been into spatial thinking and maps, and it's been a long time I dream of a tool/app that would allow me to have my own customized interactive maps. But now that I'm about to enter real estate development, it has become an absolute necessity.

I just discovered the world of GIS a couple hours ago, and it seems to allow way more than what I thought my needs would be. I'll take my time and learn as much as I can and need on my own, but to make it interesting and practical, I want to start addressing my initial requirements :

A tool that would allow me to create several layers on updated and detailed city maps, using polygons and paths, completed with my own customized pop-ups/data forms when I click or point these polygons. I would like to use/import data from google maps as well (city boundaries for instance).
Example: For a given neighbourhood, I identify every plot/building with a polygon, and to each polygon I assign a customized data form: status, owner, price, type etc. It would allow me to spatialize the data gathered on the said neighbourhood and help me know better where to invest. The goal is to obtain an interactive map on which I can click and obtain plot by plot data (using colored polygons and popups/forms).

What would be the most user-friendly tool/software/app to address this need ? I'd help me target accurately my learning process.


texts: Adding Origin to Geotiff using Python GDAL?

texts:
Alright, so I'm getting closer to figuring this out but I'm stuck at one (that I know of) thing... how do I use GDAL to apply a known origin of a GeoTiff to another that has been processed in Photoshop?

The snippet below is what I'm trying to get working. It's a pretty rough 'sketch' of what I want to accomplish and will work it out once I get going in the right direction. I'm just hoping that you can tell what I mean by looking.

So, I take the data from a good gtiff image, and want to copy it to the same image that has been broken. Is there any way to do that?

def applyRefs(input, output):
print "processing %s" %input
in_file = gdal.open(input)
gt = in_file.GetGeoTransform()
if not geotransform is None:
print 'Origin = (',geotransform[0], ',',geotransform[3],')'
else:
print" this is where the function to copy world file info to processed tif `will go"`

https://github.com/minorsecond/GIS/blob/master/tifFix.py


texts: MapInfo License Server - Can It Serve Multiple Versions of MapInfo

texts:
I have a MapInfo License Server Utility that serves 30+ MI v10.5 concurrent licenses across my organisation.

v11.5 has now landed on my desk (with new serial numbers and access codes), so can I activate my v11.5 licenses within the same instance of License Server utility, whilst retaining the v10.5 licenses, in order to gradually migrate users onto v11.5 over the coming months?


texts:
I've been asked a question regarding installation of QGIS from an unpacked file on Windows, as the organization that wants to install it needs to run a thorough malware scan.

Currently, the installation of QGIS in standalone form requires downloading of the setup file which then proceeds to download and unpack all the QGIS files. This doesn't allow a malware scan beforehand.

Is there a way to install a "pre-unpacked" QGIS program?


texts:
Simply I've got a number of PostGIS layers that have a geometry type of Geometry, as a posed to point, line, polygon.

wkb_geometry geometry(Geometry,27700)

QGIS 2.0+ doesn't like this very much and wants me to set the geometry type and in general makes my PostGIS connect in QGIS slow.

How do I change the Geometry type to point for instance?

Si


texts: I have taken a number of images from a very large rock overhang that contains ancient paintings. I want to piece the images together to create a single image of the entire rockshelter. Through ARCGIS I tried to do this by georeferencing the images, which have no spatial data. It worked quite well. However, I don't know how to save it so that when I move onto the next image, the previous georeferencing remains. How can I save it?


texts: googlemutant and leaflet problem view data

texts:
I just create a new map with leaflet and GoogleMutan with this three CDN libriaries:

1) maps.googleapis.com/maps/api/js?key=MYKEY" async defer

2) unpkg.com/leaflet@1.0.1/dist/leaflet-src.js

3) unpkg.com/leaflet.gridlayer.googlemutant@latest/Leaflet.GoogleMutant.js

Why the map shows me numbers above of the tiles? there some conflict from the cdn libraries?


texts: Prism mobile phone metadata: does it include location?

texts: I've read that Prism collects metadata, but no details of what that metadata includes.  Does anyone know if location is part of the metadata?


texts: Changing geometry type from Point to Multipoint within existing table in PostGIS?

texts:
Is there a PostGIS function that can change geometry type for an existing table?

We need to change from POINT to MULTIPOINT.

The table will be empty when we change the geometry type and we cannot just drop/create the table.


texts: Finding coordinates for 'TRMM 3B42: 3-Hourly Precipitation Estimates' data covering region in Malawi

texts:
I have weekly precicipation data covering a certain region in Malawi.
The data is from google earth engine.

I would like to know the coordinates from the data was extracted.

The data was extracted using this script:

var dataset = ee.ImageCollection('TRMM/3B42')
                  .filter(ee.Filter.date('2018-04-01', '2018-04-10'));
var precipitation =
    dataset.select(['precipitation', 'HQprecipitation', 'IRprecipitation']);
var precipitationVis = {
  min: 0.0,
  max: 12.0,
  gamma: 5.0,
};
Map.setCenter(-79.98, 23.32, 4);
Map.addLayer(precipitation, precipitationVis, 'Precipitation');

Is there a way I can get the coordinates of this region?


texts: Why is ArcGIS Desktop labeling wrong polygons?

texts: I am working with a very large data set. When I label my polygons, the labels display like a highway of labels through the middle of my map. I'm using ArcGIS 9.3. I have tried to fix it with Maplex, but then won't even display in any fashion. When I use the standard label method, they show up, but they are all over the place, some polygons show 3 or 4 labels, only one of them being the label for that polygon. Any idea what the problem, or better yet, solution might be?


texts: Possible to base GIS layers on results of calling API

texts: Is it possible to call an API which returns spatial data as JSON, GML, XML for example and use it as as data source within QGIS or other GIS? I don't know enough about APIs to understand if I am asking too much of them and whether there needs to be some other steps to this. For example - how would multiple sets of parameters be passed and collated as single result sets? Should this be done via Python to call API multiple times, obtain data (how?) and then assemble in to a shapefile some how?


texts: Change PostGIS Column from Geometry to Geography

texts:
I have a PostGIS table with a column of type Geometry(POINT) and SRID 0 that must be changed to the type Geography(POINT) as it will be used to store lat/lon values.
How can we make this alteration to the table without dropping it?

Table was originally created using SQLAlchemy, Geoalchemy2 and Python, with the following model defination class:
class User(Base):
    __tablename__ = 'users'
    id                      = Column(Integer, primary_key=True)
    name                    = Column(String())
    coordinates             = Column(Geometry('POINT'))

which should now be changed to
class User(Base):
    __tablename__ = 'users'
    id                      = Column(Integer, primary_key=True)
    name                    = Column(String())
    coordinates             = Column(Geography('POINT'))


texts: Update field the_geom name city in CartoDB?

texts:
With "update table set the_geom = cdb_latlng (lat , lon ) " I update the_geom field through latitude and longitude, but I also interested to know if there is any chance that instead of passing a latitude and longitude , pass a string is 'Madrid' for example, and you will then update the field the_geom with the coordinates of Madrid.

Con "update table set the_geom = cdb_latlng(lat, lon)" actualizo el campo the_geom a través de latitud y longitud , pero tambien me interesaria saber si hay alguna opción de que en vez de pasarle una latitud y una longitud, pasarle un string que sea 'Madrid' por ejemplo y que te te haga update del campo the_geom con las coordenadas de Madrid.


texts: Reshaping polygon on curve shape on road in QGIS

texts: In ArcMap Editing Tool -> which is End Point Arc Segment. I need similar kind of tool in QGIS.


texts:
I'm trying to generate a graticule with user-defined range of latitude, longitude, and interval; however, my code only displays one line of latitude and one line of longitude when loading it in QGIS. How do I create a grid of lines of latitude and longitude?

multiline = ogr.Geometry(ogr.wkbMultiLineString):
for lon in range(min_lon, max_lon, interval):
    line1 = ogr.Geometry(ogr.wkbLineString)
    lat = min_long + interval
    line1.AddPoint(min_lon, lat)
    line1.AddPoint(max_lon, lat)
    multiline.AddGeometry(line1)

for lat in range(min_lat, max_lat, interval):
    line2 = ogr.Geometry(ogr.wkbLineString)
    lat = min_lat + interval
    line2.AddPoint(lon, min_lat)
    line2.AddPoint(lon, max_lat)
    multiline.AddGeometry(line2)

print(multiline.ExportToWkt())
with open ('multiline.txt', 'w') as f:
    f.write(multiline.ExportToWkt())


texts: Customizing ISO 19139 metadata profile in ArcCatalog?

texts:
I'm using the ISO 19115/19139 metadata standard, and I need to convert some optional elements into mandatory elements. Is it possible to achieve this using ArcCatalog as the metadata editor?

About metadata standards and profiles says "A profile may further restrict the optionality of a metadata element, making it mandatory where before it was optional."

In ArcCatalog I can set ISO 19139 as the metadata standard:

but how can I edit my instance of ISO 19139 (to make a custom profile) so that (for example) lineage is mandatory instead of optional?

Under ArcGIS\Desktop10.4\Metadata\Translator\Rules there is a file ISO19139_min_schema.txt.

Can I edit this file, and if so what is the syntax for making lineage mandatory?


texts: IDW interpolation of values from raster in ArcGIS

texts:
I've been tasked with obtaining z-values for a large point set from a 2 m raster, but rather than extracting values directly from the raster, or with bilinear interpolation, I've been asked to use inverse distance weighted interpolation of cell centers in a 5 m radius.

This is beyond the capabilities of the Extract Values to Points tool or Zonal Statistics. As a workaround, I can convert the raster to points, resample to a finer resolution raster using IDW, and then extract values, but this isn't quite true IDW since it's limited to the cell size of the new raster (and by my patience for processing time).

Setting aside the fact that there are better approaches to smoothing elevations from a raster --- this is the approach I've been asked to take --- I could also script it in python fairly easily, but thought I might be missing something easy?

I have access to 3D Analyst, Spatial Analyst, etc in ArcGIS Desktop 10.3.


texts: using field calculator on string field that has some empty values?

texts:
I have a CityState field that has the city and states of landowners. Bellingham WA, for example. I am making two new fields, OwnerCity and OwnerState in which I strip the last two characters of CityState for OwnerCity (RTrim([CityState],2) and capture the last two characters for OwnerState (Right([CityState],2). Some of the CityState cells are empty, so I can't run the calculator.

What can I add so that it runs the calculation where there is text, and keeps the value empty if there is not any text?  I'm sure there's an easy answer to this, but I have't figured it out yet.


texts: Using Google Earth Engine to exclude upper 1% of data?

texts:
I would like to exclude the upper 1% of my Data ('VV' Band) because of noise reduction

var s1_data = ee.ImageCollection('COPERNICUS/S1_GRD')

var s1_grd_iw_vv_des = s1_data
  .filterDate(begin,end)
  .filter(ee.Filter.eq('relativeOrbitNumber_stop',orbnr))
  .filter(ee.Filter.eq('sliceNumber',slinr))
  .filter(ee.Filter.eq('instrumentMode','IW'))
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VV'))
  .filter(ee.Filter.eq("orbitProperties_pass", "DESCENDING"))
  .sort('system:time_start');

I guess i have to write a function like

var upper_percent = function(image) {
   var one_percent = ????
   var maskedImage = image.mask().and(one_percent.not());
   return image.updateMask(maskedImage);}

I dont know the command to call? Is there a pre-defined one?

EDIT:
I try this

var filter2 = function(image) {
  var nn = image.select('VV').reduce(ee.Reducer.percentile([99]));
  var out = image.select('VV').gt(nn.select('p99');
  var maskedImage = image.mask().and(out.not());
  return image.updateMask(maskedImage);}

but if i look up the histogram there were no pixel excluded (masked). Dont matter if i try 99, or 90 or even 50 ?


texts:
In ArcMap, I can easily enforce my geodatabase relationships rules (ex. one-to-one relationships) by selecting features and selecting, in the editor toolbar, the Validate Features button. I can then have a list of features that does not meet my Relationship Class rules

In ArcGIS Pro, is there a similar way to do so, without requiring a Data Reviewer Licence (ArcMap's Validate Features doesn't)?


texts: Box Select using leaflet by modifying box zoom with turf intersect

texts:
I have been trying to select multiple polygons using box select. I have modified box zoom for this purpose and included turf.js for intersect function. The problem is if i select a particular area and use interset function it shows polygons from different area. If i select two polygons it show four polygons from another area.below is my code

Tried postgis intersection, bounds intersect in leaflet, turf intersect all results in same.

Note - geojson data contains 7000 polygons.

 $(function() {
        L.Map.BoxZoom.prototype._onMouseUp = function (e) {
            this._finish();

            if (!this._moved) { return; }

             var nw = this._map.layerPointToLatLng(this._startLayerPoint);
            var se = this._map.layerPointToLatLng(this._map.mouseEventToContainerPoint(e));
            var bounds = new L.LatLngBounds(nw,se);
            getparcels(bounds);
        }
    })

function getparcels(bounds) {
       var geodata = L.rectangle(bounds).toGeoJSON();
       geojsonlayer = L.geoJson(Data, {
          onEachFeature: function (feature, layer) {
              try{
                if(turf.intersect(geodata,feature)){
                   layer.setStyle({
                      fillColor: 'black',
                      fillOpacity: 0.3
                   })
                   map.addLayer(layer)
                 }
                }catch (e) {
                   console.log(e)
                 }
          }
     })
}

need a function to select multiple polygons by box-selections


texts:
I opened a georeferenced NAIP photo to a new project in QGIS.  The coordinates were expressed as: 2199240.6, 616858.7.  The same point in Google Earth was expressed in degrees, minutes, seconds, as 41°51'25.33"N,  76°54'12.14"W.  In what form is QGIS displaying the coordinates?  Are these correct?

When I opened the NAIP tif in QGIS, the CRS set itself to NAD83/Pennsylvania North.


texts: Calculating surface and lines in 3D projection in R

texts:
For my internship research I am trying to compare the surfacearea and length of geometric shapes and lines in one coordinate system to another. I have a SpatialPolygonDataFrame in the European equal area projection ETRS-LAEA and I have a function that can calculate the surfacearea and the length of all polygons. I can transform it to another projection and calculate the surfaceareas and the lengths of all polygons in that projection.

My next step would be to calculate the surfacearea and the length of the polygons in ETRS89(x,y,z). Since this is a 3D projection the functions I am using for the 2D projections are not working. Are there functions for 3D projections? I can not find them online.

I am hoping there has been some smart guy with a similar problem and already fixed it. I am expecting this is not the case and that it is (for me) impossible in R.


texts: Open source, Python option for Arcpy's Con tool?

texts: I'm wondering if there is a method to perform a function similar to Arc's Con tool, using open source options in python. I have two rasters, one that shows locations of clouds and the other is just a landsat band. If a value is >x in the cloud raster, I want to take the pixel from the landsat band, and if it's


texts: How test whether Catalog window of ArcMap is open (and docked)?

texts: From ArcMap 10 when the Catalog window is open, how can I test whether  ArcCatalog isOpen (and docked) using ArcObjects and C#.


texts:
I am using QGIS 1.8.0 and it runs slowly (with huge layers the size of NSW, Australia) and it appears that QGIS is limited in its RAM allocation to 128.9MB. How can I increase the speed with which it deals with such big files by allocating a larger amount of my computers power to QGIS?

I am running Windows 8 Pro (64-bit) on:
 - Intel Pentium P6200 @ 2.13GHz (dual core);
 - 2.00GB RAM
 - AMD Radeon HD 6300M Series

Engine clock speed: 500-750 MHz
Processing power (single precision): 80-120 GigaFLOPS
Polygon throughput: 125-187.5M polygons/sec
Data fetch rate (32-bit): 16-24 billion fetches/sec Texel fill rate (bilinear filtered): 4-6 Gigatexels/sec Pixel fill rate: 2-3
Gigapixels/sec Anti-aliased pixel fill rate: 8-12 Gigasamples/sec
Memory clock speed: 800-900 MHz Memory data rate: 1.6-1.8 Gbps Memory bandwidth: 12.8-14.4 GB/sec


texts:
I have several raster layers with costs associated to distances from certain cells. Some distances are positive (i.e. the further from the cell, the better) and others are the opposite (i.e. the closer, the better).

How can I use GRASS's commands to convert some layers to find the negative distance-costs? Could r.mapcalc work?


texts: ArcSDE Problems: Microsoft SQL Server 2008 R2 on Windows 2008 R2 Standard (64-Bit)

texts:
I´m trying to set up ArcSDE 10.2 on a  Microsoft SQL Server 2008 R2 running on a Windows Server 2008 R2 Standard 64-Bit machine.
I can create the database connection, but when I want to enable it for spatial data I keep on getting the following error message:

Error: The version of your database server is not supported
ERROR installing/upgrading ArcSDE, Error = -494

This is strange because the ESRI page on System Requirements states, that this combination should work.

Any idea what I could do, to get it fixed?
Thanks

EDIT:

I´m using Microsoft SQL Server 2008 R2 64-Bit, or at least that is what Windows tells me (see image). I´ve tried without SP, than with SP1 and then SP2. Neither option worked.


texts: How to make polygons from interpolated points?

texts:
My goal is to make a shapefile with polygons for variable rate fertilizer application within a farm field.

I have a point shapefile with soil sample values (pH, P, K and so on) and a line/polygon file with field borders. IDW is recognizes as a suitable interpolation method for soil samples. I have gotten as far as having the interpolated points as a raster file and I have experimented with the "Contour" function to make a suitable classification.

I'm looking for a solution in QGIS to convert/make/bake the classified raster in to polygons, and clip the results to the field border. Does any one have a suitable solution for this using QGIS?

Edit:

I don't know if my question is clear enough. To illustrate my goal the attached picture shows how I need my final shapefile to be.

The picture illustrates the final result I'm after. The points are interpolated, clipped to the field border and then (some how) the classified raster is converted to a polygon shapefile.


texts: Set restricted segments in a network with QGIS

texts:
I'm trying to set up a network with QGIS. When I set the direction field there are three options for the values in that field:

Forward direction
Backward direction
Both directions

But I need to set some segments of the network as restricted. I have try with the QGIS Processing network analysis tools and with QNEAT3 plugin, but I can't find an option to restrict some segments.
Is there any way to do this?


texts: ArcGIS Online: lines are not displaying correctly on Esri Map Viewer mobile app

texts: I have a polyline shapefile uploaded to my ArcGIS Online Web Map that displays correctly via the web map and web app, but incorrectly (random shards and missing segments) in the Esri "Map Viewer" for mobile apps. Any idea what could be causing this?


texts: HTML file not loading for CartoDB.js using Leaflet

texts:
I am trying to create a web map that uses a layer hosted on CartoDB and a basemap created in Mapbox studio. Never mind adding the basemap to the L.map, I cannot get the cartoDB layer to draw at all. My code is below (without the Mapbox portion):

    <title>RCC Global Map</title>
    <meta name="viewport" content="initial-scale=1.0, user-scalable=no" />
    <meta http-equiv="content-type" content="text/html; charset=UTF-8"/>
    <link rel="shortcut icon" href="http://cartodb.com/assets/favicon.ico" />

    <link rel="stylesheet" href="http://libs.cartocdn.com/cartodb.js/v3/3.15/themes/css/cartodb.css" />


texts: Odd characters showing up in symbology view

texts: I am trying to symbolize based on company name but there are weird Greek symbols showing up in my Symbology view of the actual field values. Has anyone ever had this problem? The weird symbols do not show up in the Table View.


texts: Loop for in arcpy: for several rasters, clip with the same shapefile and create an output

texts:
I am building an arcpy script, where I would like to create a for loop, that for each raster, it clips it from a shapefile, and creates a new output, with name based on the file itself. The code:

import arcpy
arcpy.env.workspace= r"C:\Angela\LC"

#import the raster files for land cover
#Local variables:
n1992_USA = "C:\\Angela\\OK\\nlcd_1992_30meter_whole (1)\\nlcd_1992_30meter_whole.img"
n2001_USA = "C:\\Angela\\OK\\nlcd_2001_landcover_2011_edition_2014_10_10\\nlcd_2001_landcover_2011_edition_2014_10_10\\nlcd_2001_landcover_2011_edition_2014_10_10.img"
n2006_USA = "C:\\Angela\\OK\\nlcd_2006_landcover_2011_edition_2014_10_10\\nlcd_2006_landcover_2011_edition_2014_10_10\\nlcd_2006_landcover_2011_edition_2014_10_10.img"
n2011_USA = "C:\\Angela\\OK\\nlcd_2011_landcover_2011_edition_2014_10_10\\nlcd_2011_landcover_2011_edition_2014_10_10\\nlcd_2011_landcover_2011_edition_2014_10_10.img"

#import rasters
arcpy.MakeRasterLayer_management(n1992_USA, "1992t", "", "-2361915 267885 2266005 3177735", "")
arcpy.MakeRasterLayer_management(n2001_USA, "2001t", "", "-2361915 267885 2266005 3177735", "")
arcpy.MakeRasterLayer_management(n2006_USA, "2006t", "", "-2361915 267885 2266005 3177735", "")
arcpy.MakeRasterLayer_management(n2011_USA, "2011t", "", "-2361915 267885 2266005 3177735", "")

#import admin boundaries
admin_boundaries_usa = "C:\\Angela\\tmp\\cb_2013_us_county_5m.shp"
arcpy.MakeFeatureLayer_management(admin_boundaries_usa, "admin_boundaries_usa", "",
                                  "", "FID FID VISIBLE NONE;Shape Shape VISIBLE NONE;STATEFP STATEFP VISIBLE NONE;COUNTYFP COUNTYFP VISIBLE NONE;COUNTYNS COUNTYNS VISIBLE NONE;AFFGEOID AFFGEOID VISIBLE NONE;GEOID GEOID VISIBLE NONE;NAME NAME VISIBLE NONE;LSAD LSAD VISIBLE NONE;ALAND ALAND VISIBLE NONE;AWATER AWATER VISIBLE NONE")

#crop the latter to texas
# Process: Select Layer By Attribute
arcpy.SelectLayerByAttribute_management("admin_boundaries_usa", "NEW_SELECTION", "\"STATEFP\" = '48'")

# Process: Copy Features
arcpy.CopyFeatures_management("admin_boundaries_usa", "tx_shp", "", "0", "0", "0")
arcpy.SelectLayerByAttribute_management("admin_boundaries_usa", "CLEAR_SELECTION")

#dissolve to have the whole state
arcpy.Dissolve_management("tx_shp", "tx_shp_Dissolve", "STATEFP", "", "MULTI_PART", "DISSOLVE_LINES")

#make loop for all rasters to be clipped w/ shapefile
rasterlist=[n1992_USA, n2001_USA, n2006_USA, n2011_USA]
OUTFOLDER = "C:\\Angela\\LC"

Everything works fine, until when I input the loop in ArcGIS:

 for inraster in rasterlist:
    print inraster
    tx_name = os.path.join(OUTFOLDER, str(inraster))
    arcpy.Clip_management(inraster, "-1075276.88510708 309884.007197296 251355.2459483 1547680.1003418", tx_name, "tx_shp_Dissolve", "255", "NONE", "NO_MAINTAIN_EXTENT")

and it retrieves the following:

Runtime error  Traceback (most recent call last):   File "<string>", line 4, in <module>   File "c:\program files (x86)\arcgis\desktop10.4\arcpy\arcpy\management.py", line 13771, in Clip     raise e ExecuteError: ERROR 000670: output Output Raster Dataset is same as input Input Raster

I know I must have a rookie mistake on the loop part, but I can't see what else I could do to make it work.


texts:
I'm trying to import a file database into PostGIS. Running ogrinfo on my database I get:

..$ ogrinfo BRP_Gewaspercelen_2009.gdb

Had to open data source read-only.
INFO: Open of `BRP_Gewaspercelen_2009.gdb'
      using driver `OpenFileGDB' successful.
1: BRP_gewaspercelen_2009_definitief (Multi Polygon)

So I guess I have the required driver to load the file. But ogr2ogr throws me this error:

 ..$ ogr2ogr -f "PostgreSQL" PG:"host=localhost port=5433 dbname=crops
 user=postgres" BRP_Gewaspercelen_2009.gdb
-overwrite -progress --config PG_USE_COPY YES

FAILURE:
Unable to open datasource `BRP_Gewaspercelen_2009.gdb' with the following drivers.
  -> ESRI Shapefile
  -> MapInfo File
  -> UK .NTF
  -> SDTS
  -> TIGER
  -> S57
  -> DGN
  -> VRT
  -> REC
  -> Memory
  -> BNA
  -> CSV
  -> GML
  -> GPX
  -> KML
  -> GeoJSON
  -> GMT
  -> GPKG
  -> SQLite
  -> WAsP
  -> PCIDSK
  -> OpenFileGDB
  -> XPlane
  -> AVCBin
  -> AVCE00
  -> DXF
  -> Geoconcept
  -> GeoRSS
  -> GPSTrackMaker
  -> VFK
  -> PGDump
  -> OSM
  -> GPSBabel
  -> SUA
  -> OpenAir
  -> PDS
  -> WFS
  -> HTF
  -> AeronavFAA
  -> EDIGEO
  -> GFT
  -> SVG
  -> CouchDB
  -> Idrisi
  -> ARCGEN
  -> SEGUKOOA
  -> SEGY
  -> XLS
  -> ODS
  -> XLSX
  -> ElasticSearch
  -> PDF
  -> CartoDB
  -> SXF

What am I missing?


texts: no output when using subprocess.call in gdal?

texts:
I'm trying to clip an image in gdal and when I run this code:

gdalwarp -cutline path/to/polygon.shp \
-crop_to_cutline -dstnodata "-9999.0" \
in_raster.tif \
out_raster.tif

it runs ok and output is working ok too.
However, when running this code:

import ogr
import subprocess

inraster = 'in_raster.tif'
inshape = 'polygon.shp'
outraster = 'out_raster.tif'

subprocess.call(['gdalwarp', '-cutline', inshape, '-crop_to_cutline', '-dstnodata -9999.0', inraster, outraster])

it doesn't create any output. I reckon there is something wrong in  -dstnodata <nodata_value>.


texts: Geoserver 2.9 Rest API cURL Set SRS for Directory of Shapefiles

texts:
I'm automating the publishing of layers in Geoserver 2.9. I've had good luck with getting my GeoTIFFs to publish nicely, despite having garbage for incoming SRS info.

With the GeoTIFFs, I POST a coverageStore definition and then POST a coverage definition. It works great.

export t="American_Redstart_100_1_abund_aut.tif"
export sp=American_Redstart_100_1
export dataStoreName=${t%.*}

curl -v -u admin:geoserver -XPOST -H 'Content-Type: text/xml'
-d "<coverageStore>
      <name>${dataStoreName}</name>
      <enabled>true</enabled>
      <type>GeoTIFF</type>
      <url>file:///my_directory/data/species/${sp}/${t}</url>
      <workspace>stemlt</workspace>
    </coverageStore>"
"http://localhost:8080/geoserver/rest/workspaces/stemlt/coveragestores?configure=all"

curl -u admin:geoserver -XPOST -H 'Content-type: text/xml'
-d "<coverage>
      <name>${dataStoreName}</name>
      <title>${dataStoreName}</title>
      <nativeCRS>UNKNOWN</nativeCRS>
      <srs>EPSG:3857</srs>
      <projectionPolicy>FORCE_DECLARED</projectionPolicy>
      <latLonBoundingBox>
        <minx>-20051271</minx>
        <maxx>20055048</maxx>
        <miny>-2722745</miny>
        <maxy>20268684</maxy>
        <crs>EPSG:3857</crs>
      </latLonBoundingBox></coverage>"
"http://localhost:8080/geoserver/rest/workspaces/stemlt/coveragestores/${dataStoreName}/coverages"

However, I want to do something similar with a directory of shapefiles, but can't find the API parameters to do it. I'd like to be able to point to a directory of shapefiles and set the nativeCRS, srs, projectionPolicy, and latLonBoundingBox as they go in.

Currently, I am using PUT and external.shp?configure=all to load the directory of shapefiles, but then I have to manually update the projection information.

sp=American_Redstart_100_1

curl -v -u admin:geoserver -XPUT -H "Content-type: text/plain"
-d "file:///my_directory/data/species/${sp}/"
"http://localhost:8080/geoserver/rest/workspaces/stemlt/datastores/${sp}_grids/external.shp?configure=all"

How do I pass along projection information as I PUT or POST the directory of Shapefiles? I'm also okay with looping through loaded shapefiles, but can't find the right function for operating on those either.


texts: PyQGIS problem with QgsRasterCalculator CRS definition

texts:
from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry import processing
from osgeo import osr, gdal from gdalconst import *

dataset1 = gdal.Open(input1, GA_Update) band1 = dataset1.GetRasterBand(1) dataset1.SetProjection( srs.ExportToWkt() ) dataset1 = None

layer1=QgsRasterLayer(dataset1,'test1')
#layer1.setCrs( QgsCoordinateReferenceSystem(3003, QgsCoordinateReferenceSystem.EpsgCrsId) )

print('load first layer')

input2='D:\\ascii2grids\\dati_ascii\\PR_2018091500UTCplus1.txt'

dataset2 = gdal.Open(input2, GA_Update) band2 = dataset2.GetRasterBand(1) dataset2.SetProjection( srs.ExportToWkt() ) dataset2 = None

layer2=QgsRasterLayer(dataset2,'test2')

print('load II layer')

exp = "%s - %s" % (layer1,layer2)

output='C:\\Users\\assis\\test.tif'

entries=[layer1,layer2]

calc = QgsRasterCalculator(exp,
                    output,
                    'GTiff',
                    layer1.extent,
                    QgsCoordinateReferenceSystem("EPSG:3003"),
                    layer1.width,
                    layer1.height,
                    entries
                    )

res = calc.processCalculation(feedback) if res == QgsRasterCalculator.ParserError:
    raise QgsProcessingException(self.tr("Error parsing formula"))

This pyqgis script give me the following error:

Traceback (most recent call last):   File "C:\OSGEO4~1\apps\Python37\lib\code.py", line 90, in runcode
    exec(code, self.locals)   File "<input>", line 1, in <module>   File "<string>", line 57, in <module> TypeError: QgsRasterCalculator(): arguments did not match any overloaded call:   overload 1: argument 4 has unexpected type 'builtin_function_or_method'   overload 2: argument 4 has unexpected type 'builtin_function_or_method'   overload 3: argument 1 has unexpected type 'str'

I do not understand the error. Perhaps the CRS definition in the QgsRasterCalculator function. Unfortunately I do not find any example of how to define CRS


texts: How to propagate attributes down several line features?

texts:
I have a PostGIS database of road features which are split at intersections.  Some of the features on an individual road have certain attribute info (i.e. speed limits) and others do not.  I want to fill in the gaps, so that if one feature is missing a speed limit, it will automatically be assigned based on the fact that maybe four blocks south there is another feature which has the same name, roadtype, and is connected through other features to the original.

Right now all I can think of is a clunky and slow python script which will calculate predecessors and successor features, then use a bunch of loops to propagate down the road.  It feels like that would kind of be reinventing the wheel a little bit, and I am wondering if there might be any tools to help with this?  Performance is a bit of a concern because I would like to do this for cities containing over a million features, although it is an offline process so it doesnt necessarily need to be ultra fast.

I am using Python 2 and have access to ArcGIS tools.


texts:
I am fairly new in creating a web app builder, and I am trying to find a way to create a DataGrid using dojo.
In one of the in panel customized widgets, I have  a combobox that has a list of companies and I would like to style that list in a grid form that contains four columns.

the js file:

    populateList: function (results) {

        var zone;
        var values = [];
        var testVals = {};
        var xno;
        var zoneLowerCase;
        var EboID;
        //Loop through the QueryTask results and populate an array
        //console.log(results.features);

        var features = results.features;
        array.forEach(features, function (feature) {
            zone = feature.attributes.gesbez_kurz;
            zoneLowerCase = zone.toLowerCase();
            xno = feature.attributes.strXNo;
            tt = feature.attributes.strTreatyType;
            cntry = feature.attributes.UWTCountry;
            Renewal = feature.attributes.Renewal;
            EboID = feature.attributes.intEboID;

            if (((cntry == CountrySelection) || (CountrySelection == "All Countries")) &&  (Renewal == ToggleYear)) {
                values.push({
                name: zone + ' | ' + xno + ' | ' + tt + ' ['+ EboID +']',  //labelling,
                expid: feature.attributes.intExpID,
                sortName: zoneLowerCase
                });

            }
        });

    var dataItems = {
            identifier: "name",
            label: "name",
            items: values
            };

        var store = new Memory({data: dataItems, idProperty: "items"});
        var storeSorted = store.query({}, {sort: [{attribute: "sortName"},{attribute: "name"}]})
        var store = new Memory({data: storeSorted, idProperty: "items"});
        console.log(store);

        this.mySelect.set("store", store);

    }


texts:
I have built a plugin to try to read Excel files and save the chosen layer to CSV and I want to use openpyxl in QGIS, but there is a error:

ImportError: No module named openpyxl


texts: Which classification method to use when comparing two maps?

texts:
I'm making two maps of world population density. The two are using different projections, so each country's population density is a bit different in the two maps.

Which classification method should I use for these choropleth maps since my aim is to be able to compare them?

I would instinctively choose 'natural breaks' if I weren't aiming to compare the two, but I've read that using natural breaks makes it hard to compare multiple maps.


texts: Why doesn't the item properties and command history window not open up in QGIS print composer?

texts: When I open a print composer in QGIS, all I see is my map. I do not see the item properties or command history. I can add a legend, scale bar, etc., but I can't edit them. I have not had this problem in the past.


texts: Zonal Statistics of Unique Patches (Clusters)

texts:
How do you calculate the number of occurrences of a category that occur in each cluster using Google Earth Engine? My goal is to figure out the percentage (or # of pixels) of each landcover type that occur within each patch.

The script below shows a working example using a modified version of the developers guide (https://developers.google.com/earth-engine/image_objects and https://developers.google.com/earth-engine/reducers_grouping) to create clusters with an unique id and its corresponding landcover classification. However when I group a count reducer on the image the patch id is not correct.

The ultimate goal is to figure out at what point if one adds the percentage of landcover types sequentially (e.g. landcover 0 + 1 + 2...) does the proportion exceed 50% (i.e. the proportion >50% at landcover 9 for patch #19). Is there a better method to approach this problem?

//Select Region
var roi = /* color: #98ff00 */ee.Geometry.Polygon(
        [[[-122.29786986020576, 37.4658267849554],
          [-122.07007521298897, 37.46092160855244],
          [-122.0692169061042, 37.539093572188996],
          [-122.29889982846748, 37.53882133810604]]]);

// Load MODIS land cover categories in 2001.
var landcover = ee.Image('MODIS/051/MCD12Q1/2001_01_01')
  // Select the IGBP classification band.
  .select('Land_Cover_Type_1');

Map.addLayer(landcover, {}, 'landcover');

// Load a Landsat 8 image
var image = ee.Image('LANDSAT/LC08/C01/T1_TOA/LC08_044034_20140318');
Map.setCenter(-122.1899, 37.5010, 13); // SF Bay

// Threshold the thermal band to find "hot" objects.
var hotspots = image.select('B10').gt(300);

// Mask "cold" pixels.
hotspots = hotspots.updateMask(hotspots);

// Uniquely label the patches
var patchid = hotspots.connectedComponents(ee.Kernel.plus(1), 256);

//Add landcover band to patches and display
var img = patchid.addBands(landcover);
Map.addLayer(img.randomVisualizer(), {}, 'patch data');

// Grouped a count reducer: count number of land cover category pixels by patches
var count = img.select(['labels','Land_Cover_Type_1','B10']).reduceRegion({
  reducer: ee.Reducer.count().group({
    groupField: 0,
    groupName: 'patch id',
  }).group({
    groupField: 1,
    groupName: 'landcover',
  }),
  geometry: roi,
  scale: 30,
  maxPixels: 1e8
});

print(count);


texts: ArcGIS JavaScript API: SEC7127: Redirect was blocked for CORS request after submitting geoprocessing job

texts:
Using the ArcGIS for JavaScript API 3.19, I'm submitting a job to the ExtractData geoprocessing service at:
https://analysis7.arcgis.com/arcgis/rest/services/tasks/GPServer/ExtractData/submitJob
The job is submitted and run successfully, and the output is created in my ArcGIS Online root folder.
However, I can't figure out how to get the response back from the server. When I submit the job, the JavaScript console reports 2 errors:

SEC7127: Redirect was blocked for CORS request.
SCRIPT7002: XMLHttpRequest: Network Error 0x800c0014, A redirection problem occurred.

This sort of makes sense to me, because the server issues a redirect to the job page. I can see that the server responds with a status of 303, and includes the URL to the job in the response. However, the code that I have to handle the response (see below) is never called.
var url = _extractDataURL +
    "?inputLayers=" + inputLayers +
    "&extent=" + extent +
    "&clip=" + clip +
     "&dataFormat=" + dataFormat +
     "&outputName=" + outputName +
     "&context=" + context +
     "&f=" + f +
     "&token=" + IdentityManager.credentials[0].token;

var xmlHTTPRequest = new XMLHttpRequest();
xmlHTTPRequest.onreadystatechange = function () {
    if (xmlHTTPRequest.readyState == 4 && xmlHTTPRequest.status == 200) {
        console.log(this.responseText);
    }
}

xmlHTTPRequest.onload = function () {
    console.log(this.responseText);
}

xmlHTTPRequest.open("GET", url, true);
xmlHTTPRequest.send();

How can I configure things so that the redirect (the correct response) doesn't cause the CORS error, and so that I can handle the response and get the URL to the job?


texts: QGIS - QTableView or QTableWidget in custom form

texts:
I've made a custom ui form for entering data in a GIS database built with PostgreSQL and QGIS. This is a relational database and I use the Relations properties in QGIS to implement this relations (one-to-many or many-to-many relations).
Currently, my many to many relations are shown in the custom form via a QWidget with a dynamic property set to qgisRelation. This display the form of the child table in my parent table.
However, I would like to have, for certain relations, a much simpler way of displaying the child table, i.e. without the sidebar and the standard buttons (add/delete child, link/unlink child feature etc.) that comes from QGIS. I'd rather have a simple table that shows just one field (the ID) of my child table with the possibility to populate the table.
Is it possible to this with QTableView or QTableWidget? Do I have to set dynamic properties and maybe add some Python code in init function to make it work?


texts:
I would like to remove isolated pixels (or clumps of pixels <9) from my raster image.

library(raster)
# create some raster data
r <- raster(ncols=12, nrows=12)
set.seed(0)
r[] <- round(runif(ncell(r))*0.7 )
rc <- clump(r)

Before removing clumps of pixels <9

After removing clumps of pixels <9 :

In Erdas, there is the Sieve tool to do this, but how to replace it in R?


texts: QGIS display does not reflect projection

texts:
I have slight experience with ArcGIS, but am new to QGIS.

I'm trying to set the projection of my map, and I understand that this is done by altering the CRS (coordinate reference system) of each layer and/or the entire project.

When I right click a layer and select "Change Layer CRS", however, the visual representation of the map (the distortion) does not reflect the selected projection. My world map looks like it has a Robinson (or similar) projection, even if I select a Mercator CRS. The EPSG code in the lower right corner changes, as do the coordinates.

Setting the project CRS to the same projection does not help, nor does enabling on-the-fly projection.

My goal is to export the map to SVG, so the visual distortions are important.

Am I doing something wrong?


texts: Creating vector tiles from OpenStreetMap?

texts:
I need to understand the process of creating Mapbox Vector Tiles from OSM data in some detail -- in particular I need to fiddle with cutting out a lot of data (esp. buildings) to keep the final vector tile sizes down.

So far I have found two repositories that allow the creation of Mapbox Vector Tiles from OSM data:

Mapzen: https://github.com/tilezen/vector-datasource
osm2vectortiles: https://github.com/osm2vectortiles/osm2vectortiles

Is there  any information on the differences between these which is likely to be more tractable?


texts: I am new to GIS application world. I would like to ask what are the other alternatives to UMN MapServer  in the market preferably open source?


texts: How to determine if a cluster is at its maxZoom level?

texts:
I'm using the Leaflet MarkerCluster plugin, and I want to know whether a particular cluster represents points at the same location, or whether the cluster would break into smaller components if I zoomed in.

By default, the cluster will "spiderfy" if it represents multiple points at the same location - I'd like to detect this scenario and handle the individual points myself, rather than spiderfying (using spiderfyOnMaxZoom: false).

There is an event fired when you click on a cluster icon. From the documentation:

markers.on('clusterclick', function (a) {
    // a.layer is actually a cluster
    console.log('cluster ' + a.layer.getAllChildMarkers().length);
});

Is it possible to determine whether the cluster is at its maximum zoom level after the user has clicked on it? That is, if the user clicked on it again, would it continue to zoom in further, or does it represent points at the same location?


texts: What is the difference between Intersect & Overlap in ArcGIS Server?

texts:
When you run a Query Task on the ArcGIS Server's REST Endpoint, you need to give a spatialRel

There are two spatial relationships which are quite similar: esriSpatialRelIntersects & esriSpatialRelOverlaps

The JavaScript API documentation has this to say about them:

  SPATIAL_REL_INTERSECTS    Part of a feature from feature class 1 is
  contained in a feature from feature class 2.

  SPATIAL_REL_OVERLAPS  Features from feature class 1 overlap features in
  feature class 2.

If one is querying a Polygon featureclass by using a Polygon as an Input, what is the difference between the two?


texts: How to open QGIS project without being asked for PostGIS credentials every time

texts:
Problem:
We have several layers stored in a PostgreSQL database with PostGIS extension (Stored on AWS RDS Instance). The connection has been working perfectly to date.
I get asked for login credentials every time I open the QGIS 3.16.5 project which is stored locally.
What I'd like to happen:
I'd like to open the QGIS project and have the layers loaded from PostgreSQL+ PostGIS database without needing to authenticate every time.
Why it's important:

Saves time every time we open a project
Multiple database connections means entering credentials for every database connection
I'd like to know how to improve the UX for our users

What I have tried:

Stored credentials as plain text in the project file Edit PostGIS Connection > Authentication > Basic > Check Store checkbox
Converted to configuration. Edit PostGIS Connection > Authentication > Configurations > Create Authentication Configuration > Fill in Authentication Dialog inputs

Resources I've looked at already:

https://issues.qgis.org/issues/21582

From QGIS 3.4.5 and 3.6 asks password repeadly to connect to PostgreSQL 10

https://docs.qgis.org/3.16/en/docs/user_manual/auth_system/index.html?highlight=authentication

Change in PostGIS database password prompts for user credentials on QGIS start

QGIS asks for PostgreSQL credentials multiple times

Questions:

Is this a known bug with the 3.16-LTR, if so where can I find out more?
Am I managing the connections incorrectly or overlooking something?

Environment:
OS - Windows 10 Professional
QGIS Version - 3.16.5-Hannover


texts: Accessing QGIS Interface in standalone PyQGIS script?

texts: I try to set the extent of my QgsProject to the extent of one of the layers in it. I know how to do that with qgis.utils.iface but it works only when the script runs inside QGIS' Python console. What I need is a way to modify the interface without being in QGIS.


texts: Geometry-generated lines not showing on QGIS Cloud

texts:
I have a points layer (tested as both a delimted text file and as a point shapefile) which is styled as shown below where the lines were created using the Geometry Generator. Please ignore the elongated legend :)

When I upload the data onto QGIS Cloud, the generated lines disappear in both QGIS and from the online webmap even though they are clearly shown in the legend and in the symbology properties:

Symbology properties:

Is the symbology for the generated geometry not supported in QGIS Cloud? Is there a possible workaround?

Using QGIS 2.18.17 for Win 7 64-bit.


texts: Count points belonging to a country but outside of country polygon (coastal)

texts:
I'm still getting my head around GIS and am using QGIS 3.6.2.

I'd like to do a simple count of points within countries. The problem is that many of my points are located on coastlines or islands, so my current counts come up short because they do not fall within the defined polygons of the country world map. I guess this may be a resolution issue of the countries shape file?

Is there a way to add a buffer zone of a couple of kilometers to coastal polygon borders - but not to land borders? Or is there a world countries shape file that reflects accuracy of google maps or similar?

I've attached a screenshot of an example of two points in the Floria keys as seen by QGIS and by Google Maps.


texts: Handling nodata values in raster in Python

texts:
I would like to set value of raster -999.0 as NaN in crop_image function,  while reading/cropping file:

for band_path in range_bands:
    # open raster as numpy array by rater.io in content manager(open/close file)
    with rio.open(band_path) as band_src:
        if "band4" in band_path:
            print("band4")
            ## use crop.image() to read and crop
            band4, band4_meta = es.crop_image(band_src, bound)
            print(band4_meta)
            # Update with NaN data values
            band4_meta.update({'nodata': -999.0})
            print(band4_meta)

How is possible to do?


texts: Using existing styling rules in provider feature filter for 'pre-filtering' data from PostGIS DB

texts:
I have a PostGIS DB populated with OSM data (osm2pgsql). I defined styles for several features. Some are styled, others are not.
But  every time I update the canvas or load my project, all of the features from the DB are loaded  - also the features I don't styled (because they are not important for my task). I'm wondering, If I could save CPU and RAM by pre-filter the features I don't need. I already do that in my line-layer for rivers. The feature filter is simply:

"waterway" IN ('river', 'stream', 'canal')

Just the waterways are loaded and all others, like railways or streets, are filtered out.
It's easy in that waterway-case, because the work is done with that small query. But for my polygons and points I have got a way much more complex filter:

Now I would like to use my styling-rules to pre-filter the DB-Import. I know, I could copy and paste... but I need a more comfortable and flexible way (if I change something in the styling and add or delete a feature-class, I have to do the same in the feature filter).


texts: Adding Button as Icon and Text to toolbar in ArcMap

texts:
How do you add a BaseCommand to a custom BaseToolbar as an icon and text?

This is not an Add-In, it's the old way of doing buttons and toolbars even though it's version 10 of ArcMap.


texts: Porting Avenue code for Producing Building Shadows to ArcPy/Python for ArcGIS Desktop?

texts:
Whuber provided an answer at Producing building shadows using ArcGIS Desktop? which required using Avenue code.

Any idea how to make it work in ArcGIS Desktop 10?


texts: Library for parsing scientific PDS files?

texts:
Does a library exist for parsing the PDS files that scientific planetary data gets shipped as?

For example, here's a link to to the lunar topographic data from the Lunar Reconnaisance spacecraft (it's huge):

http://wms.lroc.asu.edu/lroc/global_product/256_ppd_DEM

I've seen the same format used elsewhere (for example, similar data from the Japanese Kaguya spacecraft).

Previously I've written custom parsers to massage the data into something I can actually use, but this is getting a bit complicated --- the LROC dataset has tiles of two different types of projection. What I'd really like is a library that I can just throw the PDS files at and say 'what is the data point for lat Y lon X' and have it give me an answer. Does such a thing exist?


texts: Unable to use gdal_merge.py from Python Script

texts:
If I execute this script in windows with x86 Python, works just fine. Also from command line directly.

import os, subprocess

gm = os.path.join('C:\\', 'Program Files (x86)', 'GDAL', 'gdal_merge.py')

merge_command = ["python", gm, "-o", "C:\\temp\\merge2.asc", "C:\\temp\\507.asc", "C:\\temp\\508.asc"]

subprocess.call(merge_command, shell=True)

However, in another machine, with windows and python x64, command line tool works, but the python script does not (it only differs from the previous one in the program files directory where gdal_merge.py is located).

import os, subprocess

gm = os.path.join('C:\\', 'Program Files', 'GDAL', 'gdal_merge.py')

merge_command = ["python", gm, "-o", "C:\\temp\\merge2.asc", "C:\\temp\\507.asc", "C:\\temp\\508.asc"]

subprocess.call(merge_command, shell=True)

What could be happening? Can anyone help?


texts:
I have a list of 10 000 points, and I wish consolidate 15 points each in a block in the vicinity.
I use arcgis, or mapinfo qgis, I just need the method.


texts: Count the area of rasterdata in specific shapefile polygons (residential density in different districts)

texts:
I'd like to count the exact area of three different residential classes (low, medium, high) by district of a shapefile which contains 225 districts. like: district 1 has: XXX sqm low, XXX sqm medium, XXX sqm high. distict 2 has: XXX....

the residential classes are a file system raster, the districs are one shapefile.

I can count the total area of the residential classes of the whole study area. but I can't count the area for every unique district I have in the shapefile. (Maybe I could count the area for 1-3 specific districts I pick but for 225 it would be too time-consuming) After all I'd like to see the area of each district in the attribute table of the district-shapefile and be able to export the data to Excel. I hope you can help, thanks.


texts: How to update Portal for ArcGIS system properties using ArcGIS API for Python

texts:
I'd like to update the Portal system properties using the Python API. The documentation seems to indicate that it can be done. However, making an assignment to System.properties doesn't appear to udpate the system properties.

from arcgis.gis import GIS

gis = GIS(gis_base_url, user, password)
system_properties = gis.admin.system.properties
system_properties['privatePortalURL'] = system_properties['WebContextURL']
gis.admin.system.properties = system_properties

So, is there a way to update the Portal's system properties in a way other than what's documented?


texts: Using raster calculator to show multiples of particular value?

texts:
How do I use raster calculator so that the new raster only shows multiples of 128?

I'm not sure what the expression should be.


texts: How to add fields to postgis layer loaded in qgis using python

texts:
I have loaded a layer from Postgis in qgis and I am adding the fields to the layer using the follo code snippet:

canvas = qgis.utils.iface.mapCanvas()
aLayer = canvas.currentLayer()
field_name = "test"

aLayer.startEditing()
if not field_name in fieldnames:
      res = aLayer.dataProvider().addAttributes([QgsField(field_name, QVariant.String)])
aLayer.updateFields()
aLayer.commitChanges()

However, it is not adding the field to the layer. What must I do if I want to add the fields to the layer?


texts: Converting list of stack images to time series files in R?

texts:
I am using bfast package for time series analysis. But it requires the files to be read as time series.

How to convert list of stacked images from 2000 to 2013 into a time series and given as input for analysis?


texts:
I've explored this forum fairly extensively trying to get answers as to why my SQL isn't working. I need help with the timestamp on a torque map. I've used the following formula, obviously changing the three variables, and it always comes up with an error message.

  UPDATE tablename set yearColumn = to_timestamp(yr_compl,'YYYY’)

My old column is a string and my new column is in date.

I'm a complete novice.