texts: Getting a reference to the current basemap

texts:
In my application I'm using basemap gallery widget for the user to switch between basemaps. However, when the basemap has been switched the map.getBasemap() keeps returning the original basemap assigned to the map in the constructor instead of the current basemap.

Also map.getLayer(map.basemapLayerIds) returns undefined after the basemap has been changed by the basemapGallery widget.

How can I return the reference of actual current basemap which is displayed in the map?


texts: Data Interoperability or FME ESRI Edition?

texts:
I've been debating keeping our Data Interoperability licence or going with the full FME ESRI Edition.

I'm looking for advice.


texts: How do I save my map on QGIS as a shapefile?

texts:
The project for my class involves downloading a zip (with sharefile documents) + csv file which provide the data and a basic map of Manhattan. As per instructions, I computed a new field and created a map based on it. Whenever I try to save this file as a shapefile, it only shows the information from the original zip file (and none of my work).

Hopefully this makes sense, I'm new to GIS.


texts:
I'm working on setting up the attributions for a map, and I'm wondering how to access more the styling options for the attribution control.

I've modified the .ol-attribution ul {} css property to display the items in a vertical list (display: table), but the default MapQuest attribution shows up on 3 lines instead of 1. Is this something I can fix? The line below is another layer I have.  Also, is there a way to configure where the logos will be placed?  I've read documentation but didn't find anything that gets so fine-grained.  I also tried overwriting the ol.source.MapQuest({ attribution: ...}) with all the same words/links in a single html line, but the stock attribution still shows up.

Also, when the MapQuest layer is removed, it text continues to wrap in a funny way:


texts: Remove polygons that do not fall completely inside overlay layer within the clip function

texts: I have two polygons: 1) successfulheights and 2) clipArea. I have used the clip function to remove the unwanted successfulheights. Though, there are features that have been clipped in half and those are now invalid for my project. Is there a parameter to set, "if a feature does not completely lay within the overlay layer, remove it"?


texts:
I want to edit the code I have below to save vector files as vectors and also save vector files as KMZs.

The code I have does a mass save to my temporary folder but only for rasters.

myDir = r'C:\temp'
layers = layers = iface.mapCanvas().layers()
pipe = QgsRasterPipe()
for layer in layers:
   extent = layer.extent()
   width, height = layer.width(), layer.height()
   renderer = layer.renderer()
   provider=layer.dataProvider()
   crs = layer.crs().toWkt()
   pipe.set(provider.clone())
   pipe.set(renderer.clone())
   file_writer = QgsRasterFileWriter("%s\\%s.tif" %(myDir, layer.name() ))
   file_writer.writeRaster(pipe,
                       width,
                       height,
                       extent,
                       layer.crs())
   print "%s\\%s.tif" %(myDir, layer.name() )


texts: Create GeoJSON layer based on remote server request (ArcGIS JS API)

texts:
I'm starting with ArcGIS JS API. I'm doing the following POST request to my remote server:
esriRequest("php/load_gisdata.php", {
    method: "post",
    query: {
        table: "participants_loc",
        fields: "*",
        geometry_column_name: "geom"
    },
    responseType: "json"
}).then(function(response) {
    // The requested data
    console.log(response);
    var geojson = response.data;

    const blob = new Blob([JSON.stringify(geojson)], {type: "application/json"});
    const url  = URL.createObjectURL(blob);
    const lyr = new GeoJSONLayer({ url });
    map.add(lyr);

});

the geoJson variable is correctly returning the GeoJSON (point type, srid: 4326) I want to plot.
It's giving the following message:
"[esri.layers.graphics.sources.GeoJSONSource] Some fields types couldn't be inferred from the features and were dropped"

Someone could help me out how to plot this geojson features from the remote server request?!


texts: How to import PostGIS layer with geography type to QGIS?

texts:
I tried to import in QGIS a layer from a PostgreSQL table with a geography type column created through the extension PostGIS.

I noticed that QGIS doesn't seem to like the geography column...

Here the step to create the column...

ALTER TABLE traffico ADD COLUMN geog geography(POINT,4326);

then to update the column...

UPDATE traffico SET geog = ST_GeogFromText('SRID=4326;POINT(' || lon || ' ' || lat || ')');

but then when I connect QGIS to my server I see this:

There's no way I can find to tell QGIS to treat my column as a geography type instead of geometry type. Indeed if I try to upload the layer like this the point coordinates do not correspond to reality.


texts: Marker-width determined by population in animated bubble map. CARTO Builder

texts:
I'm trying to put together an animated map of population data for two different species on five different islands. I would like the population of each to be displayed as a standard marker whose width is proportional to the population of the species for the given year, animated over twenty years. When played, the animation should show 10 stationary bubbles, in concentric pairs located on each island expanding and contracting to reflect their annual populations.

I am having some difficulty though setting this up in Carto Builder CSS, such that marker-width is defined by different thresholds of the 'value' as produced by the torque-aggregation-function: "round(avg(population))".
Perhaps I do not understand exactly what the 'value' is and how it changes, but I have followed examples such as this one: http://andrewxhill.com/blog/2015/04/17/torque-unknown/
but do not understand how the threshold values for marker-width graduations were determined.

The following code seems to work for a single species but the sizing thresholds don't seem to respond exactly to the population data. When I try to use the dataset for both species, using different colors with
"marker-fill: ramp([value], (#5B3F95, #1D6996), (1, 2), "=");"
the markers do not expand and contract in the same way as they do individually (generally much less responsive), and the markers show up unfilled with outlines only.

    Map {
  -torque-frame-count: 20;
  -torque-animation-duration: 10;
  -torque-time-attribute: "year";
  -torque-aggregation-function:"round(avg(population))";
  -torque-resolution: 1;
  -torque-data-aggregation: linear;
}
#layer {
  marker-fill: #d71500;
  marker-fill-opacity: 1;
  marker-line-width: 1;
  marker-line-color: #FFF;
  marker-line-opacity: 1;
  comp-op: lighter;
}
#layer{
[value>0]{marker-width: 7;}
[value>3]{marker-width: 7.5;}
[value>6]{marker-width: 8;}
[value>9]{marker-width: 8.5;}
[value>12]{marker-width: 9;}
[value>15]{marker-width: 9.5;}
[value>20]{marker-width: 10;}
....
[value>1650]{marker-width: 47;}
[value>1700]{marker-width: 47.5;}
}


texts: How to generate all the WMS layers getGetFeatureInfoUrl details into iframe?

texts:
With this code I am able to generate only the last layer GetFeatureInfo response only, I am able to display all the layer URLs but how am I able to show all the layers details in the iframes (number of WMS layers = number of iframes).

Here I have one base layer and 3 layers so I am expecting 3 individual iframes as overlays on the map.

  <link rel="stylesheet" href="https://openlayers.org/en/v3.19.1/css/ol.css">
  <link rel="stylesheet" href="https://code.jquery.com/ui/1.12.1/themes/base/jquery-ui.css">
  <link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/ol3/3.5.0/ol.css">

X


texts: indentation problem or function problem?

texts:
I am trying to create in Field Calculator (arcmap) a function to assing a numeric value in a column (heighs) based on a text value from another column(building).

'Pre Logic' input box

def a(heighs, building):
 if (building=  '111'):
  return 50

 elif (building=  '112'):
  return 100

 else:
  return 1

'field =' box

a( !heighs! , !building! )

I am always getting this error:
"Failed to execute. Parameters are not valid."
"ERROR 000989:Python syntax error: Parsing error SyntaxError: invalid syntax (line2)

Is this a indentation error or is an error in the function that I am creating?
How can I make the indentation on the Field Calculator window (tab don't work!!)? One space for each line?


texts: How to export the legend (Layers Panel) as an image in QGIS using the Python console?

texts: I want to export the Legend (the Layers Panel) as an image using the python console.


texts: We have to generate hundreds of mxd's and i am currently working on a template. In the template mxd, template legend items will be standard and while changing the source layer path, legend items are aligning to the bottom instead of top. As in the image, "Source A" has all the legend items as that of the template legend and it's getting aligned properly. But as one can see, "Source B" has only one legend and it's getting aligned to the bottom of the legend frame. How to fix this ? Tried with "Legend Properties" and can't able to find the specific option to achieve this


texts: Distance between current and previous GPS Point BigQuery

texts:
I have a table with 150.000 rows containing DateTime, GpsLongitude and GpsLatitude columns. Timestamp difference between rows is 10 seconds. I want to calculate the distance between every current and previous GPS point (so the distance in meters between GPS points every 10 seconds).
I tried with the query below:
select DateTime,
ST_DISTANCE(PointWKT1, LAG(PointWKT1) OVER (PARTITION BY DateTime ORDER BY DateTime)) as distance
from (select DateTime, ST_GEOGPOINT(GpsLongitude, GpsLatitude) as PointWKT1
      FROM `xxx.yyy` )
      order by DateTime

The problem is, that the returned distance is always null! Any suggections what could be wrong?


texts:
I am trying to highlight features in a tileset-layer based on this tutorial: https://jsfiddle.net/oyjLzfq4/1/
The grey layer is stored as a tileset in Mapbox studio. The feature under the mouse cursor is transferred to a dynamic layer in the frontend code. This works rather well, but some polygons appear to be split by Mapbox. No such split is visible when opening the same layer in QGIS.
The screenshot below shows the problem. It might be a coincidence, but the label (seperate layer in Mapbox studio) appears to be also split along this line.
It seems that this splitting could be processing or rendering artifact. Maybe the feature is not copied completly, but just the part of the feature that is on the current rendering tile. Is there any way this can be prevented?

I use this method to create the hover effect, because filtering the polygons is a lot slower and in this case unusably slow. See the discussion here: https://github.com/mapbox/mapbox-gl-js/issues/2492
Edit
Currently my best solution to the problem is to load the part of the feature on the current tile and let the filter on the hover-layer run at the same time. That means that part of the feature is highlighted right away and the full feature after about 3 seconds.
map.addSource("dynamic_source", {
    "type": "geojson",
    "data": {"type": "FeatureCollection", "features": []}
});

map.addLayer({
    "id": "dynamic_layer",
    "type": "fill",
    "source": "dynamic_source",
    "layout": {},
    "paint": {
        "fill-color": "#f00",
        "fill-opacity": 1
    }
});

map.on("mousemove", function(event) {
    var features = this.map.queryRenderedFeatures(event.point, {
        layers: ["MyLayer"]
    });

    // Instant hover, but only part of feature on current tile
    map.getSource("dynamic_source").setData({
        "type": "FeatureCollection",
        "features": features.map(function (f) {
            return {"type": "Feature", "geometry": f.geometry};
        })
    });

    // 3s delayed hover, but complete feature
    if (features.length) {
        map.setFilter("MyLayer", ["==", "name", features[0].properties.name]);
    } else {
        map.setFilter("MyLayer", ["==", "name", ""]);
    }
});

This forked example from the Github issue shows that this example also sometimes leads to incomplete hover higlighting at tile boundaries: https://jsfiddle.net/t2sks4pz/


texts:
I've tried the whole weekend to get the Orfeo Tool Box running in QGIS 2.18.19 (stand-alone version) and today I finally managed it by setting the paths to the bin-folder of the stand-alone version on Orfeo Toolbox 6.4.0. Now it finally shows the menu in the processing tool box and I thought, I finally made it work

BUT

I tried to execute the Unsupervised K-Means Image Classification and it failed, because the parameter for convergence threshold (-ct) doesn't exist in the application. So I ran the otbgui_KMeansClassification.bat in the OSGEO4W Shell and saw, that there really is no parameter for convergence threshold like it used to be in OTB 5.0.

Im not sure if I already made a mistake with how I got the OTB running in QGIS (my way feels kind of awkward and there might be an easier solution). But now Im completely clueless on how I can execute the Image classification.

Maybe somebody knows how to deal with this.

This is the log:

Algorithm Unsupervised KMeans image classification starting...
C:\Users\nasta\Downloads\OTB-6.4.0-Win64\bin\otbcli_KMeansClassification
-in
"C:/Users/nasta/Documents/Uni/Boku/Master/00_data/Sentinel_ARC/IMAGE_2016-09-13_33UXP_122_B02_10m_Id54562101_L2A.jp2"
-ram
128
-ts
100
-nc
5
-maxit
1000
-ct
0.0001
-out
"C:\Users\nasta\AppData\Local\Temp\processing1070d9dbfe9e434b8a708b490a3ff54a\541e075ec4dd496a9ae0f42076b9391e\out.tif"
-outmeans
"C:\Users\nasta\AppData\Local\Temp\processing1070d9dbfe9e434b8a708b490a3ff54a\ffb59a65eb5b424fb22276c3a7ef7046\outmeans.file"
ERROR: Parameter -ct does not exist in the application.
Converting outputs
Loading resulting layers
The following layers were not correctly generated.
Output Image
You can check the log messages to find more information about the execution of the algorithm


texts: Using PostGIS to find the overall difference between two (large) polygon datasets

texts:
I need to get the Difference of an input layer with ~30k polygons and an overlay layer with ~500k polygons (of which ~200k intersect the input layer).
I am trying to learn PostGIS to do this, as it will (hopefully) be much faster than e.g. QGIS. Newbie at PostGIS so probably making stupid mistakes (sorry).
I have found it easy (and fast) to use ST_Difference across the datasets as set out here. E.g.
SELECT COALESCE(ST_Difference(input.geom,overlay.geom),input.geom)
FROM input LEFT JOIN overlay ON ST_Intersects(input.geom,overlay.geom)

However, this returns every pairwise Difference. So rather than returning only those elements of the input layer which are not intersected by any element of the overlay layer, it instead gives multiple copies of the input polygons, stacked on top of each other, each with a different bite out of them corresponding to a single intersecting polygon from the overlay layer.
EDIT: Don't bother to read the rest of my original question, which is all about the bad way I tried to solve this. Skip to the excellent answer from @JGH and @geozelot.
I have tried to solve this using a double-loop approach building on the the second answer here. Essentially:
For each input polygon (first loop):
   Find the subset of the overlay polygons which intersect that input polygon (using ST_Intersects)
      (Second, nested loop) Run ST_Difference recursively with each intersecting overlay polygon
   Replace the input polygon with the result of the recursive ST_Difference

My SQL code attempting to implement this is at the bottom of this post and I know it's broken (I am a complete SQL newbie and have made a stupid mistake I don't understand - EXPLAIN output also given below).
Is this a viable approach (if I can fix the code), or is there a simpler/better way to do this? (It seems like a common enough GIS requirement...?)
My current (broken) code:
create function ST_Layer_Difference(input_geom geometry,overlay_geom geometry) returns geometry as $$
declare
   i integer;
   j integer;
   input_array geometry[];
   output_array geometry[];
   local_overlays_array geometry[];
   temp_geom geometry;
   output_geom geometry;
begin
    input_array := array(select input_geom);
    FOR i IN 1..array_length(input_array,1) LOOP
      temp_geom := input_array[i];
      local_overlays_array := array(select ST_Intersects(input_array[i],overlay_geom));
          FOR j IN 1..array_length(local_overlays_array,1) LOOP
            temp_geom := ST_Difference(temp_geom,local_overlays_array[j]);
          END LOOP;
      output_array[i] := temp_geom;
    END LOOP;
    output_geom := unnest(output_array);
    return output_geom;
end;
$$
LANGUAGE plpgsql;

create table public."output"
as
select ST_Layer_Difference(input.geom,overlay.geom)
from input, overlay;

And the output of explain (showing I've gone wrong on line 1!):
ERROR:  syntax error at or near "function"
LINE 1: ...FERS false, SUMMARY false, SETTINGS false) create function S...
                                                             ^
SQL state: 42601
Character: 132


texts:
I did an image classification using random forest algorithm in R. I got the result but I do not understand what does the confusion matrix shows ? I mean what does the rows and column show?

These are my classes, and the number of pixels for each class is the number of pixels of the training area of that class, so I do not understand what what are the two elements being compared in confusion matrix?


texts: How to manually edit 16 bit GDEM elevation data on a Geotiff file

texts:
I'm building a virtual 3D terrain scenery using both satellite imagery (google maps satellite) and elevation data from GDEM Aster elevation data (geoTIFF format files).

My problem is that the elevation data does not perfectly match the satellite images. For example, on some small coastal zones, the sea actually reaches the top of some cliffs...

So I'm desperately looking for a way to manually correct the GDEM Aster elevation data (geoTIFF format files) for these zones.
I tried using different photo editing software (trial versions and free versions) but none of them lets you set signed or unsigned 16 bit values to individual cells. Even though they are perfectly capable of processing 16 bit greyscale images, they only allow you to set 8 bit values (0-255) in their respective color pickers...
I tried :
Photoshop CS5,
Photoline,,
Picture window Pro
Paint shop pro 10,
Gimp 2

Thank you for any info on any software (Windows) that allows to do this.

P.S. I already use Geotiff tools to save and restore the geographic metadata.


texts: How to provide a default value in an arcpy script?

texts:
How do I provide a default value like recommended in every Arcpy script? This script is exported from ModelBuilder where every output except the last one is set to temporary and managed. I want this path to show to the default geodatabase every user has set on his system regardless of the exact path.
I am new to programming so I want to apologize if this is a very basic task.

    Input__DEM_ohne_Hecke = arcpy.GetParameterAsText(0)
    if Input__DEM_ohne_Hecke == '#' or not Input__DEM_ohne_Hecke:
    Input__DEM_ohne_Hecke = "H:\\Zwischenlager Patriot\\Projektarbeit_6\\Test_Verschattung.gdb\\Hoehen_Untersuchungsgebiet" # provide a default value if unspecified


texts: Time delta between consecutive points of GPX file?

texts:
I'm in a course learning about GIS systems, and we have recently recorded a GPS track. I've exported this track as .gpx and imported the points into QGIS.

Each point has a timestamp. One of the things I'm interested in is finding the temporal resolution of my file, thereby needing some sort of 'delta time' between each pair of points.

So, what I want to know: How to do calculations between different points in a vector layer? I'm not experienced in Python (but MatLab would be possible I guess), so I'd like to do it with the field calculator.


texts:
I'm facing a problem using the v.what.rast comand.

I'm trying to update the vector "canali" which represents some streams, with the slope value for every single point of the lines. The slope values are contained in the raster map "Pendenza".

I then created a column by entering:

v.db.addcol canali column='slope integer'

and that works fine: it creates an empty column named "slope".

Then I run the command line:

v.what.rast vector=canali@PERMANENT raster=Pendenza@PERMANENT column=slope

but it tells me the following:

No record for category 0 in table
1250 categories loaded from table
1 categories loaded from vector
1 categories from vector missing in table
0 duplicate categories in vector
0 records updated
0 update errors

Anyone out there can help me? I'm new to Grass and I don't want to give up!


texts: How do I make lines look "continues" (images inside)?

texts:
Possible Duplicate:
  How can I merge segmented roads into one line?

A TAB file in MapInfo looks like this :

but when I import it to QGIS (and convert to shp) it looks like this :

How can I make the QGIS layer look like the MapInfo layer (I don't know the term of this functionality)


texts: download Landsat data from Google Earth Engine for the Continental United States

texts:
I need to download one image of the Landsat 7 NDVI 32 day from Google earth engine for the continental United States. Just this specific one image from 2005-08-13 to 2005-09-14. Please see my code below:

//Import Landsat Image
var L7200508NDVI = ee.ImageCollection('LANDSAT/LE7_L1T_32DAY_NDVI')
.filterDate('2005-08-13', '2005-09-14')
.filterBounds(region);

print('L7200508NDVI', L7200508NDVI);
Map.addLayer(L7200508NDVI);

// Load a region representing the United States //this doesn't work as my boundary
var region = ee.FeatureCollection('ft:1tdSwUL7MVpOauSgRzqVTOwdfy17KDbw-1d9omPw')
.filter(ee.Filter.eq('Country', 'United States'));
Map.addLayer(region);

// Create a small area to test: geometry representing an export region.
var geometry = ee.Geometry.Rectangle([116.2621, 39.8412, 116.4849, 40.01236]);

// Export the image, specifying scale and region.
Export.image.toDrive({
  image: L7200508NDVI,
  description: 'imageToDriveExample',
  scale: 30,
  region: geometry,
  maxPixels:1e11
});

if I use the test geometry I got the error:
Error: Invalid argument: 'image' must be of type Image.
i am not sure how to create a polygon just right for the Continental United States. or do I need import a polygon from my ArcGIS layer?


texts: How import a shapefile on-the-fly and have GeoServer render it instantly?

texts:
I am tasked with adding one feature to an existing app. It was built with OpenLayers and GeoServer. The task is to import a shapefile on-the-fly and have GeoServer render it instantly. However, I can't find any doc or tutorial.

Since a user can import any shapefile, the geometry type or the data structure can't be determined ahead of time.

Thanks,
Fred


texts: How to save a georeferenced image in jpeg2000 in qgis?

texts:
I'm georeferencing an jpeg2000 image in qgis, but qgis is only capable to save it in geotif, but gdal would be able to save it in jpeg2000. How should i modify this script:
"gdal_translate -of GTiff -gcp 27915.9 20968.3......(ect...) "D:/HTI_NEMETORSZAG_jpeg/HTI_NEMETORSZAG/BIVa_0854/Elsass_lothringen.jp2" "C:/Users/~1/AppData/Local/Temp/Elsass_lothringen.jp2"
gdalwarp -r near -order 2 -co COMPRESS=NONE  "C:/Users/~1/AppData/Local/Temp/Elsass_lothringen.jp2" "D:/HTI_NEMETORSZAG_jpeg/HTI_NEMETORSZAG/BIVa_0854/Elsass_lothringen_ref.tif""

It would be 100gbs in gtif...


texts: Service published but cache generation failed on ArcGIS Server

texts:
I am trying to publish a map service containing a cached raster layer on our ArcGIS Server Enterpise 10.1.
I analyze the map service before publishing and no error.

The error I get when publishing is:

  The service has been published successfully. The cache for the service has been enabled but the cache
  generation failed."

See attached image

Further details:

the GIS Server and the folder are on the same machine, with Windows OS.
we're talking about a virtual machine in the cloud.
the raster layer is created from a mosaic raster in .tif format (size= 1,5 GB).
I use ArcGIS Desktop 10.1
I am trying to create a 5 scale cache.

I found this error on the Esri forum, but it's for ArcGIS Online servers:
http://forums.arcgis.com/threads/60488-Error-Publishing-a-Map-Service-for-ArcMap-to-ArcGIS-Online


texts: Using QGIS Expression String Builder to display selected multiple layers in a map in print composer

texts:
In QGIS, in a Map's item properties, inside print composer, there is an option to "lock layers" (display layers) based on an expression, by using QGIS "expression string builder".

I can use as expression the name of a layer in my QGIS project, and that layer will display in the map in the print composer. I'm sure I am missing some basic operator, but how can I display more than one layer? Which function or operator do I need to use?


texts: Using OSGB (ESPG: 27700) projection in Openlayers with MousePosition

texts:
I am trying to display Ordnance Survey Grid references (ESPG:27700) on a set of BaseMaps with the MousePosition control.

From looking through the OpenLayers samples, and this forum, samples I have managed to figure out most of the code that is needed. I understand that you need to reference the projection and include the proj4.js reference library.  But I think I am missing something as I can get coordinates in ESPG: 4326 and 900913, but not 27700, if I set the "displayProjection:" option accordingly.

Please see my code:

    <meta http-equiv="Content-Type" content="text/html; charset=utf-8">
    <title>Sample Map</title>


texts:
Does there exist any other software than business analyst that can do store-level trade area analysis without much programming overhead? I am particularly interested in:

Creating customer-derived trade areas that represent x% of sales
Location analysis using a Huff-style model with competition and site attributes (like size)

The data that comes with the software is not really an added benefit for me since I already have access to similar products.


texts: Selecting percentage of points from multipoint layer (ArcGis)

texts:
I am currently working with las dataset of small area, and I want to run some quick analysis on this area. Since my PC is not the best, I am looking for tool like "extract 50% of multipoints to new layer" so I could run some analysis on this PC. Putting in another way, I am looking for tool which will decrease the density of multipoint layer.

I know I was using something similar in the past, but I can't remember nor find any clues.

It would rapidly decrease the processing time for analysis.

Using ArcMap 10.2.


texts:
I have a table with many points that need to be associated with street segments to calculate an average value.
should I do a relationship for all points to a street segment that is closest.
I have dealt with several queries but most is for a specific point with coordinates, I do take care for the whole set of points (61000), similar to the "spatial join" for proximity of ArcGIS

Attributes of the streets
id integer NOT NULL DEFAULT nextval ('calles_osm_gid_seq' :: regclass),
objectid integer,     osm_id numeric,     name character varying       (48),     Inventory character varying (16),     type character       varying (16),     oneway smallint,     bridge smallint,     maxspeed     smallint,     shape_leng numeric,     maxvel integer,     MinSpd       smallint,     azimuth numeric,     geom geometry (MultiLineString,       4326),

Point attributes
    "Latitude"numeric (50,10),
   "Longitude" numeric (50,10),
   "Azimuth" integer,
   "Speed" integer,
   "Date" timestamp without time zone,
   "Geometry" geometry (Point, 4326),
   id_calles integer

Any idea?


texts: How do you view the current spatial index grid size for a file geodatabase feature class?

texts:
I am receiving an "The spatial index grid size is invalid." error during a copy/paste of polygon features in ArcMap and I successfully corrected it with AddSpatialIndex. But in an effort to figure out how and why I ended up with invalid grid size, I thought I would start with the actual grid size.

I would like to know what the current grid size is and I don't see a way to view that in ArcCatalog.  Is it possible to view the current grid size of a file geodatabase featureclass?


texts: "ERROR 000840: The value is not a Data Element" while scripting Rename_management

texts:
I am trying to use a script to go through a folder full of geodatabases, open each geodatabase, and rename the feature classes inside the same name as the geodatabase it resides in.

For example, for a folder like this....

MyFolder
   Geodatabase1.gdb
       Points, Polylines
   Geodatabase2.gdb
       Polylines
   Geodatabase3.gdb
       Points

I want...

MyFolder
   Geodatabase1.gdb
       Geodatabase1_Points, Geodatabase1_Polylines
   Geodatabase2.gdb
       Geodatabase2_Polylines
   Geodatabase3.gdb
       Geodatabase3_Points

My code to do this is giving me the error arcgisscripting.ExecuteError:Failed to execute. Parameters are not valid. ERROR 000840: The value is not a Data Element. Failed to execute (Rename). Can anyone see what might be the issue?

Here is my code:

outLocation = r"C:\Path\to\MyFolder"
arcpy.env.workspace = outLocation
MasterGDB = outlocation + 'excludedGDB.gdb'

for fgdb in arcpy.ListFiles():

    if fgdb.endswith('.gdb') and fgdb != MasterGDB:

        arcpy.env.workspace = os.path.join(outLocation, fgdb)
        featureDataset = arcpy.ListDatasets()[0]
        featureClasses = arcpy.ListFeatureClasses(feature_dataset= featureDataset)

        for fc in featureClasses:
                in_data = fc
                out_data = fgdb +'_'  + fc
                data_type= "FeatureClass"
                arcpy.Rename_management(in_data, out_data, data_type)

    else:
        print fgdb + " is not a gdb or has been otherwise excluded."


texts: Installing PostGIS Layer Viewer for pgAdmin III

texts:
I have difficulties installing pgAdmin's plugin PostGIS Layer Viewer. The official repository on GitHub gives this Web site for install procedure, then the former site gives this Web site for installing prerequisites.

I have well installed QGIS 2.8 through OSGeo4W. I have problems executing the python script that says from qgis.core import *, etc. See error message below:

I get the 'No module named...' message for either qgis.core, qgis.gui, PyQt4.QtCore and PyQt4.QtGui.

Could anybody help me?


texts:
what type of license are subject to the following QGIS libraries?

libqgis_analysis.so

libqgis_core.so

libqgis_gui.so

I would like to use the following libraries released in QGIS package in a "close" (not open source) application.
Are there any derived-license isusses ??


texts:
I am developing a plugin for QGIS that loads the layers of a WFS.

For security reasons I would like the user could not know the layer source.
How can I remove the path of the layer in Vector Properties Dialog-> Layer Info -> Layer Source? and in the tooltip of the layer in the Legend?

I have been searching but I can not find a way to modify the Vector Layer Dialog to disable functions or make other modifications.


texts: How to export image collection data with different projection from Google Earth Engine

texts:
I am exporting data from Google Earth Engine using Python for a specific region in a given temporal window.
This what I am doing:
from geetools import batch
data = ee.ImageCollection('MODIS/006/MYD11A1').filter(ee.Filter.date('2019-06-01', '2019-08-01'))\
  .select('LST_Day_1km')
tasks = batch.Export.imagecollection.toDrive(data, 'JunJuly2019', region=geometry, scale=1000)

This allowsme to save .tiff images in my Google Drive in a sinusoidal projection. However I would like to change projection. For instance I would like to implement something like this before exporting.
myData.reproject('EPSG:4326')


texts: Create a seamless raster from multiple polygons?

texts:
I have been attempting to create a total cover map of a 4x4 km target area.

I clipped the target area to my desired size and used the create polygon feature of ArcGIS to define all of the areas as 'meadow', 'roadside', 'uninhabitable' and 'wooded' ( I was very careful to make sure that all of my vertices on my multiple polygon layers connected so that there was no undefined space) .

I am trying to convert this entire 4x4 section into a seamless ASCII grid file so that it can be run by the program Circuitscape.

I have tried converting all of the polygons into individual rasters and running the mosaic rasters tool provided by ArcGIS, but this results in a huge mess.

any advice that can be offered?


texts: Selecting subset of separate bands in GRASS GIS?

texts:
As you know landsat data has separate files for different bands. Creating a
subset is a wise option to work on region of our interest. Therefore, I used
this method to subset my image -

d.zoom || g.region
g.region -p
r.mapcalc *subset*=original

But with this I can subset one band data at a time. Is it somehow possible
to use i.group and region definition to obtain subsets of all bands? Scripting might be an answer but I have no much idea about it!


texts: How to tell if map service is fully vector or raster?

texts:
I'm wondering how to tell if a map service (in web browser) is vector based or raster based.

What I think to understood so far is that Google Maps for example is vector based, because when you zoom in, street names do not get bigger but stay the same, therefore their data is not saved in raster format. Is that correct?
And if that is the case, is Here Maps then definitely raster based? Because when I zoom in on their maps, labels get bigger at first just like when you zoom in on a picture. Only afterwards the new material is loaded and the labels appear smaller and more crisp again.


texts: Trustworthy of decimal place of converted DD from DMS?

texts:
I want -122d 54m 8.41s W (DMS) to be converted as decimal (DD). I tried two different available online converters to do the calculation. Surprisingly outputs are different. It looks like they match up to 7th decimal places.

-122.902336120571
-122.90233611

Out of curiosity, I used open-source dms-conversion JavaScript library that handles conversion to test what value would be printed out. It came out as:

-122.90233611111111

Only common ground from this observation I can find is up to 7th decimal places matches. I am new to this subject area. I am not sure which decimal value I have to trust if I want it to be converted back as original DMS value.

I wonder each on-line converter and external JS library I tried out has different logic when it comes to truncate seconds?

If I need to implement conversion programmatically either using an external library or math formula available on-line such as this: Converting latitude and longitude to decimal values, how far I have to trust in decimal places for DD value?


texts: Calling ArcGIS 10 geoprocessing method in FME 2011

texts:
I am trying to create an attribute index in a feature class in FME 2011. I couldn't find any trasnformer that does it. So I'm thinking may be I can use "PythonCaller" transformer to call "Add Attribute Index" geoprocessing method from ArcGIS. Is it possible?

Here's the syntax I copied from ArcGIS Help menu:

arcpy.AddIndex_management ("counties.shp", "NAME;STATE_FIPS;CNTY_FIPS", "#", "NON_UNIQUE", "NON_ASCENDING")

Thank you in advance.

Sam

Thank you Mark and blah238. So, I intalled FME2012 on my machine and used python26.dll as my custon python interpreter and copied Mark's code. For Entry point, I entered the name of the class "MyIndexingScript". I'm getting the following error:

  Python Exception : Object: Error in executing tool
  Error encountered while calling method `input'
  PythonFactory failed to process feature

import fmeobjects
import arcpy

class MyIndexingScript(object):
   def __init__(self):
       self.featureList = []
   def input(self,feature):
       self.featureList.append(feature)
       arcpy.AddIndex_management (self, "YORK_PARCEL_ID;PARCEL_GIS_ID", "PIndex", "UNIQUE", "ASCENDING")
   def close(self):
       self.pyoutput(feature)


texts: Installing QGIS 3.4 and 2.18 on Ubuntu 18.04

texts:
I see from this older post that it was not possible to install two versions of QGIS on the same Ubuntu box with out using a virtual environment.  Now that we are a bit farther along I am wondering if this is now possible.

I would like to have my current QGIS version (3.4) but I would also like to have QGIS 2.18 installed.  Can I do this using regular installation methods or do I need to create a virtual box for the second QGIS installation?


texts:
I use the python module "shapely" to make Polygons (like below). Can I convert the sharp edges of a polygon easily to round edges? Such that when I have unified some polygons I could issue a command like

make_round_edges(polygon, edge_rounding_radius)

?


texts: ArcGIS Calculate flooding area based on input amount of water

texts:
I have created a DEM based on historical maps from 1816. I'm supposed to calculate the flooded area for certain amounts of water (in m³/s).

I intended to do it with ArcHydro - however, as far as I can see, ArcHydro takes the water gauge at certain points and based on this height it floods the area. The height is what I'm interested in though, to compare it with modern flooding data.

Since there exist no water gauge data for this date I'm at a loss on how to proceed. Is it possible with ArcHydro to set something like an input amount of water, or are there maybe different tools to solve this problem?


texts: GDAL missing functions [QGIS LTR Linux Ubuntu]

texts:
I recently installed QGIS, the 3.10.x A Coruña LTR version following those instructions. I'm running it on Linux Ubuntu 18.04 Bionic Beaver.

GDAL functions (notably in the plugin Processing), like sieve, merge, polygonize... couldn't be used as I got this error each time in the Log (here for the sieve function):

QGIS version: 3.10.4-A Coruña
QGIS code revision: 10bc9f2ff2
Qt version: 5.9.5
GDAL version: 2.2.3
GEOS version: 3.7.1-CAPI-1.11.1 27a5e771
PROJ version: 493
Processing algorithm…
Algorithm 'Sieve' starting…
Input parameters:
{ 'EIGHT_CONNECTEDNESS' : False, 'EXTRA' : '', 'INPUT' : '/home/usr/Desktop/folder/qgis/QGIS-Training-Data-2.0/exercise_data/raster_analysis/all_conditions.tif', 'MASK_LAYER' : None, 'NO_MASK' : False, 'OUTPUT' : 'TEMPORARY_OUTPUT', 'THRESHOLD' : 10 }

GDAL command:
gdal_sieve.py -st 10 -4 -of GTiff /home/usr/Desktop/folder/qgis/QGIS-Training-Data-2.0/exercise_data/raster_analysis/all_conditions.tif /tmp/processing_454f018dd02843ff825e6f2494820e75/e60f78775c294aafb296c078daaf2f48/OUTPUT.tif
GDAL command output:
/bin/sh: 1: gdal_sieve.py: not found

Execution completed in 0.08 seconds
Results:
{'OUTPUT': '/tmp/processing_454f018dd02843ff825e6f2494820e75/e60f78775c294aafb296c078daaf2f48/OUTPUT.tif'}

Loading resulting layers
The following layers were not correctly generated.

/tmp/processing_454f018dd02843ff825e6f2494820e75/e60f78775c294aafb296c078daaf2f48/OUTPUT.tif

You can check the 'Log Messages Panel' in QGIS main window to find more information about the execution of the algorithm.

So, I focused on the /bin/sh: 1: gdal_sieve.py: not found line to try to solve this problem.

I can't find any way to retrieve these functions from the QGIS UI;
I tried installing the libgdal package from the terminal, but it removes essential packages of QGIS. I had to downgrade the GDAL packages in order to reinstall QGIS...;
I managed to find the missing python scripts and tried to a) adding a new folder path in the folder scripts option to the folder containing the scripts b) adding directly .py scripts to the already existing folder (path: /home/usr/.local/share/QGIS/QGIS3/profiles/default/python).

Both attempts 3. a) and b) made QGIS crash and deactivate the Processing plugin when the software was rebooted.

Are those information sufficient to tell what I'm doing wrong and how I could get those functions?


texts: Statistics on grid cells over GRID format file layer

texts:
I have state wide data here that I would like to overlay a grid layer of 1kmx 1km onto. I used the processing toolbox in QGIS 3 (Create Grid) in order to create a grid layer of 1km x 1km lines under the Vector Creation category. After this, I imported the fmod data and changed the symbology so that unique values are mapped to the different colors in the layer. I tried to use the Zonal Statistics tool under the processing toolbox to gather a majority (mode) statistic on each grid cell by using the following settings:

Raster Layer: modified fmod layer
Raster band: band 1
Vector Layer containing zones: a gpkg layer I made via the previously mentioned Create Grid processing tool.
Output Column Prefix: _ (default)
Statistics to calculate: Majority(mode)

However, I believe (looking at the log) the tool only calculates some aggregate stat that is not what I intended on. I'm not very sure how to go about this in QGIS 3 -- currently downloading QGIS 2, but could anyone tell me what I am doing wrong?

I would like to have a grid layer where each cell in the grid contains a value determined by majority vote from the data in another layer (specifically the fmod data I am using from the link above).


texts:
I'm quite new to PyQGIS and have what is probably a pretty basic question, but I haven't been able to pick it up from the documentation. I have a QgsVectorLayer object from a PostGIS server, and I want to render it in the Map Canvas. How do I go from vector layer object to rendered map?

Thanks!


texts: Strange coordinate shift between small-scale and large-scale Ordnance Survey vector data sets (Strategi vs. OpenMap Local/VectorMap District)

texts:
I'm trying to combine data from Ordnance Survey's small-scale Strategi data set with some of their large-scale OpenData data sets (in particular OS Open Map Local and OS VectorMap District) to fill in local detail. But even though the small-scale data sets are clearly derived from the same underlying Ordnance data, there is a systemic coordinate shift between the data sets, as you can see in the screenshot below. The area shown is a closeup of the West coast of Shetland, British National Grid reference square HP, the layers with the HP_ prefix are from the OS VectorMap District data set, the remaining layers whose data is shifted around 120 meters to the West are from Strategi.

Apart from the obvious simplifications and omissions of the small-scale data set I have not been able to find any information on how or why the data from the two data sets should differ in terms of projection or position, as all Ordnance Survey data sets come with the identical EPSG:27700 projection.

Does anyone have any idea where the coordinate shift might come from? Ideally I would like to be able to correct for it in order to combine both types of data sets in my map. It seems highly unlikely that geometry simplification alone could cause an accidental 120 meter shift of the data.


texts:
I am new in this area. I was practicing the function GetParameterAsText(). I know that when you  run the code (for example in PythonWin) it works like the function raw_input(), therefore it should show a small window where the user can input the required fields.

My problem is that I tried yesterday to use the function GetParameterAsText(), but I have this error, and I think it is about the libraries:

I used the example in the ArcGIS Help:
  http://resources.arcgis.com/en/help/main/10.1/index.html#//018v00000047000000

  *File "C:\Python27\ArcGIS10.2\Lib\site-packages\pythonwin\pywin\framework\scriptutils.py",
  line 322, in RunScript debugger.run(codeObject, main.dict,
  start_stepping=0) File
  "C:\Python27\ArcGIS10.2\Lib\site-packages\pythonwin\pywin\debugger__init__.py",
  line 60, in run GetCurrentDebugger().run(cmd, globals,locals,
  start_stepping)    File
  "C:\Python27\ArcGIS10.2\Lib\site-packages\pythonwin\pywin\debugger\debugger.py",
  line 655, in run exec cmd in globals, locals

  File "C:\WCGIS\Geog485\Lesson1\scripts\practice_4.py", line 1, in
   import os

  File "C:\Program Files
  (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\management.py", line 2429, in
  CopyFeatures raise e ExecuteError: Failed to execute. Parameters are
  not valid. ERROR 000735: Input Features: Value is required ERROR
  000735: Output Feature Class: Value is required Failed to execute
  (CopyFeatures).*


texts:
I am using the function:

traverse(Graph graph, int sourceId, int[] targetIds, float maxCost,                     java.util.Properties params)

and I found out that if I feed more than 144 targetIds it will throw an ArrayOutOfBounds exception, anything below 144 works fine.

Is this an intrinsic limitation of the function or am I doing something wrong?


texts: Legend labels for Vector layer symbols - how to name them (GeoEXT)

texts:
I'm coding with GeoExt and I have some problems trying to display the names of diferent features in a vector layer.

This is my code:

mapPanel = new GeoExt.MapPanel({
    map: map
});

legendPanel = new GeoExt.LegendPanel({
    defaults: {
            labelCls: 'mylabel',
            style: 'padding:5px'
        },
    bodyStyle: 'padding:5px; border:1px solid #8db2e3; border-top-width:0px;border-left-width:0px;border-right-width:0px',
    width: 348,
    autoScroll: true
});

In my legend I have 6 sublayers and all of them are titled with "Untitled" + 1,2,3,4,5, ... I would like to change the name of these different features, removing Untitled 1, 2, 3, ...

Does anybody can help me to do that? I think it is easy but I don't know how.

Thanks in advance!
Rafael


texts:
I have two rasters Ras1 and Ras2.

I want to change the Ras1.

In Ras1, there are few areas, lines, points I am interested having nodata cells and wrong data cells.(limited to interested areas only not talking about the entire Ras1 image)

For those interested areas of Ras1 having nodata cells and wrong data cells, I want to replace cell values with Ras2.

I have seen few tools replace values, but those are value specific ( ie. Replacing all "10" value cells with "56" in entire image) and conditon tools in arcgis, deals with null/no data values.

Is there any way, I can do Replacing Raster Values in Ras1 with Ras2 of my chosen area in the image??
Is there any way I can use arcgis / envi / erdas for this?
Please guide me.


texts:
I have have a line feature which is split in equal pieces (v.split) of 20m. Now I try to Buffer each
line segment only on the sides with 5m to get equal sized polygons in QGIS 3.1, like shown in the picture below.

Running the buffer function in QGIS with:

{ 'DISSOLVE' : False,
  'DISTANCE' : 5,
  'END_CAP_STYLE' : 2,
  'INPUT' : '[...mypath...]',
  'JOIN_STYLE' : 1,
  'MITER_LIMIT' : 3,
  'OUTPUT' : 'TEMPORARY_OUTPUT',
  'SEGMENTS' : 5 }

creates similar looking result but buffers also the line end which mix up the length defined by the line segment and leads to overlapping polygons (see I)

Using this input:

{ 'DISSOLVE' : False,
  'DISTANCE' : 5,
  'END_CAP_STYLE' : 1,
  'INPUT' : '[...mypath...]',
  'JOIN_STYLE' : 1,
  'MITER_LIMIT' : 2,
  'OUTPUT' : 'TEMPORARY_OUTPUT',
  'SEGMENTS' : 5 }

also leads to overlap and this edgy transition zones in the corner. (see II)

How can I create a result looks like described in the first picture?

Update:

I tried to implement xunilk's answer, but I'm stucking at some point. I have a line feature separated in 20m pieces (except the last is 18m), with an individual id.

but I'm getting this error message:

File "<input>", line 1, in <module>
  File "<string>", line 113, in <module>
TypeError: LineString geometry cannot be converted to a point.
Only Point types are permitted.


texts: Why are my forest areas orange and non-forest areas green in Google Earth Engine?

texts:
I'm doing GEARS' 3rd lab tutorial on Remote Sensing in Google Earth Engine.
I calculated the NDVI of an area in Batchellor, Australia. The result is a black and white image so I added a color palette. I used GEE's palette to add color to the NDVI layer and visualize vegetation.
// Define a palette for the 18 distinct land cover classes.
var igbpPalette = [
  'aec3d4', // water
  '152106', '225129', '369b47', '30eb5b', '387242', // forest
  '6a2325', 'c3aa69', 'b76031', 'd9903d', '91af40',  // shrub, grass
  '111149', // wetlands
  'cdb33b', // croplands
  'cc0013', // urban
  '33280d', // crop mosaic
  'd7cdcc', // snow and ice
  'f7e084', // barren
  '6f6f6f'  // tundra
];

Result is this. Notice the \-shaped chunk of land on the left side, in the colored image, it is green while the surrounding forests are orange.
True color image

With color pallete

Here's my whole code in GEE
var image = ee.Image(sent2
.filterDate("2018-05-01", "2018-05-31")
.filterBounds(park)
.sort("CLOUD_COVERAGE_ASSESSMENT")
.first());
print("Sentinel-2 scene (May):", image);

var trueColor = {
  bands: ["B4", "B3", "B2"],
  min: 0,
  max: 3000
  };

Map.addLayer(image, trueColor, "true color image")

var falseColor = {
  bands: ["B8", "B3", "B2"],
  min: 0,
  max: 3000
  };

Map.addLayer(image, falseColor, "false color image")

var NDVI = image.expression(
  "(NIR - RED) / (NIR + RED)",
  {
    RED: image.select("B4"),
    NIR: image.select("B8"),
    BLUE: image.select("B2")
  });

var igbpPalette = [
  'aec3d4', // water
  '152106', '225129', '369b47', '30eb5b', '387242', // forest
  '6a2325', 'c3aa69', 'b76031', 'd9903d', '91af40',  // shrub, grass
  '111149', // wetlands
  'cdb33b', // croplands
  'cc0013', // urban
  '33280d', // crop mosaic
  'd7cdcc', // snow and ice
  'f7e084', // barren
  '6f6f6f'  // tundra
];

Map.addLayer(NDVI, {min: 0, max: 1, palette: igbpPalette}, "NDVI")


texts:
I wan't to filter features of a Shapefile layer from a user's selection with PyQGIS. When I select less than 30 features with a method like fid=1 or fid=2 or fid=3 ..., everything works fine but when I select over 30 features filter doesn't work and return me all the features of the layer. If I try a filter like fid >= 1 and fid <=150 (for example), it works. But it doesn't help for lists of unfollowing ids. I'm using QGIS 3.10.2.
Is there something wrong with this ? Does filters have size limits ? Is there a better way to filter features ?
# Layer
layer = iface.activeLayer()
features = layer.selectedFeatures()

# List ID of selected features
l = []
for feature in features:
    l.append(feature.id())

# Build arguments for the query
index = 0
query = []
prefix = 'fid='
while index < len(l)-1:
    query.append(prefix)
    query.append(l[index])
    query.append(' or ')
    index+=1

# Last element of the list
query.append(prefix)
query.append(l[index])

# Query
fullQuery = str(query).strip('[]').replace(',', '').replace("'", "")

# Filter
layer.setSubsetString(fullQuery)


texts: GeoTIFF layer display issue on GeoServer

texts: I have published a bioclim GeoTIFF layer on GeoServer but it looks very different from how it renders on QGIS. I am not able to figure what is the issue. Have attached a screenshot of the layer on both GeoServer and QGIS.


texts: Why does ArcGIS for Desktop legend not always show symbol as desired?

texts:
I have a problem with the legends in the ArcMap because the patch does not show the marker symbol inside it well.

Whenever, I move the legend it changes but never to what I want for every patch.

I have tried legend properties, and so on but nothing helped me.

You can see my legend here.


texts:
I have a series of raster data (GeoTIFF) with 5 bands (B / G / R / rededge / NIR).
I would like to calculate different indices by using the gdal raster calculator. I would like to use the batch process at the end.
I don't know what to do, because I can't set the GDAL raster calculator to produce reasonable results.
What's wrong with the following syntax settings:


texts: Error in map(ID=0) -- Cannot Interpret  as a date

texts:
I am refactoring code from the Javascript API to the Python API to aid in exporting. The following code works in the Javascript API:
var banded_images = ee.ImageCollection(ee_dates.map(build_annual_cube))

var sample_from_collection = function(d){
  var startDate = (ee.Date(d).advance(1.0,'year').millis())
  var endDate = ee.Date(d).advance(2.0,'year').millis()

  var pointsInThatYear = carp.filterDate(startDate,endDate)

  var img_in_that_year = ee.Image(banded_images.filterDate(startDate,endDate).first())

  var presence_samp = img_in_that_year.sampleRegions({
                              'collection':pointsInThatYear,
                              'scale': 30,
                              'properties':['Present'],
                              'tileScale':16
                            })

    return presence_samp
}

var presence_samp = ee_dates.map(sample_from_collection)
print(presence_samp.get(1))

Where var banded_images is an image collection containing a series of images with various bands, each embedded with a system:time_start property. var carp is a feature collection filled with location information. What I am trying to do is extract temporally explicit covariates from the banded image.
However, when translating to the python API, I get the following cryptic error message:
Error in map(ID=0): Date: Cannot interpret <Image<[List of bands]>> as a Date.
I have searched far and wide and can't find this specific error code anywhere. Here is my python-translated script:
banded_images = ee.ImageCollection(ee_dates.map(build_annual_cube))

def sample_from_collection(d):
    #Set start and end dates for filtering time dependent predictors (SR, NDVI, Phenology)
      #Advance startDate by 1 to begin with to account for water year (below)
    startDate = (ee.Date(d).advance(1.0,'year').millis())
    endDate = ee.Date(d).advance(2.0,'year').millis()

     # Filter points by dates
    pointsInThatYear = carp.filterDate(startDate,endDate) # <---- CHANGE ACCORDING TO SPECIES OF INTEREST

    # Filter collection by that year
    img_in_that_year = ee.Image(banded_images.filterDate(startDate,endDate).first())

    # Sample that image
    presence_samp = img_in_that_year.sampleRegions(**{
                                'collection':pointsInThatYear,
                                'scale': 30,
                                'properties':['Present'],
                                'tileScale':16
                              })

    return presence_samp

presence_sample = banded_images.map(sample_from_collection)

A link to the Javascript script can be found below, with all relevant assets shared
https://code.earthengine.google.com/094e0b2aa185bcb892d7c0958eea57c6
Edit: I underwent a rather tedious debugging process where I ran the following code for each date in the ee_dates list:
d = ee.Date("2017-01-01")
startDate = (ee.Date(d).advance(1.0,'year').millis())
endDate = ee.Date(d).advance(2.0,'year').millis()

pointsInThatYear = carp.filterDate(startDate,endDate) # <---- CHANGE ACCORDING TO SPECIES OF INTEREST

data = pointsInThatYear.getInfo()

img_in_that_year = ee.Image(banded_images.filterDate(startDate,endDate).first())

# Sample that image
presence_samp = img_in_that_year.sampleRegions(**{
                            'collection':pointsInThatYear,
                            'scale': 30,
                            'properties':['Present'],
                            'tileScale':16
                          })
print(presence_samp)

I was able to run this code with every single year of interest. What is going on?


texts: How to add field attribute values from line features to polygon features they pass through using ArcPy?

texts: I have a set of polygons and I want to add a field with the name of the streets (from another feature class of streets) that are within each polygon.  Using Select Features by Location, I select streets that are within the polygons, but I can't think of a way (using python) to automate the process of populating each polygon with the name of the street that falls within it.  I don't have a matching field for either feature class to base a join on.  Is this possible or is it something that would have to be done manually?


texts: Where do I find the code for gd_calc.py?

texts: I need to create a raster that is the average of a large set of other rasters (many thousands). I've learned from another question that the script gd_calc.py can do a limited version of this, limited to 26 input files. I'd like to modify the script to take an arbitrarily large number of input files, but I don't know where to find it on my system, being new to GRASS, Python, and just about everything else. So, can anyone tell me where to find it?


texts: Change visibility of attribute field/group of field in one layer based on checkbox or value map - QGIS

texts:
i'm using QGIS 2.18.3, in Win 7.

i want to know how to make attribute field or group of fields become visible or invisible based on value of other field

on picture above, i want to make my grup field "info" (red rectang) appear and disappear depend on value in my "status" field (blue rectang), either checkbox or value map.


texts: Trimble RTX GPS NTRIP, DNS failed to resolve IP address

texts:
I'm trying to get the RTX correction service going with my Trimble BX992. It is a CenterPoint RTX correction, for which I can see that I am correctly subscribed. I'm getting down to ~0.2m 1-sigma errors, which is pretty nice, but in the I/O Configuration->Port Configuration page I can see that the NTRIP client's status is: DNS failed to resolve IP address (including a screenshot). I was told that in the I/O Configuration->Port Summary page I should see a green row for the NTRIP Client, with "CMRxe" in the input column; however, I get a yellow row, with nothing in the input column (attached a screenshot).

My device is not connected to the internet at all times, and has a fixed IP address for a local network, which I don't know if it is part of the problem.


texts:
I have a postGIS table:

id    value    geometry
0     10       'g_0_10'
0     20       'g_0_20'
1     10       'g_1_10'
1     20       'g_1_20'
2     10       'g_2_10'
2     20       'g_2_20'
...

Now I have a web application where I want the user to see the intersection of the above geometries for any list of id, grouped by value.

For example, the user wants to see the intersection of geometries for the id list = [x,y,z].

So I want a result that looks like this:

value    geoemtry
10       st_intersection(g_x_10, g_y_10, g_z_10)
20       st_intersection(g_x_20, g_y_20, g_z_20)

Can you write such a query if you don't know in advance the size of the list?

More generaly when trying to build dynamic sql queries from Python, do you have to concatenate strings like :

request = "select * from t where id=" + id + ";"

Or is there a cleaner way ?


texts: Multiple markers with specific infowindow content

texts:
When I'm looping through a list of addresses, I'm able to plot the markers for all the addresses fine. But when I click on any marker, the infowindow content shows data of only the last marker. How do I solve this?
Javscript
var map ={};
map.markers = [];
map.addresses = [
    {
        'line': '2101 K St',
        'ref_no': '160621-000005'
    },
    {
        'line': '2131 K St',
        'ref_no': '170708-000015'
    },
    {
        'line': '2321 K St',
        'ref_no': '170707-000028'
    }
];

.
.
.
map.map_object = new Map("esri_map", {
   basemap: "topo",
   center: [<lat>, <lng>],
   zoom: 12
});

var locator = new Locator("http://geocode.arcgis.com/arcgis/rest/services/World/GeocodeServer");
for(var i = 0; i < map.addresses.length; i++)
{
    var addr = map.addresses[i];
    var params = {
        countryCode: "US",
        maxLocations: 1,
        address: {"SingleLine": addr.line}
    };
    locator.addressToLocations(params, function(candidates){
        locatorDone(candidates, addr);
    });
}

function locatorDone(candidates, addr)
{
    if(candidates.length > 0)
    {
        var candidate = candidates[0];
        var pt = new esri.geometry.Point(candidate.location.x, candidate.location.y, new esri.SpatialReference({ 'wkid': candidate.location.spatialReference.wkid }));
        pt = esri.geometry.geographicToWebMercator(pt);
        var symbol = new esri.symbol.PictureMarkerSymbol("/euf/assets/nebulogic/images/alert-icon.png", 21, 39);
        var html = "
"+addr.line+"
";
        html += "
Ref#: "+addr.ref_no+"
";
        var infoTemplate = new esri.InfoTemplate(addr.ref_no, html);
        var graphic = new esri.Graphic(pt, symbol,'',infoTemplate);
        map.map_object.graphics.add(graphic);
        map.markers.push(graphic);
    }
    else
    {
        console.log('not found');
    }
}

P.S: I've solved similar problems (in case of Google Maps API) by using closures. But I'm not sure how to use that in this case.


texts:
I have a leaflet map that I need to fitBound on a given rectangle. The coordinate for the fit bound will be transfered in JSON via a REST call on a server. We are using GEOJson to encode our geometry

So my question is : what would be the correct way to encode this using geojson? The RC candidate seemed to include a BOX type that would be perfect, but this type is not present in the final specification.

So far I see :

There is a bbox attribute, but it seem to apply to a specific geometry. I couldn't use it on its own.
I could just send a dumb JSON array [-90, -180, 90, 180], but then I lose the advantages of using a standardised format.
I could use the BOX type but it appear to have been dropped.

What should I use to send this information?


texts: Accessing GeoServer WMS service in ArcGIS Runtime SDK for QT application?

texts:
I have been looking finding the documentation on how to add GeoServer WMS services in ArcGIS Runtime SDK for QT applications.

However, I found the sample below:
https://developers.arcgis.com/qt/latest/qml/sample-code/sample-qt-wmslayerurl.htm

Which only works for services coming from ArcGIS Server's WMS service.

Is there any other way of adding the GeoServer WMS service in ArcGIS Runtime SDK for QT application?


texts:
After I updated MySQL from 5.7 to 8, spatial functions have become very slow. Specifically, this query:

SELECT points3_.track_id AS track_id
FROM waypoint points3_
WHERE ST_CONTAINS(ST_POLYFROMTEXT('POLYGON ((34.82041 32.041814, 34.82041 32.158149, 35.47959 32.158149, 35.47959 32.041814, 34.82041 32.041814))'), points3_.point)

now uses a full table scan and takes 6 seconds on 500,000 records, while before it took milliseconds.

I have a spatial index on point column, but it's not activated even if I enter force index (point)


texts: Appropriate pixel size when converting vector to raster

texts:
I am converting vector polygons into raster datasets.

Is there a standard calculation for converting the scale of the vector into an appropriate pixel size (e.g. 25 metres)?

My vector has been derived from Aerial Photo Interpretation (API) at a scale of 1:50000.


texts: least cost path between multiple points: obtain a travelling time matrix based on cost raster

texts:
I'm new to GIS and would like to obtain a matrix with travelling times between multiple facilities. I already have a friction surface with travel times.

The plugin Least Cost Path seems to be very interesting, but unfortunately, you can only enter one start point and I want to find routes between around 150 facilities. Moreover, the algorithm will only return a route between the start point and closest end point. I tried to figure out the code on GitHub, but because of my very limited programming skills, the code is too complex for me.

Another post asks this question: Least cost path between various points QGIS.
The piece of code there might help me further as well, but I don't know what the code does exactly and because I already have my cost raster, I think I should adapt the code. However, I can't ask anything as a comment since I'm a new user.

Does anyone know other ways to obtain such a travelling time matrix or how I can use the two options above?
I don't necessarily need the routes itself, only the related travel times.

kind regards,
Dana


texts: Viewing bands of landsat 8 images in QGIS for NDVI calculation?

texts:
I am trying to calculate NDVI values for sites I sampled with landsat 8 imagery.
The equation is simple : NDVI = (NIR band - Red band)/(NIR band + Red Band)

Now I understand that landsat 8 images have NIR on band 5 and Red on band 4.

I downloaded images from USGS and from there I do not know what to do. I can't find band 5 and 4.

What I have been doing is this :
 - Load a new raster layer with the jpeg files
 - Set the right CRS
 - Look in the properties...
 - Only see band 1, 2 and 3 which are red, green and blue

I figures that the file of which the name ends with .TIR (right grey-ish picture in my image) could be the infrared image?

I do not find any explanation on the web about the bands. I am probably not using the files right.

Can someone help me figure it out?


texts:
I'd like to create a web service where the user can edit some basic information of vector tiles. For example change the category of a street category from 5 to other value (4,3,2,1).

I'm not so interested in editing the geometry.

Would this be possible with mapnik or other tile server? Can you give some advice where to start?

How to include metadata (i.e. street category) in the tiles?
Which javascript libraries can be used to edit such data?


texts: Read coordinates with variable Gauß-Krüger zone (meridian) in QGIS

texts:
Is there a sensible way to visualize coordinates like this as point layer in QGIS?

"X" values starting with 2 are EPSG:31466, 3 is EPSG:31467 and 4 is EPSG:31468.

I didn't try, but i think the same problem will arise with UTM zones.

hauptstadt.txt:
---------------
ID  Name    x   y
1   Kiel    3574000 6021000
2   Hamburg 3567000 5936000
3   Hannover    3551000 5805000
4   Bremen  3487000 5884000
5   Düsseldorf  2556000 5676000
6   Wiesbaden   3446000 5548000
7   Mainz   3447000 5541000
8   Stuttgart   3513000 5405000
9   München 4467000 5334000
10  Saarbrücken 2572000 5457000
11  Berlin  4593000 5822000
12  Schwerin    4461000 5945000
13  Potsdam 4573000 5807000
14  Dresden 5411000 5659000
15  Magdeburg   4475000 5777000
16  Erfurt  4433000 5649000


texts: PostGIS Best Practices with thousands of entries

texts:
I'm creating a PostGIS database which will hold thousands of different shapes (imported from Shapefile, GeoJSON, etc), and they'll stay there and be queried over time.

Do I need to create a different table for each file that's imported?

It seems a bit weird that each "file" should have it's own separate table, but this is perhaps natural?

The data is thousands of different geojson- or Shapefiles, with different properties/columns. We plan to do all kinds of queries on the data, but mostly on one file (which translates to a few rows) at a time.


texts: R: Download a large DEM, change projection, and adjust to smaller scale

texts:
This is a process that takes just a few seconds in GIS software. My attempt to do it in R uses a large amount of memory then fails. Is there something wrong in my code, or is this just something R cannot do? I have read R can work inside Grass, can I use a Grass function from inside R?

library(raster)

# I have many environmental rasters in this format
new_r <- raster(ncol=615, nrow=626, xmn=-156.2, xmx=-154.8, ymn=18.89, ymx=20.30)
res(new_r) <- 0.00225
projection(new_r) <- "+proj=longlat +ellps=GRS80 +datum=NAD83 +no_defs +towgs84=0,0,0"

R> new_r ### not too big with a few hundred cells per side
class       : RasterLayer
dimensions  : 627, 622, 1  (nrow, ncol, nlayers)
ncell       : 389994
resolution  : 0.00225, 0.00225  (x, y)
projection  : +proj=longlat +ellps=GRS80 +datum=NAD83 +no_defs +towgs84=0,0,0
extent      : -156.2, -154.8, 18.89, 20.3  (xmin, xmax, ymin, ymax)
values      : none

# I get the DEM at much higher resolution (zipfile is 182Mb)
zipurl <- "ftp://soest.hawaii.edu/coastal/webftp/Hawaii/dem/Hawaii_DEM.zip"
DEMzip <- download.file(zipurl, destfile = "DEMzip")
unzip("DEMzip", exdir = "HIDEM")
HIDEM <- raster("HIDEM/hawaii_dem")

R> HIDEM ### 10m resolution, file is way too big
class       : RasterLayer
dimensions  : 15067, 13136, 1  (nrow, ncol, nlayers)
ncell       : 197920112
resolution  : 10, 10  (x, y)
projection  : +proj=utm +zone=5 +ellps=GRS80 +datum=NAD83 +units=m +no_defs +towgs84=0,0,0
extent      : 179066, 310426, 2093087, 2243757  (xmin, xmax, ymin, ymax)
values      : HIDEM/hawaii_dem
min value   : 0
max value   : 4200

# the following line fails (after a long time)
new_HIDEM <- projectRaster(HIDEM, new_r)


texts:
I just wonder wich version of QtCreator to install, current is Qt 5.11.1. Here I've read about building QGIS master with Qt 5.9.3 debug build: https://www.itopen.it/building-qgis-master-with-qt-5-9-3-debug-build/

So what would be the correct Qt version to build QGIS-master, or does any Qt > 5.9 work?


texts: Extract LAS tool is not working inside a loop?

texts:
In order to extract points from a .lasd (here, I call it "x.lasd") by 30 individual polygons, I am using the following code for the Extract LAS tool in ArcGIS 10.4:

inCuadriculas = r'C:\Users\...\cuadriculas' #This folder contaings 30 polygons
inLas = r'C:\Users\...\x.lasd' #This is the las dataset that contains 4 las files
outLas = r'C:\Users\...\outLas' #This folder will contain all the 33 .lasd outputs

for las in inCuadriculas: #loop in the folder
    j=2 # variable to rename the future .lasd outputs
    arcpy.ExtractLas_3d(inLas,outLas,"",inCuadriculas,"PROCESS_EXTENT","","MAINTAIN_VLR","REARRANGE_POINTS","NO_COMPUTE_STATS",os.path.join(outLas,"las{0}".format(j)))
    print inCuadriculas # Check the name of the polygons.

I get this error:

  Message   File Name   Line    Position     Traceback
        C:\Users...\clip_las.py    49
      ExtractLas    c:\program files (x86)\arcgis\desktop10.4\ArcPy\arcpy\ddd.py    2003         ExecuteError:
  Failed to execute. Parameters are not valid. ERROR 000732: Extraction
  Boundary: Dataset C:\Users...\cuadriculas does not exist or is not
  supported Failed to execute (ExtractLas).

According to the tool guide in:
http://desktop.arcgis.com/en/arcmap/latest/tools/3d-analyst-toolbox/extract-las.htm, the third parameter (extent) is optional.


texts:
I've been reading here and there that snapping while measuring should work in QGIS, but i can not make it work for the life of me.
https://hub.qgis.org/issues/8537

All the snapping options are set correctly in the snapping option dialog box (and they snap just fine while digitizing), but i don't get a snap while measuring.

Any idea's where i need to look?


texts: How do I perform a raster overlay analysis in ArcObjects?

texts:
Using the "IMultivariateOp.MLClassify" Method, it will create a field call "Count",which can be used for calculating area, as shown in the figure below:

But now I want to calculate the increase area of two rasters, for example a raster in "2000" year and "2013" year.

So I have to analyze the overlay, as shown in the figure below:

How do I execute this analysis in ArcObjects using C#?

Thanks very much.


texts: How do you set subtypes for Collector to create multiple choice dropdowns

texts:
I am trying to set up a Collector map for a wetland area to follow the format the Georgia Adopt-A-Wetland Visual Survey form has. I have been following the Collector setup guide here, but I am having trouble creating the subtypes and knowing what codes to use for the survey format.

For example, I want there to be a dropdown for Weather Conditions, and the options would be: Clear, Cloudy, Rain, & Rain within the last 24 to 48 hours. I understand that this would not really be a location point as far as feature classes go, but it would be a start.

The rest of the domain and subtypes would continue to follow the Visual Survey form.


texts: Filling nodata in a DEM using values from another DEM (in QGIS)

texts:
I have a 50cm Digital Terrain Model (tiff) with some nodata values. Instead of interpolating them, I would like to fill them in with values from a lower resolution (1m) DTM of the same area. I think this can be achieved using the raster calculator, but can't for the life of me work it out.

I'm working in QGIS 2.18.


texts:
I have the ordnance survey postcode polygons, i can get the centre point on the postcode (Using polygon centroids) and all the corner coordinates (using extract nodes)

Does anyone know how to export the coordinates of each postcode boundary to csv so it can be drawn on a google map?


texts: I have a pointlayer and some points overlay and i want to remove them, so that only one point remains. the points that overlay have the same coordinates


texts: Create moving map application from topographic maps

texts:
I'm trying to create a moving map application for an embedded ATOM device.  Data source is a 400MB tif file (500dpi) plus projection and worldfile.  My naive approach was to load that tiff into Global Mapper and export it to "Bing", using 4 zoomlevels. It works fine, but the tilecache uses up way too much diskspace.

Can anyone point me to a better way of achieveing something like this? Would it be better to use some kind of runtime decompression?


texts:
I'm trying to print a layout in QGIS that contain the WMS cartography of IGM, but each time I try there's a different hole in the result.
I discovered that MapProxy can do the trick, downloading the cache before the printing.
So, I ask you how to use MapProxy in QGIS to obtain the result and/or if there is a better way to solve.


texts: Using GRASS i.segment in Python 3 stand alone script on Windows10

texts:
I installed GRASS GIS 7.8.1 on Windows 10 using a OSGeo4W.
I would like to use in a stand alone Python 3 script the i.segment GRASS tool.
I saw some related questions (e.g., Using GRASS GIS in Python standalone scripts and call grass modules via python from outside) but they all relevant to Python 2.7, or Linux, or some old versions of GRASS or asking you to install GRASS (which is already installed on my computer via OSGeo4W). I don't mind using subprocess to exceute the i.segment function or to use it "fully" in my script.
Any easy and straightforward way to achieve that?


texts:
I have the following issue:
I am trying to do a regression analysis between the (z) values of two rasters. What I would ideally like to have is the raster output not as grid but rather a list of x1y1z1 x2y2z2, x3y3z3,...
Is there a way to do this easily?


texts: What does "Error when writing Schema: A unique destination featuretype name could not be generated for" mean?

texts:
I have a bunch of GML datasets which I need to import to Esri featureclasses. I am using the Data Interoperability extension for ArcGIS, which should support this. However when trying to import my features I get the following error:

"Error when writing Schema: A unique destination featuretype name could not be generated for """

Thanks

Here is a snippet of the code:

-
-
- http://www.opengis.net/gml/3.2 http://schemas.opengis.net/gml/3.2.1/gml.xsd" gml:id="zones">
-
-
-
-
  inspireId
  UK_zone

-
-
  eng

-
-
  Teesside Urban Area

-
-
-
-


texts: Connecting spatiallines or SpatialPoints based on minimum distance using R

texts:
Does anybody know of a way to create a single polygon (that will eventually be converted to an Owin) from either a spatialpoints object or a spatiallines object?

Here is an example of each:

Here's some very cut down code for this:

test.verts <- list()
test.verts$x <- c(518, 594, 421, 286, 104, 129, 315)
test.verts$y <- c(273, 608, 1202, 1201, 620, 285, 135)

test.owin <- spatstat::owin(poly = test.verts)

cltest <- create.lines(test.owin, nlines=20, width = 5,spacing=NULL,angle=0)

x0list <- c(cltest$x0)
meanx0 <- lapply(x0list, mean)
y0list <- c(cltest$y0)
meany0 <- lapply(y0list, mean)

x1list <- c(cltest$x1)
meanx1 <- lapply(x1list, mean)
y1list <- c(cltest$y1)
meany1 <- lapply(y1list, mean)

pspowin <- owin(xrange=c(min(unlist(meanx0)),max(unlist(meanx1))),
                yrange=c(min(unlist(meany0)),max(unlist(meany0))))

psplines1 <- as.psp(cltest, window = pspowin)

coords <- psplines1$ends
nc <- nrow(coords)
out <- vector(mode="list", length=nc)

for (i in seq(along=out)) {
  lLine <-Line(cbind(c(coords[i,1], coords[i,3]), c(coords[i,2],
                                                    coords[i,4])))
  out[[i]] <- Lines(list(lLine), ID=i)
}

aix.lines <- SpatialLines(out)

plot(aix.lines) #lines needing to be connected

sflines <- st_as_sf(aix.lines)

pt <- st_cast(sflines, "POINT") #points needing to be connected; can also convert to SpatialPoints from sp

plot(pt)


texts: Problem while using gdal_translate command through command prompt

texts:
I am trying to geo-reference a raw satellite image using the following the code through command prompt:
C:\Users\Desktop\PRODUCT1>gdal_translate -of GTiff -a_srs EPSG:4326 -gcp 2458.0 2464.0. 81.2273 21.2776 -gcp 1.0 4928.0 81.0854 21.1940 -gcp 4917.0 4928.0 81.3161 21.1443 -gcp 4917.0 1.0 81.3694 21.3612 -gcp 1.0 1.0 81.1383 21.4109 Imagery.L-4 tmp.tif

After running the above code I am getting the following:
Input file size is 4917, 4928
0...10...20...30...40...50...60...70...80...90...ERROR 3: Read of 4917 bytes for scanline 4928 failed.

ERROR 1: Imagery.L-4, band 3: IReadBlock failed at X offset 0, Y offset 4927
ERROR 1: GetBlockRef failed at X block offset 0, Y block offset 4927

Can somebody tell what could be the issue & how it can be rectified?


texts: ArcMap - Spatial Join changed the values

texts:
I just run spatial join on two big tables I have.
for some reason, one of the columns has changed its values-
 for example, on the original input layer I had this value. :

72145367277699000

the datatype of this is 'real' , column ID_Agri

and on the output join layer it changed to:

72145400000000000

Does anyone know what makes this difference.


texts: Projection Error with Google Maps Basemap

texts:
I have a shapefile of building lots in Manaus, Brazil, that I have re-projected into the WGS84 (EPSG:4326) CRS using the "project" tool in ArcGIS. I am now working in QGIS, and am trying to add the Google Maps XYZ Tile to this layer. However, it appears that these two layers are not projecting properly:

How can I transform the CRS to make the projections match?

My friend was able to fix the problem in ArcGIS by changing the properties of the data frame - is there an analogous operation in QGIS?


texts: How can I create random concave polygons for test purposes?

texts:
In the context of generating random test data using QGIS (in a shapefile):

I can generate random, continuous (adjacent) polygons by first creating n random points using ftools and then generating a Voronoi diagram from those.  Using this technique I get n convex polygons within some extent.  Additionally, I know that all the polygons are contiguous and it is possible to build an adjacency matrix.

Is it possible to generate random convex and concave polygons that are (1) bounded by some defined extent and (2) continuous?

My initial thought is to generate n x 3 polygons and then begin randomly merging them (using a Python script) until only n remain.  This should generate at least a few convex polygons.

Does QGIS (or GRASS) offer a means to do this?  How about an external package?

Edited for @whuber's comment:

It is possible to trace a path from any polygon to any other polygon without exiting a polygon in the set, e.g. the contiguous U.S.  Adding the Hawaii multipolygon would fail the contiguity test because one has to exit California (for example), enter unbounded planar space, and then re-enter the polygon set.
The probability distribution of the random points is not important, in that, the density of polygons generated from the randomly distributed points does not matter.  The goal is simply to vary the number of neighbors each polygon has - the distribution of the neighbor cardinality is not important, but it should not be equal across all polygons as with a regular lattice.
This is a test of regionalization algorithms using randomly generated data  These tests have used lattices and Voronoi diagrams.  I want to explore a particular algorithm on test data using convex and concave polygons.  The generation of convex polygons that meet my constraints is handled using the technique I describe above.

Thanks!


texts: ArcObjects (VB.Net) Python Integration with object passing?

texts:
This link illustrates one method for calling python scripts from within an ArcObjects environment.  The results of these tasks are returned as strings.  I wonder if it is possible to call a python script and pass objects back and forth.

For example, assume I have defined a 2-dimensional array (in VB.Net) of x,y coordinate pairs.  I want to pass this object to a .py script that will perform some computation (doesn't matter what really) and return a dictionary of lists, i.e. {"1":[2,3,4], "2":[1,3,4], "8":[9,11,5]}.


texts: How to add custom Icon to buttons on buttonpads in mapinfo?

texts:
I have made buttons with icons availabe in mapinfo but now i want to add custom icon. Is it possible my code is as follows:

  Include "icons.def"
  Declare Sub Main
  Declare Sub button_prompt
  Sub Main
 Create ButtonPad "Custom" As
PushButton
Icon MI_ICON_ZOOM_QUESTION
Calling button_prompt
HelpMsg "Item1"
fixed

Alter ButtonPad "Custom"
Add Separator
Add PushButton
Icon MI_ICON_GRABBER
Calling button_prompt
HelpMsg "Item2"
fixed
 End Sub
 Sub button_prompt
 End Sub


texts: ogr2ogr: extract text elements from PDF into attributes

texts:
I have PDF created by the local geodetic authority which is kind of geospatial and contains layers of cadastre underground services (water, sewer, etc) as lines, points, and textual attributes describing elevations.

Here are the OGRINFO details:

INFO: Open of `report.pdf'
      using driver `PDF' successful.
Metadata:
  CREATION_DATE=D:20180607081101Z
  CREATOR=Esri ArcSOC 10.2.0.3370
  NEATLINE=POLYGON ((573296.811062859 4797735.35436032,573297.174743531 4797926.39625182,573472.311327169 4797926.0737274,573472.76331659 4797735.04056676,573296.811062859 4797735.35436032))
...
Layer SRS WKT:
PROJCS["ETRS_1989_UTM_Zone_34N",
    GEOGCS["GCS_ETRS_1989",
        DATUM["European_Terrestrial_Reference_System_1989",
            SPHEROID["GRS_1980",6378137.0,298.257222101]],
        PRIMEM["Greenwich",0.0],
        UNIT["Degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["False_Easting",500000.0],
    PARAMETER["False_Northing",0.0],
    PARAMETER["Central_Meridian",21.0],
    PARAMETER["Scale_Factor",0.9996],
    PARAMETER["Latitude_Of_Origin",0.0],
    UNIT["Meter",1.0]]
...

I am using OSGeo4W Shell GDAL 2.2.3, released 2017/11/20 and ogr2ogr command to convert PDF into SQLite:

ogr2ogr -a_srs EPSG:25834 -f sqlite sqlite\out.sqlite report.pdf

And the result in QGIS is OK except that textual elements are not transferred at all.

I even tried to use explicitly sql query with a hope that OGR - Feature Style Specification could extract those missing text as LABEL (this works with AutoCAD DXF files and textual objects in it):

-sql "select geometry, ogr_style from Layers_Points_WN_Anno"

but with no success. ogr_style is only filled with PEN definitions (even if this layer is full with text annotations).

Question is:

Is it possible to transform somehow those text labels in PDF into some objects with OGR2OGR tool so they will show up in converted format ?


texts: Visibility analyses of two objects with ArcMap?

texts:
I'm struggling with a visibility analyses task in ArcMap.

I have two objects: a spire and a WTG (wind turbine generator). Both are point features. The object's heights are deposited in the attributes table (WTG ~200 m, spire ~15m). Surface data is provided by a raster DEM25.

What I want to find out is: Which WTGs are not visible from a certain observation point (spire) (due to the mountain effects of the DEM).

Is the only/best way to do this using the viewshed tool or is there a smarter tool available?


texts: Is it possible to embed python in code for HTML pop-ups in ArcGIS

texts:
When clicking a point feature I want the html pop-up to not only contain the data of the feature but also an image of the area around the feature (eg. 2 km^2). The image would have the base map 'World Imagery and Place' and also another raster layer 'MultiBeam Bathymetry' (which isn't shown in the workspace), the original feature layer wouldn't be shown.

This is far beyond the scope of what html could do and I wanted to use the tools in the arcpy catalog so I was wondering if it is possible to incorporate the python into the html code.


texts: Why is GeoNetwork not displaying any maps?

texts:
I've just installed GeoNetwork 2.6.4 on CentOS 6.2 using MySQL as the database backend.  Now, when I go to the web interface, I see the main web site, but no maps are displayed.  On the geonetwork machine (actually it's a VM) I see them, but over the network I don't.

Can anyone think of anything I should check or change? I set the jetty.xml file to listen on all interfaces and turned iptables off (it was on by default) so that I could access the web site at all.


texts: I'm working with NetLogo now and I have to divide my 'world' into two parts. As you can see in the picture, my world is already divided into two parts by the river. However, my turtles still cross the river as if they can walk on water. Does someone know what I have to do to make sure the turtles can't cross the river anymore and can only stay on their own side?


texts:
I have downloaded GME from spatialecology.com but I can not find a toolbar for it in ArcGIS 10. Should I be able to use GME within ArcGIS? If so, how do I add the extension. If not, how do I run GME?

Thanks, James


texts: QGIS virtual layer with SQL GROUP_CONCAT performance

texts:
I have a polygon shapefile layer, Property, with property ownership parcels and a SQLite database with projects. My goal is to highlight the parcels where an entry in the Projects database and show a label or map tip with the list of project ids for the parcel.

I have created a virtual layer using an SQL query:

SELECT p.objectid, p.pid, GROUP_CONCAT(j.ProjectID) AS ProjectIDs
FROM Property p, Projects j
WHERE j.PID IS NOT NULL AND j.PID = p.pid
GROUP BY p.objectid;

I tested with a small region from the Property layer (including the geometry), and it appears to work well, however, when I use the full Property shapefile it does not. It takes several minutes to add the virtual layer with the listed query and every operation with it requires similar or more time.

The Property layer contains approx 210000 polygons and 11 fields.
The Projects database has approx 7500 records with 14 fields.
The query will result in about 1640 records.
I have a spatial index for the Property layer.

I have tried the above with the Property layer as a spatialite database with similar results.

With sqlite I am able to attach the spatialite Property database to the Projects one and create a table in a fraction of a second using:

CREATE TABLE proj_by_parcel AS
 SELECT p.objectid, p.pid, GROUP_CONCAT(j.ProjectID) AS ProjectIDs
 FROM nsprd.Property p, Projects j
 WHERE j.PID IS NOT NULL AND j.PID = p.pid
 GROUP BY p.objectid;

I have general database and QGIS knowledge.

Is there a different approach that I should be taking in QGIS?

Maybe it will just be best to recreate a static data set periodically.


texts: Raster legend in R: how to colour specific values?

texts:
I have an elevation model plotted in R

r <- raster("example.dem")
plot(r, col = topo.colors(20))

Some of the values are below sea level (0), and I'd like to have those display in red. How can I assign specific ranges to specific colours in the plot()?


texts: Loading a specific county dataset in PostGIS

texts:
I want to load a specific county dataset in PostGIS. I can load a state script using:

SELECT loader_generate_script(ARRAY['CA'], 'linux');

However, I want to load dataset only for Los Angeles County instead of whole California state for performance reason. Geocoding with whole of California data is slow. My assumption is geocoding would be faster if I have only 1 county to work with (I need only LA data).

How can I do that?


texts: QGIS Composer: Fixed Coordinates Distance

texts:
I am working on QGIS 3 project defined to UTM 34S EPSG: 32734 with some layers on WGS84 EPSG 4326.

I need to create a printable map on print composer with grids with DMS coordinates. The grid values (DMS values) must be shown on every exact degree and seconds interval with one minute jump. For instance: 20º34'30.000", 20º35'30.000", 20º36'30.000", 20º37'30.000"... and so on.

Since I need DMS coordinates and the project is using UTM 34S which uses meters, I have setup the grid to use the SRC EPSG: 4326 (WGS84). This way I can choose "Map Units" as Units Interval, and the grid Format as Degrees, Minutes, Seconds with 3 as the coordinate precision.

The problem is I cannot find the exact value to set on the field "Interval" (X/Y), in order to get a continuous value based on the one jump on the minute. The best I could make was to fit the degree and minute: 20º34', 20º35', 20º36'... not the seconds.

As you can see on the above image, degrees and minutes are as desired (with one minute "jump"). The seconds are not fixed and the problem resides here.

To be able to do it, I had found a value "0,016500000" as per above image, but cannot find the exact value that will allow me to fix the seconds issue.


texts:
https://www.worldclim.org/data/monthlywth.html
On this page, I want to download data.
1st issue, when I click on a link, for example tmin_1990-1999, nothing happens.
I solved this using a piece of R code.
url <- "http://biogeo.ucdavis.edu/data/worldclim/v2.1/fut/2.5m/wc2.1_2.5m_tmin_CanESM5_ssp126_2041-2060.zip"
destfile <- "D:/wc2.1_2.5m_tmin_CanESM5_ssp126_2041-2060.zip"
download.file(url, destfile)


texts: Setting fixed scale of map in MapInfo layout?

texts:
How do I set the scale of a map in MapInfo v16 64bit in a simple way?

I can zoom in and out using my mouse but I cannot find how to set the scale to say 1:1000 or 1:2000.


texts: Combine two image collections into one Image collection Earth Engine

texts:
I am very new to GEE and I'm even having difficulties deciphering other people's potential solutions.

After processing each data source into one imageCollection each, I am trying to make a simple calculation between these two. One of the ideas that I had is to combine both imageCollections into one image collection where each image have bands from imageCollection1 and imageCollection2.
I tried the solution in the link, and while the join is successful, it doesn't seem to be able to add the second band. Is it because the metadata of the two imageCollections are different ?

https://code.earthengine.google.com/c0bb0bbdccf975d144b52031321795b6

var GRACE = ee.ImageCollection('NASA/GRACE/MASS_GRIDS/LAND')
                  .filter(ee.Filter.date('2003-08-01', '2003-12-01')).select('lwe_thickness_csr');

print(GRACE, 'GRACE')

var SWS = ee.ImageCollection('NASA/FLDAS/NOAH01/C/GL/M/V001');
SWS = SWS.filterDate('2003-08-01', '2003-12-01').select('SoilMoi100_200cm_tavg');

/// Get information about the GRACE projection.
var graceProjection = GRACE.first().projection();

///Resampling to meet GRACE resolution
var resampleSWS = SWS.map(function(image) {
    var mean = image
      // Force the next reprojection to aggregate instead of resampling.
     .reduceResolution({
        reducer: ee.Reducer.mean(),
       maxPixels: 1024
     }).rename('resampled SWS')
     // Request the data at the scale and projection of the GRACE image.
      .reproject({
       crs: graceProjection
     });
      return mean.set('system:time_start', image.get('system:time_start'))
})

print(resampleSWS, 'resampled SWS')

///Combine two image collections
var mod1 = GRACE.select('lwe_thickness_csr')
var mod2 = resampleSWS.select('resampled SWS')

var filter = ee.Filter.equals({
  leftField: 'system:time_start',
  rightField: 'system:time_start'
});

// Create the join.
var simpleJoin = ee.Join.simple();

// Applt join
var mod1join = ee.ImageCollection(simpleJoin.apply(mod1, mod2, filter))
var mod2join = ee.ImageCollection(simpleJoin.apply(mod2, mod1, filter))

print('Joined', mod1join, mod2join)

I don't see to be able to get this running right. Any suggestions?


texts: QGIS QFieldSync project cannot copy layers in config settings

texts: In this one particular project which I tried to recreate from another project that will not sync with QField but just crashes QGIS. When I try to configure the project for QField packaging I only see the options for offline editing, no action and remove. I have reference layers that I just want to copy but that is not an option. The other (corrupted) project has all the options as expected. What is going on?


texts: salt marsh recession mapping with landsat/WV, IKONOS)

texts:
I am trying to detect timing and pattern of change in a particular salt marsh (low marsh sitting close to a high marsh). Basically, I'm interested in detecting changes in the leading edge for different years. Just wondering if anyone knows how to map salt marsh recession using leading edge.

The leading edge is the front line - like the boundary line that fronts another adjoining vegetation/class.
My final product will be a map showing shift/change in area between different years and a graph showing the changes (after calculating size of polygon in e.g. arcgis). The timing of change will let me identify when exactly the change has occurred/is occurring (using a time slice). The pattern will enable me to determine if change in marsh is gradual and if marsh have periods of both recession and extension over time.


texts: How does a quad-tree cope with very large objects?

texts:
If I use a  quad-tree to index my map data, how does it cope with a  road that goes from one end of the map to the other?

What other sort of indexes should I consider to get round this problem?


texts:
Today I updated my geonetwork to version 4.0, but now I can't login.
I get error 403 - Forbidden.
I'm using tomcat9 and Apache2.
Configuration:
<IfModule mod_ssl.c>
 <VirtualHost *:443>
 ServerName teste@teste.com

 SSLCertificateFile fullchain.pem
 SSLCertificateKeyFile privkey.pem
 Include options-ssl-apache.conf

 ProxyRequests Off

 <Proxy *>
    Order deny,allow
    Allow from all
 </Proxy>

 ProxyPass /geonetwork http://localhost:8080/geonetwork
 ProxyPassReverse /geonetwork http://localhost:8080/geonetwork

 ProxyPreserveHost On


texts:
I have some difficulties to change "NULL" values to zero in all "DOUBLE" fields in feature class. I have tried to make a script tool with the following Python code in ArcGIS 10.5 but it didn't work. No error messages, but nothing changes in table.

# Import system modules
import arcpy

# Set the parameters
InputFeatureClass = arcpy.GetParameterAsText(0)

fieldList = arcpy.ListFields(InputFeatureClass, "DOUBLE")
for field in fieldList:
    expression = "replaceNull(field)"
    codeblock = """def replaceNull(x):
      if x is None:
        return 0
      else:
        return x"""
    arcpy.CalculateField_management(InputFeatureClass, field, expression, "PYTHON", codeblock)

What is wrong? I have looked through many topics here, but couldn't find the solution.


texts: I'm trying to georeference a table by administrative region, but none of my rows are referencing. I have the country name (listed as my 'jurisdiction' column) but no dice. I can't figure out what CartoDB would be interpreting it as. Any ideas?


texts:
I'm trying to add a part of my filename into a field that i'd created, but it comes with the following error:

  Traceback (most recent call last):
        File "Z:\04 - Temp\Charles\NDVI_tempo\NDVI_py3.tbx#Script_NDVIpy3.py", line 50, in

        File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\management.py", line 3019, in
  CalculateField
          raise e
        File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\management.py", line 3016, in
  CalculateField
          retval = convertArcObjectToPythonObject(gp.CalculateField_management(*gp_fixargs((in_table,
  field, expression, expression_type, code_block), True)))
        File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\geoprocessing_base.py", line
  506, in
          return lambda *args: val(*gp_fixargs(args, True))
      arcgisscripting.ExecuteError:  ERROR 000539: Error running expression: test
      Traceback (most recent call last):
        File "", line 1, in
      NameError: name 'test' is not defined

Failed to execute (CalculateField).
Failed to execute (Script).

This is the script that i wrote:

intershp = 'Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL\\intershp.shp'

arcpy.AddField_management (intershp, "Data", "TEXT", "", "", "", "", "", "", "")

tiff_name = tiff.split("_")[0]
field_name = "{}".format(tiff_name[0:4])

arcpy.CalculateField_management (intershp, "Data", field_name, "", "")


texts: DBMS table not found after building geometric network?

texts:
I have created a geometric network within the Water dataset, accessing as GISadmin in SQLServer, in ArcGIS 10.4.1. It builds successfully and creates the desired line-and-point network for use in one version. However, the MXD for a different version with privileges to the same dataset shows red !'s in the ToC for each of the networked features. An attempt to set data source returns the error DBMS table not found.

That same version has properly functioning networks in other datasets with the same privileges. When I remove the water network, the layers fill properly in all versions.

What tables are being referenced in the DBMS error?


texts: Mosaic to New Raster not allowing correct overlap

texts:
I am not sure what exactly is causing hiccups for me - how I processed the data (removing the white edges in 3rd party software), the bit depth (8bit), the file format (png to support transparency), the number of bands (4 for the alpha) or something I am not thinking off.

This is what I want, although it can't be 2 layers.

This is the closest I have come to combining them properly. This is a Mosaic operator of First image added.

I thought this would work since the transparent pixels have values of 255 255 255 0, however the Mosaic Operator of Minimum was way off.

Then I thought - why the heck not and did a Maximum operator as well.

All of these are not usable to me. Can someone please help me understand how to get pixels to not overwrite legitimate data with NoData values while mosaicing? Also I could use any pointers to information about how to remove the collars around imagery like this.

Limitations --- I cannot use the Display Background Value because that removes white from the background of the maps too. I also can't use Mosaic Dataset functionality (oh how I wish!) because of licensing limitations.