texts: Any known problems after Tomcat-upgrade using oracle locator?

texts:
After trying to upgrade from Tomcat 5.5.28 to 5.5.29 or above I ran into a strange problem.
We use a linux based gis system using an Oracle EE DB (without spatial license).

Problem

After the upgrade we get for some actions the following java-error:

java.io.IOException: org.hibernatespatial.helper.FinderException: Couldn't get at the OracleSpatial Connection object from the PreparedStatement.

I was able to pin the cause down to the DBCP library. Up until Tomcat 5.5.28 this was 1.2.2 after that it was upgraded to 1.3 (which is technically the same as 1.4).

The DBCP library is the database abstraction layer between Tomcat and the ojdbc6.jar (as far as I understood).

If I do a downgrade to dbcp 1.2.2 after the Tomcat-Upgrade everything works again without problems.

Has anyone encountered this out there, too? The DBCP-developers can't help, since the java exception is not directly related to dbcp (which is strange, too).


texts: Get Folder path in Processing scripts (QGIS 3)

texts:
maybe a simple question but I cannot solve it.

In a simple Processing script I have defined a parameter as QgsProcessingParameterFolderDestination where I want to save some files (vector or raster) that will be created by the script.

The problem is that I cannot get the path of the chosen folder.

Here a simplified version of the script that actually does not do anything, just printing in the log tab the chosen path.

I'm defining a QgsProcessingParameterFolderDestination and I'm trying to get the path of the chosen folder (if chosen, else Processing will take care of the temporary path) with self.parameterAsString. The

# -*- coding: utf-8 -*-

from PyQt5.QtCore import QCoreApplication
from qgis.core import (QgsProcessing,
                       QgsFeatureSink,
                       QgsProcessingException,
                       QgsProcessingAlgorithm,
                       QgsProcessingParameterFeatureSource,
                       QgsProcessingParameterFolderDestination)
import processing

class ExampleProcessingAlgorithm(QgsProcessingAlgorithm):

    INPUT = 'INPUT'
    OUTPUT = 'OUTPUT'

    def tr(self, string):
        return QCoreApplication.translate('Processing', string)

    def createInstance(self):
        return ExampleProcessingAlgorithm()

    def name(self):
        return 'myscript'

    def displayName(self):
        return self.tr('My Script')

    def group(self):
        return self.tr('Example scripts')

    def groupId(self):
        return 'examplescripts'

    def initAlgorithm(self, config=None):

        self.addParameter(
            QgsProcessingParameterFeatureSource(
                self.INPUT,
                self.tr('Input layer'),
                [QgsProcessing.TypeVectorAnyGeometry]
            )
        )

        self.addParameter(
            QgsProcessingParameterFolderDestination(
                self.OUTPUT,
                self.tr('Output layer')
            )
        )

    def processAlgorithm(self, parameters, context, feedback):

        source = self.parameterAsSource(
            parameters,
            self.INPUT,
            context
        )

        output = self.parameterAsString(
            parameters,
            self.OUTPUT,
            context
        )

        feedback.pushDebugInfo(output)

        return {self.OUTPUT: output}

UPDATE: fixed now: https://github.com/qgis/QGIS/commit/d157b6e66173


texts: Split contents of dbase text field into words at " ", "-" or "_", capitalise first letters, and replace " " with "_"

texts:
A shape file field has text with words separated by blanks " ", dashes "-" and underscores "_". I want to standardise the words to capitalise each word, and then replace spacese (" ") with underscores.
In can do this with one type of word separator (underscores) with
arcpy.CalculateField_management(in_table="xxx.dbf", field="temp", expression="' '.join([i.capitalize() for i in !inField!.split('_', '-')]).replace(' ','_')", expression_type="PYTHON_9.3", code_block="")

but this causes words separated by spaces or dashes to not be capitalised.
I have looked at using the re package but don't know how to call re.split, which should allow splitting on multiple separators.


texts: Add feature to postgis layer with sequenced id using pyqgis

texts:
Using python, I want to add a feature into a postgis layer that has a PK generated by a sequence.

layer (id [PK], field1, field2)

with edit(layer):
            fields = layer.pendingFields()
            new_feat = QgsFeature(fields)
            new_feat[1] = 'text1'
            new_feat[2] = 'text2'

            layer.addFeatures([new_feat])

But I get an error saying the value is too long for the type of the id attribute

Question: How can I assign the next value of the postgis sequence to the id?


texts:
I have three geopandas GeoDataFrame representing a road network, telecommunication wires and a set of network devices to be connected to the network. They are respectively depicted as yellow lines, blue lines and big colored dots in
.

All data are in the same geographic projection (EPSG 3857).

As you can see the device to connect can be near a road, a wire or both, but the networks doesn't precisely overlaps.

I suppose I have to build some sort of graph to compute a path between "Start point" and both "Endpoints". I have to "walk" the road to reach Endpoint 1 , while I can reach Endpoint 2 "by wire".

Howewer road network and wire network are two different set of geometries. How can I geometrically "merge" both networks? Is it necessary to build a graph?


texts: ArcGIS Server Geoprocessing Service fail to run a WeightedOverlay process

texts:
I have a simple python script that runs a Weighted Overlay process.  This script runs fine locally but when I publish it to ArcGIS Server and submit the job through the REST service page, it gave me this error message:

ERROR 999999: Error executing function.
The name of the Field is invalid: valid names may contain letters, numbers or underscores.
Failed to execute (WeightedOverlay).

Then I found out that somehow ArcGIS Server changes the raster field name from "VALUE" (as specified in my script) to "Value".  This is the message from ArcGIS Server when running the process:

Executing: WeightedOverlay "('D:\raster1' 50 'Value' (1 1;NODATA NODATA); 'D:\raster2' 50 'Value' (1 1; 2 1; 3 1; 4 1; 5 1; 6 1; 7 1;NODATA NODATA));1 9 1" ....

But the same script when run from ArcToolbox locally will return:

Executing: WeightedOverlay "('D:\raster1' 50 'VALUE' (1 1;NODATA NODATA); 'D:\raster2' 50 'VALUE' (1 1; 2 1; 3 1; 4 1; 5 1; 6 1; 7 1;NODATA NODATA));1 9 1" ....

I've restarted the ArcGIS Server service and also restarted the server as well, but no luck.

UPDATE:
Here is the code:

sRasterGLA = arcpy.GetParameterAsText(0)

lRasters = []
lRemapObjects = []
lWOTable = []
lOverlayMapping_11NN = [[1,1],["NODATA","NODATA"]]
lOverlayMapping_1NN1 = [[1,"NODATA"],["NODATA",1]]

remap_11NN = arcpy.sa.RemapValue(lOverlayMapping_11NN)
remap_1NN1 = arcpy.sa.RemapValue(lOverlayMapping_1NN1)

iTemp = 1
lRasters.append(sRasterMandatory)  # that's the directory path to the raster data
lRemapObjects.append(remap_11NN)

iTemp += 1
lRasters.append(sRasterGLA)
lRemapObjects.append(remap_1NN1)

iScale1 = int(100/iTemp)
iScale2 = 100 - iScale1 * iTemp
iParam1 = iScale1 + iScale2
iParam2 = iScale1

i = 0
while i < iTemp:
    if i == 0:
        lWOTable.append([lRasters[i], iParam1, "VALUE", lRemapObjects[i]])
    else:
        lWOTable.append([lRasters[i], iParam2, "VALUE", lRemapObjects[i]])
    i += 1

WOTableObj = arcpy.sa.WOTable(lWOTable, lOverlayScaleStep_191)
outWeightedOverlay = arcpy.sa.WeightedOverlay(WOTableObj)
sOutput = os.path.join(arcpy.env.scratchFolder, "WOTEST1")
outWeightedOverlay.save(sOutput)
arcpy.SetParameterAsText(1, sOutput)


texts: Definition Query Current Date SQL Server Native Syntax Error 10.2?

texts:
I have this definition query that worked for a layer from a file geodatabase. Now the same query will not work from the same layer in an SDE. I am getting a 'SQL Server Native Client Syntax' error. Any thoughts?

ExprDate <= (CURRENT_DATE + 30) AND ExprDate >= CURRENT_DATE


texts: Opening/moving focus to layout in ArcGIS Pro using ArcPy

texts: I have an ArcGIS Pro project which contains multiple layouts.  I've got a Python toolbox with a script which, among other things, chooses which of the layouts is appropriate for a specific map.  At the end of the script, I would like the selected layout to open (if it isn't already), and have ArcGIS Pro move the 'focus' to the chosen layout. (i.e. if the layout is open but in a tab which is not in the front, I would like it to move to the front. 'Focus' might not be the right word to describe this.)  I've looked through the ArcGIS Project object, and the layout object, without finding anything which will do this.


texts: Loading Polygon From Postgres (w/o PostGIS) directly into QGIS

texts:
I have a postgres 9.2 install that doesn't have postgis installed. There's a table which contains field name "Zone" of type polygon that I'd like to visualize using QGIS.

When I try to load the table as a PostGIS layer using the Also list tables with no geometry option, the table show up in the legend. When I open up attribute editor, I can see all of the columns except the polygon. In the table view of the DB Manager, I can see the cooridnate list in the Zone column, but am still unable to load a layer from there that has zone as part of the attribute list.

Is it possible to do what I'm trying to do without dumping to an intermediate file type? Maybe some special string in the Layer Source to let me say column "Zone" is the column containing the geometry to plot?

This is another vendor's product, so schema changes and extension installations on the database side are unlikely.


texts: Postgres ST_Distance_Sphere Runtime Improvements

texts:
I need to select all entries from a table, that are within a distance from a point.
Each entry is a point geometry (I set a gist index on geom).

Current runtime is anywhere between 3-15 seconds, depends on how many entries are in a given area (I think), but I'm trying to improve that to under 3sec.

I've seen a few posts on different ways to speed this up, including using ST_DWithin but it does not seem to improve performance, neither does casting to geography.

Any ideas how to proceed?

Here's a sample query:

SELECT  *
FROM    T1
WHERE   entry_date >= current_date - interval '1' month
  AND   ST_DistanceSphere(geom, ST_MakePoint(-74.00072343, 40.7209235)) <= 1000;


texts: how to get (number of) clustered features before map completely rendered?

texts:
When you add features into ol.source.Cluster, you can access original ol.Feature object but clustered ol.Feature which is currently displayed on the map.

check this codepen with console log. it doesn't return clustered features.
I guess calculation of cluster is kind of math so there's a way to get the number of clustered feature at least but couldn't find a way without using event such as postrender or postcompose.


texts:
I have installed ArcGIS Server 10.4 on Windows Server 2008 R2. I am intending to run standalone Python scripts. However when I tried to run the Python scripts, it crashed after importing arcpy.

I have tried uninstalling and reinstall multiple times but it still crashes.
I have only 1 installation (ArcGIS Server 10.4 (64-bit Python), no Desktop).
I have added environment path for PYTHONPATH and .pth file as well.

I would like to understand in my case why it crashes despite only having one version of Python and is there a way to solve the crashing?


texts: How can you align two shapefile layers in QGIS

texts:
I am a very beginner with QGIS.
I would like to know how to align two shapefiles in QGIS. The image below shows the mismatch:


texts: Where can I find historical political affiliation data for the United States?

texts:
I am collaborating on a research paper where one of our variables is political affiliation across the United States. I will be creating a series of historical maps. I'm looking for the percentage of the population that is Republican, Democrat, and Independent in each zip code in the US. If zip code isn't available, county would work too, or statewide if nothing else. I've searched online and have only found 2008 and 2012 election data.

Does anyone know where I could find this historical political affiliation data? I need data from 1980 to present day.


texts: ArcPy how to get field names dynamically from shapefile in customized plugin?

texts:
On selection of a feature layer ArcMap tools display all the fields, Zonal Statistics is one of the example. I have entered my feature as fishnet and it has shown me all the fields in a list below.

I know how to give static drop down list, is there any way to make this list dynamic?

My tool looks like below without dynamic field display


texts: Dividing shapefile into multiple shapefiles using an attribute value in QGIS Graphical Modeler

texts:
I have a shapefile with a column that have project references, I want to create a different shapefile for every project reference there is.
Note that the project references are not unique in the sense that there are multiple lines with the same reference.
Input example:

ID
ref
clmnC

1
001
abcd

2
002
addd

3
001
dsdd

4
003
sdsd

5
003
rzer

6
002
qdqd

7
001
ssqz

Output wanted (3 tables in this example):
  ref clmnC          ref clmnC           ref clmnC
  001 abcd           002 addd            003 sdsd
  001 dsdd           002 qdqd            003 rzer
  001 ssqz

I tried listing unique values but that only return a table with only unique project references values, any solution to solve this problem using the graphical modeler in QGIS 2.18?


texts: Map tiles rendering from local map server to disk

texts:
I am looking for a solution for following use-case:

I need to run a local map server using OpenStreetMap (OSM) data
Map server needs to be able to render a set of map tiles to local directory on disk(TMS scheme)
Input parameters for tile rendering should be GPS bounds of certain area and zoom levels.

Is there a way to achieve something like this using MapServer?


texts: Clipping Multiple Polygons in Single Layer using QGIS?

texts:
I have about 10K polygons for a single layer. Many of these polygons are overlapping. I would like to clip the polygons and either have the overlap merge to one of the polygons or even split in the middle if possible.

There is an example below:

Overlap

Clipped and merged

I've used the Clipper Plugin tool which works great, but is very manual, and would like to see if anyone has found a batch process?


texts: Can Viewshed produce area results per point feature?

texts:
I want to quantify how much land individual point features can "see" using the Viewshed tool. I can generate a complete/universal raster surface, but not one that itemizes the number of cells visible per each individual point feature. I know that I can do this by selecting individual points and then run Viewshed, but given the number of points that I wish to calculate, this would be extremely time consuming.

Does anyone know a way that I can calculate how much land can be seen (cell counts are fine) for individual point features without having to run the calculation individually on every point in a database. I was hoping for an output table with each point feature listed, followed by the number of visible cells (for the point feature). Perhaps a different tool??

Thanks in advance


texts: SQL in CartoDB to select using OR more efficiently

texts:
I made this map a CartoDB map to to show the countries that have banned their citizens working as maids in the Gulf states, I also want to show the Gulf states as one colour.

I created this SQL query for the world borders layer, the one you make available for use to all CartoDB users, in my map:

SELECT * FROM world_borders where name in ('Nepal', 'Kenya','Syrian Arab Republic','Pakistan','Bangladesh')

However, whenever I exit the map to the dashboard and then enter it again for editing, it defaults back to this query:

SELECT * FROM world_borders WHERE name ILIKE '%syria%'

This is a query which came from a filter I had applied in the past. Is the filter over riding the query somehow? I checked to make sure the filter was cleared. I have a free account. SQL queries seem to still be allowed for this level of account. What am I doing wrong? Additionally, how do I add/delete fields from the world borders table I've imported into my map?

Wil


texts: Appending one shapefile to another using ArcGIS Desktop?

texts:
I am using ArcGIS 10.6.1.  I have one shapefile with ~5,600 points, and another with ~4,200,000 points.  I would like to combine these as one, then export as a flat text file/CSV.  I've tried using the "Append" tool , and made sure all fields were named / formatted exactly the same, but it is not working.  It says the schema does not match, even though all fields are exactly the same between the 2 shapefiles.

Is the issue that the Object ID fields added by default share duplicative values between the 2 shapefiles?


texts: Layer's data source is not registered with the server and data will be copied to the server

texts:
I came across this warning message when publishing a service.

When double-click the message above I was asking to register the database on the server but then got an error message saying the database is already registered.

The feature class in my mxd was a new one I imported to the database. Since I do not want to copy data to server, is there anyway I can fix this warning and why this is happening?


texts: Find and list all polygons that overlap with another polygon

texts:
I am working with two polygon layers which overlap but are unrelated (one is congressional district boundaries, the other is MPO area boundaries).

I would like to generate a list of all the polygons in one layer which overlap with each polygon in the second layer.  I.e. I want to generate a list of all the Congressional districts which are at least partially within each MPO area.

Which vector tools or plugins can I use to get there?

I'm using QGIS 1.8 on OSX 10.6.


texts: Leaflet markers do not open popup on click (with markerclusterer)

texts:
I just started using Leaflet and Marker Clusterer to organize the markers.

Problem #1: When an unclustered marker is clicked, no popup appears.

Problem #2: When a cluster is clicked several times, all the markers within that cluster appears, and when one of this marker is clicked, its popup appears! However, after closing the popup by clicking on the map, clicking on any of these clustered markers again do not open any popups!

If I only have 3 unclustered markers, the popup works fine. However, as more markers are added, once a cluster forms, clicking on the marker within any cluster will not cause the popup to open!

Initializing markerclusterer

markers = new L.MarkerClusterGroup();
map.addLayer(markers);

All markers added to markercluster markers

A loop calls on the render function to create the marker and add it to the markerclusterer's array markers. (ignore the backbone.js code)

ListingMarkerView = Backbone.View.extend({
    template: _.template( $('#tpl_ListingMarkerView').html() ),

    render: function() {

        // Create marker
        var content = this.template( this.model.toJSON() );
        var marker = new L.marker(
            [this.model.get('lat'), this.model.get('lng')],
            {content: content});
        marker.bindPopup(content);

        // Add to markerclusterer
        markers.addLayer(marker);
    }
});

Without markerclusterer

If I add the marker directly to map instead of the markerclusterer array markers, the popups work fine, so I guess the problem has something to do with markerclusterer.

Did I do something wrong that resulted in such behavior of the popups? All help appreciated, thanks!


texts: Loading AHN files correctly into MicroDEM?

texts:
In the netherlands we have the ahn elevation data files(to be found on: http://geodata.nationaalgeoregister.nl/ahn2/atom/ahn2_05m_ruw.xml) (based on rd projection).

If i imort these files into microdem by clicking:

file-->open-->open dem

acknowledge the projection/date (amersfoort / rd new | amersfoort)

next, at the pick projection parameters, selecting

-utm zone 31
-northern hemisphere
-wgs 84
- clicking ok

although the coordinate should be in the 31u ...... it is 30nzk.

What am i doing wrong?

I use microdem in order to create Dted files. So if someone has an better option....


texts: How to replicate trunc function on QGIS?

texts:
I've search a bit before asking this question, but I haven't found the answer.
Is there a Trunc function on QGis? I've tried using floor, but it doesn't work (like I want) if the number is negative, for instance:

If

$x=-8.639

Floor will give me

floor($x)=-9

And what I'm looking to achieve is -8 (it's the trunc function).

(I'm passing angles to DMS notation)

I think this may be resolved with a conditional statement for negative values and positive values I also found this expression that makes the job regexp_replace($x, '\\..*', '')... but I guess there should be another (cleaner) way to do it.

I'm using QGis 2.0 on a 64bit Win7 machine.


texts:
I'm trying to create 745 zoning cards (one card per zone) in which different information is displayed to help citizens understand what is allowed or not in their zone. I was able to generate these cards using MS Excel and VBA code but without a map of the zone. I'm trying to find a way to generate these cards with a map of the zone using QGIS or anything you can advise. Here's an example of a zoning card (cards are in French, I hope you still understand what I'm trying to do):

I am trying to get this result:

In summary, I want to find a way to integrate a map of the right zone for the right zoning card. I have to generate 745 of these cards so it has to be automated.


texts: Deleting selected feature in OpenLayers 3?

texts:
I have problems with deleting a selected feature in openlayers 3. Last week I've spend reading this forum and ol API and trying to find an answer to my problem.
I want to draw some polygons, be able to select them (or at least one) and then delete the selected feature. With my code I can draw polygons and select them but I've not been able to delete them. What goes wrong?
I started playing with openlayers 3 two weeks ago. I hope my problem is not a missing semicolon somewhere...

var features = new ol.Collection();
var featureOverlay = new ol.layer.Vector({
        source: new ol.source.Vector({features: features}),
        style: new ol.style.Style({
            fill: new ol.style.Fill({
                color: 'rgba(255, 255, 255, 0.2)'
            }),
            stroke: new ol.style.Stroke({
                color: '#ffcc33',
                width: 2
            }),
            image: new ol.style.Circle({
                radius: 7,
                fill: new ol.style.Fill({
                    color: '#ffcc33'
                })
            })
        })
  });
featureOverlay.setMap(map);
var draw;
var typeSelect = document.getElementById('type');
function addInteraction(){
    if ((typeSelect.value) !== 'None') {
        draw = new ol.interaction.Draw({
            features: features,
            type: (typeSelect.value)
            });
        map.addInteraction(draw);
    }
    else {
        draw.finishDrawing();
    }
}
var select = null;
var selectSingleClick = new ol.interaction.Select({
    condition: ol.events.condition.singleClick
});
var selectElement = document.getElementById('type2');
var selectedFeature =selectSingleClick.getFeatures();
var changeInteraction = function() {
    if (select !== null) {
        map.removeInteraction(select);
        }
    var value = selectElement.value;
    if (value == 'singleclick') {
        select = selectSingleClick;
        }
        else {
            select = null;
        }
    if (select !== null) {
        map.addInteraction(select);
        }
    };
selectElement.onchange = changeInteraction;
changeInteraction();
typeSelect.onchange = function() {
    map.removeInteraction(draw);
    addInteraction();
    };
addInteraction();
function deleteHandler(){featureOverlay.removeFeature(selectedFeature)};
var deleteFeature = function(evt){
    if(evt.keyCode == 46){
        deleteHandler()
        }
    };
document.addEventListener('keydown', deleteFeature, false);
      }


texts: Intersection point of polygon and the line created by two points ArcPy

texts:
I have polygon layer. I selected the most distant vertex from centroid (picture: point A). (How to do that is in the answer here: Select the most distant vertex from polygon's centroid using ArcPy)
What I need to do next is that I have to find intersection point (point B) of the outline of polygon and line created by centroid and the most distant point (line OA).


texts: confused when iterating shp point to extract multi raster data

texts:
I want to use a points shp file to extract some raster files values.
It contain two loops: first is raster files loop, second is shp file points loop.
The main loop code showed as follows:

allyear={}
for filename in os.listdir(inpath):
    if os.path.splitext(filename)[1]=='.tif':
        ras_file=os.path.join(inpath,filename)

        ras_ds=gdal.Open(ras_file,0)
        gt=ras_ds.GetGeoTransform()
        rb=ras_ds.GetRasterBand(1)
        arr=rb.ReadAsArray()
        arr=nodata2nan(arr,-32768)   #this is custom function

        print(filename+' is begining processing')

        oneyearvalues={}

        for feat in lyr:
            geom = feat.GetGeometryRef()
            mx,my=geom.GetX(), geom.GetY()  #coord in map units

            px = int((mx - gt[0]) / gt[1]) #x pixel
            py = int((my - gt[3]) / gt[5]) #y pixel

            pointmeanvalue=bufferpoint(px,py,arr)  #this is custom function
            print(feat['site'],pointmeanvalue,px,py)

            # every feature site is key name and point mean value is key value
            oneyearvalues[str(int(feat['site']))]=pointmeanvalue

         #every raster filename is key name and one extraction is key value
        allyear[os.path.splitext(filename)[0]]=oneyearvalues

total_df=pd.DataFrame.from_dict(allyear)

I checked the printed result,and found that this code just operated the first .tif file. The second for loop just carried out the first raster file,and the others could not be operated. Such as follows:

rasterfile1.tif is begining processing

 ......#processing a short time

rasterfile2.tif is begining processing
rasterfile3.tif is begining processing
rasterfile4.tif is begining processing
rasterfile5.tif is begining processing

What is the problem?


texts:
I want to clip a polyline feature class by a polygon feature class (buffers created from points). Several polygons overlap, so I want to get each polygon feature clipped with each polyline feature.

Using the Clip_analysis tool, I get the orange lines:

whereas I want the following red lines:

In a python script, I tried the following code to clip every polyline with each polygon:

        temp = arcpy.ValidateTableName("clip_temp")
    bufferCursor = arcpy.SearchCursor(buffer)
    n = 0
    for pRow in bufferCursor:
        arcpy.AddMessage(n)
        feat = pRow.Shape
        arcpy.Clip_analysis(inputPolylineLayer, feat, temp, "#")
        if n == 0:
            arcpy.Copy_management(temp, clip)
        else:
            arcpy.Append_management(temp, clip)
        n=n+1

But the execution is very very slow and it crashes before finishing.

I also tried the Intersect_analysis tool but I didn't get the right geometries:

Do you have any idea of a better methods to do it?


texts: C# Arcobjects 10.3: Change the source version on a layer file

texts:
I'm using arcobjects (C#) to load my layers to the table of contents, and I need to update the data source to change the version.

The layer file is set to the correct source at the DEFAULT version. When a 'user' opens the MXD, i need to update the version just before it loads.

I'm not sure where to start with this.
Here is the loading code:

 var groupL = new GroupLayer();
 groupL.Name = ALayerName;

 var fileEntries = Directory.GetFiles(alayerDir, FileFilter).ToList();

 fileEntries.ForEach(delegate(string path)
 {
    var layer = ActiveView.OpenLyrFromFile(path);
    groupL.Add(layer);
 });
 return groupL;

So what to i need to cast 'layer' to in order to give the new version name and then restart the connection?


texts: Line Layer seems to be floating and shifted after adding terrain effect in Cesium

texts:
I had made add a line layer using cesium. In which i had added a terrain effect using the standard cesium terrain provider as below

var dataSource1 = new Cesium.GeoJsonDataSource();
var roads = dataSource1.load('../../SampleData/roads.geojson');
roads.then(function(dataSource1) {
    viewer.dataSources.add(dataSource1);
    viewer.zoomTo(dataSource1);

    //Get the array of entities
    var entities1 = dataSource1.entities.values;

    var colorHash = {};
    for (var i = 0; i < entities1.length; i++) {
        //For each entity, create a random color based on the state name.
        //Some states have multiple entities, so we store the color in a
        //hash so that we use the same color for the entire state.
        var entity = entities1[i];
        var name = entity.name;
        var color = colorHash[name];
        if (!color) {
            color = Cesium.Color.BLACK;
            colorHash[name] = color;
        }
        //Set the polygon material to our random color.
        entity.polyline.material = color;
        //Remove the outlines.
        entity.polyline.outline = false;

    }
}).otherwise(function(error){
    //Display any errrors encountered while loading.
    window.alert(error);
});

After adding the terrain, my geojson line layer seems to be floated and shifted as below

I would like to know how to avoid this floating effect of line layer alone after adding the terrain.


texts:
I have created an enterprise geodatabase in my local instance of SQL Server and I want to create one select statement of multiple select statements within multiple tables and feed this into a view.

I currently have a point feature class (point1) with (name, address, and feature_id, SHAPE).
An example of name (attribute) is Jim.

Related table (table1) that is related to the point feature class.
(feature_id, category, open)  I have managed to create a select statement showing the data I need between both of these. Which is as below:

SELECT POINT1.FEATURE_ID, NAME, TABLE1.CATEGORY, TABLE1.OPEN
FROM POINT1
INNER JOIN TABLE1 ON POINT1.FEATURE_ID=TABLE1.FEATURE_ID

   FEATURE_ID | NAME     | CATEGORY |  OPEN
    ---------------------------------------------
    1         | JIM      | FUNNY    | NO
    2         | DWIGHT   | KNIGHT   | YES
    3         | MICHAEL  | BOSS     | YES

I have two other polygon feature classes.

1.) CITIES: (NAME, NUMBER, SHAPE)
2.) ZIP: (ZIP, SHAPE)

I have also managed to create a select statement as follows. Points that intersect with the polygon.

SELECT P1.NAME, C.NUMBER Number
FROM CITIES C, NAME P1
WHERE C.Shape.STIntersects (P1.SHAPE) = 1

NAME     |  Number
----------------------
JIM      | 212
DWIGHT   | 200
MICHAEL  | 266

The other statement is as follows:

SELECT P1.NAME, Z.ZIP Zip
FROM ZIP Z, NAME P1
WHERE Z.Shape.STIntersects (P1.SHAPE) = 1

NAME     |  ZIP
----------------------
JIM      | 94203
DWIGHT   | 80201
MICHAEL  | 72201

I'm trying to merge these SELECT statements together to feed them in a view to display points with different fields from the related table, cities, and zip.

I am looking for something that looks similar to this ...

    FEATURE_ID | NAME     | CATEGORY |  OPEN |  Number  | ZIP
    -----------------------------------------------------------------
    1          | JIM      | FUNNY    | NO    |  212     | 94203
    2          | DWIGHT   | KNIGHT   | YES   |  200     | 80201
    3          | MICHAEL  | BOSS     | YES   |  266     | 72201


texts: Point from FeatureSet as Input to Python script tool?

texts:
I'm looking to use a point on the map provided through a FeatureSet as an input to both a buffer and Select by Location. I have created the geoprocessing tool that works, but am having difficulty with getting this part to not provide errors.

inPts = arcpy.GetParameterAsText(0)
distance = arcpy.GetParameterAsText(1)
arcpy.SetParameter(0, "sourcePt")
arcpy.MakeFeatureLayer_management("sourcePt", "sourcePt_lyr")

"inPts" is my FeatureSet. Distance is the field is the distance value for my buffer. When I try and make the feature layer, it is telling me that "sourcePt" has no data source.

Is there a step that I'm missing like needing to define geometry?


texts: How to find image resolution data for Open Aerial Map images through the API?

texts:
I'm prepping a data set and would like to download all images from Open Aerial Map that have a ground resolution of say greater than 5cm/pixel.

The web UI displays this data for each image:

But I can't see a way to access it through the API (Note: Link to API docs from main site is currently broken):
https://hotosm.github.io/oam-api/

I can access the metadata for any image, but can't find any kind of resolution data:

{"uuid":"http://oin-hotosm.s3.amazonaws.com/573a17d4cd0663bb003c32b6/0/575fcf142b67227a79b4fbfd.tif","title":"Carretera _ Puente Rio Chico","projection":"GEOGCS[\"WGS 84\",DATUM[\"WGS_1984\",SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]],AUTHORITY[\"EPSG\",\"6326\"]],PRIMEM[\"Greenwich\",0],UNIT[\"degree\",0.0174532925199433],AUTHORITY[\"EPSG\",\"4326\"]]","bbox":[-80.4248105,-0.9801461999990004,-80.418725308,-0.9748216999999999],"footprint":"POLYGON((-80.4248105 -0.9748216999999999,-80.418725308 -0.9748216999999999,-80.418725308 -0.9801461999990004,-80.4248105 -0.9801461999990004,-80.4248105 -0.9748216999999999))","gsd":1.8570457238770968e-7,"file_size":47649982,"acquisition_start":"2016-04-27T05:00:00.000Z","acquisition_end":"2016-04-28T17:00:03.000Z","platform":"uav","provider":"Santiago Pastor","contact":"fruizc,francisco.ruiz@engeotec.net","properties":{"sensor":"Unknow","thumbnail":"http://oin-hotosm.s3.amazonaws.com/573a17d4cd0663bb003c32b6/0/575fcf142b67227a79b4fbfd_thumb.png","tms":"http://tiles.openaerialmap.org/573a17d4cd0663bb003c32b6/0/575fcf142b67227a79b4fbfd/{z}/{x}/{y}.png","wmts":"http://tiles.openaerialmap.org/573a17d4cd0663bb003c32b6/0/575fcf142b67227a79b4fbfd/wmts"},"uploaded_at":"2016-05-16T00:00:00.000Z"}

How do I find image resolution data for Open Aerial Map images through the API?


texts: Problem with creating a model in Qgis -- select by location

texts:
I have two shapefiles (polygones) - referentiebestand and oudejaren

-> the first step is to calculate the centroids of each shapefile (gives a good result),
 then a create a buffer around one of this shapefiles, referentiebestand (gives a good result).

But then i receive a problem: The goal is find the points that lies outside the buffer,
I tried it with the tool 'Select by location' (see figure 1). But it gives nothing...

The complete model (figure 2)

If i use the tool 'select by location' outside the model, in the intergraph
And choose 'new selection' it select the points, so that works well.

But when i choose 'remove selection' instead of new selection the result is empty, nothing changed...

When i try 'new selection' 'add to current selection' or 'remove selection' , nothing changed...


texts:
with reference to Trilateration using 3 latitude and longitude points, and 3 distances.

The accepted answer was votes as shown in wikipedia  in conjunction with a conversion to radians before and after.

Is there a way to incorporate Vincenty's formulae  to this process to bring a more accurate answer?


texts: Distance computation for geospatial clustering of lat/lon points

texts:
I am working with a pretty large dataset (~500K data points) in Python (GeoPandas) and I would like to perform geospatial clustering on some subsets (~60K points) of the data. My question has to do with computing the pairwise distances for the clustering algorithm. I am having difficulty understanding how to compute geospatial distance because I can use haversine distance, assuming a constant radius for the earth; I can compute Vicenty distance; or I can simply use the lat/lon values to compute distance. My preference is the latter, for the sake of speed and computational resources, but I'm not sure how to apply this. Should I compute Euclidean distance between two lat/lon pairs? How does that translate into a surface distance? That is, given a [(lat1 - lat2), (lon1 - lon2)] pair, how does this translate into a "distance of X meters on the surface"?
Sorry, if I'm not really saying any of this correctly; I'm pretty new to geospatial computations.


texts:
I need to make theme map using diagrams like in MapInfo, can I do it in QGIS?
Or can I create a surface using attributes?


texts: How to save area value returned by r.stats in a numpy array?

texts:
I want to write a python script to get and save basin area value of interested points on a drainage line from outside GRASS. I think the workflow could be:

run r.watershed and get drainage file
run r.water.outlet to get drainage area file of interested points
run r.stats to get area value of each interested point and save them in a CSV file

The python script for the first 2 stages ran well, and thus I wrote a python script r_stats trying to run GRASS GIS module r.stats from outside,

import os
import sys
from GRASSGIS_conn import GRASSGIS_conn

def r_stats(flag, map_name, sep):

    g.run_command('r.stats', flags = flag, input = map_name, separator = sep, output = 'dd.csv')

if __name__ == "__main__":

    sys.path.append(os.path.join(os.environ['GISBASE'], 'etc', 'python'))
    import grass.script as g
    gisdb = 'C:\Users\Heinz\Documents\grassdata'
    location = 'nl'
    mapset = 'nl'
    GRASSGIS_conn(gisdb, location, mapset)
    map_name = 'b1'
    flag = 'na'
    sep = 'comma'
    r_stats(flag, map_name, sep)

The resulting CSV file contains a category column and another area value column although I just want area value. Now the problem is that because I need to get area values of more than one raster files loaded in GRASS GIS, and the code above would output too many CSV files and it's inefficient to re-organize them.

Is there a way to store these area value in a numpy array? Or are there better way to get and save area value than r.stats?


texts: Construction of R-tree data structure for tiling polygon vector data

texts:
Because it is better for irregular distributed data, I want to create tiles using r-tree for my dataset. My dataset consists of 1214 polygons which represents building foorprints. After doing some research about r-tree implementations, I could not find answers for some questions. So I don't know where to start. Here are questions which I want to find their answers:

1- Every article that I read mentions minimum and maximum entries (m and M) for nodes of the tree but did not mention how to calculate them. Or if we assume some predefined values, for my data set which values must be m and M?

2- How can we determine the level of the tree?

3- I want to construct rtree and return all child nodes and their bounds such as following example ;

mbr=Rectangle [x1=10.0, y1=4.0, x2=62.0, y2=85.0]
  mbr=Rectangle [x1=28.0, y1=4.0, x2=34.0, y2=85.0]
    entry=Entry [value=2, geometry=Point [x=29.0, y=4.0]]
    entry=Entry [value=1, geometry=Point [x=28.0, y=19.0]]
    entry=Entry [value=4, geometry=Point [x=34.0, y=85.0]]
  mbr=Rectangle [x1=10.0, y1=45.0, x2=62.0, y2=63.0]
    entry=Entry [value=5, geometry=Point [x=62.0, y=45.0]]
    entry=Entry [value=3, geometry=Point [x=10.0, y=63.0]]

Are there any python library for this purpose ??

here is my data set


texts: How to aggregate/resample a higher resolution image to coarser resolution using defined threshold?

texts:
I've got a 30 m classification image and i'd like to convert it to 500 m resolution using a user-defined threshold value.

Basically what i'd like to accomplish is similar to using the Resample tool in ArcGIS based on a Majority filter, with the difference being i'd like to define the majority filter myself so that only when the total area of a particular class in a moving window exceeds 80% of the area of the moving window will that class be used to represent the aggregated cell.

I've tried multiple approaches but none gave me satisfactory results.

I wonder if anyone has encountered similar situation before?


texts:
After adding a layer to my QGIS project with a Python script, I want to show the number of objects beside the layer in the layer window.

I have tried this code, which is from

Show feature count of layer via Python console / PyQGIS

to see if it works in QGIS 3.

## create the memory layer and add to the registry
myLayer = QgsVectorLayer("Point", "myLayer", "memory")
QgsMapLayerRegistry.instance().addMapLayer(myLayer, False)

## reference to the layer tree
root = QgsProject.instance().layerTreeRoot()

## adds the memory layer to the layer node at index 0
myLayerNode = QgsLayerTreeLayer(myLayer)
root.insertChildNode(0, myLayerNode)

## set custom property
myLayerNode.setCustomProperty("showFeatureCount", True)

The module QgsMapLayerRegistry seems to have changed from QGIS 2 to QGIS 3, because when I enter >>>help(QgsMapLayerRegistry), I am getting this error message:

help(QgsMapLayerRegistry)
Traceback (most recent call last):
  File "C:\PROGRA~1\QGIS3~1.4\apps\Python37\lib\code.py", line 90, in runcode
    exec(code, self.locals)
  File "<input>", line 1, in <module>
NameError: name 'QgsMapLayerRegistry' is not defined

I also tried to use the method showFeatureCount() on my freshly created layer, but that didn't work either.

mem_layer = QgsVectorLayer(uri, 'my_file', 'delimitedtext')
        QgsProject.instance().addMapLayer(mem_layer)
        mem_layer.loadNamedStyle(r"P:\my_style.qml")
        mem_layer.showFeatureCount()
        mem_layer.triggerRepaint()

What do I need to change to make the code work?


texts: How do I define fixed class intervalls for raster data with continuous values?

texts:
I have a raster with quantitative data (continuous values from a given minimum to a maximum). In this case its elevation, but it could also be temperature, precipitation, etc. Whenever I want to visualize these kind of raster data in QGIS I run into problems. I know I should use "singleband pseudocolour" and I can use the Interpolation / Mode / Classes options and I can also manually change the class breaks under "Value".
But how can I define an interval to automatically create classes? For example I want class intervals of e.g. 250 (0-250, 25-500, 500-750, 1000-1250, etc.). In ArcGIS it's quite easy to define this interval for ALL classes. There, the amount of classes is automatically calculated from the min / max values and the defined intervall. In QGIS I have to do all the work by myself. Especially changing all the class breaks manually can be a lot of work. In case I want to change the intervall (for example from 250 to 200), I have to do the manual work all over again.
Am I overlooking this option or is it just not possible to define a fixed intervall for all classes in QGIS?


texts: Checking if polygons are properly adjacent to each other?

texts:
I have a whole stack of polygons that were supposed to be properly adjacent (vertex to vertex) to each other, but they are not.

Through the process of cleaning this mess, I would like to visualize/run diagnostics to see which edges I haven't fixed yet.

So I wonder what is the simplest way to do this?

I was joining the polygons and checking which failed to merge, but it is a poor practice (one can overlap the border, not match it).

Probably if I create a set of "duplicated" vertices or edges that would be fine, but I didn't find a straight forward way to do that either?


texts: QGIS: subtypes and domains like in ArcGIS?

texts:
Is there an opportunity to use subtypes and domains in QGIS (Desktop 1.8.0) like the way it is used in ArcGIS?

I am specifically interested in the possibility of dependencies between subtypes and domains: not only a range or coded values for one field is required, but corresponding values.

For example: a subtype defines valid classes of land use categories like

01 - wood
02 - grassland
03 - urbanization
04 - ...

Associated to the subtype there are several domains describing the condition/worth of each category like

wood1 - deforestation
wood2 - conifer forest

grass1 - pasture
grass2 - hayfield
grass3 - grass
grass4 - ...


texts: AET Calculation with MOD16A2 Monthly and Annual?

texts: How can I calculate Actual Evapotranspiration by using Raster Calculator in ArcGIS with MOD16A2 monthly and Annual images?


texts: Can gradients be applied to lines in QGIS?

texts:
Is it possible in QGIS to apply a gradient to linear features? I'd like for some roads to fade from opaque to transparent as they get further away from some set of points. I can measure the distance to any node on the line and thus for each line could define a start and end distance to base the gradient on.

I've approximated the effect I want by splitting my roads into shorter two-node segments and applying color based on the distance to the second node. The effect however is quite discontinuous and there are many more points where semi-transparent lines overlap slightly than there would be if I was using the original linestrings.

Is there any way to make this look a bit better?

(Related questions which I hope are outdated now: Gradient colour for lines, How to visualize measure values along linestrings in QGIS?)


texts: Cannot assign output field name to FieldMap

texts:
In attempting to map input fields to output fields, I am finding that I cannot correctly assign the output field name. Running the script produces no errors, but for a given FieldMap object, the outputField property is assigned the name of the input field. I have mapped fields many times (though never with these specific tables) and never before encountered this problem.

Here is the code:

field_map = arcpy.FieldMappings()
field_map.addTable(in_select)
fm = arcpy.FieldMap()
fm.addInputField(in_select, "CHANGE_STA")
REPORT_STATUS = fm.outputField
fm.outputField = REPORT_STATUS
REPORT_STATUS.name = "REPORT_STATUS"
REPORT_STATUS.type = "TEXT"
field_map.addFieldMap(fm)

When I then call fm.outputField.name, it returns the input field name 'CHANGE_STA' (and the append I run with this field map shows that it has not mapped the fields correctly)

Any ideas? What incredibly obvious thing am I missing? (Input table is a dbf table. Is this a problem with the table type?)


texts: Polygon visibility with overlapping polygons QGIS 3.4

texts: I am marking forestry boundaries over google satellite tiles.  The forest in question  has areas of farmland within the outer boundary.  One area of farmland has a patch of forestry within it's boundary.  I need to see the g sat tiles through rh hatching for forestry and lh hatching for farmland.  I have tried making farm and forestry separate layers but where farmland is I get crossed hatching.  I think I'm asking if it's possible to draw a polygon with a hole in the middle.


texts: How to use postgreSQL to store user credentials in GeoServer?

texts:
Is it possible in GeoServer that instead of using .xml to save all user’s credentials can I use PostgreSQL, and then use those credentials to login?
PS: I tried JDBC driver, but not working. I can see new users getting created in the database, but I'm not able to use those credentials to log in

#geoserver


texts:
I have a vector map that I converted from a raster, representing the rivers in the map. However, because the colours in the map were very diffuse, the resulting vector map has many interrupted lines. Is there a way to extend the lines until they touch another, and then join the touching ones into a single line?
I am using QGIS with GRASS.

Thanks


texts: Create polygon from TIF raster (between 2 elevations)

texts:
I have a raster, and I'm trying to create a polygon and export to a shapefile that is showing one solid layer between lets say 1000m and 3000m in elevation.

So it seems like I need to filter out all the data outside of those elevations, but I'm not sure how to go about that.

How do I do this using GDAL in Python?


texts: Using Snapping Manager with ArcGIS Dynamic Map Service Layer?

texts:
Using ArcGIS API for JavaScript 3.22 I am trying to add Snapping Function on adding Point to a Graphic Layer but bases on an existing ArcGISDynamicMapServiceLayer (snapping to Edges and Vertex in the ArcGISDynamicMapServiceLayer). but looks like the snapping is not functioning.

<button id="addP2P">Add Point To Points </button>
<button id="addP2L">Add Point To Lines </button>

require([
    "esri/map",
    "esri/layers/FeatureLayer",
    "esri/layers/ArcGISDynamicMapServiceLayer",
    "esri/layers/GraphicsLayer",
    "esri/graphic",
    "esri/symbols/SimpleMarkerSymbol",
    "esri/symbols/SimpleLineSymbol",
    "esri/Color",
    "esri/SnappingManager",
    "dojo/keys",
    "dojo/domReady!"
  ],
  function(
    Map,
    FeatureLayer,
    ArcGISDynamicMapServiceLayer,
    GraphicsLayer,
    Graphic,
    SimpleMarkerSymbol,
    SimpleLineSymbol,
    Color,
    SnappingManager,
    dojoKeys
  ) {
    var map = new Map("map", {
      basemap: "topo",
      center: [-88.158805, 41.786075],
      zoom: 18
    });

    var P2PLayer = new GraphicsLayer();
    var P2LLayer = new GraphicsLayer();
    map.addLayer(P2PLayer);
    map.addLayer(P2LLayer);

    var DLayer = new ArcGISDynamicMapServiceLayer('https://sampleserver6.arcgisonline.com/arcgis/rest/services/Water_Network/MapServer');
    map.addLayer(DLayer);

    var symbolPointJFlag = new SimpleMarkerSymbol(SimpleMarkerSymbol.STYLE_SQUARE, 10,
      new SimpleLineSymbol(SimpleLineSymbol.STYLE_SOLID, new Color([0, 255, 0, 1]), 1), new Color([0, 255, 0, 1]));

    $("#addP2P").one("click", function() {
      map.on("click", addGraphic);

      function addGraphic(evt) {
        var sm = new SnappingManager({
          alwaysSnap: false,
          map: map,
          snapKey: dojoKeys.CTRL,
          tolerance: 10,
          layerInfo: [{
            layer: DLayer,
            snapToEdge: true,
            snapToVertex: true
          }]
        });
        var pgraphic = new Graphic(evt.mapPoint, symbolPointJFlag);
        P2PLayer.add(pgraphic);
      }
    });
  });


texts: Python plotly choropleth does not work with GeoJSONs

texts:
I am trying to use plotly choropleth to draw the map, lets say for a random variable of num for each of the feature regions of the map in Italy. However, it does not work. below is the code that I use:

I have downloaded the GeoJSON files for Italy from here.

import random
import pandas as pd
import plotly.express as px
import plotly.io as pio
import json
pio.renderers.default='browser'

with open('it-all.geo.json') as f:
     geojson = json.load(f)

n_provinces = len(geojson['features'])

province_names = [geojson['features'][k]['properties']['name'] for k in range(n_provinces)]

randomlist = []
for i in range(0,110):
    n = random.randint(1,30)
    randomlist.append(n)

datadata = pd.DataFrame({'province':province_names, 'num':randomlist})

fig = px.choropleth(datadata, geojson=geojson, color="num",
                    locations="province", featureidkey="properties.name",
                    color_continuous_scale="Viridis")

fig.show()

What I am getting is a mixed shape map as below. Can anyone let me know what I am doing wrong?


texts: I need to do a viewshed analysis using GRASS. I've got the DEMs loaded; and I simply want to input several (2-4) points (like 41.12345N, 122.34567W from GPS). What is the simplest, most straightforward, way of accomplishing this?


texts: See CRS of all layers at once in QGIS 3.2

texts:
Is there a way to see the CRS of all layers in the layers panel at once in QGIS 3.2.0?
QGIS is running on Windows 7


texts: Adding field that counts occurences of feature layer using ArcPy?

texts:
I have seen :

arcpy.Frequency_analysis(inTable, outTable, frequencyFields, summaryFields)

But This creates a new table with the results. I just want to add another field that will do that to an existing feature layer.

Example:

Desired output

For the NTX layer, counts each row how many times occurs in the dataset.
if SSA occurs 21 times the new column 'count' will be like this:

NTX
name COUNT
SSA   21
SSC   1
SDD   98
...

In pandas this can be:

df['count_']=[df['col1'].value_counts().loc[x] for x in df['col1']]

How to do it in ArcPy?


texts:
I have created future temperature raster (tiff format) in ArcGIS. However I want all pixel's value to be multiply with 10 (constant value). I tried to use con tool and raster calculator in spatial analysis, but it didn't work or may be I set incorrect parameter.

Can anyone help me to solve this problem?


texts: Docking QGIS Attribute Table always on the right instead of bottom

texts: Is it possible to set the default dock position for the attribute table? I have a pretty long list of fields and 1-N relations to be displayed. Now I need to always drag and drop the window and I was wondering if it is possible to avoid it and have it always on the right.


texts: TIN gets distorted using hard clip to edit it?

texts: I have two cross sections from which I have created a TIN as can be seen in the following image. For removing the extra outer portion of the TIN, I created a polygon from the cross-sections and edited the TIN by hard-clipping the polygon, but the resulting TIN is always distorted and some areas are surprisingly missing.


texts:
How can I use the franchise application connected to carto. Using your app does not come access to the carto dataset.
When I want to connect to the web: https: //franchise.carto.io/logged or https://franchise.cloud/app/ I do not see the option of carto


texts:
I am using Mapbox, Leaflet.
Via Mapbox Walking directions I create directions and save origin and destination markers coordinates; Then I want to show this created direction. So I do following:

<script>
        var start = JSON.parse($('#start_field').val());
        var finish = JSON.parse($('#finish_field').val());
        L.mapbox.accessToken = 'qwerty';
        var map = L.mapbox.map('map', 'mapbox.streets', {
            zoomControl: false
        }).setView([42.8580536, 74.6224754], 12);
        L.Routing.control({
            waypoints: [
                L.latLng(start.lat, start.lng),
                L.latLng(finish.lat, finish.lng)
            ]
        }).addTo(map);

        var start_marker = L.marker([start.lat, start.lng], {
            draggable: false
        }).addTo(map);

        var finish_marker = L.marker([finish.lat, finish.lng], {
            draggable: false
        }).addTo(map);
    </script>

But I have an error:

  Uncaught TypeError: Cannot read property 'control' of undefined


texts: Single tile in Leaflet from GeoServer to avoid duplication of labels

texts:
I have struggled with duplication of labels in polygons using WMS in Leaflet.
Is there a workaround that does not involve publishing another point layer based on the polygon with a label?
I am looking for a similar approach as OpenLayers with an SingleTile option. But as far as I have seen there is not.
Using param tiled:false is not enough
https://docs.geoserver.org/latest/en/user/services/wms/vendor.html#tiled
The link mention tilesorigin, but I am not sure what this is or if it is relevant to how the tile/s are rendered?
This is how I have instantiated the layer, but that does not solve the issue with duplicated labels inside polygons.
L.tileLayer.wms(GEOSERVER_BASE_URL, {
    transparent: true,
    format: "image/png",
    layers: "layer",
    tiled: false
});

I have been searching a fair bit to find a solution, but I have not came across a solution that worked/I fully understood to implement.


texts: How can I pass geotransform to gdalwarp to handle ungeoreferenced image as georeferenced

texts:
I'm trying to read ROIs from Landsat-8 overview using GDAL /vsicurl/ , but the image is ungeoreferenced. I have geotransform from full size image, which I can copy. Is there a way to pass geotransform to gdalwarp to handle ungeoreferenced image as georeferenced?

Here is the command I want to execute from Python:

ds_source = 'http://landsat-pds.s3.amazonaws.com/L8/044/034/LC80440342015224LGN00/LC80440342015224LGN00_B3.TIF.ovr'
# here I use -te_srs,-te to crop ROI; data should be taken from 2nd OVR level
warp_options = ' -of MEM -r cubic -te_srs "+proj=longlat +datum=WGS84 +no_defs " -te -123.75 31.952162238 -112.5 40.9798980696 -dstalpha -ovr 2'
output_dataset = gdal.Warp('', '/vsicurl/%s' % ds_source, options=warp_options)

Lines above produce error:

ERROR 1: -te_srs ignored since none of -t_srs, -s_srs is specified or the input dataset has no projection.


texts: Geoserver 2.13.2 cascading WMTS source has wrong source CRS

texts:
I'm trying to cascade this WMTS Service in geoserver 2.13.2.

Adding the source causes no additional problems but when creating layers from this source they get imported with the a source CRS of EPSG:4326. Looking at the WMTS service this is not a supported CRS.

When viewing this layer through the Geoserver's WMS the tiles are misaligned and I think this is caused by a faulty reprojection.

Is there any way to fix this? I've already set the target CRS to a supported one and forced that one to be used.


texts:
Is there a command or interface where I can restart GeoNode or reload GeoNode settings?
Is it even necessary or are the settings (in geonode/settings.py) loaded instantly?


texts:
Let's suppose you have a raster dataset (in a folder or in a personal database). How can you detect what map scale it represents (assuming of course that all rasters have same dimensions/resolution)?

Mind you, I want to have this info WITHOUT having to show the rasters on the map before.

EDIT:

Based on Matthew Snape's post, I have come to the following code that works satisfactory. The only thing that I have still hard time to get is the raster's dpi. Is there any hope to find this information embedded somewhere in the db?

public void GetDatasetMapScale(string datasetPath, out double scale)
{
     var _gp = new Geoprocessor();
    _gp.SetEnvironmentValue("workspace", datasetPath);
     var rasters = _gp.ListRasters("", "");
     var raster = rasters.Next(); // get first raster from dataset

     object dt = "";
     var rasterElement = (IDERasterBand) _gp.GetDataElement(raster, ref dt);
     var meanCellHeight = rasterElement.MeanCellHeight;
     var meanCellWidth = rasterElement.MeanCellWidth;

     var myCellSize = (meanCellHeight + meanCellWidth) / 2;
     var dpi = 300; **// HOW CAN DETECT THIS?????**
     var METERS_PER_INCH = 0.0254;

     scale = (dpi * myCellSize) / METERS_PER_INCH;
}


texts: Leaflet request non-existing tiles despite maxNativeZoom being set

texts:
I have this leaflet map with a tileLayer where maxNativeZoom is not working for some reason. If I zoom in further than the zoom level it has tiles for, and then try to pan, it requests non-existing tiles from the server. I can see it in the network console.

Not only is it requesting tiles at zoom level 5 when I set maxNativeZoom: 4, it also request tiles of zoom level 40.
And errorTileUrl isn't working either.

edit
I managed to make the map publicly accessible, so you can see it here
/edit

I copied the code from another leaflet map of mine (where it is working) and only changed things (like paths) that should not interfere with maxNativeZoom in any way.

Here is the relevant code, hopefully one of you can spot what I'm missing:

    w = 4096;
    h = 4096;
    tileDepth = 4;

    /*
    =============================================================================
    // Conversion from (x, y) raster image coordinates to equivalent of latLng
    =============================================================================
    */

    var yx = L.latLng;
    var xy = function(x, y) {
      if (L.Util.isArray(x)) { // When doing xy([x, y]);
        return yx(x[1], x[0]);
      }
      return yx(y, x); // When doing xy(x, y);
    };

    /**
    =============================================================================
    Define basic map parameters and custom CRS by extending CRS.Simple
    =============================================================================
     */

    var minZoom = 2;
    var maxZoom = 12;
    var img = [
      w,
      h
    ];
    var transformationCoefficient = Math.pow(2, tileDepth);

    L.CRS.MySimple = L.extend({}, L.CRS.Simple, {
      //                      coefficients: a      b    c     d
      transformation: new L.Transformation(1/transformationCoefficient, 0, 1/transformationCoefficient, 0) // Compute a and c coefficients so that  tile 0/0/0 is from [0, 0] to [img]
    });

    var bounds = [xy(0, 0), xy(img)];

    /*
    =============================================================================
    Create the map itself
    =============================================================================
    */

    var map = L.map("map", {
      crs: L.CRS.MySimple, // http://leafletjs.com/reference-1.0.3.html#map-crs
      maxBounds: bounds, // http://leafletjs.com/reference-1.0.3.html#map-maxbounds
      minZoom: minZoom,
      maxZoom: maxZoom,
      attributionControl: false,
      zoomControl: false    //disable zoom control at standard position - we'll add a new one later
    }).setView([w/2, h/2], 0);; //this should center the map in our 1024x1024 div

    //add zoom control at new position
    L.control.zoom({
        position:'topright'
    }).addTo(map);

    /*
    =============================================================================
    Import tiles
    =============================================================================
    */

    L.tileLayer(mapPath + 'tiles/{z}/{x}/{y}.png', {
        bounds: bounds, // http://leafletjs.com/reference-1.0.3.html#gridlayer-bounds
        noWrap: true,
        minZoom: minZoom,
        maxZoom: maxZoom,
        maxNativeZoom: tileDepth,   //how deep the actual map tiles go
        errorTileUrl: imgDir + 'error_tile.png',
        attribution: mapAttribution
      }).addTo(map);


texts: Create raster distance around a single point feature class

texts:
Okay, so I'm relatively new to ArcMap and I'm trying to use the Fuzzy membership for suitability test. I have a single point feature class that I want to use for this analysis. I know that in order to be able to apply the Fuzzy membership you need to have a raster input. My question is how do I create a raster distance around my point on order to apply the fuzzy membership after? I tried all the Distance (Spatial Analyst) tools, but none seems to work for this situation.
Help is very much appreciated :)


texts: Obtain coordinates with specific NDVI value using QGIS?

texts:
I made new layer for NDVI.  Now, I have to get coordinates for some range of NDVI.
Let's say I have to get coordinates for NDVI from 0.1 to 0.3. How I can do it?
I saw several YouTube tutorials - but there was nothing about NDVI coordinates.
Should I create new layer for coordinates or some grid?


texts: ”No module named ogr” when installing QGIS plugins

texts:
In QGIS 2.6.1, when I try to install the FlowMapper and QGIS Cloud plugins, I get an error - "The plugin is broken. Python said: no module named ogr" - giving me the only option of removing them. Other plugins work fine.

Because of this, when I open 2.6.1, this message comes up:

Nu s-a putut încărca pluginul FlowMapper datorită unei erori în timpul apelării metodei classFactory()

Traceback (most recent call last):
  File "C:/PROGRA~1/QGISBR~1/apps/qgis/./python\qgis\utils.py", line 208, in startPlugin
    plugins[packageName] = package.classFactory(iface)
  File "C:/Users/Cezar/.qgis2/python/plugins\FlowMapper\__init__.py", line 24, in classFactory
    from flowmapper import FlowMapper
  File "C:/PROGRA~1/QGISBR~1/apps/qgis/./python\qgis\utils.py", line 460, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
  File "C:/Users/Cezar/.qgis2/python/plugins\FlowMapper\flowmapper.py", line 24, in
    import ogr
  File "C:/PROGRA~1/QGISBR~1/apps/qgis/./python\qgis\utils.py", line 460, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
ImportError: No module named ogr

Versiune Python:
2.7.4 (default, Apr  6 2013, 19:54:46) [MSC v.1500 32 bit (Intel)]

Versiune QGIS:
2.6.1-Brighton Brighton, e2a51df

Cale Python: ['C:/PROGRA~1/QGISBR~1/apps/qgis/./python', u'C:/Users/Cezar/.qgis2
/python', u'C:/Users/Cezar/.qgis2/python/plugins', 'C:/PROGRA~1/QGISBR~1/apps/qgis
/./python/plugins', 'C:\\PROGRA~1\\QGISBR~1\\bin\\python27.zip', 'C:\\PROGRA~1
\\QGISBR~1\\apps\\Python27\\DLLs', 'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib',
'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\plat-win', 'C:\\PROGRA~1\\QGISBR~1\\apps
\\Python27\\lib\\lib-tk', 'C:\\PROGRA~1\\QGISBR~1\\bin', 'C:\\PROGRA~1\\QGISBR~1\\apps
\\Python27', 'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages',
'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\PIL', 'C:\\PROGRA~1
\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\jinja2-2.7.2-py2.7.egg', 'C:\\PROGRA~1
\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\markupsafe-0.23-py2.7-win32.egg',
'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\python_dateutil-
2.2-py2.7.egg', 'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\pytz-
2014.2-py2.7.egg', 'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages
\\win32', 'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\win32\\lib',
'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\Pythonwin', 'C:\\PROGRA~1
\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\Shapely-1.2.18-py2.7-win32.egg',
'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\six-1.6.1-py2.7.egg',
'C:\\PROGRA~1\\QGISBR~1\\apps\\Python27\\lib\\site-packages\\wx-2.8-msw-unicode',
'C:\\Users\\Cezar\\.qgis2\\python\\plugins\\DigitizingTools\\tools']

Before this version, I used QGIS 2.4, which I forgot to delete, so I still have it. I tested the plugins there and they work very well.
Are the plugins incompatible with version 2.6.1, or is there another problem?


texts: How do I draw North-South, East-West, Northeast-Southwest, and Northwest-Southeast lines?

texts: Is there a way to draw four line segments through a point that extend in 8 directions based on my coordinate system?


texts: PostGIS strange difference in distance calc from different point defenitions

texts:
I have a working POSTGIS database, loaded Dutch CBS data using shp2pgsql and now I want to map lat,lon coordinates to multipolygons.

All examples I see use lon lat coordinates but if I do so I get huge distance to my geom objects... If I use coordinates from http://epsg.io/28992/map# it works...

SELECT gm_naam,
st_distance(geom, ST_SetSRID(ST_MakePoint(5.112472, 52.008417),28992) ) as dist,
ST_Distance(geom, ST_SetSRID(ST_Point(134989.82, 448191.66),28992) ) as dist2
FROM nlstaging.gem_2013
order by dist2 Asc

The result is like this:

gm_naam dist    dist2
Nieuwegein  465085.92145570647  0.0
IJsselstein 464133.17518803617  2018.4612972214297
Houten  463571.9453450817   2037.1431117615386
Vianen  459472.14358348626  2102.6589540861714
Utrecht 467845.0647365362   2308.3854923630483

So Nieuwegein is the right answer, dist2 is coorect calculation but the coordinates are NOT lon,lat ...

Question:
How can I specify a point using lon lat and have correct distance?


texts: How to display BLOB value image at specific point with maptips

texts:
I am trying to display a photo loaded as a BLOB value inside a SQLite database using maptips. I have searched all over the forums, but no question seems to quite answer my problem.
Basically, I have loaded in a SQLite database into qgis (geo_ref.sqlite), that have various fields in its attributes, and then joined that to a geopackage (ref_geo.gpkg)
It looks as illustrated below:

The ref_geo: skitse layer is a duplicate layer of ref_geo, that has then joined layers with the geo_ref like so:

So you can see to the right in the first figure I showed, the respective image is stored as a BLOB value at each point.
I want to get that image to popup when i hover the mouse cursor over the point using map tips, for example.
So far, I have not even been able to display the image inside QGIS at all. I have tried using actions, map tips, symbology with raster image markers, but none of these seems to work so far.
I am a bit new to QGIS and not acquainted with coding inside QGIS. Any way to get it to display would work for me.


texts: Improve Performance of ST_DWithin - H2GIS

texts:
I have a table seg_point, that I create like this:

CREATE TABLE seg_point (id INT AUTO_INCREMENT PRIMARY KEY, segmentid TEXT, geom GEOMETRY);

And I'm trying to improve the performance of this query that looks for the nearest match within 5 km (5000 m):

SELECT SEGMENTID, ST_AsText(geom) AS geom, ST_Distance(ST_Transform(seg_point.geom, 3857), ST_Transform(ST_GeomFromText('POINT(' || ? || ' ' || ? || ')', 4326), 3857)) AS DISTANCE
FROM seg_point
WHERE ST_DWithin(ST_Transform(seg_point.geom, 3857), ST_Transform(ST_GeomFromText('POINT(' || ? || ' ' || ? || ')', 4326), 3857), 1000)
ORDER BY DISTANCE
limit 1

I saw this post and tried creating a spatial index:

Improve performance on a st_dwithin query (in PostGIS)

I can create the index, but it doesn't seem that h2gis supports CLUSTER.

Any other suggestions for improving the query?


texts: Remove polygon overlap based on id in Python

texts:
I have a Shapefile consists of multiple polygons, where I want to erase the overlap  based on the id.
An example to visualize the problem:

import geopandas as gpd
from shapely.geometry import Polygon

polygon1 = Polygon([[7, 0], [8, 1], [8, 4], [6, 4], [6, 1]])
polygon2 = Polygon([[7, 3], [8, 4], [8, 7], [6, 7], [6, 4]])
polygon3 = Polygon([[7, 6], [8, 7], [8, 9], [6, 9], [6, 7]])

polygondata = gpd.GeoDataFrame()
polygondata['geometry'] = None

polygondata.loc[0, 'geometry'] = polygon1
polygondata.loc[1, 'geometry'] = polygon2
polygondata.loc[2, 'geometry'] = polygon3

This codes gives following output visible on the left side (Exported to QGIS)

Because the overlay is created systematically in one direction i want to erase every polygon from each other based on the id. The code would be:

erase i from i+1 --> erase 0 from 1, erase 1 from 2 ...

and should give result look like in picture 2 (right). In Python I have the data loaded into geopandas geodataframe, but I only know how to perform geometrical operations between geodataframes and not single features in it. I would split the data in single geometries put it in a list and perform the calculation based on a loop, but I'm not sure how it works. Is there maybe a solution in geopandas how to perform this process without splitting the data?

I'm using Windows 10 and Python 3.7


texts:
I am trying to create a temporary point layer in memory using the following python code.

# create layer
vl = QgsVectorLayer("Point", "temporary_points", "memory")
pr = vl.dataProvider()

# add fields
pr.addAttributes([QgsField("Id", QVariant.Int),
                    QgsField("Type",  QVariant.String),
                    QgsField("size", QVariant.Double)])

# add a feature
fet = QgsFeature()
fet.setGeometry(QgsGeometry.fromPoint(QgsPoint(10,10)))
fet.setAttributes(?)
pr.addFeatures([fet])

I would like to set attributes to the label "size" and not for all. Is there any way to achieve this?


texts: With ArcObjects, distinguish between Save Edits with and without closing the current edit session

texts:
Is it possible with ArcObjects to distinguish between the user selecting the Save Edits button and the users selecting the Stop Editing button and then saving edits?

The reason for this is that in the IEditEvents.OnStartEditing event I am acquiring references to certain tables based on the current editor workspace.  These tables are used to acquire configuration parameters for an ArcMap extension and for other purposes throughout the edit session.  When the edit session is closed, I clear the references to those tables and arrays generated from those tables.

The issue this has created is that when the user selects "Save Edits", IEditEvents.OnStopEditing is called followed by IEditEvents.OnStartEditing.  Essentially, ArcMap is starting a new edit session.  Thus, the work done during the OnStartEditing event is rerun, even though the user hasn't ended the current edit session.  There is some overhead to that work that I would like to avoid if possible.

I have a workaround involving setting a shared property to the current edit workspace on the StopEditing event, and comparing that property to the new edit workspace on the StartEdting event. If they are the same then I know the edit workspace hasn't changed and I don't need to rerun the work to get table references and so forth.  It seems a little messy though.  I recall reading something about this at one point but can find no discussion regarding it on this forum or ESRI's.


texts: Configuring MapServer on Ubuntu with Apache

texts:
I am using Ubuntu 18.10 /x64 / Apache2.

I wish to install MapServer on my system. The packages can be obtained from this place: https://launchpad.net/ubuntu/+source/mapserver/7.2.0-5/+build/15292173 or by activating the proposed repository.

When I went to the MapServer website, it asks about compiling the whole thing. I found that it requires PHP5, whereas the 18.10 distribution provided PHP7. There are several inconsistencies in the documentation provided therein (i.e. MapServer website), and the actual packages. For example, the CMAKE utility asks for one libsvgcairo-dev, whereas I cannot locate any repository like that.

Overall the compilation fails.

If I install from the proposed distribution, I cannot find the php_mapscript.so, which is required for the Apache WebServer. Nor can I find any documentation related to launching the MapServer as a CGI.

I do not wish to use OSGEO Live.


texts: Unable to execute any processing tool in QGIS

texts:
I am new to QGIS and want to learn some basic tools. However, after software  installation I cannot execute any tool from the processing toolbox. All the tools is giving me the same error:

The following layers were not correctly generated. You can check the 'Log Messages Panel' in QGIS main window to find more information about the execution of the algorithm.

The only tools which I could successfully perform so far are Contour tool and Clip raster by extent tool from Raster - Extraction tab.

Below is the example log for Line Smoothing tool:

Processing algorithm…
Algorithm 'Line smoothing' starting…
Input parameters:
{ 'ITERATIONS' : 10, 'LINES_IN' : 'D:/Tanaka/Ilandag_Contours.gpkg', 'LINES_OUT' : 'C:/Users/Acer/AppData/Local/Temp/processing_9da93462bf9f4ec8af00f5e24a02cd1e/efcdce475006430596a69e483d14fb24/LINES_OUT.shp', 'METHOD' : 2, 'PRESERVATION' : 10, 'SENSITIVITY' : 3, 'SIGMA' : 2 }

shapes_lines "Line Smoothing" -LINES_IN "C:/Users/Acer/AppData/Local/Temp/processing_9da93462bf9f4ec8af00f5e24a02cd1e/148495a26d744cb1957234b468b4f7c5/LINES_IN.shp" -METHOD 2 -SENSITIVITY 3 -ITERATIONS 10 -PRESERVATION 10.0 -SIGMA 2.0 -LINES_OUT "C:/Users/Acer/AppData/Local/Temp/processing_9da93462bf9f4ec8af00f5e24a02cd1e/efcdce475006430596a69e483d14fb24/LINES_OUT.shp"

C:\Users\Acer\Documents>set SAGA=D:/QGIS/apps\saga-ltr

C:\Users\Acer\Documents>set SAGA_MLB=D:/QGIS/apps\saga-ltr\modules

C:\Users\Acer\Documents>PATH=D:\QGIS\apps\qgis-ltr\bin;D:\QGIS\apps\grass\grass76\lib;D:\QGIS\apps\grass\grass76\bin;D:\QGIS\apps\Python37;D:\QGIS\apps\Python37\Scripts;D:\QGIS\apps\qt5\bin;D:\QGIS\apps\Python27\Scripts;D:\QGIS\bin;C:\Windows\system32;C:\Windows;C:\Windows\system32\WBem;D:\QGIS\apps\Python37\lib\site-packages\pywin32_system32;D:\QGIS\apps\Python37\lib\site-packages\numpy\.libs;;D:/QGIS/apps\saga-ltr;D:/QGIS/apps\saga-ltr\modules

C:\Users\Acer\Documents>call saga_cmd shapes_lines "Line Smoothing" -LINES_IN "C:/Users/Acer/AppData/Local/Temp/processing_9da93462bf9f4ec8af00f5e24a02cd1e/148495a26d744cb1957234b468b4f7c5/LINES_IN.shp" -METHOD 2 -SENSITIVITY 3 -ITERATIONS 10 -PRESERVATION 10.0 -SIGMA 2.0 -LINES_OUT "C:/Users/Acer/AppData/Local/Temp/processing_9da93462bf9f4ec8af00f5e24a02cd1e/efcdce475006430596a69e483d14fb24/LINES_OUT.shp"
____________________________

##### ## ##### ##
### ### ## ###
### # ## ## #### # ##
### ##### ## # #####
##### # ## ##### # ##
____________________________

SAGA Version: 2.3.2 (32 bit)

____________________________
library path: D:\QGIS\apps\saga-ltr\modules\
library name: shapes_lines
library : Lines
tool : Line Smoothing
author : V. Wichmann, J. Engels (c) 2014
processors : 2 [2]
____________________________

Load shapes: C:/Users/Acer/AppData/Local/Temp/processing_9da93462bf9f4ec8af00f5e24a02cd1e/148495a26d744cb1957234b468b4f7c5/LINES_IN.shp...

Parameters

Lines: LINES_IN
Smoothed Lines: Smoothed Lines
Method: Gaussian Filtering
Sigma: 2.000000

Length of the Gaussian window: 15

Error: executing tool [Line Smoothing]

C:\Users\Acer\Documents>exit
Execution completed in 15.05 seconds
Results:
{'LINES_OUT': 'C:/Users/Acer/AppData/Local/Temp/processing_9da93462bf9f4ec8af00f5e24a02cd1e/efcdce475006430596a69e483d14fb24/LINES_OUT.shp'}

Loading resulting layers
The following layers were not correctly generated.

C:/Users/Acer/AppData/Local/Temp/processing_9da93462bf9f4ec8af00f5e24a02cd1e/efcdce475006430596a69e483d14fb24/LINES_OUT.shp

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

I am using QGIS version 3.4.10-Madeira long term release.
Looking forward for your help!


texts: Exporting to GeoTIFF from data view with Data Driven Pages and ArcPy

texts:
I'm trying to work with ArcGIS Desktop with this problem. I need to export a data to geotiff with specific boundaries which are polygons, so i decided to use Data Driven Pages (DDP). Because there is hundreds of them I found a script (batch export to GeoTIFF in dataview) which should make this what i want much quicker:

mxd = arcpy.mapping.MapDocument("CURRENT")
res = arcpy.GetParameterAsText(0)
export = arcpy.GetParameterAsText(1)
preFix = arcpy.GetParameterAsText(2)
df = arcpy.mapping.ListDataFrames(mxd)[0]

for pageNum in range(1, mxd.dataDrivenPages.pageCount + 1):
      mxd.dataDrivenPages.currentPageID = pageNum
      arcpy.AddMessage("Exporting Map " + str(pageNum) + " of " + str(mxd.dataDrivenPages.pageCount))
      pageName = mxd.dataDrivenPages.pageRow.godlo
      arcpy.mapping.ExportToTIFF(mxd, r"C:\xxx\ark_" + preFix + pageName + ".tif",df,
                     df_export_width=10369,
                     df_export_height=6375,
                     resolution = res,
                     world_file=True,
                     geoTIFF_tags=True)
del mxd

But at the end there is a problem with viewing them in arcMap, because some of rasters cover each other. When code is executing the scale in map window is jumping between two values, so I got like on the picture small rasters inside bigger ones.

Maybe someone has the same problem, with getting raster (geotiff) for each polygon?


texts: Dynamically add features to map in ArcGIS Pro 2.7

texts:
I am trying to dynamically add the contents of a feature dataset to a map in Pro v 2.7.
To achieve this, my understanding of the workflow is that I need to list the features in the feature dataset ('ResultsPath' in code snippet below), then convert to layer, save as a layer file, then add layer. The code below runs without errors, but nothing happens--nothing is added to the map and no layer files are generated. I would try addDataFromPath, but since the lyrx is not generating in the save to layer file process I am stuck for the moment. Any pointers?
#Lists the feature classes in the Results Feature Dataset

#Dynamically convert to layer and add to map.

#Code working, but nothing happening?

        m = aprx.listMaps()[0]
        for fc in arcpy.ListFeatureClasses(ResultsPath):
            OutFile = ScratchWS + os.sep + str(fc)
            lf = arcpy.MakeFeatureLayer_management(fc, OutFile, None, None).getOutput(0)
            lyf = arcpy.management.SaveToLayerFile(lf, OutFile)
            m.addLayer(lyf, "AUTO_ARRANGE")

Update: Problem solved, ListFeatureClasses was not referencing the feature dataset properly so was generating an empty list so nothing was happening. The updated code below does list the contents of the feature dataset.
    #Lists features from Results feature dataset
    #Makes feature layer from features, adds feature layers to map

    m = aprx.listMaps()[0]
    datasets = arcpy.ListDatasets(feature_type='feature')
    for ds in datasets:
        ResultList= arcpy.ListFeatureClasses("*", "ALL", ds)
        for fc in ResultList:
            OutFile = os.path.join(wkspc, str(fc))
            lf = arcpy.MakeFeatureLayer_management(fc, fc, None, None)


texts: Adding a layer with specific color to map in ArcPy

texts:
I have a blank map. And I have a shapefile. I would like to add the shapefile to the map document, change its colour to grey and then save the map document.

My plan is to make a layer from the shapefile (MakeFeatureLayer_management), add the resulting layer to the map document (mapping.AddLayer) and then change its colour using ApplySymbologyFromLayer_management.

This is my code so far:

import arcpy, os

### Set directory
arcpy.env.workspace = ...
arcpy.env.overwriteOutput = 1

### Define inputs
yellow = "symbology/yellow.lyr"

# Set map document
mxd_city_year = arcpy.mapping.MapDocument(r"...\blank_map.mxd")
DF = arcpy.mapping.ListDataFrames(mxd_city_year)[0]

# Add layers
arcpy.MakeFeatureLayer_management("states/continental_US.shp", "us")
basis = arcpy.mapping.Layer("us")
arcpy.mapping.AddLayer(DF, basis, "AUTO_ARRANGE")
arcpy.ApplySymbologyFromLayer_management(basis, yellow)

# Save map
mxd_city_year.saveACopy("thresh_" + city + "_" + year + ".mxd")

The code runs through (no error). However, in the resulting map document the colour still seems to be random.

What is my mistake?


texts: Why different SPOT5 histograms from two versions of same SPOT5 raster?

texts:
I have encountered a, to me, strange behaviour in the histogram function of QGIS (2.2.0 Valmiera on Windows 7). Basically the same image produces totally different histograms.

I might have missed something fairly obvious, but I'm not able to find out what, and I have not been able to find any previous questions on the matter.

Background
I aim to do a "Black object subtraction" in order to correct for atmospheric noise in satellite images. Without getting too technical about that, I'm looking for the lowest pixel values with a frequency of at least 1000 pixels (a suggested method according to the literature).

Problem
When I calculate a histogram for my SPOT5 image, the maximum frequency of any pixel value is below 1000 (accept for 0, that indicates no data). I find that suspiciously low, as the image is 7900*7900 pixels.

Comparison of histograms

However, if I save the picture to with QGIS, just as it is, without any modification (see attached screen dump of save dialogue), I get a totally different histogram, that looks way more plausible to me (see image above).

Save dialogue

Question
How come that the basically the same raster produce totally different histograms? Is this a bug in QGIS? Faulty handling of images by me? Two different ways of rendering histograms?


texts: UpdateCursor to insert DDM values from a layer in WGS84?

texts:
I have a layer in WGS84 (Geographic, DD) that I need to label in DDM (Degrees Decimal Minutes).

I need to calculate within an arcpy UpDateCursor, a label field with DDM in it.

I saw that you can run the Convert Coordinate Notation, but I need to be able to calculate these labels on the fly within my cursor (this script is also moving points around).

After trapsing through the Esri forums I have not seen any snippets with how to achieve this.

This is likely a noob question, but this is very early steps for me any arcpy.


texts: How to change owner of a version in an enterprise geodatabase using ArcObjects?

texts:
I have ArcObjects/.NET code used to create few versions for workspace. Versions are created successfuly but their owner is who created them (the user connected to the workspace). I want to change the owners names, but I dont see an API to do that, maybe someone knows what interface I should use etc?

Everytime connecting as different user is not an option, because somehow unmanaged memory starts to leak when inserting/deleting/updating versioned data from many different connections.

Also I found table sde.VERSIONS it has column 'OWNER' but I am not sure if simply updating this column for specific versions is the correct solution, maybe it has some side effects.


texts: Using Map Tips and Ruled Expression in QGIS?

texts:
QGIS 2.18.5

I'm trying to add a condition to my expression in my html/css map tips where if a value is 'null', the result will show 'no data available'.

I've got a screenshot of my map tip so far for a feature with populated values (the date and the number of years). Not all the features have a date and integer associated with them.

What do I add to the code below?

 [%"Local Auth"%]

 Date LAP Adopted

 [% "Local Pl_1" %]

 Five Year Housing Supply

 [% "5 YHLS (#)" %]

I tried

 Date LAP Adopted

 [%CASE WHEN  "Local Pl_1" = ' ' THEN 'No Data Available' END%]

but no luck.


texts: Add/Remove 1000 separator (Comma) From Integer Field In ArcGIS Collector?

texts:
I have a team of inspectors using ArcGIS collector to collect and edit data.  When I set up the web map in AGOL I obviously wasn't thinking clearly and checked some fields to Use 1000 Separator (Years, ID numbers) and others I left unchecked (Lengths, $ amounts).  I inverted the selection I needed.  I had the team download the map and inspectors are using it in the field.  I corrected the 1000 separator mistake in AGOL, but when I sync in collector it still will not correct the problem.  This may sound a bit OCD, but the 1000 separators are actually causing a bit of disturbance in the way we collect data.

Does anybody know a way to correct this without downloading the map over again?


texts: How to find address from the Census BlockCode?

texts: I have Census BlockCode information (15 digits) and need to get corresponding street name and address in NYC. What is the best source of such information?


texts: Set geometry generator polygon to have flat end caps

texts:
I'm trying to buffer a line using Geometry Generator in QGIS (buffer($geometry,3)) However although I've set the Simple Fill join style to Miter I can't get the ends of the buffer to be flat.

I wanted to use a polygon for styling (rather than having the line width set to 6 map units) because it's much easier to have a semi-transparent fill with a consistent border and a line fill effect.

Is there any way to set the end cap for polygon styles or using expressions? I notice single_sided_buffer() has a parameter for join style but buffer() does not.


texts: Aggregating census tracts and reducing Coefficients of Variation?

texts: I am working on a project that maps out poverty rates (ACS data) by census tract. Many of the tracts have CV's that are just too high for our study and I am hoping for find a sound methodology in merging neighborhood tracts in hopes of lowering the CV's. At first, I arbitrarily chose which tracts to merge (neighboring tracts that also had high CV's) but I am wondering now if there is a more appropriate method for this. I know how to merge- I am wondering about how to choose which to merge.


texts: How to fill the void/empty space between two polygons

texts:
I am using QGIS 3.01

I am wondering how I can get get rid of the empty space that is splitting the polygon in the picture into two pieces. I know I can edit the polygons together, but is there a better tool for doing this? Like some kind of buffer tool that merge the two polygon?


texts: Generate and add L.control.layers from GeoJson properties to Leaflet

texts:
I have a problem with L.control.layer on this code, the control is added but he put another marker on top of the correct one.

Naturally, I'd like to generate the layers from the results, not adding new ones.


texts: GeoPandas read_file only reading first part of KML

texts:
I am trying to read the following file with GeoPandas:

import geopandas as gp
gp.io.file.fiona.drvsupport.supported_drivers['KML'] = 'rw'

df = gp.read_file('http://watchthemed.net/index.php/json/layer/9', driver='KML')

This returns one geometry:

gp.read_file('http://watchthemed.net/index.php/json/layer/9', driver='KML').shape

>> 1,3

Which is a multilinestring of length 1.

print(gp.read_file('http://watchthemed.net/index.php/json/layer/9',
       driver='KML').geometry.loc[0])

>> MULTILINESTRING Z ((11.526082 33.171136 0, 11.531247 33.196851 0,
11.537371 33.213983 0, 11.566443 33.298784 0, 11.585108 33.364755 0, 12
33.916667 0, 12.218734 34.175 0, 12.40997642454452 34.37477266876127 0,
23.583333 34.333333 0, 24.1667 34 0, 25.150691 31.669654 0))

However, if you visit the link, you can see that the file has many lines. How can I get the rest of them?


texts: Mapserver - label at beginning of line string (WMS)

texts:
I am creating a WMS service made up of linestrings pulled from a POSTGIS database

I want to label the start of each line with a sequence number. At present I can label the centre of each with:-

 PROCESSING 'LABEL_NO_CLIP=1'  (to prevent repetition across tiles)
 ....
 LABEL
      COLOR 132 31 31
      TYPE TRUETYPE
      FONT "arialbd"
      SIZE 12
      POSITION AUTO
 END

This places the label in the middle of the line string. Any idea how to place it at the start?

Note that GEOMTRANSFORM "start" does not work with labels


texts: Displaying Kriging Interpolation 3D in ArcScene?

texts:
I am trying to display my kriging.tiff interpolation raster in 3d within ArcScene. The reason to do this would be to visually see where the high values and low values are for my analysis. This would be similar to how DEM data is displayed in 3D for easier viewing and analysis.

I have tried using the 'raster to TIN' tool on the kriging raster but my image is still flat even though it says the base heights should be floating on a custom surface.


texts: Change automatic column st_length into st_perimeter in ArcSDE database (PostGIS/PostgreSQL)

texts:
When I import a feature class into my ArcSDE enabled PostGIS/PostgreSQL database, it automatically adds the st_area and st_length columns into ArcGIS. The problem is that st_length gives 0 when the features are polygons. According to PostGIS documentation, this should be normal (http://postgis.net/docs/ST_Length.html). I have to use st_perimeter instead. How can I change this automatic function (of sde?) so that it creates a st_perimeter column instead of the st_length?

I use PostgreSQL version 9.4, PostGIS 2.1.8 and ArcGIS 10.2.


texts: Split polygons with lines that do not cross features

texts:
I'm working on dividing a polygon layer with another line layer. My goal is to create subzones by cutting each polygon feature along these lines. At the moment I cannot make geometry adjustments, I have to use them as they are.
When using "Split with lines" these lines must (?) cross the polygon and the problem would be only in those that are inside or do not cross the polygon.
I noticed that in some cases (blue arrow, image below) some line that starts outside the polygon and crosses with another line inside the polygon makes the cut. Is there a processing order for this calculation? Am I missingre any topological rule?

transparent red thick line = cut lines
green polygons = result after "Split with Lines"


texts: Is this a Google Earth Engine bug when zooming in the Map.Layer?

texts:
I'm starting to construct some algorithms from scratch to calculate terrain parameters. My first algorithm is about the gradient slope (radians) and it apparently looks fine when zoomed out:

The result is also quite similar to the ee.Terrain.slope algorithm. The differences in intensity are only due to layer styles.

However, when zooming in into the image to check for small changes related to resolution, I found the following strange situation:

Which is very different from the ee.Terrain.slope, over the same region:

I don't know why this is happening and I think it can be a bug.
The full code is attached below, where its possible to compare the results from both the algorithms. The style for ee.Terrain.slope is custom 0-90, while my slope product is stretched by 2 sigmas.

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


texts: Applying partial opacity in WMS layer using OpenLayers?

texts:
I am working on a GIS API, and I can't share the code, so I will try to explain what I want.

I want to apply a style map on a single layer, where I have some elements with opacity 1 and others with opacity 0.5.

Please, look at the image below (where I show the desired feature using a "vector" layer):

I have tried using the opacity property from OpenLayers.Layer, but it sets the whole layer to 0.5 opacity.

I know WMS is an image layer, so it is maybe not possible to set multiple values for opacity. Even though, it thought it would be worth asking.

Is it possible to use a style map with multiple opacities in a WMS layer?

With my current SLD Body, I get the object with only opacity 1.

<sld:StyledLayerDescriptor xmlns:sld="http://www.opengis.net/sld" version="1.0.0" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:schemaLocation="http://www.opengis.net/sld http://schemas.opengis.net/sld/1.0.0/StyledLayerDescriptor.xsd" xmlns:ogc="http://www.opengis.net/ogc" xmlns:gml="http://www.opengis.net/gml"><sld:NamedLayer><sld:Name>topp:tasmania_water_bodies</sld:Name><sld:UserStyle><sld:FeatureTypeStyle><sld:Rule><ogc:Filter><ogc:Or><ogc:PropertyIsGreaterThan><ogc:PropertyName>AREA</ogc:PropertyName><ogc:Literal>1067500000</ogc:Literal></ogc:PropertyIsGreaterThan></ogc:Or></ogc:Filter><sld:PolygonSymbolizer><sld:Fill><sld:CssParameter name="fill">#0000ff</sld:CssParameter><sld:CssParameter name="fill-opacity">1</sld:CssParameter></sld:Fill><sld:Stroke><sld:CssParameter name="stroke">#000000</sld:CssParameter></sld:Stroke></sld:PolygonSymbolizer></sld:Rule><sld:Rule><ogc:Filter><ogc:Or><ogc:PropertyIsLessThanOrEqualTo><ogc:PropertyName>AREA</ogc:PropertyName><ogc:Literal>1067500000</ogc:Literal></ogc:PropertyIsLessThanOrEqualTo></ogc:Or></ogc:Filter><sld:PolygonSymbolizer><sld:Fill><sld:CssParameter name="fill">#0000ff</sld:CssParameter><sld:CssParameter name="fill-opacity">0.5</sld:CssParameter></sld:Fill><sld:Stroke><sld:CssParameter name="stroke">#000000</sld:CssParameter></sld:Stroke></sld:PolygonSymbolizer></sld:Rule></sld:FeatureTypeStyle></sld:UserStyle><sld:UserStyle><sld:FeatureTypeStyle><sld:Rule><ogc:Filter><ogc:Or><ogc:PropertyIsGreaterThan><ogc:PropertyName>AREA</ogc:PropertyName><ogc:Literal>1067500000</ogc:Literal></ogc:PropertyIsGreaterThan></ogc:Or></ogc:Filter><sld:PolygonSymbolizer><sld:Fill><sld:CssParameter name="fill">#0000ff</sld:CssParameter><sld:CssParameter name="fill-opacity">1</sld:CssParameter></sld:Fill><sld:Stroke><sld:CssParameter name="stroke">#000000</sld:CssParameter></sld:Stroke></sld:PolygonSymbolizer></sld:Rule><sld:Rule><ogc:Filter><ogc:Or><ogc:PropertyIsLessThanOrEqualTo><ogc:PropertyName>AREA</ogc:PropertyName><ogc:Literal>1067500000</ogc:Literal></ogc:PropertyIsLessThanOrEqualTo></ogc:Or></ogc:Filter><sld:PolygonSymbolizer><sld:Fill><sld:CssParameter name="fill">#0000ff</sld:CssParameter><sld:CssParameter name="fill-opacity">0.5</sld:CssParameter></sld:Fill><sld:Stroke><sld:CssParameter name="stroke">#000000</sld:CssParameter></sld:Stroke></sld:PolygonSymbolizer></sld:Rule></sld:FeatureTypeStyle></sld:UserStyle></sld:NamedLayer></sld:StyledLayerDescriptor>


texts: How to export raster of Mean Annual Land Surface Temp from Landsat Imagery?

texts:
I am using google earth engine to estimate Mean Annual land surface temp of a local scale from 1984 to 2019 from landsat Imagery.
I am having trouble determining a script for my goal, there are many example but the problem is that all of them are calculating the LST for one year.
I have tried this script but I think there is some thing wrong with it and I cannot export the raster of calculated Mean annual LST for my study.
I will be helpful, if anybody could help me to rectify the script. The GEE Code for Export raster of mean annual lsthere


texts: Closest raster VALUE to a given x,y point

texts: Is there a way in ArcGIS to find the closest VALUE in a raster to a given x,y point? There might only be NoData at a certain point, in which case I want to find the closest pixel value.  I can't seem to find a function or combination of functions to do that.


texts: Error message: Commit errors when creating a new shapefile

texts:
I want to create a shapefile (ESRI shapefile) layer in QGIS version for 3.6.2 for polygons however I keep getting error messages at first it was:

  Error: could not commit changes to layer error 1 feature(s) not added

I then went onto the table changed id column to Text(string) format as I looked up what the error message meant in the hope that it would work, but now I get another error message that says

  Error: could not commit changes to layer error 1 attribute(s) not added

Does anyone know how to avoid these messages so I can go in and out of edit mode in the shapefile?


texts: Extracting layers from gdb and merging them?

texts:
I have a gdb folder which contains ~3000 multipolygon layers (with very few features on each layer). I would like to take around 100 of those layers and merge them together while adding a field to the features to indicate what layer they originate from.

I think that ogrmerge.py can help me with the merge and with feature labeling. But I haven't found a way to either select specific layers to merge (ogrmerge seems to take an entire dataset as its input) or to extract the relevant layers to create a sub-dataset to pass to ogrmerge.

Any idea how this could be achieved (preferably with command line open source tools)?


texts: Retaining georeferencing information from CAD when opening in another mxd?

texts:
I georeferenced about 150 CADs, updating the georeferencing as I work. I created a new mxd using the same base layers and added one of the CADs that was georeferenced but it did not retain the location, rotation, or scaling information.

I am using ArcGIS Desktop 10.6.

What am I doing wrong?

I tried saving the georeferencing in the project folder, and then tried saving it in the root CAD folder. Both options gave me the same result. It opens up, but not where I georeferenced it.