texts: Dynamic MAP parameter for MapCache + MapServer

texts:
The MapServer and MapCache are installed on the Ubuntu machine, all works fine. However, the specific of developed system is that the displayed through WMS data often have different map files created. The question is - can the MAP parameter be dynamic, i.e. fetched from the original request to MapCache?

So, basically I want to make request like

  http://31.227.56.57/mapcache/?SERVICE=WMS&REQUEST=GetMap&VERSION=1.1.1&LAYERS=localmapserver&HEIGHT=1024&WIDTH=1024&SRS=EPSG%3A3857&BBOX=11271098.442818949,6261721.35712164,11897270.578531113,6887893.492833804&map=/s/some/file.map

The example of working config

  <source name="vmap1" type="wms">
      <getmap>
         <params>
            <FORMAT>image/png</FORMAT>
            <LAYERS>result_0</LAYERS>
            <TRANSPARENT>true</TRANSPARENT>
            <MAP>/s/temp/mapfiles/RSCnuQp7iH.map</MAP>
            <!--<MAP>{MAP}</MAP>-->
         </params>
      </getmap>
      <http>
         <url>http://31.227.56.55/cgi-bin/mapserv</url>
      </http>
   </source>

   <tileset name="localmapserver">
      <source>vmap1</source>
      <cache>disk</cache>
      <grid>WGS84</grid>
      <grid>GoogleMapsCompatible</grid>
      <format>PNG</format>
      <metatile>5 5</metatile>
      <metabuffer>10</metabuffer>
      <expires>3600</expires>
   </tileset>


texts: Openlayers 3.10.1 bind Overlay to moving feature

texts:
In my map I'm tracking vehicle locations. If a vehicles is clicked a Overlay with additional informations is shown. How can I "bind" the Overlay to the vehicle location. I tried

feature.getGeometry().on(
            'propertychange',
            function(propEvt) {
                console.log('I'm here');
                var extent = feature.getGeometry().getExtent();
                var centerOfFeature = getCenterOfExtent(extent);
                popup.setPosition(centerOfFeature);
            },
            feature.getGeometry());

But log message is not shown but the feature is still moving.

-- Rüdiger


texts: QGIS custom CRS for map of manufacturing shop floor

texts:
I'm working on a web app that shows a map of a manufacturing plant's shop floor.  I've built it in Inkscape and QGIS, but when I try to load it into Tilemill, I get a CRS error.  Specifically, it says:

Detected out of bounds geographic extent (14.944551, 93.956594, 180, 85.051) for layer 'aisles'.  Please ensure that the SRS for this layer is correct.  It's native extent is '14.944551, 93.956594, 195.298512, 254.314634'

How do I set a custom CRS in QGIS, OR, how do I specify the map's location within the current CRS so that it will work in tilemill?

Background: The webapp will be for a map of a factory shop floor.  Different zoom levels will have different information/ have different actions.  I'm planning on styling the map in Tilemill, then using Leaflet to show the map.  This will run on the company's internal intranet.

Research: I've read this question on stackexchange.  It's helpful, but I'm not sure it fits my needs, and I don't know what the proj4 parameters are.  Also, I've looked at this documentation, but haven't been able to understand it and apply it to my situation.

Construction details: I created the map first as SVG in Inkscape.  Then I put each layer (one at a time) into QGIS.  After that, I'm planning on loading each layer (one at a time) into Tilemill.  Obviously, the geographic scope of this is much smaller than most map projects.

This map does NOT need to correspond to real-world lat/long, but I do need to be able to show where machines are (and their status) on the map.

Finally, if I'm thinking about this wrong, and need to go another direction or use different tools, don't hesitate to let me know that, too.

Thanks.

EDIT: SOLVED.  Ended up going with the answer below.  I used the qgsAffine plugin to move/scale everything to fit between the 0 and 180 coordinates, on both the x and y axes.

One hiccup was that QGIS/qgsAffine can only edit vectors that are in the shapefile format (in my experience.  There might be some other format out there it can do).  QGIS easily converted my dxf files to shp, though, so it was an easy fix.


texts: Processing modeler: how to dissolve a vector layer and sum the attributes

texts:
I have a complex workflow prepared in the processing modeler, but I'm missing one step to link all the model.

I need to dissolve a polyline vector in a single line vector that has as attribute the SUM of a specified attribute of the polyline vector.
I can do it via plugin (dissolve with stats) but just outside the modeler, I tried to find a way to link the "dissolve" tool and the "basic statistics for numeric field" tool but I still can't figure out how to put the html output inside the vector layer.

I saw some comments about an sql query that looks easy, but I don't know if I can use it outside spatialite or sql, inside the modeller.

I tried the command SELECT SUM(column_name) FROM table_name;

but as out put i got GDAL command output:
FAILURE:
Unable to open datasource

Is it because the layer has to be inside an sql database to run this kind of query?


texts: How to start L.Control.Geocoder.nominatim() with user's current location

texts:
I'm trying to develop a routing system and I want to know how can start with the user's current location from A point?

var control = L.Routing.control({
    //   waypoints: [null],
    waypoints: [
        L.latLng( user's current location),
        L.latLng(" string address")
    ],
    routeWhileDragging: true,
    show: true,
    language: 'it',
    geocoder: L.Control.Geocoder.nominatim(),
    autoRoute: true
}).addTo(map);


texts: QGIS line style with alternating dashes - and x

texts:
I found the latest QGIS has a great tool to create complex styles, however looks like one small but very useful thing is missing, or I can't find it.
I'm trying to create a line style something like "- x - x -". I can do it in various ways, but none of them gives right results. The problem is that I need to set starting point (shift) either for "-" or "x", so to draw it not from 0. If I set X coordinate shift for the marker it works ok, except line bends, where marker can go outside the line (vertex).

QGIS style:

ArcMap style:

Is there is something to do with it?


texts:
I have little experience with pyproj, but I want to achieve the following task:

I have a trajectory of a robot, given in lat/long coordinates which I want to convert to x/y Coordinates. In aerospace sciences, we used the NED coordinate system for this, therefore I want to transfer a bunch of lat/long coordinates to local x/y coordinates whereby the first x/y coordinate shall be at 0/0,(which is typically just the start of the trajectory).

So I have been trying to use the PyProj library, which sucessfully converts lat/long to x/y, however the coordinate system in which the result is given is that of england, however ideally I want the 0/0-point to be my first lat/long point. So I though perhaps it would be possible to define a custom coordinate system originating at a specific coordinate system, but couldn't find anything on google.

Therefore: does anyone know how I can define a NED coordinate System at a specific lat/long combination?

Here is my code so far:

import pyproj as proj

crs_wgs = proj.Proj(init='epsg:4326')  # assuming you're using WGS84 geographic
    crs_bng = proj.Proj(init='epsg:27700')  # use a locally appropriate projected CRS
    # then cast your geographic coordinate pair to the projected system
    x, y = proj.transform(crs_wgs, crs_bng, s['gps_lng'], s['gps_lat'])


texts:
Is there any method or tool/model for ArcGIS to draw automatic line based upon raster data values.
E.g i have a DEM and i want to automatically generate a line feature passing through all those pixels whose values is less than 1.


texts: Iterating through list of fields using field calculator?

texts:
import arcpy, sys, itertools

arcpy.env.workspace = r'path to geodatabase'

fields1 = ['gearqty', 'cnemarea', 'fzone', 'qdsq', 'tensmsq', 'crew', 'nanglers', 'port', 'totalquant', 'trip_days', 'day_fracti', 'gear_days', 'qtykept', 'qtydisc', 'mean_hp', 'sd_hp', 'mean_len', 'sd_len', 'mean_gtons', 'sd_gtons', 'size_fleet', 'anon_id', 'fdays']
fields2 = ['gearqty1', 'cnemare1', 'fzone1', 'qdsq1', 'tensmsq1', 'crew1', 'nangler1', 'port1', 'tquant1', 'tripdays1', 'dayfract1', 'geardays1', 'qtykept1', 'qtydisc1', 'mean_hp1', 'sd_hp1', 'mean_len1', 'sd_len1', 'mean_gton1', 'sd_gtons1', 'sizefleet1', 'anonid1', 'fdays1']

# add new fields of type LONG
for fc in arcpy.ListFeatureClasses():
    for b in fields2:
        arcpy.AddField_management(fc, b, "LONG")

# calculate LONG fields from STR fields
for fc in arcpy.ListFeatureClasses():
    for a, b in itertools.izip(fields1, fields2):
        arcpy.CalculateField_management(fc, b, !a!)

# delete all unnecessary string attributes
for fc in arcpy.ListFeatureClasses():
    for a in fields1:
        arcpy.DeleteField_management(fc, a)

The list of fields, fields1, is all string data, but it should be numeric.  The first code block successfully created new LONG attributes using the names provided by fields2.

The second part is where I'm having trouble.  I just need the new fields2 attributes to be calculated from the old fields1 attributes.

The third section should delete all of the attributes from fields1, as they are useless to me as strings.  I haven't been able to test that part yet.

I've been testing this on one shapefile, but there will be 20+ shapefiles in the geodatabase.


texts: How to match a raster NetCDF data with a vector layer in QGIS?

texts:
When I load a raster layer with humidity data from NOAA database (for instance rhum.2m.2014.nc at ftp://ftp.cdc.noaa.gov/Datasets/NARR/Dailies/monolevel/) I am not able to match it with a simple map of the US counties (https://www.census.gov/geo/maps-data/data/cbf/cbf_counties.html).

I tried to give them both the same CSR but it does not work, it seems that the raster layer is a lot smaller than the US map.

Does anyone have any idea? I have been working on this issue for several days now...

NB: The raster file is a netcdf

NB2: With arcGIS the two maps match automatically but unfortunately I have to work with QGIS.


texts: Applying Google maps styles to OpenStreetMap?

texts:
I have these styles of  google maps:

[
{
    "featureType": "landscape",
    "elementType": "all",
    "stylers": [
        {
            "saturation": "-16"
        },
        {
            "lightness": "33"
        }
    ]
},
{
    "featureType": "landscape.natural",
    "elementType": "all",
    "stylers": [
        {
            "saturation": "-21"
        },
        {
            "lightness": "15"
        }
    ]
},
{
    "featureType": "poi.business",
    "elementType": "labels",
    "stylers": [
        {
            "visibility": "off"
        }
    ]
},
{
    "featureType": "road.highway",
    "elementType": "all",
    "stylers": [
        {
            "hue": "#53FF00"
        },
        {
            "saturation": -73
        },
        {
            "lightness": 40
        },
        {
            "gamma": 1
        }
    ]
},
{
    "featureType": "road.arterial",
    "elementType": "all",
    "stylers": [
        {
            "hue": "#FBFF00"
        },
        {
            "gamma": 1
        }
    ]
},
{
    "featureType": "road.local",
    "elementType": "all",
    "stylers": [
        {
            "hue": "#00FFFD"
        },
        {
            "lightness": 30
        },
        {
            "gamma": 1
        }
    ]
},
{
    "featureType": "transit",
    "elementType": "all",
    "stylers": [
        {
            "hue": "#00cdff"
        }
    ]
},
{
    "featureType": "transit.line",
    "elementType": "all",
    "stylers": [
        {
            "hue": "#ffbf00"
        }
    ]
},
{
    "featureType": "transit.line",
    "elementType": "geometry",
    "stylers": [
        {
            "visibility": "on"
        },
        {
            "color": "#9b9b9b"
        }
    ]
},
{
    "featureType": "transit.station.rail",
    "elementType": "all",
    "stylers": [
        {
            "saturation": "-37"
        },
        {
            "lightness": "-6"
        }
    ]
},
{
    "featureType": "water",
    "elementType": "all",
    "stylers": [
        {
            "saturation": "-69"
        },
        {
            "lightness": "-13"
        },
        {
            "gamma": 1
        },
        {
            "color": "#9fbed4"
        }
    ]
}
]

I would like to apply it on a OSM map in OpenLayers.


texts: Line density for mulitple commuters in QGIS

texts:
I have a layer with multiple lines, one city area commuting to another.
This means that on some roads there are multiple lines lying on top of eachother, which at some point will split.
I also have some roads going next to eachother.

I would like to create a sort of line density map based on these lines in which the value of the line segment is based on the amount of lines on top of eachother or within 10 meters.

I have tried a couple of solutions but none of them seem to work.

I have tried to join the lines to a more simple line layer, but it seemed to join on the wrong lines.

I have tried to split the lines into points and then connect to points to a more simple line layer, yet again some of them join on the wrong lines.


texts: Does the Trimble GeoXT 2008 (3000 series) work with ArcPad 10?

texts:
I've been looking around and I've found some info on the 6000 series...however, nothing much on the 3000 series.

My main concern is the firmware on the GeoXT 3000 is not compatible with the ArcPad 10.x and that it won't be able to run. Before I try going to knock down Trimble's customer service door, does anyone here have experience with using ArcPad 10.x on a GeoXT 3000?


texts: Open Source alternatives to ArcGIS Locator Report?

texts:
I'm a pretty new GIS programmer and I'm trying to run a Locator Report in ArcPy with no success.

I'm wondering if there are any other software packages that can accomplish the same thing but with open source code?


texts: Create table into PostgreSql from Geoserver SQL View query

texts:
I know how to connect Geoserver to PostgreSql and create SQL View queries with parameters from Geoserver, but ;it is possible to create a new table into PostgreSql from a Geoserver Sql View with parameters like this?:

Geoserver:

SELECT source FROM pgr_nearest_node('net',%lat%,%long%)


texts:
Is there a way to build topology and work with topology relations in PyQGIS?

I'm asking because I need to run some custom generalization algorithms on my data, and I need to preserve topology. Presumably this would be easier if there is already some functionality for handling topology in QGIS.

The data is basically shapefiles with multipart-polygons.


texts:
I would like to modify a raster (pixel values) in QGIS.

The plugins "value tools" can easily provide me with the value of the pixel of my raster.
However I don't know how to modify it.

I have found the Q&A titled Manually correcting raster DEM (OpenSource)? which says that it is possible to do it in GRASS.
I'v spent several hours on trying to understand how grass works without success.

I can't find the d.rast.edit tool in the sextante plugin.

Is there another way to acess it?


texts: Reducing tree cover area from Hansen in GEE yields me areas larger than the feature itself

texts:
I am trying to calculate the tree cover in Maine (by county) using the reduceRegions() tool. I based my code on this thread. I am interested in areas in square km. The output of my code yields, though, areas of counties' tree cover that are larger than their areas (and/or the whole Maine itself!). The calculated tree cover in the 16 counties of Maine sums up to about 6.5 million sq. km. The area of Maine, instead, is of about 95 thousand sq. km.

I set the scale to 30 meters so to fit Hansen's data pixel resolution. Here is my code:

//Load and filter the Hansen data
var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2017_v1_5')
              .select(['treecover2000','loss','gain','lossyear']);

// turn your scale into a var in case you want to change it
var scale = 30;

//load maine:
var maine = ee.FeatureCollection('ft:1S4EB6319wWW2sWQDPhDvmSBIVrD3iEmCLYB7nMM')
              .filter(ee.Filter.eq('StateName', 'Maine')).select(['CntyFips']);

//find the area of tree cover in sq km:
var treeCover = gfc2014.select(['treecover2000']);
var areaCover = treeCover.multiply(ee.Image.pixelArea())
                .divide(1000000).select([0],["areacover"]);

// Sum tree cover area per feature (county):
var mainetree = areaCover.reduceRegions({
  collection: maine,
  reducer: ee.Reducer.sum(),
  scale: scale,
});

print(mainetree);    // just for visual inspection, the sum is done in the .csv

Map.addLayer(maine, {}, 'counties');    // visualise area of interest
Map.centerObject(maine,6);

// Export the FeatureCollection.
Export.table.toDrive({
  collection: mainetree,
  description: 'maine',
  fileFormat: 'CSV'
});


texts: How to display data within a drawn polygon

texts:
I am using mapbox as my base layer and would like to implement a functionality as shown in this gif below, where if i draw a polygon , it reads from my geotiff file and display NDVI data. I would like to color code pixels appropriately ranging from red to green based on NDVI value

Can someone point me to right direction please.


texts: Make rotate button always visible in OpenLayers

texts:
The documentation mentions that the rotate button is always hidden, when the rotation of the view is 0:
https://openlayers.org/en/latest/apidoc/module-ol_control_Rotate-Rotate.html
A .ol-hidden css selector is added to the button when the rotation is 0.
After searching for hours, I wonder if OpenLayers offers any way to make this button always visible?
I know how to override CSS and I also know how I could use vanilla JS. But before doing so, I would like to make sure that I haven't overseen anything in the docu.
UPDATE 1
According to the first answer, the configuration is done on the map object on its control param.
// Vue
import { defaults as defaultControls } from 'ol/control'
export default {
  mounted() {
    this.olMap = new Map({
      controls: defaultControls({
        rotateOptions: { autoHide: false }
      })
    })
  }
}

That works. But how is "autoHide: false" as a configuration option for the control implemented, as the first part of the answer suggests?
UPDATE 2
The first half of the first answer refers to a situation, where ONLY a rotate button is supposed to be included. In that case the code would look like this:
// No need to import "default" controls...
import Rotate from 'ol/control/Rotate'
export default {
  mounted() {
    this.olMap = new Map({
      // "controls" expects an array:
      controls: [new Rotate({ autoHide: false })]
    })
  }
}

The second solution as shown, hover, does not display any other buttons (e.g. zoom) but the rotate button. They would have to be implemented separately in the same way as the rotate button.


texts: How to use filter on the Search widget in ArcGIS JS API 4.9?

texts:
var searchWidget = new Search({
  view: view,
  includeDefaultSources:false
  sources: [{
    locator: new Locator({ url: "https://geocode.arcgis.com/arcgis/rest/services/World/GeocodeServer"}),
    countryCode:"GB",
    categories:["Address", "Postal"],
    maxResults: 10
    minSuggestCharacters: 4,
    placeholder:"Find address",
    filter:{
        ...
    }
  }]
});

I'm trying to limit the results to a certain geographical area. I've been told I can use the filter property, but I can't find how or any examples.


texts: I have a .DWG (CAD dataset file) that when I bring into ArcMap 10.1, it will not display anything. I've tried converting to a shapefile and importing to geodatabase and I am not able to progress from there. Am I doing something wrong in particular?


texts:
Testing my VB.NET 2010 app in ArcMap an exception occurs when ICommandBars is instantiated.

' zoom to selected
Dim mxApp As ESRI.ArcGIS.Framework.IApplication
Dim numUID As UID = New ESRI.ArcGIS.esriSystem.UID
numUID.Value = "{AB073B49-DE5E-11D1-AA80-00C04FA37860}"
Dim commandBars As ICommandBars = mxApp.Document.CommandBars   '>>> error here!!!!!!!!
Dim cmdItem As ICommandItem = CommandBars.Find(numUID, False, False)
If Not (cmdItem Is Nothing) Then
    cmdItem.Execute()
End If

Message: Exception first-chance type 'System.NullReferenceException' in LUCAS.dll.

P.S. ArcGIS 10 SP1 on Windows 7, VS2010.


texts: Problems installing QGIS on Windows Vista

texts:
I have tried several times on several windows vista machines to install QGIS (both versions 2.10.1 and 2.8) but with no luck.
It seems to install fine without any errors but when the installation is finished QGIS is nowhere to be found. There is no short cut created on the desktop and it does not appear when I search in the start menu.

I have very basic computer knowledge and have simply followed the onscreen instructions and not altered anything else.

It works perfectly on my windows 7 and 8 machines.


texts: I use topology checker and add a rule that polylines must not have dangles. It find errors for all 36 polylines but I think some topologies are correct at the right side. Can you help me how can I fix this. Because I tried same rule in ArcGIS and results seems better. Pictures are below.


texts: Is the output of the Dijkstra one-to-many routing algorithm always sorted, cheapest to most expensive?

texts:
I have implemented the pgrouting Dijkstra one-to-many algorithm, and have noticed that the result is a table of every route together, with and edge value of -1 indicating the end of an individual route.

I am interested in selecting the cheapest route of all the routing results.  By manually parsing and comparing the results via the agg_cost values, I have found that every single time, the first group of results (i.e. the rows up until the first -1 value is the) is the shortest path.

I'm wondering now if this is by design, if it is a consequence of Dijkstra's algorithm (or at least this implementation of it), or if it is just a coincidence given my relatively simply sample dataset.

Is it necessary to manually compare the routing results to find the shortest path, or will the first route in the group always be the shortest one?


texts:
How can I put a line break in the popup window in the following statement? I am trying \n and  but its not working.

for lt, ln, el, plc in zip(lat, lon, elev, pl):
    fga.add_child(folium.Marker(location=[lt, ln], popup= "Place is " + plc + \n +"Magnitude of " + str(el),
    icon=folium.Icon(color=color_producer(float(el)))))


texts: Hotspot Analysis Output Color Disappears When I Reload the File in Another Project

texts:
I am using ArcGIS desktop 10.x. I used hot spot analysis on a couple of points feature classes, and saved the output to a gdb in shapefile format.

I then reload the shapefiles in a newly created project, but the coloration of hot/cold spot no longer exists.

Should I export the output into another format in order to reload, or other kinds of gdb? Or is it that I can not reload the analysis output in another project whatsoever?


texts: gstat's krige.cv() error: "Error: dimensions do not match:"

texts:
I am using krige.cv() for cross validation of a data set with 1394 locations. In my code, empirical variogram estimation, model fitting, and kriging (by krige()) everything works fine. But when I use krige.cv() I get the following error.

  Error: dimensions do not match: locations 2786 and data 1394

One can notice that 1394 * 2 =  2786. What could I be missing?
Please note that there are no NA or NaN or missing values in the data, variogram or kringing results. Everything works fine, and it's just krige.cv() that does not work.

Update: Here is the source code.

library(gstat)
library(sp)
library(geoR)

data <- read.table("data.txt", header=T)
names(data)
attach(data)

#variable of interest
var = z
coordinates(data) = ~x + y
data.var <- variogram(var ~  1, data=data, cutoff=1, width=.05)
model <- fit.variogram(data.var, vgm(0.20, "Sph", nugget=0.06, range=1, fit.method=2))
poly <- chull(coordinates(data))
poly.in <- polygrid(seq(69.0, 74.1, 0.1), seq(31.5, 37.5, 0.1), coordinates(data)[poly,])
coordinates(poly.in) <- ~ x + y
krige.out <- krige(var ~ 1, data, poly.in, model=model)
krige.out.cv= krige.cv(var ~ 1, data, model=model)


texts: Defining custom projection for moving all objects on X and Y coordinates

texts:
I have a map project wih layers in WGS1984 UTM projection. I also have a bunch of layers without .prj files. When I add it on top of my layers, it's placed with a huge offset on X and Y coordinates but other than that everything is fine, features are not distorted etc. So basically I just need to "move" all the features along X and Y axes and everything will be alright.

Can I define such projection using ArcGIS "Define projection" tool or something else? I need a projection that just adds a linear offset to X and Y coordinates. I know I can fix one layer by using "Spatial Adjustment" tool and it works fine, but my goal is to have a .prj-file that I can use for any such layer.

I'm kinda new to ArcGIS and will give clarifications if needed.
Any help would be appreciated, thanks in advance.


texts: Changing Excel point coordinates in different coordinate system to points in ArcGIS Desktop

texts:
I have coordinates in a format e.g. 22°42'39.065”E, 49°4'59.864”N. I tried to add the points using "XY Table to Point" but when I add my Excel file as input, it doesn't recognize the fields as X and Y values. However, when I convert points in that format to DD format it works perfectly.
Is there a way to add those points from an Excel directly, or I need to have my points always in DD format in my Excel database?


texts:
Been using Esri JavaScript for a while but just started trying out the silverlight API but things are not going smoothly!

Tasks that I used to accomplish easily with Javascript is now proving unreasonably long and even unsolvable when using silverlight.

I am looking to create a little tool to add to the silverlight viewer.  The tool interface has a listBox a Button and a ComboBox.

The logic is that The listBox will have hard-coded serviceURL's. User selects a serviceURL, clicks on the button and the tool display the layers/map information in the comboBox.

JavaScript did this easily with GetResourceInfo etc but I just don't know where or how to carry this out with SilverLight.

Any help will be very much appreciated.


texts: Clipping/exporting multiple selected vector layers in map to extent of another raster/vector layer in ArcGIS Desktop?

texts:
Is there a way to clip/export selected vector layers (more than one) in a map to the extent of another raster/vector layer at a time?? in Arcgis or any other software.

I am using Arcgis10. I am having many vector layers as shown in figure. I want to clip for example "Australia" region of all vector layers. I can do it each and every layer separately. But I want to clip all layers/few layers in the map to Australia region at a time as a bundle. So-that I can avoid re-work of clipping each layer and save time.


texts: How to make 3d stack bar and overlap in QGIS

texts: This is a stack bar plot made by QGIS 3.18 (1st image). All I want is to flip the image and overlap this plot as 3D in QGIS (Like 2nd image). How can I do that?


texts: Could bridges prevent ArcMap from creating correct watersheds?

texts:
I have a very detailed DEM of a region in Northern Germany. For hours I tried creating correct watersheds but my results all look illogical. Is it possible that bridges disturb the ArcMap Tools? I have a shapefile including the streams as vector lines but I don't know how to include them in one of the tools.

Here is one example cutout of my DEM, I marked the mentioned bridges (I'm not sure if "bridge" is the correct word to use: pipes which are covered by earth and allow ditches to flow underground).

My aim is to analyse the path which rains takes after hitting the ground,  similar to the WMS service ofered by our local government (http://www.umweltkarten-niedersachsen.de/arcgis/services/Hydro_wms/MapServer/WMSServer?). An example from this WMS:

You can see that water from the green circuited area flows this way: "Hekelner Kanal" -> "Ollen West" -> "Motzener Kanal" -> "Weser" (which is subdivided into several sections).


texts: How to use shp2pgsql to import a very large shapefile into PostGIS?

texts:
I need to import very large datasets into the Postgis server so I need to use the -D option in the prompt. Problem is that I'm getting an error no matter what I do.

The command lines I'm using are:

shp2pgsql -D "C:\Documents\Work\Data\AddressBase\AB4.shp" AB | psql -h serverip -p 5432 -U name

The DBF has over 14GB, 1,000,000 rows and 17 columns

I'm getting the following

I have the latest pgadmin and postgis. Is my command line wrong in some way?

Also am I missing a password input?


texts: WMTS map with different origins for different zoom levels

texts:
I have WMTS map and I'm using proj4leaflet for nonstadard crs. This map has different origins for different zoom levels. First two levels have the same origin, then the next two levels have also the same origin but different from the first two and so on. As consequence every second zoom map not only zooms but also moves.

Here is GetCapabilites request for my wmts service:
https://prostor4.gov.si/ows2-gwc-pub/service/wmts?service=wmts&request=GetCapabilities

My layer is SI.GURS.DK:DTK5_D48 and tilematrixset EPSG:3912_DTK5.

The whole code:

var myUrl = 'http://prostor4.gov.si:80/ows2-gwc-pub/service/wmts?&layer=SI.GURS.DK:DTK5_D48&style=default&tilematrixset=EPSG:3912_DTK5&Service=WMTS&Request=GetTile&Version=1.0.0&Format=image/png&TileMatrix=EPSG:3912_DTK5:{z}&TileCol={x}&TileRow={y}';

var myLayer = new L.TileLayer(myUrl, {
    maxZoom: 11,
    minZoom: 0
});

var crsResol = [420, 210, 140, 70, 42, 28, 14, 7, 2.8, 1.4, 0.7];

var crs3912 = new L.Proj.CRS('EPSG:3912',
             '+proj=tmerc +lat_0=0 +lon_0=15 +k=0.9999 +x_0=500000 +y_0=-5000000 +ellps=bessel +towgs84=682,-203,480,0,0,0,0 +units=m +no_defs',
             {
               origin: [370210.0, 235895.0],
               resolutions: crsResol
             });

var myMap = new L.Map('map', {
  crs: crs3912
});

myMap.addLayer(myLayer);
myMap.setView([46.047711, 14.507136], 0);

I tried to dynamically change crs attributes (origin and transform) but it didn't help. Is there any solution for this?


texts: Folium draw polygons with distinct colours

texts:
I have trouble with my assignment. I have to draw over 100 polygons from the xml file and every of them has to have other color. I don't know how to accomplish it.
First I tried to draw polygons in loop, but eventually all polygons have the same color.
colors = []
for x in range(len(gminy)):
    color = np.random.randint(16, 256, size=3)
    color = [str(hex(i))[2:] for i in color]
    color = '#'+''.join(color).upper()
    colors.append(color)
cnt = 0
print(colors)
for x in gminy.iterrows():
    styl = {'fillColor':colors[cnt], 'color':'#000000', 'weight':0.1}
    folium.GeoJson(data=x[1]['geometry'], style_function = lambda y: styl).add_to(m)
    cnt += 1

Then I tried to create folium.Choropleth object with bins parameter set to 116 (that's the count of polygons) like this:
folium.Choropleth(
    geo_data=gminy,
    name='idTerytTerc',
    data=gminy_nr,
    columns=['idTerytTerc', 'Liczba'],
    key_on='feature.properties.idTerytTerc',
    fill_color='OrRd',
    fill_opacity=0.7,
    line_opacity=0.1,
    legend_name='Liczba mieszkańców',
    highlight=True,
    bins=116
).add_to(m)
m

But when I did this I get an error message:
---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
<ipython-input-53-42d7422910e7> in <module>
----> 1 folium.Choropleth(
      2     geo_data=gminy,
      3     name='idTerytTerc',
      4     data=gminy_nr,
      5     columns=['idTerytTerc', 'Liczba'],

C:\Programy\Anaconda\envs\geo\lib\site-packages\folium\features.py in __init__(self, geo_data, data, columns, key_on, bins, fill_color, nan_fill_color, fill_opacity, nan_fill_opacity, line_color, line_weight, line_opacity, name, legend_name, overlay, control, show, topojson, smooth_factor, highlight, **kwargs)
   1223             # We add the colorscale
   1224             nb_bins = len(bin_edges) - 1
-> 1225             color_range = color_brewer(fill_color, n=nb_bins)
   1226             self.color_scale = StepColormap(
   1227                 color_range,

C:\Programy\Anaconda\envs\geo\lib\site-packages\branca\utilities.py in color_brewer(color_code, n)
    150     if not explicit_scheme:
    151         # Check to make sure that it is not a qualitative scheme.
--> 152         if scheme_info[base_code] == 'Qualitative':
    153             matching_quals = []
    154             for key in schemes:

KeyError: 'OrRd'

So, how can I fix this?


texts: Check if a point coordinate is between linestring coordinates

texts:
I know that, linestring is a vector and if i want to check that,  if point coordinates are between linestring beginning and ending coordinates then i need to use vector comparison
(x2 > x && x2 < x1 && y2 > y && y2 < y1) that comparison show me if x2(Point) coordinate is between a linestring(vector) codrinates.

But I don't know how do compare linestring coordinates and point coordinates,

I wonder if somebody can give me some tips on how do compare if a point corodinate is between linestring coordinates?


texts: Arcpy MakeFeatureLayer_management issue using SQL condition

texts:
I have a lines shapefile. I made a small code to select some of its elements and turn them into a standalone layer. The code is like this:

arcpy.env.workspace = r"C:\MyDirectory\Mygdb.gdb"
network = r"Linesshapefile"
arcpy.MakeFeatureLayer_management(network,'Winners_lyr', "OBJECTID_12" in Winners)
arcpy.SavetoLayerFile_management('Winners_lyr','Winners_lyr')

The last line is because, according to the documentation, the layer is not saved but only locally created, since I want to be able to access it also after the code.

Winners is an list containing integer numbers, like this:

Winners = [1, 4, 59, 330]

ObjectID_12is the code that univocally identifies my elements ( I think it received the appendix _12 because this code is the result of several geoprecessing actions, including intersections)

when I run my code, anyway, I get this error:

 File "C:\Program Files (x86)\ArcGIS\Desktop10.4\ArcPy\arcpy\management.py", line 6477, in MakeFeatureLayer
raise e
arcgisscripting.ExecuteError: ERROR 000358: Invalid expression 0
Failed to execute (MakeFeatureLayer).

I guess the issue is the SQL condition I wrote, but I am not sure how to fix it since the MakeFeatureLayer_management is not very generous of information concerning how to use this function.
How can I solve the problem?


texts: Incorporating known GPS/INS data to automatically produce geo-registered 3D models

texts:
I am trying to reconstruct a 3D building model from 2D video frames using an open-source Structure-From-Motion(SFM) software (Bundler). The final reconstructed 3D model will be in Cartesian coordinates but ultimately it needs to be geo-registered so that I can overlay it on a virtual globe or in 3D analyst.

As far as I know, the typical SFM approach solves for relative camera poses and the 3D points. However, since I already have GPS/INS metadata for each video frame, is it possible to use those as an initial estimate for the extrinsic camera parameters (camera pose) to input into bundler, and get my 3D point cloud directly in the same world coordinate system directly?

I know a common approach (correct me if I am wrong) after the SFM process is to manually try to fit the model or find a rigid transformation to the world coordinate system, by using some marked ground control points on a reference imagery, but I just want to know whether this geo-registration can be done automatically instead by specifying known rough initial camera poses as input to Bundler (or any other similar sFm software) for subsequent refinement through bundle adjustment.

Thanks!


texts:
I'm relatively new to PostGIS and postgresql . I'm trying to build a simple query containing two relations as follows :

I've two relations R1 and R2 both of which contains geometries(Points) . I need to find out the closest point of every point in R1 which is in R2 . Say R1 is a relation of all hospitals and R2 is a relation of all medical shops , i need to find a medical shop closest to each hospital . More over i've a geometry types of srid 4326 in case you want to use ST_Distance or similar functions . Assume i need the medical shops within say 500 meters,say .


texts: Visibility analysis with mask for non visible areas?

texts:
I have LiDAR data and I've made a Digital Surface Model (DSM):

I've made a visibility analysis for 3 buildings with this model:

My goal is to know from where the buildings are visible, not what do you see from them. But as you can see in the orthophoto, there are a few forest areas where the trees block the view of the buildings from the ground level and yet, they appear as visible (from the top of the trees I guess). The same happens with the roof of another building:

Is there any way to use a layer (vector or raster) as a mask for non-visible areas when you run a visibility analysis?

I'm currently working with ArcMap 10.3 but I'm open to suggestions with QGIS 3 as well.


texts:
I need to download image patches from a WMTS server. However, WMTS only allows to request tiles specified by row and column of the tile matrix.

What I need to achieve is some kind of "WMS-like" query which lets me specify a spatial bounding box. I am aware that there is no server-side processing on WMTS servers. Therefore I am looking for a workaround.

I achieve the desired result by adding the WMTS as a layer to QGIS and use the "Save as..." function of the layer, which lets me specify the spatial extent and resolution. Since I need to do this for ten-thousands of image patches I need to find a way to automate this. Unfortunately, FME does not have a reader for WMTS.

Any ideas?


texts: Change Title Of Geoserver Legend Decoration

texts:
My colorramp legend's title gets set to be the layer's name. I'm wondering if I can override this and set the title explicitly. For example let's say I want to change the colorramp's title from "Spring Index, Daily Anomaly - Leaf Anomaly" to just "Leaf Anomaly".

Currently This is the WMS result I'm getting:

As you can see it's getting "Spring Index, Daily Anomaly - Leaf Anomaly" from the layer title:

Here is my layout:

<layout>
  <decoration type="text" affinity="bottom,left" offset="26,20">
    <option name="message" value="Difference from average in the timing of the Start of Spring, as estimated through an index of early-leafing plants"/>
    <option name="font-size" value="15"/>
    <option name="font-bold" value="true"/>
    <option name="font-family" value="Arial"/>
    <option name="halo-radius" value="4"/>
  </decoration>
  <decoration type="text" affinity="bottom,left" offset="26,0">
    <option name="message" value="Results Provisional  (www.usanpn.org/data/spring)"/>
    <option name="font-size" value="15"/>
    <option name="font-italic" value="true"/>
    <option name="font-family" value="Arial"/>
    <option name="halo-radius" value="4"/>
  </decoration>
  <decoration type="legend" affinity="bottom,right" offset="4,-55" size="auto"/>
</layout>

Here is my colorramp sld:

<?xml version="1.0" encoding="ISO-8859-1"?>
<StyledLayerDescriptor version="1.0.0" xmlns="http://www.opengis.net/sld" xmlns:ogc="http://www.opengis.net/ogc"
  xmlns:xlink="http://www.w3.org/1999/xlink" 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">
  <NamedLayer>
    <Name>days</Name>
    <UserStyle>
      <Name>days</Name>
      <Title>Days Difference</Title>
      <FeatureTypeStyle>
        <Rule>
          <RasterSymbolizer>
            <Opacity>1.0</Opacity>
            <ColorMap>
              <ColorMapEntry color="#FFFFFF" quantity="-9999" label="Day Difference" opacity="0.0"/>
              <ColorMapEntry color="#67001f" quantity="-20"   label="20 Days Early"/>
              <ColorMapEntry color="#b2182b" quantity="-8"    label="8 Days Early"/>
              <ColorMapEntry color="#d6604d" quantity="-6"    label="6 Days Early"/>
              <ColorMapEntry color="#f4a582" quantity="-4"    label="4 Days Early"/>
              <ColorMapEntry color="#fddbc7" quantity="-2"    label="2 Days Early"/>
              <ColorMapEntry color="#f7f7f7" quantity="0"     label="No Difference"/>
              <ColorMapEntry color="#d1e5f0" quantity="2"     label="2 Days Later"/>
              <ColorMapEntry color="#92c5de" quantity="4"     label="4 Days Late"/>
              <ColorMapEntry color="#4393c3" quantity="6"     label="6 Days Late"/>
              <ColorMapEntry color="#2166ac" quantity="8"     label="8 Days Late"/>
              <ColorMapEntry color="#053061" quantity="20"    label="20 Days Later"/>
            </ColorMap>
          </RasterSymbolizer>
        </Rule>
      </FeatureTypeStyle>
    </UserStyle>
  </NamedLayer>
</StyledLayerDescriptor>


texts:
I'm calculating areas of a KML in QGIS, ArcGIS and GeoPandas and all of then are different, my results are this:

Geopandas
ArcGIS
QGIS

29.3498
28.500084
25,651554

The original crs of the KML is 4326, so in QGIS and ArcGIS first, I reproject the KML (in ArcGIS first I have to make the KML a shp or feature) and then I calculate the area.
I checked this post about the same problem but with QGIS, and I calculate the area in QGIS like they said, in ArcGIS I did something similar, first I set the coordinate system to the dataframe and then re calculate the area, with GeoPandas I´m doing this:
import geopandas as gpd
import fiona

gpd.io.file.fiona.drvsupport.supported_drivers['KML'] = 'rw'
df = gpd.read_file('original_kml.kml', driver='KML')
gdf = df.to_crs('esri:102015')
gdf["area"] = gdf['geometry'].area/ 10**6

gdf.sum()

Just in case someone wants to check the process the files are here


texts: How to download a big area OSM by JOSM editor

texts: I download a small rectangle of data osm using JOSM, but when trying to download a big area the message error shows that the server will not allow it.


texts:
How can I add a geopandas geodataframe to a cartopy map? My objective is to create a large number of PNGs to turn into a GIF. The geopandas default plotting options are limited so I would like to leverage the abilities of cartopy.

For Shapefiles, I prefer the ease of use of geopandas over shapely (which geopandas uses under the hood) or other polygon libraries. Is there a way to combine the power of cartopy with the power of geopandas?

I found this article but it doesn't really combine geopandas with cartopy

This is what I tried so far. Also I noticed that cartopy can add Shapely Geometries to its Feature class but not geoDataframes. So far my conclusion is that using Cartopy's reader is easier than using geopandas if your goal is to add the data to cartopy.


texts: Finding rivers from GeoTIFF files with GDAL

texts:
I've grabbed a bunch of .tif-files (from http://dwtkns.com/srtm/) and I am trying to render it in 3D as terrain.
I've converted it to a 16-bit heightmap successfully but I am curious is there any way I can get a rivers-map/texture from GDAL? If there is no tool for this is there any algorithm to detect rivers?


texts: I'm using QGIS and I've been trying to use the fixed buffer on layer of a utility right-of-way(.shp file) in my map. It seems to be leaving parts of the input layer un-buffered. I've made my own shapefile from a line and I'm still unable to buffer that newly created line.


texts: Can I get a unique reference to a feature layer in ArcObjects?

texts:
I want to build an IDictionary with the feature layers of my TOC in ArcObjects. So I need key value pairs. The value will be an attribute of the feature layer but I am wondering what to use as key. I can use the layer name but there can theoretically be duplicate layer names in the TOC. Any suggestions on how I can get a unique reference to a feature layer?

Background is that I need to store the layers and their querydefinitions in a custom application that runs in desktop and ArcEngine environment. I don't know the layers before. It could be that the user defines layer groups that have sublayers with the same name. This would lead to an error when storing it in the dictionary.


texts:
I have downloaded data set from USGS.  I have images in tiff files and their data in XML format I want to upload it on ArcGIS Pro so that tiff files become labelled from XML files so that I can export them in chips and use them for deep learning.
What format does the ArcGIS Pro use for loading labelling?
How to convert XML to that format?
How to upload tiff file along with that format on ArcGIS Pro so that I can visualize the bounding boxes?
Below is an image of collection of all my dataset:


texts:
I am using the Identity tool which takes 3 parameters: in_features, identity_features, and out_features. Does the output FeatureClass need to have fields that are a combination of the fields in in_features and identity_features?

I keep getting an error when I the out_features only contains the fields in in_features and not both in_features and identity_features.

public static void PerformIdentityOverlay(string workSpace, string inputFeatureClassName, string identityFeatureClassName, string outputFeatureClassName)
{
     IWorkspaceFactory workspaceFactory = System.Activator.CreateInstance(System.Type.GetTypeFromProgID("esriDataSourcesGDB.FileGDBWorkspaceFactory")) as IWorkspaceFactory;
     IFeatureWorkspace featureWorkspace = workspaceFactory.OpenFromFile(workSpace, 0) as IFeatureWorkspace;
     IFeatureWorkspaceManage featureWorkspaceManage = featureWorkspace as IFeatureWorkspaceManage;

     IWorkspace workspace = workspaceFactory.OpenFromFile(workSpace, 0);

     IFeatureClass identityFeatureClass = featureWorkspace.OpenFeatureClass(identityFeatureClassName);

     IFeatureClass inputFeatureClass = featureWorkspace.OpenFeatureClass(inputFeatureClassName);
     IFields inputFieldsFromFeatureClass = inputFeatureClass.Fields;
     String shapeFieldName = inputFeatureClass.ShapeFieldName;

     IFeatureClassDescription fcDesc = new FeatureClassDescriptionClass();
     IObjectClassDescription ocDesc = (IObjectClassDescription)fcDesc;

     // Use IFieldChecker to create a validated fields collection.
     IFieldChecker fieldChecker = new FieldCheckerClass();
     IEnumFieldError enumFieldError = null;
     IFields validatedFields = null;
     fieldChecker.ValidateWorkspace = workspace;
     fieldChecker.Validate(inputFieldsFromFeatureClass, out enumFieldError, out
                                    validatedFields);

     IFeatureClass outputFeatureClass = featureWorkspace.CreateFeatureClass(outputFeatureClassName, inputFieldsFromFeatureClass, ocDesc.InstanceCLSID,
                                                ocDesc.ClassExtensionCLSID, esriFeatureType.esriFTSimple, shapeFieldName, "");

     Geoprocessor gProc = new Geoprocessor();
     Identity ident = new Identity(inputFeatureClass, identityFeatureClass, outputFeatureClass);
     IGeoProcessorResult result = gProc.Execute(ident, null) as
                IGeoProcessorResult;
}


texts: How can I get weighted Flow Accumulation to deal with negative values?

texts:
here's one I haven't seen anywhere, to my surprise. I am trying to do something a little different with Flow Accumulation, which I would think has been attempted before, but it's stumping me for now.

I want to weight flow accumulation by a loading map of pollutant z to get its accumulation as it's transported downstream. This is simple enough, except I'd like some areas to be net sinks -- that is, where a fixed amount of the pollutant is removed rather than accumulated. So, say our flow path contains cells with values 1 1 -1 1. -1 is a net sink, and so our accumulation here, with flow going left to right, would be 1 2 1 2. This too seems to be doable. I can make one weighting raster of 1 1 0 1 and another of 0 0 1 0, run the flow accumulation with both, then subtract the 2nd from the 1st. OK.

However, imagine this flowpath scenario -- 1 1 -5 1. Here I have a powerful "sink" at -5, which is in a position to remove everything that has flowed to it. Accumulation at that cell would be 2, and with 5 removed, this leaves us with -3 using the above method. But of course, having a quantity of -3 is nonsensical, and instead it should simply drop to 0, as everything is removed. At the end of this little 4-cell flowpath, accumulation should thus be 1.

Is there any way to weight or otherwise get Flow Accumulation to intelligently do the opposite of accumulating (de ... cumulating?) in certain specific areas, while also respecting a lower bound of 0? I'm hoping the method, if it exists, will be implementable in Python as well.

UPDATE:
If anyone is interested in this problem in the future, I may have figured out how to do it (though am still open to correction, as I'm not 100% sure). By default, the Flow Accumulation tool treats negative weightings as 0.

The steps are as follows:
1). Run Flow Accumulation with your weighting grid that includes negatives.
2). Apply the following conditional statement to the result:

Con(("FlowWeight" < 0) & ("FlowAccum" > 0),Con(-"FlowWeight" < "FlowAccum",
-"FlowWeight","FlowAccum"),0)

This assigns a non-zero value to all cells with a negative weight and a non-zero accumulation, such that if the weight represents uptake of some but not all of what has accumulated there, the cell gets the value of that amount, as a positive number, and if uptake is more than what has accumulated, the value is simply the accumulated amount.

3). Run Flow Accumulation weighted by the grid resulting from the above. This is like the "loss accumulation".

4). Subtract the result of step 3 from the result of step 1. Force any negative values to 0.


texts: Connect two fields belonging to two different shapefile

texts:
Is there a way in QGIS to connect two fields belonging to two different layers of the same type?

For example in my case I have a linear "father" layer containing all the information of the pipelines of my network, from this layer I extracted another "son" layer containing the same elements but with different fields, except for some that are identical (for example the lengths of the pipes).

I would like to make sure that when I update the lengths contained in the "father" layer, the lengths in the "son" layer are updated automatically. Can you help me?

I attach a screen of the two layers, in which I would like to make sure that if I update the "L_LUN (m)" field of the "father" layer, also the "Length" field of the "son" layer is updated automatically.


texts: Change the icon to a marker in Leaflet.Editable

texts:
I'm using Leaflet.Editable to add markers with different properties to a map.
Like: https://github.com/Leaflet/Leaflet.Editable/blob/master/example/index.html
How to customize the icon of the marker or color of the path according to the used tools?


texts:
I am trying to find the nearest point on the road network closest to a given point. I imported the road network into PostGIS from QGIS and am using the SQL window in QGIS DB Manager to query the road table. The query I used is:

SELECT *
FROM public."wa_transformed"
ORDER BY ST_Distance(
ST_StartPoint(geom),
st_SetSrid(ST_MakePoint(-122.2652671, 47.30995661), 4326),
true
) ASC
LIMIT 1;

The idea is to find the distance between all the roads in the road network and find the first item in the sorted list. . The red marker is the point (47.309, -122.265) whereas the thick brown line on the left is the query result (which is definitely not the closest source from red point).

I had to transform the roads (original SRID - 4326) to match the OpenStreetMap layer on QGIS using a CRS defined by (+proj=lcc +lat_1=45.83333333333334 +lat_2=47.33333333333334 +lat_0=45.33333333333334 +lon_0=-120.5 +x_0=500000.0000000001 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs). However, I get the same result whether I use the transformed table or original table.

I am using POSTGIS="2.5.2 r17328" [EXTENSION] PGSQL="110" GEOS="3.7.0-CAPI-1.11.0 3.7.1" PROJ="Rel. 4.9.3, 15 August 2016" GDAL="GDAL 2.2.4, released 2018/03/19" LIBXML="2.7.8" LIBJSON="0.12" LIBPROTOBUF="1.2.1" RASTER

The original roads shapefile is here and the transformed shapefile is here


texts: What algorithm should I use to remove outliers in trace data?

texts:
I have GPX traces and want to remove outlying observations. There are some route points that are impossible because it would require acceleration that the vehicle can not achieve.

Are there good algorithms for filtering out such extreme noisy observations?


texts:
I would like to use this function to update the marker every time the GPS is detected, it works but always adds a new marker and circle.

function onLocationFound(e) {
        var radius = e.accuracy / 2;

        L.marker((e.latlng), {icon: IsyUserIcon}).addTo(map).bindPopup("You are within " + radius + " meters from this point").openPopup();
        L.circle(e.latlng, radius).addTo(map);
    }


texts: GeoJson to Spatial Geography in sql server 2012, fix orientation of polygon

texts:
I've download from Mapzen many geoJson files related to countries with their respective administration level. Then I've parsed all the geoJson files and saved them to my sql server 2012 as geography spatial data.

I've then done a check of the Luxemburg map and noticed some strange results.
Some of the maps where created correctly, with the corresponding inner area:

In other cases the area of the selected boundary was the outer:

I do know that there is the left-hand convention for storing inner area of a polygon in sql server as spatial data. Clearly fixing each polygon one by one will be a waste of time.

GeoJson files aren't bounded to the left-hand rule, so I need a way to fix the coordinates orientation of the polygon before storing them in the db.

I want the orientation fix logic to be done in the c# code (MVC) not by the sql server.

I've read about a possible solution for convex polygons that handels the issue with cross-product, but still this doesn't work.

How can I solve this issue?
Is there a program that fixes the geoJson orientation of the polygons?
Is there an algorithm based on a mathematical principle?
Without going to the fuss of checking the orientation and all the extra coding, is there a site where I can download world data coordinates with all the administrative levels that use the left-hand convention for polygons?

Code used to get cross product

public decimal GetCrossProductZDirectionAndLength(VectorCoordinate coordA, VectorCoordinate coordB, VectorCoordinate coordC)
{
    //***********************************************
    //This method gets the cross product
    //for vectors:
    //BA X BC
    //***********************************************
    decimal crossProductDirectionLength = 0;

    //These vector coordinates are relative to the origin (0,0,0)
    VectorCoordinate vectorCoordinateBA = new VectorCoordinate();
    VectorCoordinate vectorCoordinateBC = new VectorCoordinate();

    vectorCoordinateBA.X = coordA.X - coordB.X;
    vectorCoordinateBA.Y = coordA.Y - coordB.Y;
    vectorCoordinateBA.Z = coordA.Z - coordB.Z;

    vectorCoordinateBC.X = coordC.X - coordB.X;
    vectorCoordinateBC.Y = coordC.Y - coordB.Y;
    vectorCoordinateBC.Z = coordC.Z - coordB.Z;

    crossProductDirectionLength = GetDeterminat2x2(vectorCoordinateBA.Y, vectorCoordinateBA.Z, vectorCoordinateBC.Y, vectorCoordinateBC.Z);
    crossProductDirectionLength += GetDeterminat2x2(vectorCoordinateBA.X, vectorCoordinateBA.Z, vectorCoordinateBC.X, vectorCoordinateBC.Z);
    crossProductDirectionLength += GetDeterminat2x2(vectorCoordinateBA.X, vectorCoordinateBA.Y, vectorCoordinateBC.X, vectorCoordinateBC.Y);

    return crossProductDirectionLength;
}

private decimal GetDeterminat2x2(decimal a, decimal b, decimal c, decimal d)
{
    //get the determinat for a 2x2 matrix
    //|a b|
    //|c d|

    return (a*d)-(b*c);
}


texts: QGIS crashes when packaging QGIS Project for export to QField

texts:
QGIS is crashing when I try to package a project for export to QField, using the QFieldSync plugin. Previous projects have succesfully exported without issue. I am running QGIS 3.4.11 Madeira on a Windows 10 laptop. The error report is pasted below. Any idea why this is happening?

Crash ID: 57306b64118c1c5d75ffb118873c33fc668bdf6e

Stack Trace

QBitArray::QBitArray :
QgsRelation::referencedLayerId :
QgsCoordinateTransformContext::calculateDatumTransforms :
QgsCoordinateTransform::transformPolygon :
QgsCoordinateTransform::QgsCoordinateTransform :
QgsMapCanvasAnnotationItem::updatePosition :
QgsMapCanvas::updateCanvasItemPositions :
QMetaObject::activate :
QgsMapCanvas::setExtent :
QgsMapCanvas::readProject :
QMetaObject::activate :
QgsProject::readProject :
QgsProject::readProjectFile :
QgsProject::read :
PyInit__core :
PyMethodDef_RawFastCallKeywords :
PyMethodDef_RawFastCallKeywords :
PyEval_EvalFrameDefault :
PyMethodDef_RawFastCallKeywords :
PyEval_EvalFrameDefault :
PyFunction_FastCallDict :
PyMethodDef_RawFastCallDict :
PyObject_Call :
PyInit_QtCore :
PyInit_QtCore :
PyInit_QtCore :
PyInit_QtCore :
PyInit_QtCore :
QMetaObject::activate :
QAbstractButton::clicked :
QAbstractButton::click :
QAbstractButton::mouseReleaseEvent :
QWidget::event :
QApplicationPrivate::notify_helper :
QApplication::notify :
QgsApplication::notify :
QCoreApplication::notifyInternal2 :
QApplicationPrivate::sendMouseEvent :
QSizePolicy::QSizePolicy :
QSizePolicy::QSizePolicy :
QApplicationPrivate::notify_helper :
QApplication::notify :
QgsApplication::notify :
QCoreApplication::notifyInternal2 :
QGuiApplicationPrivate::processMouseEvent :
QWindowSystemInterface::sendWindowSystemEvents :
QEventDispatcherWin32::processEvents :
CallWindowProcW :
DispatchMessageW :
QEventDispatcherWin32::processEvents :
qt_plugin_query_metadata :
QEventLoop::exec :
QCoreApplication::exec :
main :
BaseThreadInitThunk :
RtlUserThreadStart :

QGIS Info
QGIS Version: 3.4.11-Madeira
QGIS code revision: 9a8a6d4687
Compiled against Qt: 5.11.2
Running against Qt: 5.11.2
Compiled against GDAL: 2.4.1
Running against GDAL: 2.4.1

System Info
CPU Type: x86_64
Kernel Type: winnt
Kernel Version: 10.0.18362


texts:
I have a data from terrestrial scanning (3D) projected to 2D plane but still in XYZ format available here: http://ulozto.cz/xYkfywxY/xyz-to-2drst-txt.

Finally, my data are all on same plane, but every point has its XYZ information. I would like to export the XYZ point shapefile to raster file in 2D where 1 point in .shp = 1 pixel (pixel size not important for a while).

Please, any idea how to do that? I can use ArcGIS 10, R, ERDAS or some other software?

I found some information here: http://support.esri.com/es/knowledgebase/techarticles/detail/39771 but they finish with raster interpolation and I want the raster information in format: point = pixel, no point = NoData.

Every advice is really welcome ! XYZ point data in ArcScene representing two trees


texts: Preserve holes after Contour extraction in QGIS

texts:
I'm trying to get a clean vector layer from an elevation raster layer. I'm using Raster > Extraction > Contour... to get lines from the elevation raster here:

This tool produces this layer (I kept only one elevation level to keep it simple here):

I'm then using Vector > Geometry tools > Lines to polygons to transform these lines into polygons and it produces this result:

The problem is: All the parts which are below the main elevation contour are filled. But those should be holes in the main polygon.

How can I get the polygon with the holes?


texts: geoTIFF File conversion in QGIS - creation parameters not defined

texts: I want to convert a geoTIFF file into ESRI Grid Format (Arcinfo Binary Grid)using QGIS --Raster--Conversion---translate tool but I get a error message saying "creation parameters not defined". What is this and how to proceed next?


texts: How to combine MapInfo tab files into one spatialite DB?

texts:
I have many folders, each with a collection of MapInfo tab files (each folder has the same table names). For each folder I want to merge the tab files into one spatialite DB.

My current method is to use ogr2ogr to convert the tab files into shapefiles, and then use .loadshp to load each one into the DB. This truncates my column names, which is irritating.

Is there a way to avoid the shapefile part of this process? ogr2ogr can convert to spatialite, but I do not know how to go about merging multiple spatialite files.


texts: Select all fields from table and use with arcpy.DeleteIdentical_management()?

texts:
I would like to use arcpy.DeleteIdentical_management() function without having to specify all fields in the table.

Is there a keyword used to select all field in a table? I tried "*" or "ALL" but none worked.

Here is the code :

import os
import arcpy
from arcpy import env

table = "table1"
fields = ["a", "b", "c", "d", "e"]

arcpy.DeleteIdentical_management(table, fields)


texts: GeoServer create a datastore using PHP and the REST API

texts:
First of all I hope you forgive my bad English,I'm trying to create a datastore in GeoServer using PHP so can you please give me an example that can guide me to achieve it.

The code I have is:

$app->match('/send', function() use ($app){
$curl = new Curl();
$encoder =  [new XmlEncoder('datastore'), new \Symfony\Component\Serializer\Encoder\JsonEncoder()];
$normalizer = [new \Symfony\Component\Serializer\Normalizer\GetSetMethodNormalizer()];
$serializer = new \Symfony\Component\Serializer\Serializer($normalizer, $encoder);
$service = "http://localhost:8080/geoserver/"; // replace with your URL
$reqt = "rest/workspaces.xml"; // to add a new workspace
$url = $service . $reqt;
$curl->setOpt(CURLOPT_POST, True);
$curl->setBasicAuthentication('admin', 'geoserver');
$json = json_encode($curl->get($url));
$data = json_decode($json, TRUE);
$request = Request::createFromGlobals();

if($request->isMethod('POST')){
    $params  = [
            'name' => $request->get('name'),
            'type' => 'GeoTIFF',
            'description' => $request->get('description'),
            'url' => __DIR__ . '/' .$file->getClientOriginalName(),
            'enabled' => 'on' == $request->get('enabled') ? 'true' : 'false',

    ];
    $curl->setOpt(CURLOPT_UPLOAD, TRUE);
    $json = $app['serializer']->serialize($data, 'json');
    $curl->setHeader("Content-type"," application/json");
    $target = "http://localhost:8080/geoserver/rest/workspaces/{$request->get('workspace')}/datastores";
    return new Response($curl->post($target, $json));
}

return $app['twig']->render('form.html.twig', ['data' => array_shift($data)]);
})->bind('send');


texts: Exported DXF from QGIS won't open in AutoCAD

texts:
I exported a DXF from QGIS by right-clicking on the layer. Then tried to open in AutoCAD 2000.
That did not work. AutoCAD gave a failure-message:

Unknown Value "AC1018" found in drawingversion.

Then I tried to view the file in DWG-Truview2013 - which didn't work either giving the following error:

Invalid database handle 0 on line 1548.
Invalid or incomplete DXF input -- drawing discarded.

Does anyone know a solution? Or is there anywhere in QGIS change the DXF-Version to export?


texts: How to convert CityGML (LoD2) to Digital Surface Model?

texts:
I'm very new to GIS and doing this for a university project. For evaluating Solar Potentials on Building Roofs I want to use QGIS 3 (the UMEP Plugin, with SEBE tool). Amongst other data I need a Digital Surface Model preferably as a .tif file. For that I already have elevation data but where I'm really struggling are the CityGML files. The goal is to get moderately detailed information about the solar potential on differently shaped roofs in their geographic environment, hence roof shapes, position towards the sun and shades are very important.
Here is my problem:
I want to merge the CityGML (LoD2) buildings on the elevation map and from that create a 2D Digital Surface Model, with accurate representation of heights of the roofs. The DSM should reflect the shaping of the roofs not only cubics or similar.
I though about gdal_rasterize but am convinced it won't help me too much as it doesn't work for CityGML as far as I'm concerned. I also experimented a bit with CityJSON but at the moment I'm running out of ideas on what to google.
The building data is from this page (unfortunatelly in German):
https://www.hannover.de/Leben-in-der-Region-Hannover/Verwaltungen-Kommunen/Die-Verwaltung-der-Landeshauptstadt-Hannover/Dezernate-und-Fachbereiche-der-LHH/Baudezernat/Fachbereich-Planen-und-Stadtentwicklung/Geoinformation/Open-GeoData/3D-Stadtmodell-und-Gel%C3%A4ndemodell/Digitales-3D-Stadtmodell


texts: How to display Values as Labels in Thematic Map - Mapinfo

texts:
I was able to produce the Thematic Map for Agricultural Production for the various boundaries. How can I display the production value associated with each polygon as a label (to make it look like the the pic below)

?


texts: Getting geometry of Image in Google Earth Engine

texts:
I have an image which is from type 'ee.Image' and I want to extract its geometry so that:

im_geo = im.geometry();

but the geometry function of 'ee.Image' is DEPRECATED!

Is there another way for it?


texts: Is it possible to upload data directly to QGIS Cloud from a mobile phone?

texts:
I am going to be performing some field work, but I don't have a GPS unit as of yet due to cost.

Ideally I would like to pre-upload some features onto my QGIS cloud online map first, then go into the field and take elevation data from various points from my phone and have them uploaded on the spot to my online published map (as a point file).

Is there a way to do this directly without any other hardware?


texts: What kind of margin of error can I expect in figuring distance?

texts:
I have reviewed Why is law of cosines more preferable than haversine when calculating distance between two latitude-longitude points? and that does not equate to the specific question that I wish to ask because I'm looking for a simple answer about accuracy of multiple formulas.  While there is an answer there that discusses accuracy, it does it more as a comparison between two methods and depends on the reader understanding trigonometry to follow the answer.  I'm looking, mostly, for simple margins of error.

This is my first time trying to compute distance from two points on the ground.  I've been researching this and see that a lot of people use the Haversine Formula, but I also see comments, whenever I read about it, that there are accuracy issues.  Apparently, from what I see, Vincenty is more accurate, but I've also seen comments that the Law of Cosines is better.

What I don't see is any way to find out just how inaccurate I can expect a method to be.  Apparently there can be problems at extremely small distances in Haversine, but I'm not clear on that.

I'd like to just go on and use GeoPy in my program.  I'm basically working with distances of 100 miles or less.  Most of the time I'm comparing coordinates of two locations to see if they're close enough to be considered the same address.  If the accuracy is going to be within a few feet, that's close enough.

Can someone either point to a good resource or explain what kind of accuracy I can expect from the different methods I've mentioned and if there's a way to know what the conditions are that will give me more or less accurate answers?

I do see that GeoPy has a feature to tell me the confidence of a calculation, but I need to be able to know, in general, what kind of accuracy I can get from any one method overall, not for each individual calculation.

I see a lot about one method being better than another, but I don't see much on what kind of accuracy I can expect from any one method and what factors I can look at to tell me how accurate that method is.  That's what I'd like help finding.


texts: Specifying URL to ArcGIS Online organization or Portal for ArcGIS?

texts: I have downloaded Esri web appbuilder now I want run it but I'm stuck here in the URL.


texts:
I have developed a python script to copy and paste in the field calculator. Having an input of three values, it compares all of them and returns a list of unique values.
The problem is that when I use it in ArcMap, it only returns the first result for all the rows, even they should be different. Then, how can I evaluate each row?

Here is the code:

index =  0
list_yes = []
def myFunc(var_1,var_2,var_3):

        global list_yes
        global index

    myList_1 = [var_1,var_2,var_3]
    for i in myList_1:
        myList_2 = myList_1[index:]
        for j in myList_2:
            if i == j:
                if not j:
                    pass
                else:
                    list_yes.append(j)
            elif not j:
                pass
            else:
                pass

        index = index +1

    list_set = set(list_yes)
    unique_list = list(list_set)
    return str(unique_list)


texts: Is it possible to have multiple versions of the same ESRI Addin on the same machine?

texts:
I will be adding new functionality to a ESRI addin that I have created and I would like to give it to users to test.

They already have a working version.

Is it possible to install 2 versions of the same addin, one that is the current version (v2.1) and another that they can beta test(v3.0)?

If it is possible, how do I get this to work?

I am using ArcGIS 10.2.2.


texts: Building model for extraction of urban features on Google Earth Engine

texts:
I'm trying to build a model for extraction of urban features on Google Earth Engine, I've calculated NDVI & NDBI and subtracted them to get the built-up area but I'm unable to display it on the map. It is showing the error:-  Layer 3: Layer error: Parameter 'expression' is required.
var landsat = ee.ImageCollection('LANDSAT/LC08/C01/T1_TOA')
                .filterBounds(varanasi);
Map.centerObject(varanasi,9);
// Printing an image collection

var imgfunc = function(image){
  var result = image.normalizedDifference(['B5','B4'])
                    .rename('ndvi')
                    return image.addBands(result);
};
var ndviimage = landsat.map(imgfunc);
var ndviband = ndviimage.select('ndvi');
Map.addLayer(ndviband.mean().clip(varanasi), { min: -1, max: 1,palette: ['blue', 'white', 'green']},'NDVI');
//NDBI
var ndbifunc = function(ndbiimage){
  var ndbiresult = ndbiimage.normalizedDifference(['B6','B5'])
                    .rename('ndbi')
                    return ndbiimage.addBands(ndbiresult);

};
var ndbiimage = landsat.map(ndbifunc);
var ndbiband = ndbiimage.select('ndbi');
Map.addLayer (ndbiband.mean().clip(varanasi), { min: -1, max: 1,palette: [ 'yellow', 'blue']},'NDBI');
//var diff = ndviband.mean().subtract(ndviband);
//Map.addLayer(diff,'Built up Area');
var diff = landsat.mean().expression(ndbiband-ndviband);
Map.addLayer(diff, {min: -10, max: 10, palette: ['FF0000', '00FF00']});

Link to my GEE Code


texts:
I'm making a map on ArcMap and I converted the legend to graphics in order to manually manage my symbols and texts. For instance I want to increase the size of all the points in the legend to the same value. If I ungroup all elements I can then select all points one by one keeping control key pressed and then change the size for all of them in one move. Then I would need to regroup all the elements together in order to be able to move the legend around as a single piece. The problem is I couldn't find a way to select more graphics at once, for example with a rectangular selection. That would save me much time especially when regrouping.
Alternatively, is there a faster method to do this whole operation (changing legend items size)?


texts: How to find bounding coordinates (lat/long) for UTM zones?

texts:
I want to error check incoming lat/lng points by checking if the point resides within the given UTM zone. (This would work for the data im using because sometimes there is a missing sign or lat/lng is reversed and the UTM zone can still be correct)

So given a UTM zone, how would you find the bounding coordinates for that zone?(Either the lat/lng for SW and NE corners or a range of lat/lng values for a given zone)

Also: Since I am supplied with both UTM and lat/lng I can convert both and cross compare as a way to validate incoming points but I thought it would be interesting to know how to calculate a UTM zone's bounding box in lat/lng form just given a zone (and possibly a hemisphere).


texts: I've just been looking at Ext 5.0 which has some changes. I'm wondering if it's compatible with GeoExt 2.0. There seems problems with GeoExt.Popup


texts: Cleaning LiDAR points that classification algorithm misses?

texts:
There are numerous methods for classifying ground and not-ground points from raw LAS point clouds. They all do a great job, depending on local circumstances, but inevitably even when algorithm and terrain are well matched there are some exceptions. How to clean up those outliers? Are there any tools out there in which one can interactively lasso a collection of points and manually change some attributes? E.g. like in a paint program?

Horizontal 3D view of LAS point cloud with all "non-ground" points removed, however some points are mistakenly still classed as "ground".

Viewpoint top down view


texts: Running Buffer_analysis in PyScripter fails to execute

texts:
For some reason, I can't get Buffer_analysis to run successfully from PyScripter(2.5.3.0) for a particular shapefile.  It fails to complete and exits with error 999999.    I can run the exact same command and data in ArcMap(10.1 SP1) and it completes with no errors.

Can somebody suggest a reason why this may happen and possibly a solution?  I have noticed this with other commands as well.  It runs fine in ArcMap, but has problems in PyScripter.


texts: Running GRASS functions in QGIS 2.4 from python console

texts:
I am using qgis 2.4 and i am looking for the "easiest" way to execute e.g. the r.in.gdal function from the python console in qgis.

Specifically for r.in.gdal, I could not found this function in processing module (i.e. import processing; processing.alglist('r.in')), thus, processing module seems not to be the way to go ?!?

Before executing a grass command the mapset needs to be opened. Which command does this job?


texts: How to detect openlayers feature vector intersections in any zoom level?

texts:
I have to change vector feature when vectors are intersect each other,
I wrote the following code:

// check all vectors features (arrows) in the current layer for intersections
  allArrowFeatures = pointsLayer.features;
  var hideArr = new Array();
  var l = 0;
  for(var i = 0 ;i < allArrowFeatures.length;i++){

    arrowObj = allArrowFeatures[i];
    for(var j = 1 ;j < allArrowFeatures.length;j++){
      if(j===i){ // same vector
        continue;
      }
      if(allArrowFeatures[j].geometry.bounds.intersectsBounds(arrowObj.geometry.bounds)){
        hideArr[l] = pointsLayer.features[i].style = {display:'none'};
        l++;
      }

    }
  }

When you zooming in to the max zoom you wont see any intersections so I assume this code is OK but the problem is when I zooming out the vectors still intersects each other.

Its look like I don't use the correct function for my case, and I tried several ways this was the closest to my case.
Any idea what openlayers method I have to use to detect if there was intersection also in zoom out probably not base on the bounds?


texts: geocoding CSV using MMQGIS and error massage

texts: While geocoding a csv file using MMQGIS utility in QGIS, I provide an input and the follwing errror message pops out, what could be the problem?


texts: WGS84 geographic vs projected representation in GIS

texts:
I'm new at this (studying for a masters of GIS) and trying to understand the WGS84 geographic coordinate system as presented in ArcPro
WGS84 is a datum and reference ellipsoid
EPSG 4326 is WGS84 using latitute and longitude as a coordinate system
The data is have is being represented in ArcPro in WGS84 with lat/lon coordinates, so really EPSG 4326.
However, WGS84 is an ellipsoid, while the representation I see in ArcPro is flat. So the coordinates are being projected in some way. My lecturer assures me that there is no projection. When I use the "measure distance" tool the line between two points is curved since it's a geographic system. I get that, that part is fine. But the map is still flat, so I keep insisting there has to be some sort of projection. If there was no projection I'd be looking at an ellipsoid. What am I missing?


texts: Visualizing specific class from unsupervised classification in Earth Engine

texts:
I have classified an image into 20 classes using unsupervised classification for a specific region in Google Earth Engine. I want to check a specific class (say class 9) and visualize only that specific class.

How can I do that?

//loading the true image
var poi=ee.Geometry.Point([78.032188,30.316496]).buffer(35000);

var dataset = ee.ImageCollection('LANDSAT/LT05/C01/T1_TOA')
                  .filterDate('2010-01-01', '2010-12-31');
var trueColor321 = dataset.select(['B3', 'B2', 'B1']);
var tru=trueColor321.reduce(ee.Reducer.mean());
var truclip= tru.clip(poi);
var trueColor321Vis = {
  min: 0.0,
  max: 0.4,
  gamma: 1.2,
};
Map.setCenter(78.032188,30.316496, 9);
Map.addLayer(truclip, trueColor321Vis, 'True Color (321)');
//ndvi calculation
var dataset1=ee.ImageCollection('LANDSAT/LT05/C01/T1_32DAY_NDVI').filterDate('2010-01-01','2010-12-31');
var NDVI= dataset1.select('NDVI');
var ndvimean=NDVI.reduce(ee.Reducer.mean());
var ndviclip=ndvimean.clip(poi);
var visparam={
  min: 0.0,
  max: 1.0,
  palette: [
    'FFFFFF', 'CE7E45', 'DF923D', 'F1B555', 'FCD163', '99B718', '74A901',
    '66A000', '529400', '3E8601', '207401', '056201', '004C00', '023B01',
    '012E01', '011D01', '011301'
  ]
};
Map.addLayer(ndviclip,visparam,'NDVI',0);

//clustering
var training=ndviclip.sample({
  region : poi,
  scale : 30,
  numPixels : 5000
});
var clusterer=ee.Clusterer.wekaKMeans(100).train(training);

var result = ndviclip.cluster(clusterer);

Map.addLayer(result.randomVisualizer(),{},"clusters");


texts: Why outside border of service area always unreasonable?

texts:
In ArcGIS 10.1 Network Analyst, why results of service area are different if number of break points is changed?

For example, the grid is 300m*300m, I try to produce service area around a node.

settings are like this...

Then I got this result

then I only change breaks to 300,600, result:

Why the outside border is always in strange shape?

Why the service areas of 600m are different in these 2 cases?

I try to change polygon generation to "detailed", it still generate this kind of results.

Really confused, and I have problems understanding any service area result because of this.


texts: Problem inserting new features to file geodatabase feature class using ArcObjects 10 SDK

texts:
I'm trying to add new feature to file GDB feature class and keep getting this error: Objects in this class cannot be updated outside an edit session. Feature class has attachments enabled (exists relationships). Googled alot about this problem and still can't figure it out. FYI I'm developing a standalone application. My code:

Public Shared Sub CreateFeature(ByVal workspace As IWorkspace, ByVal featureClass As IFeatureClass, ByVal point As ArcGIS.Geometry.IPoint, ByVal data As String)
    Dim workspaceEdit As IWorkspaceEdit = CType(workspace, IWorkspaceEdit)
    workspaceEdit.StartEditing(True)
    workspaceEdit.StartEditOperation()
    Dim featureBuffer As IFeatureBuffer = featureClass.CreateFeatureBuffer()
    Dim featureCursor As IFeatureCursor = featureClass.Insert(True) ' <- at this point error occurs
    Dim featureOID As Object
    featureBuffer.Value(featureBuffer.Fields.FindField("Data")) = data
    featureBuffer.Shape = point
    featureOID = featureCursor.InsertFeature(featureBuffer)
    featureCursor.Flush()
    workspaceEdit.StopEditOperation()
    workspaceEdit.StopEditing(True)
End Sub


texts: How do I convert an active kmz file to a shapefile using either ArcMap or QGIS

texts:
I need to convert an active kmz file to shapefile. I have a layer file and a kmz file which feed data from a server. If I try to convert using the kml to layer tool in ArcMap it fails as there is no underlying data. I am also not able to export the layer file to a shapefile for the same reason. Both the layer and the kmz file work well as the data is being fed directly from the online server. A few posts suggested using Arcmap2SLD but this crashes as soon as I try to run it. I have tried running it on a couple of computers with no luck. I am running Windows 7.

I'm attaching an image of a couple of kmz files loaded into Google Earth Pro. The first two are ordinary kml files. The last one (in Temporary Places) is an active kmz file which link directly to an online server (in this instance a government server). You can see this by the difference in the icons.


texts: how to calculate the average temperature of each band of a .grd file?

texts: I am having a GRD file in which the temperature at different longitude and latitude are given in a 32x35 array for one day(band),i am having 365 such bands and want to calculate the average of all individual bands.


texts: Implementation of Esri ArcGIS Plugin Data Source using Monogdb

texts:
I worked with plugin data source using file data and it worked fine as given here. For MongoDB I downloaded the code given in GitHub and compiled it and registered it with categories.

How can I use this in ArcMap like will I get a new connection option in Add new database connection? I do not know how I will implement this.

For Simple file I used Java code and deployed the jar file to  ArcGIS by placing the JAR file in the <ArcGIS Install Dir>\java\lib\ext folder and when I opened ArcCatalog I saw the custom file as a feature class.

For Mongo I used the .Net code but I do not know after registering the dll with categories page how it can be accessed in ArcMap. i.e. How i can connect to MongoDB in ArcCatalog?


texts:
I'm searching for SVG of different (german) hydrant types ("Unterflurhydrant", "Überflurhydrant", "Hydrantenschacht") etc.

The symbols are currently stored as polygons in MapInfo's proprietary .TAB (the equivalent of .shp) format. You specify the tab, tell which feature you want as symbol and voila - that symbol is used for all hydrants of type xy, stored as polygons.

I have now to setup our look'n'feel in QGIS. So I need our polygon-symbols as SVG OR I get the symbols as SVG from somewhere.

I need to use them for visualization purposes in QGIS. Any possible/suggested data sources?


texts: Avoiding Legends too large to display in ArcGIS message from ArcGIS for Server WMS service?

texts:
I have created a WMS on ArcGIS for Server 10.2 with a polygon layer which contains a unique value renderer with 4 categories. The resulting legend item can be called via service call, which responds with an image of 22 width and 64 height in pixels, this legend item is not displayed in the table of contents. (other layers symbolized with a single symbol of 12x12 pixels have a legend in ToC). When I open the properties of the WMS Server and click on the layer which has no legend in Table of Contents (ToC), I can see:

  Legend Information: Available(Legends too large to display in ArcGIS)

A little googling results in this page: Using WMS service layers

It is documented that...

  ...when the image received for a layer in a WMS service exceeds this
  limit, no legend is displayed for the layer

Is there a way to alter the limits of the legend item shown in ArcMap?


texts: pyqgis - refreshing only QgsMapCanvasItems not the entire QgsMapCanvas

texts: I have a fairly complex map, consisting of many layers, and despite they're optimized, they still take a certain time to refresh. I also have a QgsMapCanvasItem representing my position, which i would like to have updated about 20 times per second. This could be easily done, if it was possible to only update the position of the QgsMapCanvasItem, and use a cached image of what is displayed in the map canvas. Is this possible?


texts: How to get started with pgRouting and OSM?

texts:
I have followed these instructions to install PostgreSQL and pgRouting on our Ubuntu server:

http://obsessivecoder.com/2010/02/01/installing-postgresql-8-4-postgis-1-4-1-and-pgrouting-1-0-3-on-ubuntu-9-10-karmic-koala/

I am working with united_kingdom.highway.osm, It all seemed to go ok.

I have pgAdminIII on my Windows machine, and can connect to my PostgreSQL db and I see the tables: geometry_columns, planet_osm_line, planet_osm_nodes, planet_osm_point, planet_osm_polygon, planet_osm_rels, planet_osm_roads, planet_osm_ways, spatial_ref_sys.

Now, I have not done anything with PostgreSQL before, and I am looking for the simplest possible way of getting something out to show it is working, like finding the distance between two UK points.


texts:
Is it possible to set env.workspace = "in_memory" in ArcPy using both ArcGIS Pro and the ArcGIS 10.2.2 (or 10.3) architectures?

What I am trying to do is get the output of a snap pour point operation written to memory as opposed to disk.

I realize I could write to disk then bring it into memory but this would not help. I am trying to optimize a series of tasks minimizing the use of writing to physical media the output of a process that is merely required as an input to the next process.


texts: How to visualize pgrouting query results in mapguide open source

texts: I would like to make a spatial query using postgis 2.0 and pgrouting 1.05 then visualize results in Mapguide open source.Any suggestions?


texts: Using values from a tool parameter defined as multiple value

texts:
I am trying to use values from a parameter that accepts multiple value and use the in the iteration but am getting an error. I have failed to understand where the error is coming from because the code was working fine. The values are obtained from the interface below

The user can select one or more values.The code is below

import arcpy,os

inputraster = arcpy.GetParameterAsText(0)
fc = arcpy.GetParameterAsText(1)
sea_level = arcpy.GetParameterAsText(2)
clipraster = arcpy.GetParameterAsText(3)
outputfolder = arcpy.GetParameterAsText(4)
valuelist = [x.strip() for x in sea_level.split(";")]
for value in valuelist:

 arcpy.Clip_management(inputraster,"#",clipraster,fc,"#","ClippingGeometry","MAINTAIN_EXTENT")
        clip_raster = arcpy.Raster(clipraster)
        floodedarea = os.path.join(outputfolder + "floods" + str(value))
        areasbelow = clip_raster <= value
        areasbelow.save(floodedarea)
        arcpy.AddMessage(areasbelow)

Edit:

The error message:

ERROR 000732: Input Raster: Dataset 1 does not exist or is not supported


texts:
I'm trying to convert Lat/Lon coordinate to decimal and vice-versa.

Now using https://stackoverflow.com/questions/1140189/converting-latitude-and-longitude-to-decimal-values.

When I convert 16° 17′ 31″ S 48° 46′ 31″ W to decimal I got -48.775277777777774, -16.291944444444447 but I can see the original coordinate from OpenLayers is -48.775058943087124, -16.291811459160023.

My steps:

Take coordinate from map (-48.775058943087124, -16.291811459160023).
Convert to DMS using hdms = ol.coordinate.toStringHDMS( coordinate ). Got 16° 17′ 31″ S 48° 46′ 31″ W.
Convert DMS back to decimal using the function from link above. Result is -48.775277777777774, -16.291944444444447.

How can I achieve the correct precision?


texts: How to let ArcGIS Server see related records in a geodatabase

texts:
I have an ArcSDE geodatabase which contains multiple relationship classes.

If I add this layer to ArcMap, the relationship is automatically applied (eg, the Identify dialog allows me to drill-down to the related features).

However, if I publish this map document to ArcGIS Server, the REST endpoints for the layer don't show the relationship classes.

Similarly, if I create a new ArcMap document and create a Relate within ArcMap (using > Joins and Relates in the TOC) this also fails to appear in ArcGIS Server's REST directory.

This sample REST directory does show the related data. How has this been set up?


texts:
I have GWC layer on my GeoServer. Now I want to  print this layer using GeoServer print module but don't know how to do this.

When I try WMS layer all works fine. But how to print the GWC layer?

MORE INFO

In client side I have a GXP and try to print using gxp_print plugin. My problem is that I can't set URL to the GWC layer.

I read this and this: so I try this URL:

layer.setUrl("/geoserver/gwc/service/wms?BBOX=7637927.0361051%2C6845332.4516735%2C7750900.9638949%2C6971453.5483265&HEIGHT=256&WIDTH=256&TRANSPARENT=true&map_resolution=305.74811309814453&VERSION=1.1.1&FORMAT=image%2Fpng&TILED=true&SERVICE=WMS&format_options=dpi%3A100&LAYERS=gisurfo:atd&REQUEST=GetMap&STYLES=&SRS=EPSG%3A900913&");

But GeoServer returns me error:

400: Requested horizontal resolution: 441.30440542890574 , best match:
305.74811309814453 exceeds 10% threshold. Perhaps the client is configured with an
incorrect set of scales (resolutions), or the DPI setting is off compared to the one in GWC ?

This same like in links above. So how to print GWC?

UPDATE

Good news. I changed scales in my config.yaml file to scales that I get from the GWC layer (but I get only left part of each scale). Now I try to print GWC layer but it not printing. In GeoServer log I found Bad Request:

java.io.IOException: Error (status=400) while reading the image from
http://geourfo.ru:9091/geoserver/gwc/service/wms?BBOX=6866019.0%2C7636631.0%2C6950767.0%2C7752196.5    &HEIGHT=625&WIDTH=458&VERSION=1.1.1&FORMAT=image%2Fpng&SERVICE=WMS&
format_options=dpi%3A75&LAYERS=gisurfo%3Aatd&REQUEST=GetMap&STYLES=&SRS=EPSG%3A900913
: Bad Request

After I change width and height to 256/256 and get my tile!

So question: why iprint plugin take tiles with HEIGHT=625 and WIDTH=458 but not with 256x256?


texts:
I have just started using QGIS for an assignment. I'm just new on that so I will try to explain what my question is.
I have added a map  of Vic Australia from near map. I used the GDA94/ MGA zone 55 for CRS. I then added an OpenStreetMap map (using the same CRS) which is aligned ok. On the next step, I'm trying to add some vectors (shapefile) I have download from a government site. However, when I add the vectors, their coordinates are wrong (although I use the same CRS). The location from these vectors is correct (you can see the features) but they are aligned.

Can anyone explain me what am I doing wrong?


texts: How to search for several different features in the attribute table?

texts:
I guess my problem is a quite simple one, but the query builder and I have our difficulties, so maybe someone can help me.
I have a long attribute table in QGIS and am trying to extract about 200 fields with a certain ID, e.g. ID 4440 and 3376 and 3156 and so on. It's too much to do it manually and I'm pretty sure that there is a neat command, but somehow everything I tried is not working and thus I'd be really grateful for some help!
Thx!


texts: How to choose a field and Apply a General Symbology in Python

texts:
I did a very simple field to apply a symbology layer:

import arcpy

mxd = arcpy.mapping.MapDocument('CURRENT')
df = arcpy.mapping.ListDataFrames(mxd, "")[0]
in_layer = arcpy.mapping.ListLayers(mxd, "", df)[0]
in_symbology_layer = arcpy.mapping.Layer(r"F:\PercACH.lyr")

arcpy.ApplySymbologyFromLayer_management(in_layer,in_symbology_layer)

This *.lyr works well, but only in one field (year 2003). I want to do the process below:

Choose a feature class in GDB -> Choose the Field -> Apply a Symbology -> Generate a layout and export in PDF or PNG.

The export process I already did, but not with the field chosen step.

EDIT: all the fields has the same range.


texts: How to access the result(s) of a Processing module?

texts:
I am experimenting with the QGis Processing toolbox. I open up the module dialogue, set it then parameters and click. If the run succeeds the dialogue is closed automatically and no result(s) are ever shown.

Under the Processing menu there is a dialogue called Results Viewer, but it is always empty.

How can I access the results of a Prossessing module after the dialogue is automatically closed?


texts: Style a feature in one layer based on an attribute in another layer in Print Composer Atlas

texts:
I have two layers; one is the foreground and the other a background. Both layers have fields containing similar values (field values in the foreground layer (mun_code) is also present in the background layer's fields (ADM3_PCODE). For instance, both fields may contain the value 012802000).

For the background layer, I am trying to use Rule-based styling to add a thicker border using the expression:

@mun_code = attribute( $currentfeature, 'ADM3_PCODE' )

but @mun_code returns NULL (most likely because mun_code as a field is only present in the foreground layer while the background layer contains the field ADM3_PCODE):

I already Generated an atlas choosing the foreground layer as the Coverage layer. My understanding is that $currentfeature (in the code above) is already automatically controlled by the Print Composer Atlas, so I am wondering why it returns NULL.

How can I fix this?


texts: Draw shortest line between a point and a line

texts:
I am trying to draw a line between a thousands of points and a line using QGIS. I then want to carry out some calculations on the new lines to see which polygons they intersect in a different layer.

I have used NNjoin to find the shortest distance between each point and the line but I want a layer displaying these lines.


texts: Using Computer Games to display GIS data?

texts:
Besides being a GIS user I'm a big computer game fan.

Games like Civilization, Minecraft, World at War and many others have a very powerful engine to display huge maps in a beautiful way.

Is it possible to use a game engine to display GIS data for a better comprehension or visualization of the environment?


texts: Seeding multiple bounding boxes using GeoWebCache

texts: We have a system with nationwide coverage but only need to seed major cities. I am using the RestAPI but cannot seem to find a way to seed in multiple bounding boxes with a single xml config file. It would be really useful to be able to specify multiple bounding boxes in a single seeding operation.


texts: Center Openlayers map based on amount of polygons

texts:
I have number of polygons on my map, how can i center the map and set zoom so all fits in the map,

I managed to get the center of each polygon and put a marker on the center using  getExtent() would it be the same for the map?

heres my code
https://codepen.io/anon/pen/bqyqog?editors=1010


texts: Percentage of polygon in one shapefile within polygon of another

texts:
I'm a newbie, apologies if this is obvious/has already been asked and answered but I couldn't find anything.

I have two shapefiles: 1. an administrative boundary layer for a county in the UK known as an LSOA boundary that has 500 little zones in it 2. a flood zone.

Ideally I want to find out which of the little LSOA zones are ≥50% within the flood zone and end up with a yes/no or a 1/0 for each of the 500 LSOA zones.

But I don't know how to do this. I figured I could Join the two shapefiles, but there's no common attribute between them. Then I thought I could use the Join Attribute by Location function, which worked and shows me which LSOA are in the flood zone, but that's nearly all of them (see image 2).

I think this is an SQL problem but I don't know. I'm new to QGIS and have never used PostgreSQL.

Any help would be greatly appreciated. I can provide whatever info you lovely people need to help me out.