texts: Weighted overlay output raster does not match the extent of my input rasters

texts:
I am having trouble understanding why the resulting extent of my raster produced by my weighted overlay doesn't not cover the extent of the inputs I used. All my inputs are in the same coordinate system and using the same datum.

The environment settings has my processing extent set to the grey input raster, as does the raster analysis cell size. All my other inputs have the same extent/outline as the grey raster shown in the image, so I am wondering why the weighted overlay output does not cover the entire area?


texts: Performing a Python scripted tool based on whether a field is populated or not

texts:
My question is if I'm running the ArcMap Dissolve tool through a Python script, how can I tell it to check for attributes in one field, and if that field is empty, to perform the tool on another field.

To be more specific, I would like the Dissolve tool to check for route numbers in a field called [rt_shrt_nm] (which are provided when the route name is a numerical value).  But if that field is empty (which is sometimes is if the field has textual route names) to perform the tool on the [rt_long_nm] field instead.

This is what I'm working with now, just to give you a better idea, even though I know its probably wildly wrong:

if [rt_shrt_nm] is null:
        arcpy.Dissolve_management(OutShapesFCname, outGDB, ["rt_long_nm"], "", "", "")
elif:
        arcpy.Dissolve_management(OutShapesFCname, outGDB, ["rt_shrt_nm"], "", "", "")


texts: ArcObjects c# .Net Merge with Geoprocessing?

texts:
I'm trying to run the Merge tool using ArcObjects c# .Net. I'm fairly new to the skill but I have been able to run other tools such as Buffer and Points to Line.
However, I have met my match with Merge because it requires more than one input (two layers). Would anyone have an idea how I'm suppose to put two layers into "MergeTool.inputs"? I think right now, my Layer1 is being replaced by Layer2 so when the tool runs, I end up with just one layer being merged into nothing and I get the same layer returned.

IMxDocument pMxDoc = ArcMap.Application.Document as IMxDocument;
        IMap map = pMxDoc.FocusMap;

        Geoprocessor GPMerge2 = new Geoprocessor();

        ESRI.ArcGIS.DataManagementTools.Merge MergeTool = new ESRI.ArcGIS.DataManagementTools.Merge();

        MergeTool.inputs = Layer1;
        MergeTool.inputs = Layer2;
        MergeTool.output = OutputName;

        GPMerge2.Execute(MergeTool, null);

Edit: I've seen some solutions to use semicolon and quotes. Something like
MergeTool.inputs = "Layer1;Layer2"
But I'm actually using some variables and I can't get the exact path to the layers to make it into a string. The difference can be seen below if I try to do that. I also want to try and keep this set of code if possible and not set up parameter array.


texts:
I am using the Decarta JS API.  I am displaying a map with layers.  I want to refresh the map and redraw new layers, but I see no way to do this in the documentation.

Is there a way to refresh the map using the API?


texts:
I am trying to split an irregular polygon into 120 equal area blocks. I also need a buffer around each new block within the polygon.

Is this possible in QGIS 2.12.1?


texts: QGIS: Vector Layer Split: Renaming of the end files based on the value of the Unique ID

texts:
I have a model that convert fields and does a vector layer split, to split the CSV files into multiple files for easy geo-reoference (these are addresses). The thing is the unique value is of course an Integer, so I would love it that my split vector layers and rename them to the value of unique ID. This column is called "Blok" in my Vector/CSV file. Because I get seriously annoyed that it goes like:

Blok 0
Blok 1
etc.

While Blok 0 has all the attributes who 26 as value. It makes it kind of hard to keep track which files contains which values.
So how do I do this using QGIS 3.14.15 by the way?


texts: QGIS WMS ver 1.3.0 Exception with GeoServer

texts:
I'm trying to load WMS layers from a local GeoServer into QGIS so that I can edit them. I can access and edit the same layers if I access them directly from a PostGIS database. However, trying to load them as WMS layers in QGIS fails with an XML Parsing Exception:

 <<?xml version="1.0" encoding="UTF-8"?><ServiceExceptionReport version="1.3.0"    xmlns="http://www.opengis.net/ogc" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"   xsi:schemaLocation="http://www.opengis.net/ogc http://localhost:8080/geoserver/schemas /wms/1.3.0/exceptions_1_3_0.xsd">   <ServiceException>
  javax.xml.transform.TransformerException: Translator error
Translator error
Error occurred trying to write out metadata for layer: Pk50095
Unable to acquire a reader for this coverage with format: WorldImage
</ServiceException></ServiceExceptionReport>
Tried URL: http://localhost:8080/geoserver/wms?SERVICE=WMS&REQUEST=GetCapabilities

In addition, when I try to access the GetCapabilites URL directly, an XML Parsing error is generated for WMS ver 1.3.0. However, version 1.1.1 GetCapabilities URL executes fine.

Any idea on what may be wrong with my server installation?


texts: Highlighting selections using ArcObjects?

texts:
I have a problem you want help:

In ArcMap 10.2.2

In my code I want to highlight selections but I can't.

My code:

void drawPoint(IFeature pfeature)
{
  IRgbColor rgbColor = new RgbColorClass();
  rgbColor.Red = 255;
  rgbColor.Blue = 0;
  rgbColor.Green = 0;

  IMarkerSymbol simpleMarkerSymbol = new SimpleMarkerSymbol();
  simpleMarkerSymbol.Color = rgbColor;// color RED
  simpleMarkerSymbol.Size = 15;
  IMarkerElement markerElement = null;

  markerElement = new MarkerElementClass();
  markerElement.Symbol = simpleMarkerSymbol;
  IGraphicsContainer graphicsContainer =
  GeneralVariableMap.mapControl.Map as IGraphicsContainer;
  IElement element = null;
  element = (IElement)markerElement;
  element.Geometry = pfeature.Shape;
  graphicsContainer.AddElement(element, 0);
}


texts: Calculating distance from point to points along roads using ArcGIS Network Analyst?

texts:
i have start poin(living city) and 10 other points and road network, I want calculate distance from my start point (my living city) to each other city along the road network. **In the result i want 10 distances for every point.  **  I dont know which analysis use for it.

I am using ArcGIS 10.4 with all extensions.


texts:
I'm trying to create KML files manually (XML format) by extracting data from a postGIS database with "geography" values in a column.

But i believe i need to transform the projection of each row in the database from EPSG:4326 to EPSG: 900913.. How can i do that?

Code:

 map = new OpenLayers.Map({
maxExtent: new OpenLayers.Bounds(-20037508.34,-20037508.34,20037508.34,20037508.34),
numZoomLevels:16,
zoom: 0,
maxResolution:156543.0339,
units:'m',
projection: "EPSG:900913",
displayProjection: new OpenLayers.Projection("EPSG:4326"),
div: "map_canvas",
controls: [
                new OpenLayers.Control.Navigation(),
                new OpenLayers.Control.PanZoomBar(),
                new OpenLayers.Control.MousePosition(),
                new OpenLayers.Control.LayerSwitcher()
        ],
});

var KML = new OpenLayers.Layer.Vector("KML", {
        strategies: [new OpenLayers.Strategy.Fixed()],
        protocol: new OpenLayers.Protocol.HTTP({
            url: "assets/text.kml",
            format: new OpenLayers.Format.KML({
                extractStyles: false,
                extractAttributes: true

            })
        })
    });

Coordinate system in database: (one "Area") / Polygon.

  -179.999995473772 -89.9999977368861,-89.9999954737723 -89.9999977368861,4.52622771263123e-    06 -89.9999977368861,90.0000045262277 -89.9999977368861,179.999995473772     -89.9999977368861,179.999995473772 89.9999977368861,90.0000045262277     89.9999977368861,4.52622771263123e-06 89.9999977368861,-89.9999954737723     89.9999977368861,-179.999995473772 89.9999977368861,-179.999995473772 -89.9999977368861

Current error:

  Uncaught Bad LineString point coordinates: -179.999995473772


texts:
I want to sum the areas of several polygons in the same shape layer by using the Aggregate function in Select features using an expression. However, I can't find the Aggregate menu nor the SUM function within it. Where do I find it? I think I have installed most plug-ins connected to this function.

2.14.15 on a mac


texts: QGIS - PYQGIS Import gpx as vector layer

texts:
I'm trying to use the python console to import some gpx files as vector layers.
When I do it manually I have to choose which layers to add, my options are:

I've tried various versions of the code from the manual:

uri = "path/to/gpx/file.gpx?type=track"
vlayer = QgsVectorLayer(uri, "layer name you like", "gpx")

I also tried adding the layer name I want like:

path/to/gpx/file.gpx|layername=track_points?type=track", "layer name you like", "gpx")

Edit:
On each try I don't get any errors, it just doesn't import anything, ie no new layers in the layers panel. I would like to import the files using python as I have quite a large number of them. I could do it manually but I would have to choose the track_points layer for each file

Edit:
From Joseph's answer I adapted his code and got what I wanted, not entirely sure how but if it helps anyone else:

import os
path = 'Path/To/GPX Files/'
names = ["track_points"]
for dirpath, subdirs, files in os.walk(path):
    for f in files:
    layername = f[:-4]
        for name in names:
            iface.addVectorLayer(os.path.join(dirpath, f)+"?type="+name, layername, "gpx")


texts: Advanced Python Field Calculation gives error: TypeError: 'in ' requires string as left operand, not QVariant

texts:
I am using the following function in QGIS 3.10 Advanced Python Field Calculator as part of a Processing Modelling file. The model processed 4 layers and generates a single layer for use in another program.  I want to rerun the Process Model frequently as the base data changes regularly. I checked the function in Jupyter Notebook and it ran fine (note the indents in the version below are out):
def ADD(zone,area):
    if zone in "R1":
        return 0.7
    elif zone in ["R2","R5","RU1","RU2","RU4","RU5"]:
        return 0.96
    elif zone=="R1":
        return 0.4
    elif zone=="IN2":
        return 11.5*area
    elif zone=="IN2":
        return 26.5*area
    elif zone in ["IN3"]:
        return 68.5*area
    elif zone in ["B1","B3","B4","B5","B6"]:
        return 11.5*area
    else:
        return "0.0"

Where "zone" and "area" are from the input attribute table and are string and float types respectively.  I am expecting the returned value to be float also.  However, I receive the following error:
Traceback (most recent call last):
File "C:/PROGRA~1/QGIS3~1.10/apps/qgis/./python/plugins\processing\algs\qgis\FieldPyculator.py", line 176, in processAlgorithm
exec(bytecode, new_ns)
File "<string>", line 1, in <module>
File "<string>", line 2, in ADD
TypeError: 'in <string>' requires string as left operand, not QVariant
Error encountered while running Determine ADD flows in KLD day
Error encountered while running Determine ADD flows in KLD day
Execution failed after 25.19 seconds

Strangely the Processing model terminates, but it looks like it completes the calculation correctly and adds in the new attribute as intended.  I actually have two functions in a row, both in similar structure, and the second one won't run after the model crashes.  I tried swapping the order and the same thing occurs, it crashes after the first one.  I need to add in quite a few more python functions so really need it to not crash!
I am not an expert in either Python or QGIS, and I am wondering if I have an obvious error anywhere or something I am missing.


texts: Web-GIS cannot display WMS data because of CORS header ‘Access-Control-Allow-Origin’ missing

texts:
I have an OpenLayers application in a secured connection (with apache2 and port 443) in a machine A and I want to display data from a GeoServer WMS which it exists in a machine B with insecured connection (with Tomcat 8 and port 8080).
But while WMS server is CORS enabled (it serves the same layers in other web-GIS applications),  browser console raises the popular error: CORS header ‘Access-Control-Allow-Origin’ missing.
The weird thing is that the corresponding getLegend service of each layer works fine.
Here is my code for layer definition:
var serverNameType = 'geoserver';
var serverLink = 'http://<domain>:8080/geoserver/wms';

var tmpInfo = new Array(
              {
                        layerName: "<workspace:<layername>",
                        altLayerName: "<workspace:<layername>",
                        hasIdentifyInfo: true,
                        crossOrigin: 'anonymous',
                        legendName: getTranslation(lang, "Layer title"),
                        legendSrc: serverLink + "?REQUEST=GetLegendGraphic&VERSION=1.0.0&FORMAT=image/png&WIDTH=20&HEIGHT=20&LAYER=<workspace:<layername>&legend_options=fontName:Times%20New%20Roman;fontAntiAliasing:true;fontColor:0x000033;fontSize:14;bgColor:0xf0f8fb;dpi:180"
              }
      );

for(var i = 0; i < LayerControl.layerControls.length; ++i)
      {
                LayerControl.layerControls[i].addLayer(
                        new Layer(
                                tmpInfo[i].layerName,
                                new ol.layer.Tile({
                                          source: new ol.source.TileWMS({
                                                    preload: Infinity,
                                                    url: serverLink,
                                                    serverType: serverNameType,
                                                    crossOrigin: "anonymous",
                                                    params: {
                                                              'LAYERS': tmpInfo[i].layerName, 'TILED': true
                                                    }
                                          }),
                                          visible: false
                                }),
                                tmpInfo[i].hasIdentifyInfo, tmpInfo[i].legendName, tmpInfo[i].legendSrc
                                )
                        );

Why is this happening ? Is it possible to fix it in front-end level ?


texts: SQL query in Python console to sum column based on values of other two columns

texts:
We have 100's of shapefiles for which we have to generate excel reports based on a few (combination of) parameters in the attribute table. I am using SQL query to achieve this and when I run the query in the "Create Virtual layer" process, it works fine and the output has been generated, which will be saved as an excel file. Instead of this process, I am trying to run this through the Python console. Though the code didn't show any error, the output is not being printed. Do I have to change the SQL query format in the code?
Image 1: For the shapefile final_bdry, Total_Area has to be calculated based on two other columns - MS_Area, Ground_Con (similar to Pivot_Table in excel, the sum of a column based on a unique count of the other two columns). No issues with the query and I am getting the desired output in the virtual_layer.
.
Image 2: But, if I run the same query within the Python console, the result is "None"

Ref: Used this as ref SQL queries in shapefiles using python/gdal return none result


texts: Where is Copy As Python Snippet of ArcGIS Pro?

texts:
When writing ArcPy code to use  with ArcGIS Desktop 10.x applications like ArcMap I will often run a Geoprocessing tool and then go to the Geoprocessing | Results window to right-click on that successful run to Copy As Python Snippet so that I can Paste that code into my script.

When I look for the same option in ArcGIS Pro there is no Results window to start from.

What is the equivalent workflow for ArcGIS Pro?


texts: Problem with my google map, place search with a text bar

texts:
I am doing a small project with google maps api, I am trying to locate all the interesting services around me. For example, i want to search "restaurant" then it will shows all the restaurant near me within my radius.

Javascript framework that I am using

Ionic framework
Angular JS
Google maps API

Here's the link, it won't shows the result since either you have to run it on localhost or have an emulator to run it.

http://jsfiddle.net/#&togetherjs=DhsmJ98GNe


texts: Exporting map from composer in QGIS - text labels have moved in output?

texts:
I am exporting a map from QGIS.  I have some text labels for the title of my map, rendered with html.

All looks good before I export but when I open my exported map (I have tried exporting to both PDF and image), the map title has completely moved to the side and is partly off the page altogether.  The picture below shows the output for the bottom of my map where the title is and you can see it is cut off. It is supposed to say ONEROA. and the other lines are missing letters as well.


texts: Using different fields for symbology and legend labeling and maintaining structure using ArcPy?

texts:
I have multiple fields in my attribute table so I want to use one general layer and apply to many fields (same range). The problem is everytime I apply the layer or it`s get the old field and apply to my map or apply to new field but do not get the settings.

I can get the map with the correct range but always do not came with the legend format (in percentage (%)).

import arcpy
import csv
import pythonaddins

TabelaP1 = arcpy.GetParameterAsText(0)
Ano = arcpy.GetParameterAsText(1)
IndexL = int(arcpy.GetParameterAsText(2))
SiglaProduto = arcpy.GetParameterAsText(3)
mxd = arcpy.mapping.MapDocument('CURRENT')
df = arcpy.mapping.ListDataFrames(mxd, "")[0]
lyr = arcpy.mapping.ListLayers(mxd, "", df)[IndexL]
Layer1 = arcpy.mapping.Layer(r'F:\Agro\Layers\PercAPH2.lyr')
arcpy.mapping.UpdateLayer(df, lyr, Layer1, True)

arcpy.AddMessage('Calculando Porcentagem de Area Plantada')

if '2002' < Ano < '2016':

    arcpy.AddMessage('Iniciando Calculo para o Ano {0}'.format(Ano))

    NovoCampo1 = 'PercAPH_'+SiglaProduto+'_'+Ano
    Campo1 = 'APH_'+SiglaProduto+'_'+Ano
    arcpy.AddMessage('Campo selecionado: {0}'.format(Campo1))
    arcpy.AddField_management(TabelaP1, NovoCampo1, 'Float')
    arcpy.CalculateField_management(TabelaP1, NovoCampo1, '(!{0}!/!AreaTotalHec!)*100'.format(Campo1), "PYTHON_9.3")
    lyr.symbologyType == "GRADUATED_COLORS"

    lyr.symbology.valueField = NovoCampo1
    arcpy.RefreshActiveView()
    arcpy.RefreshTOC()

if Ano < '2003' or Ano > '2015':

    arcpy.AddMessage('O Ano {0} informado nao possui Serie Historica disponivel'.format(Ano))

arcpy.RefreshActiveView()
arcpy.RefreshTOC()

[

The original layer give me the label in % with round(2) for the numbers. I tried many ways, including make by ArcCatolog for the GDB but also do not work. There is any possible way to solve this?


texts: Connecting GeoJSON on GeoServer to Carto

texts:
GeoServer can output GeoJSON via WFS.

Carto can read GeoJSON via URL.

Should it then be possible to connect a GeoJSON layer in Carto from a GeoServer WFS?

I'm trying this setup, but I'm getting a cryptic "unknown error" in Carto when trying to connect to GeoJSON output.

If I save the same output to file and upload it to Carto, it loads fine. So the problem seems to be in the connection rather than the data.

[Edit] Here is a link to the relevant GeoJSON: GeoJSON output from GeoServer


texts:
I have been working through the ASF's "Recipes" for the creation of a DSM from Sentinel 1 (A and B) imagery.

I am able to produce a DSM, and am relatively happy with the results in terms of the surface produced - it appears to follow the trends of my known datasets (hills where there should be hills and valleys where there are valleys).

I am however seeing discrepancies between my derived surface and a known point (derived from LiDAR). I generated a DSM using imagery acquired in the same month (a few days either side to improve the Baseline and remove temporal discrepancies) as when my LiDAR dataset was acquired.

I am not going to dispute that the LiDAR Elevation is going to give me a better Z-Value than a Satellite derived DSM. I am however looking to use Sentinel in areas where there is not a decent surface or topographic dataset I can use.

I would like to understand why there is a difference (consistently around 120m) between my DSM and LiDAR surfaces.

How can I improve the vertical accuracy of my DSM? Is the issue a vertical datum issue or perhaps a setting / attribute I should be looking for in the acquisition of the imagery?


texts: Summing attributes intersected by under a radius

texts:
Apologies - new user.

I am working with QGIS on Mac OS X and with shapefiles downloaded from the US Census Bureau that show population data by state and by county.

How do I:

select certain points (e.g. cities)
draw a defined distance radius around those points (e.g. 500km - At the moment I can only get the units to be in degrees)
Sum the populations in the counties that are intersected by the circle?

Any assistance is greatly appreciated!


texts: Leaflet.js displaying incomplete base map

texts:
I am working on a tile server with a leaflet.js UI. I have successfully created the tiles for my base map and the tiles for the overlay. The result can be seen here (zoom in further to see overlay). Unfortunately, as can be observed, the base map is incomplete and is wrapping in a weird way. This is not the result of a fault tiling process I believe as I am using the same tile hierarchy here, using openlayers.js. As you can see, in that UI, the base map is complete.

This is my current leaflet interface:

    <title>Mapea - Mapeamento Meoteorológica Inteligente</title>
    <link rel="stylesheet" href="style/leaflet.css" />
    <link rel="stylesheet" href="style/map.css" />
    <link rel="stylesheet" href="style/L.Control.MousePosition.css"  />

Although I doubt this has to anything to do with the tiling process, as I wrote the following based on this [google tutorial][3], here is my tiling script.*

Sorry, I could not post the link directly, I'm new to this stackexchange.

#!/bin/bash
##
## tiler: creates tiled hierarchy from an ungeoreferenced image.
## usage: tiler <input_image> <N> <E> <S> <W>
##

prefix=$(basename ${1%%.*})
N=$2
E=$3
S=$4
W=$5

# First, recover the image's dimensions.
IMG_CHARS=`identify "$1" | cut -f 3 -d' '`
xsize=`echo $IMG_CHARS | cut -d'x' -f 1`
ysize=`echo $IMG_CHARS | cut -d'x' -f 2`

# Insert georeferencing points, prepare for map projection.
gdal_translate -of VRT -a_srs EPSG:4326 -gcp 0 0 $W $N -gcp $xsize 0 $E $N -gcp $xsize $ysize $E $S $1 ${prefix}.vrt

# Projection itself.
gdalwarp -of VRT -t_srs EPSG:4326 ${prefix}.vrt ${prefix}.vrt

# Create the tiles according to rules above.
gdal2tiles.py -t $prefix -p geodetic ${prefix}.vrt

UPDATE - I have followed the recommendations of the solutions below, but to no avail (see comments). Still, @Sorin is correct in his analysis.

This is the tree output of the first few levels of my base map tile hierarchy.

natural
├── 0
│   ├── 0
│   │   └── 0.png
│   └── 1
│       └── 0.png
├── 1
│   ├── 0
│   │   ├── 0.png
│   │   └── 1.png
│   ├── 1
│   │   ├── 0.png
│   │   └── 1.png
│   ├── 2
│   │   ├── 0.png
│   │   └── 1.png
│   └── 3
│       ├── 0.png
│       ├── 1.kml
.
.
.

As @Sorin said, the folders seem to have something wrong. For example, zoom 0 should only have one folder, the 0 folder, which only has one file, the 0.png file. Instead it has a 1 folder, which itself has another 0.png. Inspecting these files, I have come upon the other missing half of the map (0/1/0.png).

link [3]: developers.google.com/kml/articles/raster


texts: Getting the distances from multiple points to multiple lines in QGIS

texts:
I am looking for advice on how to tackle this problem. I have many points and I want to find a list of the nearest lines to each of those points, in order, and export it to excel.

Here is a description of my inputs:

5 lines (from KMZ files) within a few meters of eachother, sometimes
crossing, many miles long
QGIS layer with many points, each point on one of the aforementioned lines

Desired output:

5 additional columns in the attribute table, each noting the nearest line, and if possible, the perpendicular distance to that line.

Things I have considered, and issues:

Drawing circles with different radii on each point. Then trying an
overlay analysis.
Nearest Neighbors method, but I am new to QGIS and Python and cant
figure out how to do this for multiple lines and points. QGIS only lets me (as far as I can tell) get the first nearest neighbor, whereas I need multiple


texts: Join multiple spatial tables by distance

texts:
I'm trying to join two tables, one which contains road names and the other just a point. After exhausting research I'm now stuck.

The idea is to join up the two tables using only distance, ie, select all road names which is within 500 meteres of a postcode

My current query does nothing except waste time and I'm wondering where I'm going wrong?

SELECT "R"."name", "P"."Postcode"

FROM "planet_osm_roads" AS "R" JOIN "gb_postcodes" AS "P"

ON st_dwithin(st_transform("P"."Geom", 900913), "R"."way" , 500)

ORDER BY "P"."Geom" <-> st_setsrid(st_makepoint(0.06703376770019531,  51.51646108300825), 900913)

LIMIT 10;


texts: Module: CSV table, points to lines, intersect, add field, calculate geometry

texts:
I am looking to take point csv data and use this to create lines (point to Lines tool) then intersect this with a polygon feature. With the resulting table I need to add a field (Length) and then calculate geometry.
I have tried building this out as a model and it works, somewhat. However, I'm unsure what to set as parameters and what to leave out. Ideally I would like to trim down the optional parameters to be only what I need them. I'm not overly familiar with how this would look to have something like this as a python script instead of a model but I believe that this may be the better way to handle things.
Within the ModelBuilder there are a number of intermediate data created that I don’t necessarily need.
Are there are insights as to how I could clean this up or tackle it in a more efficient manner?


texts: Extracting pixel values by points and converting to table in Google Earth Engine?

texts:
I realise that this question has been asked before as Extract pixel values by points and convert to a table in Google Earth Engine at Stack Overflow.  However, I can't get the code to work. It runs, but the value of the raster isn't added to the table it just prints the point attributes

var yearly = ee.ImageCollection('JRC/GSW1_0/YearlyHistory');

Map.addLayer(yearly, {}, 'yearly', false);

var fg_points = ee.FeatureCollection('ft:1quTYHcLHQNy1yrssUPaubflVBEkGGMoPFLdbc_-J')

// Empty Collection to fill
var ft = ee.FeatureCollection(ee.List([]))

var fill = function(img, ini) {
  // type cast
  var inift = ee.FeatureCollection(ini)

  // gets the values for the points in the current img
  var ft2 = img.reduceRegions(fg_points, ee.Reducer.first(),30)
  // gets the date of the img
  var date = img.date().format()

  // writes the date in each feature
  var ft3 = ft2.map(function(f){return f.set("date", date)})

  // merges the FeatureCollections
  return inift.merge(ft3)
}

// Iterates over the ImageCollection
var newft = ee.FeatureCollection(yearly.iterate(fill, ft));

// Export
Export.table.toDrive(newft,
"export_Points",
"export_Points",
"export_Points");


texts:
I have been install QGIS 2.0.1.when I use Polygonize (Raster to vector)for convert jpg to shp I get error message such as:

  C:\windows\123.exe is not recognized as an internal or external
  command,  operable program or batch file.

Please help me to fix it.
thanks


texts: I have the following yellow line layer and I would expect the v.clean tool to split my lines at each intersection. Running the tool gives me the red lines. The settings are all on default. Why is that, and how can I change it? Is there perhaps a detailed tutorial on the tool I haven't found yet?


texts: Previewing data type of vector layers in QGIS Browser?

texts:
One of the features I really like about ArcGIS is that in the Catalog pane, it will show you what vector data type (whether it be point, line, or polygon) a layer is in the little icon, which is quite helpful especially when dealing with municipal or federal data that is poorly named and I have to guess what it will be. Being able to tell what type of shapefile it is saves me time waiting for a large file to render.

I know QGIS can tell what type of vector shapefile each layer is in the properties, is there a way to make the icon update according to the data type? A screenshot is shown below as to what I mean.


texts: Rearrange plotorder of SpatialPolygonsDataframe

texts:
I have a SpatialPolygonDataframe (loaded with readOGR) I want to plot in either tmap or leaflet. A lot of the polygons are overlapping, and the smaller polygons are completely covered by larger polygons. Is there a way of rearranging poltorder based on polygonsize?

Sorry for not making a reproducible example…


texts: Mapnik not generating tiles for tilecache

texts:
Imported data by using Osm2pgsql, to get rid off all problems related to data conversion, I stored this as latlong (-l parameter). I created a stylesheet, and I tested with Python script and it was working (so both data and stylesheet should be fine).

Here comes problem... when I tried to use it with OpenLayers, I got my basic layer covered by pinkish "image not found" tile, but when I try to access tile by calling cgi script "http://localhost/tc/tilecache.cgi/1.0.0/osm/5/31/26.png" I have a nicely generated tile.

Does any one have idea what I'm doing wrong?

files:

osm.xml

<Map background-color="#00000000" srs="+proj=latlong +datum=WGS84">
  <Style name="roadStyle">
    <Rule>
      <LineSymbolizer stroke="#ff0000" stroke-width="100" />
    </Rule>
  </Style>
  <Layer name="roads" srs="+proj=latlong +datum=WGS84">
    <StyleName>roadStyle</StyleName>
      <Datasource>
       <Parameter name="table">
       (select way,highway from planet_osm_line where highway in ('motorway','motorway_link','trunk','trunk_link','primary','primary_link','secondary','secondary_link','tertiary','tertiary_link','residential','unclassified') order by z_order) as roads
       </Parameter>
      <Parameter name="type">postgis</Parameter>
      <Parameter name="user">www-data</Parameter>
      <Parameter name="password">gis</Parameter>
      <Parameter name="dbname">gis</Parameter>
      <Parameter name="estimate_extent">true</Parameter>
    </Datasource>
  </Layer>
</Map>

tilecatche.cfg

[cache]
type=Disk
base=/tmp/tilecache

[osm]
type=Mapnik
mapfile=/var/www/tc/osm/osm.xml
extension=png

[basic]
type=WMS
url=http://labs.metacarta.com/wms/vmap0
extension=png

Definition of my layers:

var map = new OpenLayers.Map('map',
      {
        numZoomLevels: 12,
        projection: new OpenLayers.Projection("EPSG:4326"),
        displayProjection: new OpenLayers.Projection("EPSG:4326")
      });
var layer = new OpenLayers.Layer.WMS( "VMap0",
        "tilecache.cgi?", {layers: 'basic', format: 'image/png' } );
var tms = new OpenLayers.Layer.WMS( "osm","tilecache/tilecache.cgi?",
      { layers:'osm',  format: 'image/png'},{ isBaseLayer: false}  );


texts: How do I know when a NOMAD GFS Grib2 file is ready for download

texts:
The NOMAD GFS system creates files on a 4 times a day cycle.

The file naming specification is:

cc is the model cycle runtime (i.e. 00, 06, 12, 18)
xxx is the forecast hour of product from 00 - 384
YYYYMMDD is the Year, Month and Day

'cc' is stated as the cycle run time, so it is reasonable to assume that that is the time the processing begins.

How do I know what time the cycle actually is finished so that the grib2 files can be downloaded?


texts: Random Forest land cover classification in ArcMap?

texts:
I am an undergraduate student working with a very large data set (to me), and who has NO experience with R, but is fairly comfortable with ArcMap.

I have about 60GB of 8band 1m satellite imagery of African rainforest, and a shapefile of about 90 GPS points taken by park rangers on the ground. The GPS shapefile contains a column with the land cover type at each point. There are 9 different land-cover types, which are not easily distinguished by eye from the imagery. My goal is to create a map showing the distribution of these cover types over the entire scene.

Can anyone provide me with near-step-by-step instructions on how to use Duke's Marine Geospatial Ecology Tools (MGET), or potentially any other free Arc add-in, to run a Random Forest classification? Note that I have absolutely no idea how to use R, which seems to be the most-referenced way to process this information.

I know that running a simpler form of classification would be simpler, but since this is not just a typical undergraduate arbitrary product for an assignment (It will actually end up being used by several NGOs working with wildlife in the area), I would like to produce the most accurate product that I can, and a Random Forest model seems to be the most appropriate for this goal.


texts: Convert simple PostGIS polygon layer to polygon topogeometry

texts:
I have created a topology schema 'parcele_topo' with the default tables edge_data, face, node and relation. I also imported a polygon shapefile to this schema using shp2pgsql-gui as table 'parcele_tsur_116'. I added a topogeometry column (named 'topo') to this table using:

select topology.addtopogeometrycolumn('parcele_topo','parcele_topo','parcele_tsur_116','topo','polygon');

This is the topology summary:

Topology parcele_topo (id 3, SRID 2462, precision 0.005)0 nodes, 0 edges, 0 faces, 0 topogeoms in 0 layers Layer 1, type Polygonal (3), 0 topogeoms Deploy: parcele_topo.parcele_tsur_116.topo

And this is the screenshot:

Now I want to convert this simple table to topogeometry layer inside this topology schema so I can use topological analyses.


texts: Parameter 'input' is required by using GEE

texts:
I am new to using google earth engine. Exploring forums I found a script that is very useful to me, but when I run it I get the following error:
Error: Error in map(ID=20201020T135119_20201020T135116_T21KVQ):
String.replace: Parameter 'input' is required.

What is the mistake?
Here is the link: https://code.earthengine.google.com/7d322c73d4efc6c12203133bf025c04d[enter link description here]1
var Paraguay = ee.FeatureCollection("users/marcelos/Poligonos_segundo_informe_1-polygon")

Map.centerObject(Paraguay)
Map.addLayer(Paraguay)

var S2 = ee.ImageCollection('COPERNICUS/S2')
  .filterDate('2020-10-20', '2020-10-30')
  .filterBounds(Paraguay)
  .filterMetadata('CLOUDY_PIXEL_PERCENTAGE', 'less_than', 20)

var maskcloud = function(image) {
  var QA60 = image.select(['QA60']);
  return image.updateMask(QA60.lt(1))
}

var NDVI = function(image) {
  return image.expression(
    '(NIR - RED) / (NIR + RED)', {
      'NIR': image.select('B8'),
      'RED': image.select('B4'),
    }).rename('NDVI').copyProperties(image, image.propertyNames());
}

var NDVIcol = S2.map(maskcloud)
  .map(NDVI)

var ts = NDVIcol.map(calculateFarmersNDVI)

function calculateFarmersNDVI(image) {
  var properties = ee.Dictionary(
    Paraguay.iterate(function(region, acc) {
      return ee.Dictionary(acc)
        .combine(
          ndviForArea(image, ee.Feature(region))
        )
    }, ee.Dictionary())
  )
  return ee.Feature(null, properties
    .set('Year', image.date().get('year'))
    .set('Month', image.date().get('month'))
    .set('Day', image.date().get('day'))
  )
}

function ndviForArea(image, region) {
  var ndvi = image.reduceRegion({
    reducer: ee.Reducer.mean(),
    geometry: region.geometry(),
    scale: 30
  }).get('NDVI')
  // Not all characters are valid property names in EE
  // dictionaries. Replace them with _
  var farmers = region.getString('Farmers')
    .replace('\\W', '_', 'g')

  // There might be images without any value for the region.
  // ee.Dictionary.set() must have a non-null value.
  // The If() statement replaces null with an empty string.
  var ndviOrBlank = ee.Algorithms.If(ndvi, ndvi, '')
  return ee.Dictionary().set(farmers, ndviOrBlank)
}

var listOfFarmers = ee.List(Paraguay.aggregate_array('Farmers'))
  .map(function (farmers) {
    return ee.String(farmers).replace('\\W', '_', 'g')
  })
ee.List(['Year', 'Month', 'Day'])
  .cat(listOfFarmers)
  // We want to provide toDrive() selectors, to
  // get the CSV columns ordered. This must be a
  // client-side object though, so we need to call
  // evaluate() on it before using it.
  .evaluate(function (selectors) {
    Export.table.toDrive({
      collection: ts,
      description: 'NDVI_test4',
      fileFormat: 'CSV',
      selectors: selectors
    })
  })


texts: Can not open shapefiles from QGIS 2.18 in MapInfo 16.0, projection was missing

texts:
I am working with QGIS 2.18 where I do have several shapefiles in EPSG 4326 and encoded in ISO 8859-1, see image below.

I can check their geometry state for instance with Geometry Validator or Check Geometries plugins. As the result, I do not have any vital issues, and it basically means that the geometry itself is not missing. See the output below.

But when I open both of them in MapInfo Pro 16.0 I can't see any of them. I did it through the following command (see the image below)

  Home > Open > Relation

Then MapInfo asks me to save my shapefiles with .tab extension. And I proceed with it.

And after I can not see my shapefiles in MapInfo. Check the following image.

In QGIS

I tried to Save As... MapInfo TAB which did not solve the issue.

I deployed ogr2ogr from the GDAL/OGR Processing Toolbox, where I used [OGR] Conversion > Convert format. However, the output still does not want to show up in MapInfo.

Moreover, I changed the encoding between System, UTF-8 and ISO 8859-1, nothing helps.

I have tried to Save As... GML file, also does not help.

In MapInfo

I tried "Universal Translator". Nevertheless, I do not have an option to chose .shp format from the "FME Reader Galarie" list, but there are some other options, for instance .gml, .kml, .gdb, .mdb etc.

I have also tried with MapInfo Professional 12.5, and MapInfo Pro 15 - nothing.

So, my Question is: What could be actually the problem? Geometry? Coordinate System? Original data (which was created from .dxf)? Wrong conversion from QGIS into MapInfo? Corrupt data? Wrong execution of the data input in MapInfo?

I have tried several similar topics but they do not transparently respond to my issue.

Shapefile data not displaying on MapInfo workspace?
Tool for converting .shp file to .tab file?
QGIS from SHP to TAB
How to import a shapefile (.shp) into MapInfo Pro


texts: Install GDAL frameworks and QGIS -- redundant?

texts: Just wondering if GDAL tools are actually already within QGIS -- was it redundant for me to install them on my OS X Yosemite independently? Recently I used GDAL tools to explode .kml data for use within QGIS -- I was able to do this since following cookie crumbs led me to install GDAL for the solution.


texts: Incomplete shapefile created by arcpy.Union_analysis()?

texts:
When I attempt to use arcpy.Union_analysis() with 4 or more input files (all shapefiles; works with 3 or less) the operation is successful (no errors, exit code 0) but the final output is only the .sbn, .shx, and .shp.xml files (a few short of a shapefile set...).

If I watch my directory while the union is executing, I can see it creating the rest of the shapefile files (including the .shp) but right before the code completes, these files are deleted.

I've check my licence (arcpy.ProductInfo() = ArcInfo). I understand union is limited to 2 files if you have the basic licence. This is not the issue (it works with 3 files but not 4 or more and I don't get the 'you don't have the right licence' error).

I've replicated the issue on two machines with different, but valid, licenses.

Any solutions?

Code snippet is as follows. Unfortunately I cannot share the actual shapefiles.
arcpy code:

env.workspace = os.getcwd()
files = ['1.shp',
        '2.shp',
        '3.shp',
        '4.shp',
        '5.shp',
        '6.shp',
        '7.shp',
        '8.shp',
        '9.shp',
        '10.shp',
        '11.shp']

arcpy.SetProduct('arcinfo')
print arcpy.ProductInfo()
arcpy.Union_analysis(files, 'TestUnion', 'ONLY_FID')

I am using ArcGIS Version: 10.4.1, arcpy 2.7.10

With respect to a comment:

  "With ArcGIS Desktop Basic and Desktop Standard licenses, the number
  of input feature classes or layers is limited to two."

We verified that we have the advanced license. This is not the issue. Also, it works with three files, proving this is not the issue since if it were we'd be limited to two.

We are trying the suggestion of converting to Geodatabases before the union since shapefiles are not known for their reliability. I'll update again if it works.

As shown in the code snippet, we tried with the JoinAttributes set to ONLY_FID. It still did not complete successfully.

How can we check where the problem is?

We also verified we are staying under the 2 GB shapefile size restriction.


texts: How to draw ellipse in openlayers? I have not problems with drawing circles, or multipolygons, but i need to draw simple ellipse and not found any functions for this action.


texts: Problem using object as output path for arcpy.SaveToLayerFile_Management

texts:
I am simply trying to create a selection of a street layer from street names that appear within a separate .dbf of street names and create individual layers of each selection with a for loop. Everything works until the last step where I use the object to serve as the path string for the resulting layer.

Is this not allowed in arcpy.SaveToLayerFile_management? I have looked everywhere for the answer. I have also tried both \ and / within the pathname and with/without an "r" in front of the pathname. Here is the latest version of the code:

import arcpy
import sys, os
arcpy.env.overwriteOutput = True
arcpy.env.workspace = os.path.join(os.path.dirname(sys.argv[0]), 'Test.gdb')
rows = arcpy.SearchCursor("COR_Length_Sum_Greater_Half_Mile")
arcpy.MakeFeatureLayer_management("COR_Maintained", "CORlyr")
for row in rows:
    street = row.getValue("Display_Street")
    expression = "Display_Street" + " " + "=" + " " + '\'' + str(street) + '\''
    lyrName = street.replace(" ", "_")
    arcpy.MakeFeatureLayer_management("COR_Maintained", lyrName, expression)
    outLayer = "r" + '\"' + "A:\\" + "StreetMaintenance\GIS\Justin Castrati" + "\\" + lyrName + '\"'
    arcpy.SaveToLayerFile_management(lyrName, outLayer)


texts: How to download Open Street Map data in QGIS?

texts:
I'm trying download data from OSM in QGIS via Openlayers plugin. Downloading is successfull, but when I add layer to project, many of the features are missing.
Is there any limitation, or does anybody know why the data is not complete?


texts: posList conversion to area in square meters


texts:
I have 102 different sites in two cities, which I have created buffer zones for.
I also have layers that have polygons for different landscape features (e.g. water surface area).
I am trying to find a way to find the area of the polygons within each individual buffer, so that I can conduct statistical analysis to determine if the area of each landscape feature differentiates two groups of sites.
I am new to GIS, and so far have been able to narrow down the extent of some polygons to those that occur within the buffer zone, but not being able to differentiate between each individual site. I used Vector Selection > Extract by location > Extract from: vegetation, by comparing: buffer zone.
I would ideally like to produce a spreadsheet with each individual site, with the area of each polygon type within that buffer zone. Is it possible to keep the labels of the sites within the results as well? I am hoping to collate this information into a spreadsheet and run it through R for statistical analysis.


texts: Is it possible to update only latitude or longitude in postgis?

texts:
I have a postgis database that contains Points created based on a couple (longitude, latitude), and I need to update only one of those two for a given point.

I found that it's possible to do so but to modify both of the coordinates using :

UPDATE poi SET place = ST_GeographyFromText('SRID=4326; POINT(' || "+lg+" ||' '|| "+lat+" || ')') where mykey='"+key+"';

Any help would be very appreciated.


texts: How can I import the export file (*.sdx) from ArcSDE9.3 to ArcSDE10.1?

texts:
We got some files (*.sdx) from our clients using sdeexport from ArcSDE 9.3. We only have ArcSDE 10.1. We tried to use sdeimport command like this:

sdeimport -o create -t NHDFCODE -f c:\temp\NHD.NHDFCODE -v -i sde:oracle11g:**\** -u sde** -p **.

When run this command, there is the only message we see:

Copyright (c) 2009 Microsoft Corporation.  All rights reserved.

ArcSDE 10.1  for Oracle11g Build 679 Thu Apr  5 11:43:21  2012
SDEX File Import         Administration Utility
-----------------------------------------------------

Since I cannot tell what was going on based on this,I am just wondering if we can import the sdx file from ArcSDE 9.3 to ArcSDE 10.1. Thanks.

Update:

sdexinfo -o describe -f c:\temp\NHD.NHDFCODE
------------------------------------------------------------------------
SDEX File:             c:\temp\NHDNHDFCODE.sdx
Version:               SDEX 10.0.0_1
Type:                  Final
Volume Type:           Single

SDE Attribute Columns: 15
Name                                 Type    Width numDecimal  NULL?   RowID
-------------------------------- ----------- ----- ---------- -------- -----
OBJECTID                         SE_INT32       10      0     NOT NULL  SDE
FCODE                            SE_INT32       10      0     NOT NULL
DESCRIPTION                      SE_NSTRING    255      0     NOT NULL
CANALDITCHTYPE                   SE_NSTRING     32      0
CONSTRUCTIONMATERIAL             SE_NSTRING     32      0
HYDROGRAPHICCATEGORY             SE_NSTRING     32      0
INUNDATIONCONTROLSTATUS          SE_NSTRING     32      0
OPERATIONALSTATUS                SE_NSTRING     32      0
PIPELINETYPE                     SE_NSTRING     32      0
POSITIONALACCURACY               SE_NSTRING     32      0
RELATIONSHIPTOSURFACE            SE_NSTRING     32      0
RESERVOIRTYPE                    SE_NSTRING     32      0
STAGE                            SE_NSTRING     32      0
SPECIALUSEZONETYPE               SE_NSTRING     32      0
KCODE                            SE_NSTRING     30      0     NOT NULL


texts: Polygonizing only pixel values of 1 from raster using QGIS?

texts:
I have a raster mask that I have created in QGIS (2.2.0) that consist of 0, 1 and nodata values.  I would like to end up with a polygon dataset that consist of only the value of one pixels from this same raster dataset.  I have tried the polygonize tool, but it doesn't seem to be producing what I am looking for as it is taking both the 0 and 1 pixel values into the polygons.

How can I accomplish this task?


texts: How to add state polygon information to a city point?

texts:
I have the following layers:

Country border
federal state borders
cities (located as points)

I would need the following result file (txt or csv..): A column with the name of the city and another with the federal-state in which the city is located.

It would be great if someone could help me (I alredy tried: vector/merge shape files - but does not work).


texts: Translate wkt reference system to proj4 format

texts:
I want to translate this wkt version of the SR-ORG:8399 code reference system, to its corresponding proj4 format, so I can use it in QGIS. The exact characteristics of the system are the following: http://spatialreference.org/ref/sr-org/8399/prettywkt/

I have tried to translate to proj4 writing the following format,

+proj=longlat +lat_0=40 +lon_0=-170 +k=1.0 +x_0=0.0 +y_0=0.0 +units=m +no_defs

but QGIS Custom CRS gives a "This projection definition is not valid" message.

Can you help me correct it?


texts: Are long lines in shapefiles to be considered geodesics or straight lines in the 2D latlong space?

texts:
Is there a definition on how lines connect vertices in the shapefile format?

In the simplest case, imagine a line with only 2 points, from 40,-118 to 40,-112 - that's somewhere random in the US - with standard WGS84 geographic coordinate system. Here is the content of the .prj file:

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

Is the point say 40.1,-116 north or south of the line?

If we consider lines are linearly interpolated in the latlong space, it follows the 40 deg parallel (small circle) and the point is north of the line.
If we consider lines are shortest paths on the Earth surface, it's a geodesic (great circle) with a maximum latitude at the middle of the line, higher than 40.1 deg. Then the point is south of the line.
Or is it simply undefined? Since shapefile format has no notion of curves, but only straight segments connecting lines. Line needs to be densified (points added along the line) to clarify this answer.

If I create such scenario in QGIS, the line follows the 40 deg parallel, and would tell me answer is 1. But I wouldn't take this as a definite answer and would like to hear a more solid one.


texts: Intentionally increase coordinates in PostGIS geometry

texts:
I have a table with a geometry column (polygon) in a projected coordinate system in meters. Let's say, I want to shift all the polygons by some constant. For example, every point in all polygons must be shifted by X + 4 and Y - 4. Is there a way to do that using PostGIS?
In QGIS, I can do that by selecting all the polygons and using the Move Feature tool. However, QGIS crashes as the amount of vector features is huge, so I try to complete it in PostgreSQL.


texts:
I'm trying to add a new layer to GeoServer programmatically, So at first i have added a new workspace and then a new datastores, so far eveything is good, but when i'm trying to add a new layer i get error 500 which is:

Cannot cast org.geoserver.catalog.impl.LayerInfoImpl to org.geoserver.catalog.FeatureTypeInfo

What i'm doing is that i'm using below code:

<layer>
<name>Shahrestan</name>
<type>VECTOR</type>
<defaultStyle>
<name>polygon</name>
<atom:link xmlns:atom="http://www.w3.org/2005/Atom" rel="alternate" href="http://localhost:1024/geoserver/rest/styles/polygon.xml" type="application/xml"/>
</defaultStyle>
<resource class="featureType">
<name>cite:Shahrestan</name>
<atom:link xmlns:atom="http://www.w3.org/2005/Atom" rel="alternate" href="http://localhost:1024/geoserver/rest/workspaces/cite/datastores/kamal/featuretypes/Shahrestan.xml" type="application/xml"/>
</resource>
<attribution>
<logoWidth>0</logoWidth>
<logoHeight>0</logoHeight>
</attribution>
</layer>

And send it to the url below:

http://localhost:1024/geoserver/rest/workspaces/cite/datastores/kamal/featuretypes

Please note that i have already created a new layer using geoserver ui, and then i went through the geoserver rest index to get the xml codes

http://localhost:1024/geoserver/rest

But still i'm getting error. What's wrong?


texts:
Does anyone know where the Processing menu moved in QGIS 2.4? Specifically, I am trying to install the LAStools box for QGIS following these instructions:

http://rapidlasso.com/2013/09/29/how-to-install-lastools-toolbox-in-qgis/


texts:
I have multiple GPS traces (polylines) from cars. I would like merge them somehow to create roads. I don't have access to maps in this solution so basically I'm trying to create a poor man's map, with roads by looking where cars go.

So I have multiple polylines created from pretty inaccurate GPS receivers and I need somehow determine that which ones are from cars traveling on the same road.

Then I would like to merge them into a road using some kind of average.

Surely this must be a common problem, but I don't even know what to google... There must be some GIS term for this that I'm not aware of...

I have the data in Postgis but could do this in memory if it is easier...


texts: How can I build GDAL without the static libraries?

texts:
Background

I'm building GDAL myself, because I need it linked against libkml, with:

wget -O gdal.tar.gz "https://github.com/OSGeo/gdal/archive/${GDAL_VERSION}.tar.gz"
tar --extract --file gdal.tar.gz --strip-components 1
cd gdal
./configure --prefix=/my/install/dir --with-libkml=/where/libkml/is
make
make install

However, the build on my hosting provider (heroku) is time-critical - on the borderline where sometimes it times out, sometimes not. If I could accelerate the hour-long build by even 10 minutes, that'd solve my issue.

The goal

I believe I don't need the static library libgdal.aas well as the .so which is what I use (for geodjango). If I don't need it, I don't have to build it. You see where this is going...

The question

What configure option, if  any, would avoid building the unneeded static version of GDAL and instead build only the dynamic verison?


texts: ArcGIS Parcel Fabric: Reset adjustments/transformations?

texts:
So I'm learning how to use parcel fabrics, and getting frustrated so far. I attempted to do an adjustment on one of my subdivision plans (600 parcels). I thought it worked, so I saved it. Then on closer inspection, I realized that ArcGIS had totally mangled my parcels. However, since I saved it, UNDO is no longer an option. Is there a way to reset the geometry for a parcel/plan/fabric back to the internal (COGO) geometry so I can start again?

It appears that the only help I can find from Esri on this is to delete the parcels and enter them again, but this seems like an idiotic "solution" (especially for 600 parcels). I have tried unjoining and rejoining them, but this keeps the warped geometries. I may be remembering things wrong, but it seems like the old survey analyst handled things like this much better.


texts: ArcGIS geoprocessor memory leak: what's the underlying cause and why is it so hard to fix?

texts:
@blah238's very nice answer to Porting Avenue code for Producing Building Shadows to ArcPy/Python for ArcGIS Desktop? reminded me of all the time I've spent over the past two years battling against the dreaded geoprocessor memory leak in ArcGIS Desktop.

I only have personal experience of problems with version 9.3.1, but there are numerous threads on the ESRI support pages (e.g. here, here) which suggest that the problem goes back to at least v9.1. blah238's answer and this thread indicate that the problem still isn't fixed at v10.

I know that workarounds are possible, usually involving the use of subprocesses to force the gp to release the memory it's hogging, but they're clunky and often very slow.

This question isn't meant as an attack on ESRI, but I'd be really interested to know:

Why has this problem proved so difficult to fix?

For my work at least, it's a major issue and one that's pushed me away from using ArcGIS for spatial modelling. This has been a great learning experience, but I'm surprised that fixing the problem isn't higher on ESRI's priority list.


texts:
I'm trying to find out what the formula would be to convert the X, Y of an image to the Easting and Northing coordinates.

I've downloaded some of the raster maps from OS UK web site and I'm looking at adding the miniscale tif file to one of my applications. I've been given a number of coordinates in the easting northing format and I need to convert this to an x y of the image so I can place a label or image on top of the mini scale tif file.

I've been searching online quite a bit but everything seems to point towards applications like ArcGIS and other GIS software. However I'd quite like to make this a simple project, so I'd like to avoid embedding 3rd party sources within the application.

Does any body have any good starter guides or examples on this, or what the forumla would be to work out the x y from the easting and northing I have?


texts: Measuring length of network reached within specific search tolerance?

texts:
The purpose of this analysis is to make a measure of how well developed the street network is in different areas of a given city. The result I'm looking for is to give different parts of the network different shades of colour, to indicate which areas have the most developed pedestrian Networks. This Id think could be further developed into a polygon layer (areas of the city With best road network) for overlay analysis.

I have a good pedestrian network, but don't know how to do this analysis. Any ideas?


texts: Is there a way to have the file dialog come up and let a user choose a file name and location for output?  I have the table to Excel script as an geoprocessing tool in the Flex Viewer and would like to have the user select the file name and location rather than use hard coded ones.


texts: Disable access to GeoNode resources per default

texts:
I am testing GeoNode for geodata storage. Now I face the problem, that per default some pages are viewable by unauthorized users, especially all registered users with additional details, all registered groups, etc.

I want to disable this and found the parameter LOCKDOWN_GEONODE in the documentation. Other parameters from this document I can set in my local_settings.py file, but this parameter won't go effective when I put it at the end of this file: LOCKDOWN_GEONODE = True.

What am I doing wrong? Is there another place to put this information?
Or is there probably another way to achieve the same?


texts: How to publish map to web using ArcGIS Desktop and Server?

texts:
Just starting out with ArcGIS server at work, I've had no prior experience but I been given the responsibility to display our map on the web.

I have a ArcGIS server running on a server and I am using ArcMap 10 on my desktop.  I created a Folder in the server where I saved the map document.  I can add a new service using this document and I'm given a url when I do so.  But when I go to that URL I get this error:

Server Error in '/ArcGIS/Services' Application.
Runtime Error

And, when going to the same url from the server I get:
Server Error in '/ArcGIS/Services' Application.
No Content

I'm using ArcMap10 and ArcGIS server 10


texts: Import grid text file to shapefile in QGIS

texts:
I'm working with the Climate Change forecasts from the Hadley Center Climate Model (HADCM3) available here: http://www.ipcc-data.org/sres/hadcm3_download.html

They have many different forecasts and I'm working with the file: HADCM3_A1F_TEMP_1980.tar.gz. It seems to be a grid text file, and includes information for several different months, so it has headers at the start and also in the middle of the text:

What is the best way to import these seemingly txt files (maybe they are some other format I'm not familiar with) to QGIS, and save them as shapefiles? I've managed to extract them using 7zip, but I don't know how to proceed from there. I'm working on Windows and QGIS 2.16.2


texts: Timestamp video/time-series in Google Earth Engine

texts:
I successfully managed to output a video time-series of processed Sentinel 1 SAR images in GEE for my areas and dates of interest by adapting an existing script. However, I've not been able to timestamp it.
I have tried working off some help provided to the following two discussions and scripts within, but with no success.
Adding timestamp to video export in Google Earth Engine
How to add a Timeline tool to Time Lapse in GEE?
I have tried many things with many errors, but my error with the script at its current state is:
Collection.iterate: Error in map(ID=1):
Date: Parameter 'value' is required.
    in <global>, line 165
    in <global>, line 175

This is my full code: https://code.earthengine.google.com/b1f180b3f823680d8d138dee1d0536bc and the portion that I use to try to timestamp is below:
var text = require('users/gena/packages:text');
var S1_comp =  S1_comp.map(function(image){
  var start = ee.Date(image.get('system:time_start'));
  var end = ee.Date(image.get('system:time_end'));
  var label = start.format('YYYY-MM-dd').cat(' - ').cat(end.format('YYYY-MM-dd'));

  return image.visualize({
    forceRgbOutput: true,
    //palette: palettes.BrBG[9], //original palette was "000000", "fdbb84"
    //min: 0,
    //max: 1
  }).clip(roi).set({label: label});
});

// annotate
var annotations = [
  {
    position: 'left', offset: '1%', margin: '1%', property: 'label', scale: Map.getScale() * 2
  }
]

S1_comp = S1_comp.map(function(image) {
  return text.annotateImage(image, {}, roi, annotations)
});


texts: Need assistance setting center coordinate with GDALWARP

texts:
I'm having some issues converting my project into an MBTiles format. I'm trying to use gdalwarp to convert my main file to geotiff format, but I can't seem to get the center coordinate to how it is with my google maps custom map.

My center coordinate is latitude:-53.722541 and longitude:-100.625000 using Google Maps. I'm not sure how to use gdalwarp to ensure my center coordinate to remain at the same position. My files dimensions are 8192x8192.

The command I've been running is:

gdalwarp -s_srs EPSG:4326 -t_srs EPSG:3785 -r bilinear -te -4096 -4096 4096 4096 input.png output.tif

The SRS projection I'm using with TileMill is:

+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0.0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs +over

I'm at a loss at which of those two I would need to change in order to keep my markers properly placed, and not need to redo each and every one of them.

Thanks in advance for any support with this,


texts:
I'm trying to import a shapefile to my PostGIS database but I can not, give a look:

Importing with configuration: 12mie250gc_sir, public, geom, / home/samuelribeiro/Desktop/MALHA_DIGITAL/ac/12MIE250GC_SIR, mode = c, dump = 1, simple = 0, geography = 0, index = 1, shape = 1, srid = 4326
Shapefile type: PolygonM
PostGIS type: MULTIPOLYGONM [3]
Failed SQL begins: "SET TO UTF8 CLIENT_ENCODING;
STANDARD_CONFORMING_STRINGS SET TO ON;
BEGIN;
CREATE TABLE "public". "12mie250gc_sir" (gid serial,
"id" float8
"nm_micro" varchar (100),
"cd_geocodu" varchar (2));
ALTER TABLE "public". "12mie250gc_sir" ADD PRIMARY KEY (gid);
SE "
Failed in pgui_exec (): ERROR: function AddGeometryColumn (unknown, unknown, unknown, unknown, unknown, integer) is not unique
LINE 9: SELECT AddGeometryColumn ('public', '12mie250gc_sir ',' geom ', '4 ...
                ^
HINT: Could not choose a best candidate function. You might need to add explicit type casts.

Shapefile import failed.


texts: Depicting multiple point features in same location using ArcGIS Desktop?

texts:
On a map that I am creating (in ArcGIS Desktop 10), I need to include symbols for point features that are sometimes right next to (or even on top of) each other. This, of course, creates potential for all sorts of clutter or even for some features to not be displayed at all when zooming out.

I am thinking of positioning the symbols in a neat row (or column) off to the side and having some kind of arrow point from them toward the features' actual location. Is there some kind of symbology setting for this in ArcGIS (I want to be sure that the distance between the symbols is always consistent every time this is done), or do I have to align the symbols and draw the arrow (or whatever the proper term is) by hand?


texts:
I am running a fairly large script where I export to a PDF. Now I added an export to Excel along with export to PDF from the same MXD. The Python window says python.exe has stopped working, ends then says restarts. I never had python crash on me before I added the export to Excel feature. can this be a memory issue?
I ran the file from the command prompt and got this error or whatever you want to call it. Not sure what it means exactly

Unhandled Exception: System.IO.FileNotFoundException: Could not load file or assembly 'GpMetadataFunctions, Version=10.2.0.0, Culture=ne
eyToken=8fc3cc631e44ad86' or one of its dependencies. The system cannot find the file specified.
File name: 'GpMetadataFunctions, Version=10.2.0.0, Culture=neutral, PublicKeyToken=8fc3cc631e44ad86'
   at System.Reflection.Assembly._nLoad(AssemblyName fileName, String codeBase, Evidence assemblySecurity, Assembly locationHint, StackC
kMark, Boolean throwOnFileNotFound, Boolean forIntrospection)
   at System.Reflection.Assembly.InternalLoad(AssemblyName assemblyRef, Evidence assemblySecurity, StackCrawlMark& stackMark, Boolean fo
)
   at System.Reflection.Assembly.InternalLoad(String assemblyString, Evidence assemblySecurity, StackCrawlMark& stackMark, Boolean forIn
   at System.Reflection.Assembly.Load(String assemblyString)
   at GetManagedType(Char* _assembly, Char* _path, Char* _className)

WRN: Assembly binding logging is turned OFF.
To enable assembly bind failure logging, set the registry value [HKLM\Software\Microsoft\Fusion!EnableLog] (DWORD) to 1.
Note: There is some performance penalty associated with assembly bind failure logging.
To turn this feature off, remove the registry value [HKLM\Software\Microsoft\Fusion!EnableLog].

*i do not save or close the mxd

*my code is very large so if anyone needs any more information I can provide it

export to excel code

        in_excel = fc
        out_excel = "C:path\\{}{}.xls".format(table['County'],table['MunName'])
        arcpy.TableToExcel_conversion(in_excel, out_excel)

this is my export to PDF code

        table['zoom'] = arcpy.mapping.Layer("UnionLayer")
        arcpy.SelectLayerByAttribute_management(table['zoom'])
        df.extent = table['zoom'].getSelectedExtent()
        df.scale = df.scale * 1.1
        arcpy.RefreshActiveView()
        PDFMap = arcpy.mapping.ExportToPDF(mxd, table['County']+" "+table['MunName']+" "+"Zoom Out")

The script always crashes on the second iteration right before the export to excel. the first iteration it works fine. the export to excel occurs before the export to PDF.

UPDATE

Instead of using the table to excel conversion I used a search cursor to extract the info i wanted and tried writing to a csv and I encountered the same problem. The answer does not seem to be in my code. Has anyone ever encountered a similar error? I have found this post https://geonet.esri.com/thread/61030 which alludes to the same message I get..the person on the post says "I had to copy GpMetadataFunctions.dll and MetadataTranslator.dll to the Global Assembly Cache (GAC) directory."


texts: Add unique id to merged attributes tables QGIS

texts:
I need to find a way to stack up numeric id fields from two tables without duplications.
To put it simply here is what I want to do schematically :
Table 1

Id

1

2

Table 2

Id

1

2

3

Expected result (Addition and reset of table 2 Id)

Id

1

2

3

4

5

I'm currently copy/pasting and editing rows manually but it's completely suboptimal and time consuming.
Any ideas, plugins, functions...?


texts: Autofill field based on location in QGIS

texts:
I want to create a project where a user drops a point in QGIS within a predefined polygon and one of the fields is automatically completed with the value from a column in that specific polygon.

I have been playing around with relation reference in the field properties but can't get anything to work.

Please help!

Thomas


texts: Finding closest segment to point which contains same attribute using NNJoin plugin of QGIS?

texts:
I'm using QGIS 2.18

I have a point layer with the 4405 address positions of a municipality, as well as a road layer with road segments.

I would like to join the attribute information of the closest road segment to each point UNDER THE CONDITION that the attribute "street name" in the road layer has the same value as the attribute "street name" in the address point layer.
In other words, the tool should look for the closest road segment which has the same street name and join the attribute values of this road segment to the point.

Using the NNJoin plugin in QGIS does the job for the first part of my question, but I'm not able to include my condition (as a consequence, for 302 points the attribute values of a road segment are joined for which the street name is not the same as the one indicated in the point).

Is there some kind of extended plugin of NNJoin which could do the job?


texts: Reading attribute from shapefile using ArcPy

texts:
Currently I am making a script that allows the user to draw a line, the line is then created into end points with lat/long available.
From this lat/long, I have to run some math equations to figure out where the next coordinates should be.
I have gotten it to create the points based on a line, calculate the bearing, and used ArcPy to place the next point into the shapefile using my math equation as a test.
I am using Arcmap 10.4 for testing. I have access to Pro and 10.7, but I am using 10.4 just to learn.
My question is how do I read the attributes of a file so that I can use the X/Y of a point to run my equations on?
I have (light) experience using Pandas, working with CSV's etc just very new to ArcPy. Below is a list of what I am doing to make it easy.

Create line. Use COGO to add bearing.
Create points from lines.
Read attribute table. Read X/Y.  (This is what I need help with).
Run equation on X/Y to get the next 4 points.
Add points to point shapefile.

When I have it figured out I can provide an update to share the code with for future use.


texts: I have aerial view images (.TIFF) which I have georeferenced using ArcMap 10.4. What I am trying to do with the images is export them with the georeferenced data attached and import them to Google Earth Pro. What would be the best way to go about this task?


texts: Dissolve polygons based on max attribute value and retain attributes for that feature?

texts:
I want to dissolve a layer containing 9 features and have remaining in the attribute table only the attributes associated with max area feature.

I have tried several dissolve tools within QGIS but none seem to allow this.

In the example I wish to dissolve the features into one and only have the max area (at top) associated with the layer.


texts:
I'm trying to overlay satellite images in a web page, images are properly georeferenced.

When I use gdal2tiles by default it throws an error:

gdal2tiles.py --zoom 0-5 --s_srs EPSG:4326 20173481715B02G16.tif dist/

ERROR 6: EPSG PCS/GCS code 900913 not found in EPSG support files.  Is this a valid
EPSG coordinate system?
ERROR 6: No translation for an empty SRS to PROJ.4 format is known.
Traceback (most recent call last):
  File "/usr/bin/gdal2tiles.py", line 2278, in <module>
    gdal2tiles.process()
  File "/usr/bin/gdal2tiles.py", line 482, in process
    self.open_input()
  File "/usr/bin/gdal2tiles.py", line 856, in open_input
    self.out_ds.SetMetadataItem('NODATA_VALUES','%i %i %i' % (self.in_nodata[0],self.in_nodata[1],self.in_nodata[2]))
IndexError: list index out of range

However if I use a profile such as raster or geodetic it works fine but when I overlay that layer over an OpenStreetMap layer, using leaflets or openlayers the result its a mess, the layer goes to Antarctica.

gdal2tiles.py --zoom 0-5 --profile geodetic --no-kml 20173481715B02G16.tif dist/
gdal2tiles.py --zoom 0-5 --profile raster --no-kml 20173481715B02G16.tif dist/

This example it's using the raster profile, with geodetic the result it's pretty similar.

I did try to force projection with --s_srs EPSG:4326 parameter, change the projection in openlayers but nothing seams to work.

What am I doing wrong or how can I fix this?


texts: MapServer not displaying data from other network locations?

texts:
I have installed MapServer and try to access raster file(.img) from other Network locations (NAS), which is displaying map background color for the mapcanvas.  I have also tried to load the same raster file from local disk and it's loading fine.

Please tell me what are the possible solutions.


texts: Obtaining range of values in zones created using Slice tool

texts: I have to reclassify rasters in my workflow into two zones using natural jenks. The Slice tool does the job perfectly, but I would like to print the range of values out in my script to be used as variables and written in csvs. Is there some sort of command I can use that will print out the threshold value? I should mention that the input raster is a 32bit float which complicates things a little. My alternative is to use Zonal Statistics and derive min and max for these ranges.


texts:
Is there any chance to edit the label size of features from OSM Basemaps in QGIS 3?
For example, the road names are displayed way too small and are barely readable when exporting maps as pdf.


texts: Find points that intersect a user-defined polygon in GDAL

texts:
I'm trying to do something that I don't know is possible in GDAL/Python.

I have a list, L, of pre-defined lat/long points that make up a uniform grid.

A user will send me a set of points that make up the vertices of a polygon. This polygon will usually be irregular, not square. Any number of points could make up the polygon vertices, but it will likely not exceed 10 or so.

I would like to:
1) Create a polygon, P, based on the user-defined vertices.
2) Find the lat/long points in L that intersect polygon P.

Is this possible? I don't normally work with vectors in GDAL, so I'm not familiar with the best way to do this.


texts: File header incorrectly reporting 0 point records in R?

texts:
I have a set of 10 LiDAR files that I'm attempting to retile and eventually use to create a Canopy Height model. I'm able to import them as a catalog, however when I go to retile my catalog, I get the following error:

Error: File not found
In addition: Warning message:
In normalizePath(path.expand(path), winslash, mustWork) :
  path[1]="NA": The system cannot find the file specified

Additionally, the progress figure looks like this (indicating there were not data points found in each file):

To investigate, I imported a single LiDAR file, and was successfully able to plot the point cloud. Proving that there definitely ARE data in the file. However, upon importing the file, I got the following warning:

Warning messages:
1: Invalid file: header states the file contains 0 points but 3408961 were found.
2: Invalid file: the header states the file contains 0 1st returns but 3322882 were found.
3: Invalid file: the header states the file contains 0 2nd returns but 81189 were found.
4: Invalid file: the header states the file contains 0 3rd returns but 4852 were found.
5: Invalid file: the header states the file contains 0 4th returns but 58 were found.

And when looking at the header, it tells me the number of point records is 0. However, using summary(), on the file, I'm told I have "3408961 points,  NA pulses."

This is all very perplexing to me. Does anyone know of any possible solutions?


texts: Creating a portable map package file folder

texts:
I am searching for an expression or some assistance in changing filepaths and links to photos within attributes to relative paths.

I have used QConsolidate to create a folder with all the files and links I want.

I would like to create a single folder so I can give the folder to a client and when they install on their system they will have all the connections retained.

I have been finding info on relative paths and thought I had cracked it with @project_folder ||'/data/'|| filename

But this expression does not change the file path when the folder is moved to a different location.


texts: How can I calculate the Kolmogorov-Smirnov Tests (ks.test) with raster data in R?

texts:
I want to test if my data (especially the temperatures of the landsurface) are normal distributed.
ks.test(rasterstack[],"pnorm",alternative="two.sided",exact=NULL)

This is the result:
One-sample Kolmogorov-Smirnov test
data:  r_vector
D = 0.086658, p-value = 0.4405
alternative hypothesis: two-sided
Is this the right code to test raster with a ks.test?


texts: How to examine length of vector lines in qgis/grass?

texts: I have vector layer of lines. How can I examine length of lines?


texts: Calculate field in a joined table in ArcPy

texts:
EDIT: I was able to solve the problem with the dictionary approach suggested by @BERA. However, I still don't get why it does not work with the CalculateField method.
I'm trying to calculate a field in a joined table in ArcPy. Code snippet below.
The Calculate Field always gives back the error message:

ERROR 000539: Invalid field MAX

When I do the same thing manually one after the other in ArcGIS everything works perfectly well.
poly_riskobj = r"E:\SRRM_Walldorf\Risikoobjekte.gdb\Ro_Mitte"
out_stat_table = r"E:\SRRM_Walldorf\Risikoobjekte.gdb\Ro_Mitte_HQ100"

arcpy.MakeTableView_management(out_stat_table, "table_view")
arcpy.MakeFeatureLayer_management(poly_riskobj, "obj_layer")
arcpy.AddJoin_management("obj_layer", "OBJECTID", "table_view", "OBJECTID_1")

arcpy.CalculateField_management(poly_riskobj, "WSP_HQ100", "!MAX!")

When I checked the fields with arcpy.ListFields I realized that the joined fields are missing in the poly_riskobj fc I'm refering to in the CalculateField. They only exist in the obj_layer. Which makes sense...
So I adapted the CalculateField to
arcpy.CalculateField_management("obj_layer", "WSP_HQ100", "!MAX!")

and saved it with CopyFeatures. Hoewever, when I check the result in ArcGIS the field "WSP_HQ100" is not updated with the MAX values.
I don't know what I'm missing. I read through these posts but could't figuere it out:

Understanding arcpy.CalculateField_management with Join?
How to Calculate Field in joined table using ArcPy?


texts: Nominating/upvoting/promoting QGIS bug-reports at QGIS-Hub for faster fixing?

texts: Since I started working with QGIS, I encountered several bugs that were already reported at the  QGIS-Hub. All of the bugs that limit my daily workflow in QGIS are just considered "normal" - not "blocking" or worse and have not been touched for 2 years now - I think the promotion of bugfixes would be more appropriate than just filing another ticket. I wonder if there was any mechanism to draw more attention to specific bugs, for example a voting system like here at gis.stackexchange.com? (The specific bug I am referring to is: https://hub.qgis.org/issues/11106)


texts: Enlarge/roughen the Vector Size of Tiles

texts:
I'm currently running a GeoServer service and my task is to make the vectors rougher than they are right now. I really don't know where and how GeoServer processes the data from the database and how it will be displayed.

I want to do this to speed up the render time for loading the map. If the server uses 20 points right now, I want it to use 10 for example.

I don't know where to start. If my information provided is not enough, just tell me and I'll tell everything you want to know.


texts: Translate Excel formula to be used in field calculator

texts:
I have a series of equations to create corner points to make rotated polygons based on angles in the attribute table, but I have created it in excel and can't figure out the syntax to get it to work properly in the field calculator.  I have 8 different trigonometric formulas to create XY coordinates for the 4 corner points.  One of them is as follows:

=COS(ATAN(1500/2250)+RADIANS(MOD(360-([Angle]+270),360)))*(SQRT((1500^2)+(2250^2)))+[X]

Once I understand the syntax for one, I should be able to use it for the other 7 formulas, even though they vary slightly.
Can anyone advise the best way to use this formula in field calculator?
Thanks.

Edit, after receiving @Paul's answer:

I figured out the correct syntax of some of the other ones, but two of them are still giving me trouble.  One is producing a different result than Excel and one is saying that the syntax is incorrect.

=COS(RADIANS(270)-ATAN(750/1500)+RADIANS(MOD(360-([A]+270),360)))*SQRT((1500^2)+(750^2))+[X]

=COS(-ATAN(1500/2250)+RADIANS(MOD(360-([A]+270),360)))*SQRT((1500^2)+(2250^2))+[X]

I think it has to do with the negative signs and me probably putting them in the wrong spot...


texts: Semi-Automatic Classification Plugin Scipy Dependency Error

texts:
I've been using QGIS with the Semi-Automatic Classification plugin (SCP) at work, and am attempting to install and use it at home. Unfortunately, I am encountering an error which causes QGIS to crash when installing, disabling, or uninstalling SCP.

Upon opening QGIS 2.18.6, I get the following messages:

Semi-Automatic Classification Plugin: Please, restart QGIS for executing the Semi-Automatic Classification Plugin. Possible missing dependecies: SciPy

Semi-Automatic Classification Plugin: Please, restart QGIS for executing the Semi-Automatic Classification Plugin

From investigating the plugin's source, it looks like it is trying the following import before throwing the SciPy dependency exception:

import scipy.stats.distributions as statdistr

I tried the same import from the OSGeo4W Shell, and get the following traceback:

Taceback (most recen call last):
  File "<stdin>", line 1, in <module>
  File "C:\OSGEO4~1\apps\Python27\lib\site-packages\scipy\stats\__init__.py", line 334,
    in <module> from .stats import *
  File "C:\OSGEO4~1\apps\Python27\lib\site-packages\scipy\stats\stats.py", line 181,
    in <module> import scipy.special as special
  File "C:\OSGEO4~1\apps\Python27\lib\site-packages\scipy\special\__init__.py", line 546,
    in <module> from ._ufuncs import *

Additionally, doing a pip list does not display scipy but does display numpy and matplotlib, the other two dependencies for the plugin.

I've attempted to reinstall QGIS and supporting apps / libraries using the advanced and express installers, but run into the same problem each time I try to install SCP. The error also causes QGIS to crash and write a minidump on close. I'm also using Windows 10.

I am able to import scipy from the OSGeo4W Shell, but any module imports do not work, i.e. import scipy.stats or import scipy.linalg

How can I resolve this specific Scipy dependency issue within QGIS?


texts: Reducing Sentinel-1 timeseries in Google Earth engine

texts:
For my research I want to have sentinel-1 layers for each month and in order to reduce the amount of border noise I used the imageCollection.median() on a Sentinel-1 timeseries within the Google Earth Engine.

The problem is that when I zoom in the pixels seem to be overlapping and not aligning properly.

Does anybody know how to fix this or is it simple not possible to reduce sentinel-1 timeseries?

See code below:

// Load the Sentinel-1 ImageCollection.
var sentinel1 = ee.ImageCollection('COPERNICUS/S1_GRD');
var geometry2 = /* color: #d63000 */ee.Geometry.Polygon(
    [[[17.12881578170004, 46.81643443671282],
      [17.101364827626412, 46.55282544717593],
      [18.4462854747203, 46.861496672885984],
      [18.094957936090168, 47.11612461249994]]]),

// Filter by metadata properties.
var vvvh = sentinel1
// Filter to get images with VV and VH dual polarization.
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VV'))
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VH'))
  .filter(ee.Filter.eq('instrumentMode', 'IW'))
  .filter(ee.Filter.eq('resolution_meters', 10))
  .filter(ee.Filter.eq('relativeOrbitNumber_start', 51))
  .filter(ee.Filter.bounds(geometry2));

var S1 = vvvh.select('VH','VV');

//Create timeseries for whole year 2015
var S1y = S1.filterDate('2015-01-01', '2015-12-31');

// Create timeseries for summer and winter
var S1s =  S1.filterDate('2015-06-01', '2015-09-01');
var S1w =  S1.filterDate('2015-10-01', '2016-03-01');

// Create timeseries for each month in the year 2015
var S1_1 =  S1.filterDate('2015-01-01', '2015-02-01').median();
var S1_2 =  S1.filterDate('2015-02-01', '2015-03-01').median();
var S1_3 =  S1.filterDate('2015-03-01', '2015-04-01').median();
var S1_4 =  S1.filterDate('2015-04-01', '2015-05-01').median();
var S1_5 =  S1.filterDate('2015-05-01', '2015-06-01').median();
var S1_6 =  S1.filterDate('2015-06-01', '2015-07-01').median();
var S1_7 =  S1.filterDate('2015-07-01', '2015-08-01').median();
var S1_8 =  S1.filterDate('2015-08-01', '2015-09-01').median();
var S1_9 =  S1.filterDate('2015-09-01', '2015-10-01').median();
var S1_10 =  S1.filterDate('2015-10-01', '2015-11-01').median();
var S1_11 =  S1.filterDate('2015-11-01', '2015-12-01').median();
var S1_12 =  S1.filterDate('2015-12-01', '2016-01-01').median();

var compall =
ee.Image.cat([S1_1,S1_2,S1_3,S1_4,S1_5,S1_6,
S1_7,S1_8,S1_9,S1_10,S1_11,S1_12]);

Map.addLayer((compall), {min: -25, max: 0}, 'RGB');

All right, I tried to use only 1 relative orbit and now they seem to allign but still overlap:


texts: Raster Copy: Output bit size does not match specified value

texts:
Using ArcGIS (Map/Catalog) 10.0.  I have a raster file (kdens) developed from point data using the kernel density tool.  It is 32-bit float file format.  I need this file converted to a polygon shapefile.  (Note I have used the Math/Times tool to preserve the decimal units in a new raster: kdens100K).  Using the copy raster tool, I have specified the output raster to 32-bit unsigned.  However, the output raster is identified as a 16-bit signed when I check the raster properties and the raster values are truncated to match the 16-bit value range.

Any insight into to why the output is being formatted to 16-bit signed when I specify 32-bit unsigned?


texts:
I am having trouble exporting the raster image via Earth Engine. The error is in the print screen below, as well as the script used. The data is from the product MCD64A1 from MODIS.
// Exporting to Google Drive
Export.image.toDrive({
  image: burnedArea_2001,
  description: '2001',
  scale: 500,
  region: Amazonia
});


texts: Identifying source of reference mosaic dataset in ArcPy?

texts:
I want to check the source mosaic dataset for a reference mosaic dataset within Python.

Using:

arcpy.Describe(reference_mosaic).referenced

I can confirm that the reference mosaic is indeed a reference mosaic, however, the describe object does not expose the path to the source mosaic dataset as far as I can tell, even though this information is easily accessible in ArcMap:

Using:

arcpy.Describe(reference_mosaic).childrenNames
[u'AMD_test_reference_BND', u'AMD_test_reference_LOG', u'AMD_test_CCA', u'AMD_test_reference_DEF']

I can see that the reference mosaic has a child (AMD_test_CCA) that refers to the source mosaic, however, I cannot interact any further with this item as arcpy says that it does not exist. Although the name of this child does confirm the name of the source dataset, I would prefer to confirm the absolute path of the source dataset.

Does anyone have any ideas or suggestions as to how I could achieve this?

As the reference mosaic dataset has an image service being published from it, I will do not want to have to recreate the reference mosaic and go through the rigmarole of stopping and/or reconfiguring the service.

Using ArcGIS Desktop 10.5.1 on Windows 7 64-Bit.


texts: Calculate intersection for two FeatureCollections containing points in GEE

texts:
I have a FeatureCollection of points that I want to sample data from. However, I only want data for the points where information is available from both the winter and summer image. The solution I've come up with is to sample all data from both and then intersect those point classes. I've tried implementing the intersection code from here, which was written for intersections of polygons, but I'm getting two incorrect results that I don't understand.
First, depending on the order that I nest the functions in, the intersection FeatureCollection will have a different number of elements. I don't understand why that is, since it seems like it shouldn't matter what order the mapping happens in. In my specific code, winterOutside produces 7 features and summerOutside produces 5 features (5 was the correct number).
Secondly, on the map, both FeatureCollections look identical, but incorrect - there are only four points on the map, and they're both missing a fifth point that should be included.
var winter = landsat8.filterDate('2017-12-01', '2017-12-10').mosaic()
var summer = landsat8.filterDate('2017-06-01', '2017-06-10').mosaic()
var visParams = {bands: ['B4', 'B3', 'B2'], max:20000};

Map.addLayer(winter, visParams, 'winter')
Map.addLayer(summer, visParams, 'summer')

// makes FeatureCollection with band values of winter values for invariant points
var winterSamples = winter.sampleRegions({
  collection: points,
  scale: 6.5,
  geometries: true
});
// makes FeatureCollection with band values of summer values for invariant points
var summerSamples = summer.sampleRegions({
  collection: points,
  scale: 6.5,
  geometries: true
});

// with winter on the outside:
var winterOutside = winterSamples.map(function(feat1){
  feat1 = ee.Feature(feat1);
  var mapped1 = summerSamples.map(function(feat2){
    feat2 = ee.Feature(feat2);
    var intersection = feat1.intersection(feat2, ee.ErrorMargin(1));
    return ee.Feature(intersection);
  });
  return mapped1;
}).flatten().distinct(['B1', 'B2', 'B3', 'B4', 'B5'])

print(winterOutside);
Map.addLayer(winterOutside, {color: 'blue'}, 'winterOutside')

// with summer on the outside:
var summerOutside = summerSamples.map(function(feat1){
  feat1 = ee.Feature(feat1);
  var mapped1 = winterSamples.map(function(feat2){
    feat2 = ee.Feature(feat2);
    var intersection = feat1.intersection(feat2, ee.ErrorMargin(1));
    return ee.Feature(intersection);
  });
  return mapped1;
}).flatten().distinct(['B1', 'B2', 'B3', 'B4', 'B5'])

print(summerOutside);
Map.addLayer(summerOutside, {color: 'red'}, 'summerOutside')

Here's the link to my code if you want to grab my imports. My question is really just "How do you intersect two FeatureCollections that only have points in them?" but also "Why isn't the above working?"


texts: WFS shows only a small proportion of its area

texts:
I'm using QGIS and when I want to add features of this WFS layer
"https://www.wfs.nrw.de/geobasis/wfs_nw_atkis-basis-dlm_aaa-modell-basiert?"
it covers only a small proportion of the area it should cover - namely the most south western part of the German federal state Northrhine-Westphalia (NRW), although it should cover areas in the whole state.
Even using a vector layer of NRW as backgrund and using the option "only request features overlapping the current view extent" does not help.


texts:
(as simply as you can) what are the differences between INSPIRE, GEOSS and OGC?

I will explain what I understood:

The OGC provides standards for web-mapping that like WMS, WFS (are these standards or services? is there any difference?) Does OGC provide data
If INSPIRE provide data to EU countries and GEOSS/OGC does that but for all the world so what is the the point of INSPIRE, then? is there and differences in SDI architecture between INSPIRE and GEOSS?


texts: GeoPandas setting geometry for a row with .loc and Multipolygons

texts:
Given the following example:
def merge_with_largest(input_gdf,poly):

    # select the polygons bordering polygon
    close_by_gdf = input_gdf[input_gdf.geometry.intersects(poly)]
    # selects the index of the row with the largest area
    selected_index = close_by_gdf['area'].argmax()
    # now we merge the geometry of the selected index with the polygon geometry
    selected_geom = input_gdf.loc[selected_index]['geometry']
    # setting the
    input_gdf.loc[selected_index,['geometry']] = gpd.GeoSeries((selected_geom,poly)).unary_union

What I want to do here is to merge the geometry of poly with the largest of the polygons that do intersect it with. I want to do this by using .unary_union between the two selected polygons and then update the geometry of the selected polygon with what .unary_union returns. The problem though is that sometimes .unary_union returns a mulitpolygon and then the .loc method throws an error because it's an iterable.

ValueError: Must have equal len keys and value when setting with an iterable

Is there a way to overwrite this and set the multipolygon as the the geometry for the selected column?


texts: How does oversampling work in WMTS layers?

texts:
(This question arises as a sub-topic of the answer provided here.)

In the symbology options for a raster layer in QGIS, there is a setting called "Oversampling" (see image below). How does this setting affect the display of a raster layer from a tile-based map service (WMTS), specifically one that displays different tilesets depending on the zoom level?

For example, the USA Topo Maps service contains tilesets for three different scale-based map series (1:250k, 1:100k and 1:24k), and QGIS loads the appropriate tileset for the current display scale and magnification. However, changing the Oversampling setting has an effect on which tileset is loaded.

In the examples below, the display scale is set to 1:125,000 at 2x magnification, and QGIS loads the following tilesets depending on the Oversampling value:

At 2.0, the 1:100k tileset is loaded
At 3.0, the 1:24k tileset is loaded
At 0.0, no tileset is loaded

The same results will occur at other equivalent scale/magnification factors (such as 1:62,500 at 1x magnification).

Why does changing the Oversampling setting have this effect? Is this a bug, or is it how the feature is designed to work?


texts: Problem with plugins installation (2.18.3)

texts:
I'm using QGIS version 2.18.3. When I try to install a plugin, a message appears after the installation progress

This problem occured not only with this plugin, but with every plugin I tried to install (I tried at least 30 different). How can i solve this problem?


texts: How to fill holes/missing values in XYZ data / a DEM with QGIS?

texts:
I have equidistant XYZ point (text) data that I used to create a point grid. However, it's not perfect - there are some holes. I then created a DSM raster from it, which obviously has the holes too. The holes range from a single point to larger areas (~50 or more points missing).

How do I patch up those holes? I tried to interpolate the raster, but that reduces the quality, and doesn't fix the zero "holes" - they're just slightly raised afterwards. Ideally I'd like to just fill those zero values with Z values from their neighbouring squares.

Another issue is that the XYZ data set isn't rectangular but like in the third attached image (red path). If I use a fill method, how do I avoid filling up the whole surrounding area?


texts: How to identify polylines with multiple end/starting points?

texts:
I'm looking for a way to programmatically identify polyline features which do not represent a 'simple' line from A to B but instead have several starting and/or end points. These are called complex network edges when working with geometric networks.

Why am I asking this? Because I have a geometric network with a feature class that is defined as only having simple network edges. But it somehow occurred that there are some 'bad' (i.e. complex) edges in there as well and I need to sort them out. Checking if the interface IComplexEdgeFeature is implemented doesn't work on the feature objects because by definition all features within that class are simple edges.


texts: Building ArcGIS Model for identifying location?

texts:
I have a point shapefile of some plant species in an area of about 25000 sq.km. I have some other layers:

unit boundary (private/public)
roads
research plots
raster (land use)

If I overlay all the shapefiles and raster files, I will know the what kind of plant species is found in what land use, proximity of road to the plants, and unit type and so on.

For example, I have plant Species (A) in location 1 which has

landtype: forest
unit boundary: private
research plot: control

Now I want to build a model that can identify other locations in the area with the same conditions as location 1.

Has anyone built such models in ModelBuilder, or have any idea what tools can be used to achieve this goal?


texts: Determining percentage of overlapping lines using QGIS?

texts:
i have a lines layer that represents the bikeways of my city, and also have another overlapping layer that represents the bikeways that need repairs.

How can i get an analysis of the percentage or distance that needs reparing compared to the total length or 100% of the bikeways of the city?


texts: Extracting NetCDF raster values to time enabled shapefiles - what's wrong with my python script?

texts:
I'm trying to make a python script tool that extracts the values from a NetCDF raster to each data point from a time enabled shapefile at the correct time. My NetCDF raster consists of monthly sea surface temperature (SST) averages and I also have a plastic debris point shapefile that has been time enabled with time step interval of one month. I have checked with Time Slider that both the raster and the shapefile values change through time correctly.

Can someone please help me figure out why I keep getting an error with this script:

mxd = arcpy.mapping.MapDocument('SST.mxd')
df = arcpy.mapping.ListDataFrames(mxd)[0]
rasterLayerName = "SST"
rasterLayer = arcpy.mapping.ListLayers(mxd, rasterLayerName, df)[0]
plasticLayerName = "Plastic"
plasticLayer = arcpy.mapping.ListLayers(mxd, plasticLayerName, df)[0]
df.time.currentTime = df.time.startTime
while df.time.currentTime <= df.time.endTime:
   uniqueName = "Plastic_" + df.time.currentTime.strftime("%Y_%m")
   outPointFeatures = "U:/gis/" + uniqueName
   arcpy.sa.ExtractValuesToPoints(plasticLayer, rasterLayer,
         outPointFeatures, "NONE", "VALUE_ONLY")
df.time.currentTime += df.time.timeStepInterval

Here is the error:

Runtime error
  Traceback (most recent call last):
   File "<string>", line 6, in <module>
   File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\sa\Functions.py", line 1311, in ExtractValuesToPoints
     add_attributes)
   File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\sa\Utils.py", line 47, in swapper
     result = wrapper(*args, **kwargs)
   File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\sa\Functions.py", line 1304, in wrapper
     add_attributes)
   File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\geoprocessing\_base.py", line 498, in <lambda>
     return lambda *args: val(*gp_fixargs(args, True))
 ExecuteError: ERROR 001143: Background server threw an exception.

I'm using ArcGIS Desktop 10.1 and Python 2.7 on Windows 7.


texts: How to turn off caching of fill patterns in TileMill?

texts:
I'm using CartoCSS in TileMill and png files on disk for filling polygon symbols with patterns.

polygon-pattern-file: url('sprites/yard.png');

The display is ok the first time, but when I change the file itself in a raster editor, polygon fills don't change. They do change only if I rename my png file to a different name and make a corresponding change in CSS.

So I guess TileMill caches files for given URLs somehow. How to turn that off?