texts:
I've a vector layer with a bunch of draw controls on it. You can add polygons, points and you can modify polygon.

Now I changed the style of the points so they have a nice pin-like image. But them I noticed when I select the modify polygon tool and select a polygon, the handles on the polygon are also replaced by my point-image. How do I fix that?

this.panel = new OpenLayers.Control.Panel({
        allowDepress: true,
        autoActivate: false
    });
    this.controls = {
        "polygon": new OpenLayers.Control.DrawFeature(this.layer, OpenLayers.Handler.Polygon, {
            "title": "New Zone"
        }),
        "POI": new OpenLayers.Control.DrawFeature(this.layer, OpenLayers.Handler.Point, {
            "title": "New POI",
            "displayClass": "olControlDrawPoi"
        }),
        "modify": new OpenLayers.Control.ModifyFeature(this.layer, {
            "standalone": false,
            "title": "Edit Zone"
        }),
        "clear": new OpenLayers.Control.Button({
            displayClass: "olControlClearLayer",
            trigger: myLocations.clearLayer,
            title: "Clearlayer"
        })
    };
    for (key in this.controls) {
        this.panel.addControls(this.controls[key]);
    }

    myMap.map.addControl(this.panel);

I tried adding the style in the Point Handler but that has no effect.
I also tried to set a feature property so I can use a symbolizer in the style, but i can't figure out how to add this property before the feature is actually drawn. (so this works, but only after a layer.redraw().

So I want to have defautl styling for the polygons and its handels when its in modify mode, but I do want a nice marker image for my newly drawn points.

Any ideas?


texts: How to cluster points based only on location

texts:
I have profile data points which lie almost on a line (s. picture below)

Is there an easy way / plugin to add an attribute to each point which states to which profile line it belongs?

It should not be very difficult from a mathematical point of view in my opinion?

In other words: I want to create clusters from my points and each cluster should contain the points that belong to one profile line. So that in a second step I can fit a line through each of the profile points


texts: How to handle User Created Project Variables with qgis API

texts:
I am developing a project and a plugin for a company working with QGIS 2.14.4.

I want to update - through the plugin - a user created project variable which is used for rule-based styling purposes by i am not able to find the way by reading the api.

I found a workaround by using QGIS' own project variables in particular QgsProject.instance().setTitle() and changing the layer styles' rules to use @project_title variable but this is far away from clean.

Can anyone help with how to access and change user created project (or even layer) variables?


texts: Opening CTF and DAT files in ArcGIS for Desktop?

texts:
How does one open CTF and DAT files as listed here using ArcGIS 10.3?

I looked for instructions on their WELCOME to the Geoscientific Database Guide website but there isn't much guidance.


texts: Specifying MINUSERPIXELVALUE when clipping raster file?

texts:
I am trying to clip a raster file. When I do so, the following message come up;

  ERROR 5: MINUSERPIXELVALUE must be specified.


texts: Set unique values to different groups programmatically

texts:
I have a polyline feature layer. After doing some computation, I update a field of that layer and want to set unique color for different value of that field. We can do this from symbology using unique values. But I want to do it programmatically.

I know there should be a very easy way to do this thing with arcObjects. But I couldn't find the way.

I go through How to define a renderer for your layer document. But I did not find the way to set "Value Field" and "Add All Values" options here. And of course I have tried to use IUniqueValueRenderer.AddReferenceValue Method. But the document is not giving me all necessary information.

Question : How to set unique color value for different value of a specific field using c#? It will be vary helpful if anyone can give me some pointer about this.

Thanks.


texts: QGIS hyperlink to server location- working in 3.0 but not 2.18.20?

texts:
This is not necessarily a pressing issue, as there is a workaround, but it has piqued my curiosity about how QGIS manages hyperlinks through its 'Action' functionality. Please feel free to correct my if I get some terminology wrong here, I'm not strong on the IT side of things.

I have a shapefile with an attribute field containing a full filepath to a server location (\\IP address\folder\folder\file.format). This is to prevent issues with trying to access these files from different machines where the drives may be mapped differently, as opposed to directly to the drive, for example (C:\etc.).

This has worked perfectly well both in ArcMap and QGIS 3.0, setting up the hyperlink and clicking the feature behaves as expected, the file opens in the correct program. Note that as the full filepath is stored in the attribute the Action script just has the field [%Hyperlink%] and nothing else.

However, using QGIS 2.18.20- which is what most people who are accessing this data are using, and this is why it is important- once the action has been set up, using the action tool on a feature elicits no response at all, not even an error message. Is this something to do with mapping the file location to the server directly? If so, how and why do the two versions of QGIS treat this differently?

Edit: Upon implementation of the recommendation by eurojam, using a direct python script to access the file, I have at least got an error message. It seems that possibly there is some issue with how the filepath data is being stored that is not reflected in the attribute? the error message reads:

<
An error occurred during execution of following code:
import os
os.startfile('\172.17.50.18\vsdcid117_vol3_archive_srrc\Assets\IMG_US_Roads_March_2018\9C2B20A700335B4_1_P\9C2B20A700335B4_1_P.gps', 'open')

Traceback (most recent call last):
  File "", line 2, in
WindowsError: [Error 2] The system cannot find the file specified: '\172.17.50.18\x0bsdcid117_vol3_archive_srrc\Assets\IMG_US_Roads_March_2018\9C2B20A700335B4_1_P\9C2B20A700335B4_1_P.gps'
>

The filepath in the first section is how it should read, and this is what shows in the attribute table, however the second iteration, which it appears is the one that Python is attempting to work with, has some interesting formatting variants- the double backslash, and 'x0bsdcid117' where it should read 'vsdcid117'. I used the ArcMap field calculator to build the filepath string, and used the double backslash as an escape character so that the string would end up reading as it does once the filepath string had been concatenated with the file name. Does this have something to do with it?


texts: Getting and setting the symbol sizes for a graduated layer in QGIS 3

texts:
I have a whole set of maps that have graduated colors based on a set ramp and I also like to change the symbol sizes as it makes anomalous areas stand out even more.
Id like to change the symbol sizes with the console (if I use the automatic classification by value the symbol sizes are reset).

In QGIS 2.18 I could do it by:

layer = qgis.utils.iface.mapCanvas().currentLayer()
symbols = layer.rendererV2().symbols()
symbols[0].setSize(1.5)
symbols[1].setSize(1.75)
symbols[2].setSize(2)
symbols[3].setSize(2.25)
symbols[4].setSize(2.5)
symbols[5].setSize(2.75)

In QGIS 3 the rendererV2 has been renamed (to renderer), but I am stuck on how to call the symbols. I know its likely a python 3 syntax I am missing, but it hangs on

symbols = layer.rendererV2().symbols()

TypeError: symbols(self, context: QgsRenderContext): not enough arguments

I did manage to get at the symbols a different way:

layer = qgis.utils.iface.mapCanvas().currentLayer()
renderer = layer.renderer()
renderer.ranges()[1].symbol().size()
ranges = renderer.ranges()
ranges[0].symbol().size()

but if I do:

ranges[0].symbol().setSize(1.5)

it does not actually change anything and a layer refresh does not help. If I query the size again it tells me the original size, and I guess it is the default value the graduated renderer is using.

In the change log for v3, QgsSymbolLayer (renamed from QgsSymbolLayerV2) says usedAttributes() now requires a QgsRenderContext argument.
I've been poking at this for a long time and I feel I am close, but I just am not familiar enough with python to get over this hurdle.

QGIS API


texts: Is it possible to render map tiles with a Nominatim PostGIS database?

texts: I've imported a OSM database using Nominatim's setup.php. Is it possible to use some software to render the tiles using that database, or do I have to re-import in a different database for Mapnik and similar?


texts:
I would like to merge 2 shapefiles to create just one:

1 - a shapefile for the 9 GOR English regions with
2 - a shapefile which has the rest of the UK: Scotland, Wales, Northern Ireland and England

The aim is to end up with a shapefile with 12 parts (9 GOR + Scotland + Wales + NI) so I can see difference on certain factors across the UK.

Does anyone know how to do this?

Just to be clear, I don't want to overlay them, I want them to sit side by side so I can see all the parts at once.

Alternatively if anyone has a .shp file for the 9 GOR English regions + Scotland + Wales (NI not necessary), that would also be great.

I have tried to follow these instructions: https://www.statsilk.com/maps/merge-multiple-map-layers-single-shapefile-using-quantum-gis#shapefile-preparation but the "Merge shapefiles to one" is not in my vector menu and I can't seem to be able to download the Merge Shape plugin due to "bad gateway".


texts: Raster map in different zone () than region ()

texts:
I am new to QGIS (just installed QGIS 2.0.1) and am having problems running any spatial analysis tools. I followed the instructions on how to call GRASS tools into QGIS, but when I try to run the r.walk tool I get an error.

￼Raster map is in different zone (14) than current region (0)

￼Unable to open raster map

￼Raster map is in different zone (14) than current region (0)

￼Unable to open raster map

￼Unable to open raster map

I used the r.region tool to find try to set the region to the raster and also checked all the layers to make sure they had the same CRS which they did. I created the slope data from my DEM in GRASS, so i know that they are the same, and the vector data already had the same CRS.  I have tried a number of different approaches from what I have found on forums, but none have solved my issue.  I have no idea how to proceed. Any suggestions?


texts:
ArcGIS 10.2 Server enforces case sensitivity on map service URLs.
This has been the case for some URLs in previous versions but was reported as bugs?

Does anyone know the reason WHY ESRI makes AGS URLs case sensitive?
I don't see the logic on making URLs case sensitive as it just adds another layer of unnecessary complexity and adds confusion..

Example of URL that will not work:
http://services.arcgisonline.com/arcgis/rest/services/world_street_map/mapserver (incorrect case on service name)

More information:
http://blogs.esri.com/esri/arcgis/2013/12/17/esri-basemap-system-upgrade-to-arcgis-10-2-for-server/


texts: Calculating difference two rasters through gdal?

texts:
I'm writing a module for Qgis I'm creating my own bitmap calculator, I need to compute the pixel values of the two rasters, I want to calculate with the help of Gdal.I create two functions Where I read all the information (how many channels, projection, coordinates, etc.) . On the output there are two variables as can now find the difference of the pixel values?

from osgeo import gdal
def inputFile():
    filepath = r"C:\Users\user2\Desktop\1\rastr calc\2012.tif"

    # Open the file:
    raster = gdal.Open(filepath)

    # Check type of the variable 'raster'
    type(raster)
    # Projection
    raster.GetProjection()

    # Dimensions
    raster.RasterXSize
    raster.RasterYSize

    # Number of bands
    raster.RasterCount

    # Metadata for the raster dataset
    raster.GetMetadata()
    # Read the raster band as separate variable
    band = raster.GetRasterBand(1)

    # Check type of the variable 'band'
    type(band)

    # Data type of the values
    gdal.GetDataTypeName(band.DataType)
    # Compute statistics if needed
    if band.GetMinimum() is None or band.GetMaximum()is None:
        band.ComputeStatistics(0)
        print("Statistics computed.")

    # Fetch metadata for the band
    band.GetMetadata()

    # Print only selected metadata:
    print ("[ NO DATA VALUE ] = ", band.GetNoDataValue()) # none
    print ("[ MIN ] = ", band.GetMinimum())
    print ("[ MAX ] = ", band.GetMaximum())

def inputFile2():
    filepath = r"C:\Users\user2\Desktop\1\rastr calc\2011.tif"

    # Open the file:
    raster = gdal.Open(filepath)

    # Check type of the variable 'raster'
    type(raster)
    # Projection
    raster.GetProjection()

    # Dimensions
    raster.RasterXSize
    raster.RasterYSize

    # Number of bands
    raster.RasterCount

    # Metadata for the raster dataset
    raster.GetMetadata()
    # Read the raster band as separate variable
    band = raster.GetRasterBand(1)

    # Check type of the variable 'band'
    type(band)

    # Data type of the values
    gdal.GetDataTypeName(band.DataType)
    # Compute statistics if needed
    if band.GetMinimum() is None or band.GetMaximum()is None:
        band.ComputeStatistics(0)
        print("Statistics computed.")

    # Fetch metadata for the band
    band.GetMetadata()

    # Print only selected metadata:
    print ("[ NO DATA VALUE ] = ", band.GetNoDataValue()) # none
    print ("[ MIN ] = ", band.GetMinimum())
    print ("[ MAX ] = ", band.GetMaximum())

How can I now take deduct the values of the 2012.tif pixel from the pixel values of 2011.tif Example of calculations as in Qgis calculator rasters

In Qgis the raster calculator, it takes pixels from the coordinates, not from the beginning of the images My task is to make the same algorithm as in the Qgis raster calculator

Example Qgis

2012.tif-2011.tif=outputfile.tif Only from gdal


texts:
I need to convert this code to be used by python in ArcGIS

if [Field_1] < [Field_2] :
pd = [Field_1]
else:
pd = [Field_2]


texts: Is there a way to use some sort of "instancing" when storing polygons in a spatial database?

texts:
I am quite new to GIS database addons, and I am still exploring the possible uses for them. In general the idea of storing objects indexed by their position which enables additional operations is pretty cool.

However, it's fit for data gathered from real world scenarios, where every park, building etc. will have its own set of points. The simplest way is to simply store them as points and execute regular queries on them.

However, imagine we're putting some imaginary objects on the map, belonging to a few categories, but in large quantities. Let's take a star shape for example. Now, I can simply add the set of polygons every time I want to place my star in my "world", but this creates a huge duplication of data and requires me to modify each and every record shall I want to update the shape somehow.

The obvious solution to that is to store only the position of said shape, and a key to some other table holding shapes with 0-based coordinates. I can perform things like intersection tests by using ST_Translate on the shape and its position; so far so good.

Will this approach scale to a bigger number of instances? Will I be able to perform every query I would be able to using regular polygon storage?

And lastly, but not leastly, doesn't PostGIS provide a built-in support for that? I was unable to find that in the documentation.


texts: Does the PostGIS ST_AsMVTGeom function perform geometry simplification on the fly?

texts: PostGIS now comes with a set of functions that allow to generate Mapbox Vector Tiles. Among them, the ST_AsMVTGeom allows to transform a geometry into the coordinate space of a tile (typically 4096x4096). Given a geometry containing consecutive points that falls onto the same coordinate after transformation, does the function perform simplication on the fly, or should the ST_Simplify function be called beforehand?


texts: Are there any MSSQL Licencing restrictions regarding ArcGIS Server Services

texts:
I use ArcGIS Server in conjuction with MSSQL Server 2012.

According to my Database Administrator there could be licencing restrictions on the side of MSSQL when publishing Data over the web with ArcGIS Server when the source datasets are stored within MSSQL. It's all to do with mssql nodes.

Would anyone know if there are any such restrictions regarding data use from MSSQL Servers?


texts: EVI and SAVI indices for Landsat 7 and Landsat 8

texts:
I'm trying to run the indices Savi and EVI at landsat 7 & 8. About SAVI I read this from ESRI and explains very good this index. About SAVI I have a question, I will choose buy my self the L factor if it is 1, 0.5 or 0?
About EVI index I will use the some value for L factor like SAVI? and what about coefitiens C1 & C2.


texts: How to read user defined project variables using Python

texts:
I can change user defined project variables by QgsExpressionContextUtils.setProjectVariable('myvar','hello world')

But how can I access/read this variable in a python script?


texts: How can I apply symbology from two layers to one feature class with common features using ArcMap 10?

texts:
Say I have two layer files. Layer file 1 has symbology for feature types A, B & C. Layer file 2 has symbology for feature types B, C & D. In these two layer files, feature types B & C are the same and have the same symbology. If I merge these two layer files, I get Feature Class 3 with feature types A, B, C & D.

How can I get all 4 feature types in Feature Class 3 with the "correct" symbology?

A: symbolised from Layer 1
B: symbolised from Layer 1 OR Layer 2
C: symbolised from Layer 1 OR Layer 2
D: symbolised from Layer 2

I hope that's clear.

In the real world I have about a dozen layer files each with several dozen feature types, with an unknown level of overlap, so I'm looking for the quickest, least manual way of doing this.

I know that using a Style, then symbolising by Matching Symbols to a Style is the "correct" way of doing this. However I don't have access to the style used to create these layers, and I really don't want to have to manually save every single symbol from every single layer to a new style. However, if it's done once, then it's done forever, so it's not out of the question. If anyone knows of an easy way of generating a Style file from one or more layers, that would be much appreciated as well.

Something in the back of my mind tells me it might also be possible to do this using Cartographic Representations. I don't know how to go about it, and I suspect it would be massive overkill, but if it's the best way, I'll take it.

I'm pretty sure there's no easy way to do this, but I just want to confirm I'm not missing something. Any suggestions appreciated! Many thanks!


texts: How to make invalid geometry valid in Tibero Spatial?

texts:
Keeping invalid geometries in databases is rather considered as a bad practice. Yes they can be stored (inserted) but further process, including Intersection (which is one of strongest reason of using spatial database), gives us only error.
In Tibero Spatial there's a function to check validity of a geometry called ST_ISVALID().
But how to do the following in Tibero DBMS?

checking reason of invalidity (equivalent to PostGIS ST_ISVALIDREASON())
making an invalid geometry valid (equivalent to PostGIS ST_MAKEVALID())


texts: Is there a way to start the route on pgrouting from a coordinate that is not a node?

texts:
After execute all steps from this pgRouting tutorial, now I have a map that creates and displays routes, but as the routing methods rely on parameters representing vertex Ids that was created using osm2pgrouting, my route always starts either after or before the place I selected on the map, and ends on the same way.. after or before the place it should end.

We all know that there are tools that are based on pgrouting, that displays the route ending and starting exactly where the user clicked, not in the nearest node.

I was wondering if there is a way to create more nodes between the current ones, so I would reduce the distance between nodes to few meters, what would fix it, or if there is another solution to my issue..

Thank you in advance!


texts: OpenLayers: constraint a map extent and fit a div

texts:
I need to display a map with a world extension -18,-90,180,90. I want the map on the initial view to fit in the div. I tried to define the div proportionally to the map extent, but it does not work: the map displays vector before lon=-180° and after lon=180°.

I tried to define a maxExtent and restrictedExtent in the map, but it does not work (see image). How can I force the map to fit to restrictedExtent AND to the div (or what is the correct div dimension)?

Source code:

<!DOCTYPE html>

     <title></title>
     <link rel="stylesheet" href="http://myip/geoserver/openlayers/theme/default/style.css" type="text/css">


texts:
Are there any tools or methods in ArcMap that can measure the area of each level of the kernel density output?

Note this is just a kernel density plot (not a shape file or polygon). Kernel Density Analysis (dark green (10%) to red (90%) (1)

I have tried a few things but i need it to be accurate.


texts: CartoDB/Javascript - Layer Selection Buttons help

texts:
I'm trying to make my own layer selector (due to the CartoDB one not appearing on mobile devices) and I'm having issues, mainly due to my lack of experience with Javascript.

The issue seems to be most of the tutorials use SQL to pull out data from a single layer, but I can't figure out how to to just pull out layers from a visualization.

I'm basing my code off of http://jsbin.com/sihayino/4/edit , because it seems to work besides that. It produces buttons at the very least.

There are two issues that confuse me,

A) What does this block of code do?

var subLayerOptions = {
sql: "SELECT * FROM estest_all_3"
};

B) And how do I get this to turn into displaying layers in the visualization, rather than a

sublayers[0].setSQL("SELECT * FROM estest_all_3");

These are the only things in the file I can't figure out how to replace with my information.

I can provide more code if needed.


texts: How find a given point is with in given distance to a route get from google map direction service

texts:
I have created a route with start point, end point and the optional stops points.
I have to check if the given point is away with in a given distance from the route.
I have to match the given points others than the saved points in my DB.
I have saved route's start, end and stops (optional) in the db.

The route is created from google map direction service.


texts: Set raster metadata in gdal and view in GIS software

texts:
Using gdal, where (and how) do I store raster metadata so that it can be viewed in ArcGIS (or similar)?

Gdal datasources and raster band objects have the methods SetMetadata and GetMetadata. These seem to work fine to store the data, and I can get/set the data using Python. However I can't seem to find or view the metadata entries in ArcGIS. An example of what I tried to do is below:

import gdal
import numpy as np

# Some data
fname = 'metadata_test.tif'
data = np.ones([100,100], dtype=int)
lats = np.arange(0, 10, 10.0/100)
lons = np.arange(0, 10, 10.0/100)

# Create data source
driver = gdal.GetDriverByName('GTiff')
ds = driver.Create(fname, len(lons), len(lats), 1, gdal.GDT_Byte)

# Set metadata in datasource
ds.SetMetadata({'description':'test description'})

# Set sample data
out_band = ds.GetRasterBand(1)
out_band.WriteArray(data)

del ds

If you load the test data, you will see that the metadata has indeed been set:

ds = gdal.Open(fname)
ds.GetMetadata()
Out[8]: {'description': 'test description'}

However, I can't find it in ArcGIS.

I also tried the revese process, i.e. creating the metadata entries manually in ArcGIS (by editing the data source item description) and loading up the raster in Python with gdal. However, I can't find the metadata when I try the ds.GetMetadata() command, and in fact am struggling to find it anywhere in the object at all. I thought that the metadata might be stored in the band, however band.GetMetadata yields nothing either.


texts: Are there any js libraries for reading OpenStreetMap History data?

texts:
I want to do time based mapping on openstreetmap data in the Kathmandu, Nepal area. I think i have two options

Using OpenStreetMap Change files (OSC)
Using OpenStreetmap History files (OSH)

I have found js libraries for reading OSC files. Are there any js libraries for reading OSH data?


texts:
I am new to GIS and trying to create a map in my UWP app. I have a working MapView and used SketchEditor to enable drawing polygons, points, polylines etc. However, can someone please explain why I get a wkt string like

PROJCS["WGS_1984_Web_Mercator_Auxiliary_Sphere",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Mercator_Auxiliary_Sphere"],PARAMETER["False_Easting",0.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",0.0],PARAMETER["Standard_Parallel_1",0.0],PARAMETER["Auxiliary_Sphere_Type",0.0],UNIT["Meter",1.0]]

instead of something like

POLYGON((571178 6337246,571178 6402217,598061 6402217,598061 6337246,571178 6337246))

My code below:

var coloradoCorners = new Esri.ArcGISRuntime.Geometry.PointCollection(SpatialReferences.Wgs84);
coloradoCorners.Add(-109.048, 40.998);
coloradoCorners.Add(-102.047, 40.998);
coloradoCorners.Add(-102.037, 36.989);
coloradoCorners.Add(-109.048, 36.998);

// create a new polygon from the point collection
var coloradoPoly = new Esri.ArcGISRuntime.Geometry.Polygon(coloradoCorners);
string wktString = coloradoPoly.SpatialReference.WkText();

Gives the string:

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

The map I have coded is:

 private void Initialize()
        {
            var myNewMap = new Map(
                            BasemapType.Imagery, // use the Imagery basemap
                            42.175131,           // latitude
                            9.192313,            // longitude
                            8                    // level of detail (lower number = smaller scale)
                            );

            // Create graphics overlay to display sketch geometry
            _sketchOverlay = new GraphicsOverlay();
            MyMapView.GraphicsOverlays.Add(_sketchOverlay);

            // Assign the map to the MapView
            MyMapView.Map = myNewMap;

            // Fill the combo box with choices for the sketch modes (shapes)
            SketchModeComboBox.ItemsSource = System.Enum.GetValues(typeof(SketchCreationMode));
            SketchModeComboBox.SelectedIndex = 0;

            // Set the sketch editor as the page's data context
            DataContext = MyMapView.SketchEditor;
        }

and to draw a polygon:

SketchCreationMode creationMode = (SketchCreationMode)SketchModeComboBox.SelectedItem;
                Esri.ArcGISRuntime.Geometry.Geometry geometry = await MyMapView.SketchEditor.StartAsync(creationMode, true);

                // Create and add a graphic from the geometry the user drew
                Graphic graphic = CreateGraphic(geometry);
                _sketchOverlay.Graphics.Add(graphic);


texts: Creating polygons as % of original area using ArcGIS Desktop?

texts:
I have a polygon shapefile and need to create new polygons, within the originals, as % of original area.

I used buffer by field (Field= new area as %) but the buffers are created in circles around the polygons.

What am I doing wrong?

I am using ArcGIS Desktop 10.


texts: Saving Maplex settings in style manager?

texts:
I am using ArcGIS 10.2 for Desktop and I can see in the style manager that there is a part for Maplex.

Does anyone know how I can save Maplex settings for different features?

I tried using layer files but as the data I am working with is completely different it would be better to just click on a Maplex style.


texts: Export table with monthly total precipitation data for each year with GsMaP in GEE

texts:
I'm trying to download a CSV table with total precipitation for a series of years expressed in monthly total. I would like this information for a specified polygon (geometry)
¿Can anyone please help me by telling me what am i doing wrong? Thanks!
var iniDate = ee.Date('2014-09-01')
var endDate = ee.Date(ee.List(GSMAP.get('date_range')).get(1))
print('lastDay available', ee.Date(ee.List(GSMAP.get('date_range')).get(1)))

var filterGSMAP = GSMAP.filterDate(iniDate, endDate).select(['hourlyPrecipRateGC']);

var difdate = endDate.advance(-1, 'month').difference(iniDate, 'month');
var createList = ee.List.sequence(0, difdate)
var listdates = createList.map(function(month){
return iniDate.advance(month, 'month')})
print(listdates)

var gsmapImageCollection  = ee.ImageCollection.fromImages(listdates.map(function(summarize_month){
var filterCol = filterGSMAP.filterDate(ee.Date(summarize_month), ee.Date(summarize_month).advance(1, 'month'))
return filterCol.sum().select(['hourlyPrecipRateGC'])
.copyProperties(filterCol.first()).setMulti({
Date: ee.Date(summarize_month), 'system:time_start': ee.Date(summarize_month).millis()});}))
print(gsmapImageCollection)
print(ui.Chart.image.series(gsmapImageCollection, geometry, ee.Reducer.mean(), 5000)
.setChartType('ColumnChart'));

var reduced = filterGSMAP.reduceRegions({
    collection: geometry,
    reducer: ee.Reducer.mean(),
    scale: 30});
Export.table.toDrive({
  collection: reduced,
  description: 'export_sampled_values',
  folder: 'EE Exports',
  fileFormat: 'CSV'});


texts: viewshed analysis: how to assess the visibility of at least the top 3m of a tower?

texts:
This is a follow up question to Viewshed analysis: how to assess the amount of a building that would be visible?

Unfortunately, I have no access to ArcGIS' Viewshed2 tool.

I want to perform a viewshed analysis with the aim to assess if the top 3m of a 6m-high tower is visible from a given location.

Does it make any sense to first calculate two (binary) viewsheds with an OFFSETB (i.e., target height) of 3 and 6 m respectively, and then calculating their LOGICAL union (via Map Algebra)?

My guess is that the union would result in a raster indicating those cells (i.e., location across the landscape) where at least any part of the target tower between 6m and 3m from the top would be visible.


texts: Rasterio error too many point when reproject from epsg:4326 to Behrmann

texts:
I am trying to reproject a raster image with the raster package. The source CRS is EPSG:4326 and the destination CRS is Behrmann obtained from another raster, but I get an error. I leave you the code below and the error:

import numpy as np
import rasterio
from rasterio.warp import calculate_default_transform, reproject, Resampling
import geopandas

route = r"I:\miguel_analysis\AMPHIBIANS_2020_04\RASTER\raster_50x50_2020_05_20_behrmann\Adenomera bokermanni.tif"
src_file = r'I:\miguel_analysis\0_utils\wc_2_1_5m_4326.tif'
dst_file = r'I:\miguel_analysis\0_utils\wc_2_1_5m_berh_rast.tif'

with rasterio.open(route) as src:
    dst_crs = src.crs

with rasterio.open(src_file) as src:
    transform, width, height = calculate_default_transform(
        src.crs, dst_crs, src.width, src.height, *src.bounds)
    kwargs = src.meta.copy()
    kwargs.update({
        'crs': dst_crs,
        'transform': transform,
        'width': width,
        'height': height
    })

    with rasterio.open(dst_file, 'w', **kwargs) as dst:
        for i in range(1, src.count + 1):
            reproject(
                source=rasterio.band(src, i),
                destination=rasterio.band(dst, i),
                src_transform=src.transform,
                src_crs=src.crs,
                dst_transform=transform,
                dst_crs=dst_crs,
                resampling=Resampling.nearest)

the output is

CPLE_AppDefinedError: Too many points (7569 out of 7569) failed to transform, unable to compute output bounds.

the defintion from crs is:

       GEOGCS["WGS 84",
              DATUM["WGS_1984",
                    SPHEROID["WGS 84",6378137,298.257223563,
                             AUTHORITY["EPSG","7030"]],
                    AUTHORITY["EPSG","6326"]],
              PRIMEM["Greenwich",0],
              UNIT["Degree",0.0174532925199433]],
       PROJECTION["Behrmann"],
       PARAMETER["False_Easting",0],
       PARAMETER["False_Northing",0],
       PARAMETER["Central_Meridian",0],
       UNIT["metre",1,
            AUTHORITY["EPSG","9001"]],
       AXIS["Easting",EAST],
       AXIS["Northing",NORTH]]


texts: Problems with OpenLayers and WFS - proxy.cgi issues

texts:
I am having problems getting WFS data to display on a map using OpenLayers/GeoExt/GeoServer.

I created a proxy script like on this page:

http://trac.osgeo.org/openlayers/browser/trunk/openlayers/examples/proxy.cgi

I have put this is file in C:\Program Files\GeoServer 2.0.2\webapps\geoserver.

I have changed the path at the top line of the file to read

#!C:/Python26/python.exe -u

I have also added the address of the machine that I want to allow to the file (10.0.0.155:8080) and added an environment variable for python (the machine is running Windows XP).

When I point to the proxy.cgi in a web browser it just displays the proxy.cgi file in the web browser, when from what I have read it should be showing the OpenLayers homepage.

Does anyone know how I can get this to work?

Regards,

James


texts:
We use an ancient Macon MapKit! component to display GIS maps in a Windows gui.

Does anyone know what has become of MapKit! ? All traces of MapKit! have been removed from the current GfK GeoMarketing (their new name) web site.

And more specifically, is there any way to convert the maps (extension .mkl and .mkb) from/to other GIS formats?


texts: How to List Fields of a Featureclass in a FileGDB using ogrinfo

texts:
I have an ArcGIS File Geodatabase, and I know that there is a featureclass within it.

If I type: ogrinfo myfilegdb.gdb it lists the featureclasses within it.

I now want to just find out what are the field names & types of the Featureclass within it.

I looked at the documentation, and it says that -so should be used:

  Summary Only: suppress listing of features, show only the summary
  information like projection, schema, feature count and extents.

But when I type: ogrinfo -so myfilegdb.gdb, It just shows me the listing of featureclasses & does not show me the schema.

How do I see the field names & types of the featureclass?


texts:
I am new to OpenLayers and GeoServer... can someone please help me to show a WFS layer... here is my code:

  const vectorSource = new ol.source.Vector({
    format: new ol.format.GeoJSON(),
    url: function(extent) {
      return 'http://localhost:8080/geoserver/Geoservertest/ows?service=WFS&' +
          'version=1.1.0&request=GetFeature&typeName=Geoservertest%3Apoi&' +
          'maxfeatures=50&outputformat=json&srsname=EPSG:4326&' +
          'bbox=' + extent.join(',') + ',EPSG:4326';
    },
    strategy: ol.loadingstrategy.bbox
  });

  const vector = new ol.layer.Vector({
    title: 'Sample Layer',
    source: vectorSource,
    style: new ol.style.Style({
    image: new ol.style.Circle( /** @type {olx.style.IconOptions} */ ({
        radius: 5,
        fill: new ol.style.Fill({
            color: '#ffff00'
        })
    }))
})
  });

  var map = new ol.Map({
    target: 'map',
    layers: [
      new ol.layer.Group({
            'title': 'Base maps',
            layers: [
                new ol.layer.Tile({
                    'title': 'Google Maps Satellite',
                    'type': 'base',
                    visible: true,
                    'opacity': 1.000000,
                    source: new ol.source.TileImage({ url: 'http://mt1.google.com/vt/lyrs=s&hl=pl&&x={x}&y={y}&z={z}' }),
                }),
                new ol.layer.Tile({
                    'title': 'Google Road without Building',
                    'type': 'base',
                    visible: false,
                    'opacity': 1.000000,
                    source: new ol.source.TileImage({ url: 'http://mt1.google.com/vt/lyrs=r&x={x}&y={y}&z={z}' }),
                }),
                new ol.layer.Tile({
                    'title': 'Google Satellite & Roads',
                    'type': 'base',
                    visible: false,
                    'opacity': 1.000000,
                    source: new ol.source.TileImage({ url: 'http://mt1.google.com/vt/lyrs=y&x={x}&y={y}&z={z}' }),
                }),
                new ol.layer.Tile({
                    title: 'Water color',
                    type: 'base',
                    visible: false,
                    source: new ol.source.Stamen({
                    layer: 'watercolor'

                    })
                }),
                new ol.layer.Tile({
                    title: 'OSM',
                    type: 'base',
                    visible: false,
                    source: new ol.source.OSM()
                }),
            ]
        }),

                    new ol.layer.Group({
                title: 'Layers',
                layers: [vector]}),
    ],
    overlays: [overlay],
   // target: 'map',
    view: new ol.View({
      center: ol.proj.transform([-73.985130, 40.705896], 'EPSG:4326', 'EPSG:3857'),
      zoom: 13,
      minZoom: 2,
      maxZoom: 25
    }),

  });


texts:
Prior to ArcGIS 10.0, I worked on several attribute tables at once and I often minimized them.

Now, the attribute tables have no icon to minimize them using ArcGIS 10.1 for Desktop on our Windows 7 computer.

We have them on PC and in Citrix as well, and the software behaviour is the same there.

Is there a way to get attribute tables to minimize?


texts: How to start using the SharpMap library?

texts:
I am a very newbie on GIS and also programming especially C#. I am trying to develop a program to:

display map from GeoTIFF, DTED (level 2 I guess), shapefiles and JPEG type of files.
display movement a vehicle with regards to velocity input
display latitude, longitude, elevation, pitch and roll
vehicle icon representing correct geographic location
calculate actual distance travelled

I have been told to use SharpMap library and I am doing it on Microsoft Visual Studio 2010.

Any tips how to start using the SharpMap library and to program? Any pro on SharpMap that might help me?


texts: Simplifying Shapefile without risking Spatial Accuracy?

texts:
I am building a feature for a web application. Basically would use the ABS (Australian Business Statistics) data to render the regions on the Google Maps.

I am able to implement the basic working prototype however the ABS dataset (65 MB Shapefile) is huge and contemplated to use the tools

Map Shaper
QGIS

for simplification.
The rendering time is poor without simplification as the data set is huge.

Map Shaper was not really helpful because the shapefile (65 MB) didn't get imported into the database after simplification because of error 333 object not found (something like this).

QGIS helped me to simplify the shapefile from 65 MB to 4 MB, 2 MB, 1 MB depending upon the tolerance level (0.0001, 0.0020,0.0080 etc). I was able to successfully import this into the database and used it to render on maps. The rendering time has been significantly improved, however the spatial accuracy is going for toss because of the level of simplification. The functionality here is to store regions into the database for a user and then notify the users when a delivery matches.

We are using Grails, Postgres (using PostGIS), GeoJSON here to implement it.

Could anyone here assist me in simplifying the shapefile (Statistical Area Level 2) of size 65 MB to at least 1-2 MB without risking the spatial accuracy when we zoom in to the lowest level?


texts: Symbolising multiple lines in the same location

texts:
I would like to produce maps showing multiple lines that are in the same location.

This poses the problem that the lines overlap each other making one of multiple lines visible when displayed in a map.

Is there a method of making all the lines visible but keeping their location.


texts: Transformation types in geo-referencing of QGIS?

texts:
I am a beginner in geomatics.

I am trying to georeference some images in QGIS.

I would like to understand what is happening behind the different transformation types (linear, polynomial, Helmert...).
What I don't understand is the red lines that appear after inserting a GCP. These lines take as a starting point the point the user put, but I don't understand their orientation.

Does this show a gap between the coordinates given by the user and the coordinates of the map?

Example :
http://hpics.li/dbe04f9


texts: Distance between points and polygons in QGIS / PostGIS

texts:
I need to work out the distance between a table of polygons and a table of points.

The table of polygons has different site boundaries, and the point table has points relating to each of the site boundaries.

I need to work out the distance between each point and the site it relates to.

What is the best way to do this? Is this best done with a PostGIS query or is there a way to do it in QGIS (I notice there is distance matrix tool in QGIS but it doesn't support polygons). Ideally I would like a way to do this in QGIS as its something my users do quite regularly.


texts: Converting format and CRS of many shapefiles to kml (or csv, json..)

texts:
I tried to use QGIS to transform shp files in kml and to translate the CRS to use in a google maps application.
My problem is that I have a lot of shp files that i can open in a bulk manner by drag and drop but i can only export one layer at a time.

I have used the MMQGIS plugin to merge layers and export one that containes all. This works fine, except that the only way i can keep track of the shp files is by their file name (the objname is truncated to 10 chars), and when i merge the layes the order gets somehow scrambled and i cannot make sense of the what is which.

So, bottom line, how can i convert a lot of shp files (and also to change the CRS) in a more friendly format (kml, csv, json, etc)?


texts: Is there a rule of thumb regarding the number of services which should run on a single instance of ArcGIS Server?

texts:
I have a single installation of ArcGIS Server serving all the GIS data for my organization, including internal and open data. It is installed on a Windows Server with 32GB RAM. We are up to well over 100 services. I currently have problems with image services failing, although map services seem to be stable. I cannot say for sure that by decreasing the amount of published services would increase performance or make the current services more stable as the network plays a role too.
How or when would I need to be thinking about running multiple instances of ArcGIS Server?


texts:
We have SDE (10.3) running on SQL Server (SQL Server 2012 on Windows Server 2012 R2)

Over time we've noticed our GIS performance degrade, and can have good days and bad days.  We run versioned editing, replication, have geometric and network datasets, and have archiving enabled.  We run regular (several times a day) compresses to keep things running smoothly.

My question is whether having archiving enabled on our heavily edited datasets would have a noticible effect on performance?  We are not in a position to test this, so wonder if anyone else has experience with testing the performance of archiving heavily edited data.  Some of our feature classes could have anywhere from 100 to 10,000 edited records in a day.

We've had our geodatabase with all our different types of datasets and functions running for a number of years, and it's only since we migrated from arcsde 9.3 to 10.3 (about a year ago) that we've noticed the significant performance issues.


texts:
How do I show only some buildings?

I only need to show 3 cabins in the forrest, but not every house in the whole village around, and the city close by. Maybe some of the symbol buildings in the city area to, but thats less important.
Actually it would be best if the built area could just be gray (with the roads layer on top)(in other words less detailed outside the forrest).

The buildings are made up of 2 layers, lines and areas (color inside the building lines).
And one of those cabins need to be a different color than the two other, or it can be a symbol (I have an icon for it if needed, and its even preferable (a (DNT) tourist hut)).

Im using QGIS Dofour 2.0.1 for this but also have QGIS Lisboa 1.8.0 (wich i used sligtly more), but I'm more used to ArcGIS (and the much older ArcView (3.x)(nope not the new one you might think of :-)).


texts:
I am getting the following OGR error when using Grass r.mapcalc in QGIS 2.18 on Windows 7:

Data source is invalid ()

The process seems to run, the window shuts but there is no output.

This problem has been diagnosed in the following link as caused by non-ASCII characters in path used:
  Qgis gdalogr:contour is not generating output.

My user also uses this Firstname.Surname format but the directory I am trying to save to is not in this path.

Is my username the issue?

Both rasters tested have the same number of pixels, extent and spatial resolution. Here are my parameters:

Similarly, SAGA raster calculator, produces no output and the following python error:

AttributeError: 'module' object has no attribute 'TYPE_RASTER'


texts: Product License not initialized error when trying to open an SDE connection in ArcObjects

texts:
I am trying to connect to create connection in using ArcObjects using the code below.

IPropertySet propertySet = new PropertySet();
IWorkspaceFactory sdeWsFactory = new SdeWorkspaceFactory();

propertySet.SetProperty("SERVER", server);
propertySet.SetProperty("INSTANCE", instance);
propertySet.SetProperty("DATABASE", database);
propertySet.SetProperty("USER", user);
propertySet.SetProperty("PASSWORD", password);
propertySet.SetProperty("VERSION", version);

IWorkspace WS = sdeWsFactory.Open(propertySet, 0);
IFeatureWorkspace fws = WS as IFeatureWorkspace;

When it tries to execute the Open, I get the error message: The Product License has not been initialized. How do I tell which "product" it is talking about?

We are runnuing ArcSDE 9.3.1 SP2.


texts: CAD Equivalent of Maplex and VBScript for Labeling

texts:
Is there a CAD equivalent of maplex functionality and VBScript usage to produce custom labels?

I'm pulling a shapefile into CAD, but not sure what support there is to intelligently label the information behind the features.


texts: Error when extracting specific field of shapefile in Google Earth Engine?

texts:
I'm trying to extract specific field data which is in shapefile. There was an error while executing the following code in GEE as argument invalid type. the screenshot shows the error message.

Can you help to solve the problem?

var startYear = 2018;
var endYear = 2019;
print(countries)
for(var i=0;i<countries.size().getInfo();i++){
  var startMonth =countries.filter(ee.Filter.eq('startMonth',i));
  for(var j=0;j<countries.size().getInfo();j++){
    var startDay = countries.filter(ee.Filter.eq('startDay',j));
    var start = ee.Date.fromYMD(startYear, startMonth,startDay);

  }
}
for(var i=0;i<countries.size().getInfo();i++){
  var endMonth =countries.filter(ee.Filter.eq('endMonth',i));
  for(var j=0;j<countries.size().getInfo();j++){
    var endDay = countries.filter(ee.Filter.eq('endDay',j));
    var end  = ee.Date.fromYMD(endYear, endMonth, endDay);
    //print(typeof endMonth)
  }
}
print(typeof startDay)
print(start)
print(end)

Here startMonth,startDay and endMonth,endDate are the fields in attribute table(in shapefile) I want to extract these field data based on the country.special thing is there is another column country as 'ADM0_NAME'.


texts: Conversion of NAM Lambert Conformal coordinate system to NAD 1983 UTM Zone 11N

texts:
I am attempting to convert the NAM wind data coordinate system using data from the thredds server:

Catalog = TDSCatalog('http://thredds.ucar.edu/thredds/catalog/grib/NCEP/GFS/'
                        'Global_0p25deg/catalog.xml?dataset=grib/NCEP/GFS/Global_0p25deg/Best')

using Siphon. I have extracted the lat/long values and have the pyproj string:

Proj('+proj=lcc +lat_1=12.190, +lon_0 = -133.459 +datum=NAD83 +no_defs')

but I am not sure if this is correct. Using the right string, how might I convert to the NAD 1983 UTM Zone 11N coordinate system and is an intermediate projection conversion needed?


texts:
I am using QGIS v 3.2. The data set I am using involves pipelines. I need to mark an exclusion zone which is running 50m either side of the pipeline for the entire length of the pipeline.

Is there a function which will allow me to do this or do I have to plot it manually?

The reason I ask is that the pipeline is not straight, portraying various bends and curves, which I would find it hard to create a 50m exclusion zone accurately by hand.

Would a possible way, be by duplicating the shapefile and somehow offsetting its coordinates by 50m?


texts: Add NDVI Date to stacked image in Google Earth Engine

texts:
I have created an Image collection in GEE, and I have computed NDVI index and mapped in each Imagen from Collection. Later, I create another collection with NDVI as a band, an finally Statacked in a Imagen each NDVI band.
BUT the question is: how add the NDVI´s DATES in each band of the Stacked Image?
function addNDVI(image) {
  var ndvi = image.normalizedDifference(['B5', 'B4']).rename('ndvi')
  return image.addBands([ndvi])
}
var collection = ee.ImageCollection ("COPERNICUS/S2_SR")
  .filterBounds(geometry)
  .filterMetadata ("CLOUDY_PIXEL_PERCENTAGE", "Less_Than", 30)
  .map(addNDVI)

print(collection)

var Coleccion_ndvi =  collection.select(["ndvi"])
print(Coleccion_ndvi)

var stackCollection = function(collection) {
  var first = ee.Image(collection.first()).select([]);
  var appendBands = function(image, previous) {
    return ee.Image(previous).addBands(image);
  };
  return ee.Image(collection.iterate(appendBands, first));
};
var stacked = stackCollection(Coleccion_ndvi);
print('stacked image', stacked);


texts: QGIS - Distance between cluster analysis not working with K means and taking to long with hierarchical method. Any advice?

texts: QGIS- Trying to carry out a cluster analysis regarding the density of villages in Kazakhstan.  I am placing latitude and longitude as my fields and not sure if this is correct (see image). Hierarchical method is taking too long and K mean is giving wrong results. Any advice on what I can do?


texts:
Someone knows about a tool that test the cpu usage, memory, IO Read/Write with MapProxy or a WMS, WMTS, TMS Server?

Like a 'Web Stress tool' that test Server with concurrent clients request and navigation through the site.


texts:
I need the shapefile of the map attached (or very very similar) - but I can't find. I attached it here.
Somewhere I found a picture where it is named "Terrestrial Major Habitat Types".


texts: Batch Processing 400 MODIS HDF4 Files Using Several Tools In ArcGIS 10

texts:
This is my first post, I hope I don't screw it up.  We have about 400 MODIS HDF4 files we are looking to batch process in the ArcGIS 10 or Python environment.  The tools we are looking to use are Extract Subdataset (HDF4 files have subdatasets and we are only interested in one of them, which its index number is 1), ProjectRaster(to UTM 11 NAD 83), convert HDF4 files to TIFF (although I believe this can be done when executing the Extract Subdataset tool as it lets you choose an output file format), Int (we need the raster file to have an attribute table with pixel values).  I've tried performing these tasks in modelbuilder in batch mode but i've had no success in properly setting it up.

We are interested in both the ModelBuilder and Python approach.

It would be great if someone can post some batch processing code that would apply to these tools.


texts: Transform Box from SR-ORG:6630 to WGS-84 in Python

texts:
I have an image file that I'm opening with gdal and using GetGeoTransform() to determine the minimum and maximum lats and lons:

Lon_Min=-2493045.0
Lon_Max=2342655.0
Lat_Max=3310005.0
Lat_Min=177285.0

I want to transform it to WGS-84.  I believe I've found how to do so using pyproj:

from pyproj import Proj, transform

# project from SR-ORG6630 to ESPG4326 (WSG84)
# copied input projection from Proj4 format in http://spatialreference.org/ref/sr-org/6630/proj4/
inProj = Proj('+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23 +lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs')
outProj = Proj(init='epsg:4326')

lon_min_coord, lat_min_coord = transform(inProj, outProj, Lon_Min, Lat_Min)
lon_max_coord, lat_max_coord  = transform(inProj, outProj, Lon_Max, Lat_Max)

This works fine without error; however, I noticed that depending on which corner of the box I use, I get a different result:

(Pdb) transform(inProj, outProj, Lon_Max, Lat_Min)
(-73.59459648893585, 22.07673063064898)
(Pdb) transform(inProj, outProj, Lon_Max, Lat_Max)
(-63.67219185074133, 49.177063191811406)
(Pdb) transform(inProj, outProj, Lon_Min, Lat_Min)
(-119.78610533599247, 21.74230777798903)
(Pdb) transform(inProj, outProj, Lon_Min, Lat_Max)
(-130.2328280158085, 48.70739591346951)

For example, in the first two commands the only thing I changed was using Lat_Min vs using Lat_Max.  However, the output longitude was -73 then -63.  How do I know which coordinates to use in my transformation?


texts: extract amenity from pbf file using osmosis

texts:
I want to extract some amenity(atm) from a pbf file.
Here is what i've found.

  osmosis --read-pbf myfile.pbf --tf accept-nodes amenity=atm --write-xml atm.osm

Myfile.pbf is a file for Paris area (France). It's 211Mo.
The atm.osm file is 1.8Go I think it's really big for only ATM of Paris.
So is it the good way to extract atm?
I also tried with rarely amenity and it's always a 1.8Go file.
Is there a way to read atm.osm to know if the data is good?
If I'm not on the right way I want to extract ATM from the data to insert it in the future in a database.


texts:
I have this raster (land-and-water.tif):

And this one (land-only.tif, black = no data):

And I need this (black = no data):

I'm trying to do it with GDAL. Here's what I tried so far:

gdaltindex clipper.shp land-only.tif
gdalwarp -cutline clipper.shp -dstalpha land-and-water.tif water-only.tif

But the clipper.shp is empty, and obviously the resulting water-only.tif is the same as land-and-water.tif . I think I need to tell gdaltindex what is the no-data value in my land-only.tif , but I didn't find how in the doc.

How can I do that?


texts: Javscript API esri/dijit/search JSON Response

texts:
Hello guys I am having some trouble trying to parse out the returned JSON data from a dijit search. I have the following code where s.set("value","48035") will be populated from another page so for the time being i manually set it to 48035.

    var s = new Search({
        enableButtonMode: true, //this enables the search widget to display as a single button
        enableLabel: false,
        enableInfoWindow: true,
        showInfoWindowOnSelect: false,
        map: map
    }, "search");

    s.set("value", "48035");

//have tried the following
    s.search().then(function (response) {
        var extra = response[0][0].x
        alert(extra)
        });
//and the following
s.search().then(function (response) {
                var extra = response.candidates[0].location
                alert(extra.x)
                });
    s.startup();

It seems like no matter what i try or do i can not figure out how to get the returned value of x. Here is what the returned JSON data looks like.

{"spatialReference":{"wkid":4326,"latestWkid":4326},"candidates":[{"address":"48035, Clinton Township
, Michigan","location":{"x":-82.9172080701955,"y":42.554260928999781},"score":100,"attributes":{"Addr_type"
:"Postal","Match_addr":"48035, Clinton Township, Michigan","StAddr":"","City":"Clinton Township"},"extent"
:{"xmin":-82.922208999999995,"ymin":42.549261000000001,"xmax":-82.912209000000004,"ymax":42.559260999999999
}}]}


texts: How to make ordinary krigging with a spherical variogram in SAGA GIS?

texts: I have a very special map to make: it's a sea bird density map. To do this, I interpolate the observation points with ordinary kriging with a variogram. As I've read, the variogram must be of "spherical" type. The default one in SAGA is linear, but one can modify the equation. My question is: which equation should I use?


texts: ArcGIS-Javascript-API use "color" property of each polygon feature in geoJSON file

texts:
I am using ArcGIS Javascript API v4.17.
My geoJSON file has "color" property defined for each feature. each feature is a polygon. I would like to color-code each polygon on the map using the color that is defined in geoJSON file. Could you help?
GeoJSON file:
{
   "type":"FeatureCollection",
   "features":[
      {
         "type":"Feature",
         "properties":{
            "color":"red"
         },
         "geometry":{
            "type":"Polygon",
            "coordinates":[
               [
                  [-77.06818685702217,38.81969344450444],
                  [-77.06761126974416,38.81970363664622],
                  [-77.06762429634536,38.82015392436274],
                  [-77.06819988724344,38.82014373205801],
                  [-77.06818685702217,38.81969344450444]
               ]
            ]
         }
      },
      {
         "type":"Feature",
         "properties":{
            "color":"yellow"
         },
         "geometry":{
            "type":"Polygon",
            "coordinates":[
               [
                  [-77.0666293687195,38.82557775397718],
                  [-77.06605373350021,38.82558793973532],
                  [-77.06606675302561,38.82603822744011],
                  [ -77.06664239186607,38.82602804151912],
                  [-77.0666293687195,38.82557775397718]
               ]
            ]
         }
      }
   ]}

Code I have so far:
const url = "localhost:8080/javascript/GeoJSON.json";

const renderer = {
        type : "simple",
        symbol : {
            type : "simple-fill",
            color : "green",
            outline : {
                color : "white",
                width : 0.7
            }
        }};

const geojsonLayer = new GeoJSONLayer({
            url : url,
            renderer : renderer,
            opacity : 0.70,
            geometryType : "polygon"
            });

this.map.add(geojsonLayer);

The question  is, how do I replace the color: "green" in the code to use the "color" property defined in my geoJSON file?


texts: I am trying to write a tutorial for simple QGIS functions. I would like to include a set of buttons, each with an explanation of what its function is. Can anyone tell me where I could find the buttons as an image, or would I need to print screen and cut up (digitally) each toolbar?


texts: Updating PostgreSQL database with PHP language

texts:
I want to update the text information contained in my database for a specific layer. This layer is called via a WFS stream.
At first I displayed my layer with an Overlay with OpenLayers, and with a php script I want to update one of the data only according to the unique identifier of the object selected by the user.
To modify the data I created a text html object and a "submit" button
The html text object adds a new variable "maj_infos" which takes the value given by the user. And it uses the variable "feature-id" created by the Overlay which is displayed in text in the table.
I put you all the parts of my script:

Informations Signalétique

             Id

             Code Insee

             Année

             Sentier

             Infos

Modifier les informations de l'objet "Signalétique"

    <form method=POST >
        <input class="commentaire" type="text" id="maj_infos" name="maj_infos" required maxlength='140'/>

        <input type="submit" value="Envoyer" />
        <input type="hidden" id="feature-id" />
        <?php
            $dbconn = pg_connect("host=localhost dbname=xxxx user=xxxx password=xxxx")
            or die('Connexion impossible : ' . pg_last_error());

            if (isset ($_POST['maj_infos']) && isset ($_GET['feature-id'])){
                $sql='UPDATE signaletique SET infos="'.$_POST['maj_infos'].'"WHERE id="'.$_GET['feature-id'].'"';

            pg_query ($sql) or die ('Erreur SQL'.$sql.pg_error());

                echo 'Mise à jour des informations:'($_POST['maj_infos']);
            }
            else {
                echo 'Les variables ne sont pas déclarées';
            }
            pg_close();
        ?>
    </form>

I also put a capture of the rendering.

I want the added information to update in the info column of my database. And since it is a WFS feed, my table contained on the page will also update.
I don't know php but unfortunately I have to use it. I don't understand where my error is in the php script. When I launch "submit" nothing happens no error message appears. Whereas until now there were.


texts:
Using the grass GUI I get this error message:

ERROR: G_calloc: unable to allocate 2004389888 * 8 bytes of memory at raster.c:71

What can i do to remedy this?


texts: How to add DockWidget above the Layer panel using PyQGIS

texts:
I'm trying to dock an attribute table on the left hand side of my window, above the Layers panel.

By enabling "Open new attribute tables as docked windows" (under Settings Options > Data Sources > Feature attributes and table and using code I adopted from Displaying QGIS attribute table on the right, I can do something like:

from qgis.PyQt.QtCore import Qt
from qgis.PyQt.QtWidgets import QApplication

# Display the attribute table (by default it appears at the bottom)
iface.showAttributeTable(some_layer)

# Find, and move attribute table to left hand dock widget area
attrTables = [d for d in QApplication.instance().allWidgets() if d.objectName() == u'AttributeTable']
iface.addDockWidget(Qt.LeftDockWidgetArea, attrTables[0])

However, the attribute table appears underneath the layers panel.  How do I get it to appear above it?


texts:
I would like to know how to export a GeoTiff image with n bands in Python.

I have a function but it only works with one band raster.

def save_raster ( output_name, dataset, raster_data, driver ,NaN_Value,arr_projection=None):
    if arr_projection is None:
        arr_projection=[]
    else:
        if str(type(arr_projection)) == "<class 'osgeo.gdal.Dataset'>":
            srs=arr_projection.GetProjectionRef ()
        else:
            srs = arr_projection
    """
    A function to save a 1-band raster using GDAL to the file indicated
    by ``output_name``. It requires a GDAL-accesible dataset to collect
    the projection and geotransform.
    """
    # Open the reference dataset
    g = ( dataset )
    # Get the Geotransform vector

    if raster_data is False:
        raster_data=g.ReadAsArray()
    if type(raster_data) == tuple:
        raster_data = np.array(raster_data[0])
    if str(type(g)) == "<class 'osgeo.gdal.Dataset'>":
        geo_transform = g.GetGeoTransform ()
        x_size = g.RasterXSize # Raster xsize
        y_size = g.RasterYSize # Raster ysize
        srs = g.GetProjectionRef () # Projection
    elif str(type(g)) == "<class 'affine.Affine'>":
        geo_transform = (g[2],g[0],g[1],g[5],g[3],g[4])
        RastArr = raster_data
        x_size = int(RastArr.shape[1])
        y_size = int(RastArr.shape[0])
        #PROCESS RASTERIO NUMPY
    else:
        geo_transform = (g[1][2],g[1][0],g[1][1],g[1][5],g[1][3],g[1][4])
        RastArr = np.array(g[0])
        x_size = int(RastArr.shape[2])
        y_size = int(RastArr.shape[1])
    if raster_data.ndim > 2:
        raster_data=raster_data[0]
    NaN_rast = NaN_Value
    # raster_data[raster_data == NaN_rast] = 'NaN'
    raster_data[raster_data == NaN_rast] = np.NaN
    # Need a driver object. By default, we use GeoTIFF
    driver = gdal.GetDriverByName ( driver )
    dataset_out = driver.Create ( output_name, x_size, y_size, 1, \
            gdal.GDT_Float32 )
    dataset_out.SetGeoTransform ( geo_transform )
    dataset_out.SetProjection ( srs )
    dataset_out.GetRasterBand ( 1 ).WriteArray ( \
            raster_data.astype(np.float32) )

Any advice?


texts: How to cut polygon or split polygon using polyline and ArcObjects?

texts:
I have ArcMap10.1 and ArcObject.net.

i have created tool and on mouse down event i write some code I tried it using topologicalOperator4.Cut2(newPolyline) but i get error "A polygon cut operation could not classify all parts of the polygon as left or right of the cutting line."
i also tried with topologicalOperator4.Cut but this method requires 3 parameter

polyline
LeftGeometry
RightGeometry

I have polyline but where should i get this remaining two param values?

this is my polygon which i have to divide into two part where user will create polyline


texts: Counting number of points in exact same location using ArcGIS Desktop?

texts:
I have a points shapefile where in any given location there may be one point, two points, or three points on the exact same spot.

Any suggestions on how to count those up so I get one point at each location, but a count of how many points were there before?


texts: how to put the data on the map from database

texts:
I wanted to display a polygon with the coordinated on my map and each time i click the button.

This is my database:

I get the data from the database with php

        $sql = "SELECT * FROM waypoint";

        $server = "localhost";
        $username = "root";
        $password = "";
        $db = "scandrone";
        $conn = mysqli_connect($server, $username, $password, $db);

        $result = mysqli_query($conn, $sql);

        $json_array = array();

        while ($row = mysqli_fetch_assoc($result)) {
            $json_array[] = $row;
        }

        echo json_encode($json_array);

?>

I am made an load function, i created an variable to stock my JSON data and i tested to see if it works by doing this "myJsArray[0].longitude, myJsArray[0].latitude" and it works but i have add it one by one but instead of that i want to add everything at once.

var myJsArray = <?php echo json_encode($json_array); ?>;

map.on('load', function() {
var calcButton = document.getElementById('anim');
calcButton.onclick = function() {

map.loadImage('./images/index.jgp',
    function(error, image) {
    if (error) throw error;
    map.addImage('marker', image);
        }
    );

map.addSource('france', {
'type': 'geojson',
'data': {
'type': 'Feature',
'geometry': {
'type': 'Polygon',
'coordinates': [
                [
                    [myJsArray[0].longitude, myJsArray[0].latitude],
                    [myJsArray[2].longitude, myJsArray[1].latitude],
                    [myJsArray[1].longitude, myJsArray[1].latitude],
                    [myJsArray[0].longitude, myJsArray[0].latitude]
                ]
            ]
        }
    }
});
map.addLayer({
        'id': 'france',
        'type': 'line',
        'source': 'france',
          "layout": {
            "line-cap": "round",
            "line-join": "round",
          },
          "paint": {
            "line-color": "#D20C0C",
            "line-dasharray": [0.2, 2],
            "line-width": 2
          }
    });
};

});


texts: GeoServer - Horizontal Banding Between tiles

texts:
Anyone know what is causing this banding between tiles:

This is a PostGIS layer, imported from SQL Server into PostgreSQL using ogrfdw as outlined here: Import SQL Server Spatial Data into PostGIS


texts: Script to upload several .xyz files into QGIS

texts:
The following script has the purpose to upload all .xyz files into QGIS that are within the folder xyz_data. That's the script:

from qgis.core import QgsVectorLayer, QgsProject
import glob, os

dir = "C:/Users/user/Desktop/xyz_data"

os.chdir(dir)

for file in glob.glob("*.xyz"):
    uri = "file://" + dir + file + "?delimiter=%s&xField=%s&yField=%s" % (" ", "field_1", "field_2")
    name = file.replace('.xyz', '')
    lyr = QgsVectorLayer(uri, name, "delimitedtext")
    QgsProject().instance().addMapLayer(lyr)

When I run the command in the python console the names of the files appear in the QGIS layerlist but there is no data and I get an attention sign next to the layer. There is no attribute table I could open nor can I zoom to the layer

Hovering over the attention sign gives me the following error message:

  No available Layer! Layer data source was not found.

Apparently it has something to do with the uri. However, I don't see any errors in the code.
The directory is correct, I went to the folder and copy & pasted the url.

Perhaps there is something wrong with writing + file +?

Perhaps I don't get the whole filename as a string?

The delimiter is correct as well.
The xField and yField names are correct as well I think. The file itself already starts with the x-coordinate and doesn't have a line called xField, however when I upload each file manually, QGIS represents the columns as "field_1" and so on.

Does anyone see the error?


texts:
I am using rio clip like this:

rio clip C:\\Users\\rit\\Documents\\Projects\\GG\\Data\\geonode.tif C:\\Users\\rit\\Documents\\Projects\\GG\\Data\\test.tif --like C:\\Users\\rit\\Documents\\Projects\\GG\\Data\\tmpHungary.tif'

The tmpHungary.tif file has a no data value of 255. However in the output test.tif file, even cells where tmpHungary.tif has nodata value are extracted.

How can I fix this?


texts: Joining nearest neighbor small polygons using R

texts:
The figure shows a plot of numerous polygons (zoom in) made from a raster file by

require(sp)
require(raster)
pol100 <- rasterToPolygons(r, fun=function(x){x==1}, dissolve=FALSE)

There polygon has 2059 features

> pol100
class       : SpatialPolygonsDataFrame
nfeatures   : 2059
extent      : 254358.6, 268808.1, 2619318, 2642087  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=46 +datum=WGS84 +units=m +no_defs
nvariables  : 1
names       : value
min values  :   1
max values  :   1

Each of the small polygons can be explored/plotted individually

#if want to plot a particular polygon
plot(pol100[[1]][[500]])
#if want to know the area of a particular polygon (say polygon#100)
#(all polygons are in fact of same size)
pol100@polygons[[100]]@area

I want to know if there is any automated way to merge the nearest neighbor polygons? Looking at the figure tells that the process will create approximately 15-20 large polygons, for instance. Thanks in advance.


texts: Can't get SRS from PostGIS using Proj4 text

texts:
I'm having a problem transforming a point from the Australian MGA UTM datum (eg MGA Zone 51 http://spatialreference.org/ref/epsg/28351/) to WGS84 using GeoDjango. I've narrowed that down to a problem where PostGIS can not select a Spatial Reference System using Proj4 text.

One of the main limitations is that I'm dealing with many geospatial locations, and I don't actually know which CRS I'm working with for each one; I only know that they're using the MGA UTM datum, they're located in Australia, and they have a northing, easting and zone which I am given.

I figure I should be able to use the Proj4 text string in order to get the correct CRS like this, because the Proj4 string is consistent for these CRSs except for the zone:

from django.contrib.gis.gdal import SpatialReference
proj = "+proj=utm +zone=%d +south +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs" % zone
srs = SpatialReference(proj)

Is this a reasonable approach? If not, that would at least answer the question.

However, if so, it's not working as I think it should, because srs.srid gives nothing.

>>> srs.srid
>>>

Doing it the other way works fine:

>>> srs = SpatialReference(28351)
>>> srs.srid
28351
>>> srs.proj4
u'+proj=utm +zone=51 +south +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs '

Digging deeper, I went into the database to check things out.

database=> SELECT proj4text FROM spatial_ref_sys WHERE srid = 928351;
                                    proj4text
----------------------------------------------------------------------------------
 +proj=utm +zone=51 +south +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
(1 row)

That works fine. But going back the other way doesn't:

database=> SELECT * FROM spatial_ref_sys WHERE proj4text = '+proj=utm +zone=51 +south +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs';
 srid | auth_name | auth_srid | srtext | proj4text
------+-----------+-----------+--------+-----------
(0 rows)

I've tried updating the spatial reference table with the insert statements from http://spatialreference.org/.

Any ideas about this?


texts: OpenLayers does not render WMS layer generated by MapServer

texts:
I need to add a WMS layer to my OpenLayers map and it should be generated by MapServer installed locally.

OpenLayers HTML file

<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/html4/loose.dtd">

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

Map File

MAP
    NAME parcels_map
    STATUS ON
    SIZE 420 300
    EXTENT 489885 5450946 490904 5451637
    UNITS METERS
    SHAPEPATH "shapefiles"
    IMAGECOLOR 255 255 255
    WEB
        TEMPLATE "template.html"
        IMAGEPATH "C:\\xampp\\htdocs\\output\\"
        IMAGEURL "/output/"
        METADATA
            "wms_title" "WMS Demo Server"
            "wms_onlineresource" "http://localhost:8080/cgi-bin/mapserv.exe?map=C:\\xampp\\htdocs\\mapfile.map"
            "wms_srs" "EPSG:4326"
            "wms_feature_info_mime_type"  "text/html"
        END
    END
    LAYER
        NAME parcel_boundaries
        DATA parcels
        STATUS DEFAULT
        TYPE POLYGON
        CLASS
            NAME "parcels_color"
            COLOR 255 215 0
        END
    END
END

Obviously MapServer can generate the map using the given shapefile. I have tested it by entering following URL to the browser.

http://localhost/cgi-bin/mapserv.exe?map=C:\\xampp\\htdocs\\mapfile.map

But, when adding this URL as WMS service URL, it just gives  broken image tiles.

Any help would be very much appreciated.


texts: Field Calculator - use a string field to determine how a calculation is applied to a numeric field

texts:
I have a String field ('Prefix') and a Double field ('Result'). I am trying to group data by a specific range, <0.002, 0.002 to 0.07, >0.07, however some values have a '<' or '>' prefix contained in the String field. This causes me issues where i have the value 0.07 in one cell and the value 0.07 in another cell with the prefix '<', which represent 0.07 and <0.07 respectively.

What I would like to do is create a new Double attribute called 'Result_Edit'(easy enough) and use Field Calculator to look for all cases where the 'Prefix' field (String) has a '<' or '>' and the 'Result' is 0.002 or 0.07, and either return 0.0019 or 0.071 to the corresponding 'Result' field (Double). The values that have no prefix will remain as they are.

Prefix | Result | Result_Edit

..>.......0.07.......0.071

...........0.07.......0.07

..<.......0.07.......0.07

..<.......0.002......0.0019

...........0.002......0.002

..>.......0.002......0.002

I only care about <0.002 and >0.07.

This is what I have, but it doesn't work.....

def EDIT ( Prefix , Result , Result_Edit ):
  if (Prefix = '<' & Result = 0.002):
    return 0.0019
  elif (Prefix = '>' & Result =0.07):
    return  0.071
  else:
    return Result

EDIT ( !Prefix! , !Result! , !Result_Edit! )

Using ArcGIS Desktop v 10.5 with Python checked.


texts: QGIS Updating composer label text via python - problem with Item Properties panel

texts:
I am writing a plugin that replaces the text of a composer label based on its itemID.

When I change the text of the label using:

item.setText("Hello World")

the label changes on the composer but the text in the Item Properties panel for the label doesn't, it stays as what it was originally.

If I call item.displayText() before and after I change the label text this returns the correct information.

How do I get the text to change in the Item Properties panel? Using qgis 1.8 and 1.9


texts: Configure MapServer 7 for both time and elevation

texts:
I am trying to use MapServer to serve up multi-dimensional raster data.
All the data is stored in separate TIFF images - 1 image stores (for example) air temperature at a single time and elevation.

Using the current MapServer docs, I have a working service that uses a time dimension.  When I try to add a second dimension for elevation, no error is thrown, but the GetCapabilities request only includes the time dimension. Any other dimension I add is ignored.

My tile index postgres request:

NAME "time_idx"
TYPE POLYGON
DATA "geom from ms.grid_metadata using unique model using srid=4326"

and the layer metadata:

  wms_dimensionlist "time, elevation"
  wms_timeitem "time"
  wms_elevation_item "elevation"

  wms_timeextent "2018-10-20/2018-10-24/PT12H"
  wms_timedefault "2018-10-23T00:00:00Z"

  wms_elevation_extent "2, 10"
  wms_elevation_unit "meter"
  wms_elevation_default "2"

DEBUG level 5 does not print any error or warning output.

I currently am using the same LAYER METADATA for both the tile index and raster layers.

I tried creating a layer that includes only the elevation dimension.
In that case, no dimensions were advertised.

It seems I have some stupid map file problem, but I'm not seeing it.

Full disclosure: I did build MapServer from source, but this has not caused any problems so far and I have no reason to think this is my problem.

output of mapserv -v:

MapServer version 7.2.1 OUTPUT=PNG OUTPUT=JPEG SUPPORTS=PROJ SUPPORTS=AGG SUPPORTS=FREETYPE SUPPORTS=ICONV SUPPORTS=WMS_SERVER SUPPORTS=WFS_SERVER SUPPORTS=WCS_SERVER SUPPORTS=FASTCGI SUPPORTS=GEOS INPUT=JPEG INPUT=POSTGIS INPUT=OGR INPUT=GDAL INPUT=SHAPEFILE

Running on RHEL 6.


texts: Arcpy module to read csv files and create polyline shapefiles

texts:
I am creating a module in arcpy that will loop through a folder with csv files in it. I have both of them written, and they work fine on my computer. However when I send them to anyone else they won't run. I do not have any hard coded paths in them and cannot figure out the problem. The error says 'IOError: [Errno 2] no such file or directory:u'COTS_234AR.csv'' in regards to the line

#Open File
vehicle = open(input, "r")

Am I overlooking something?

Code that calls module:

#Import Modules
import arcpy
from arcpy import env
import sys
import csv
import Script1

#environment settings
arcpy.env.overwriteOutput = True

#set workspace parameter: ONLY USER DEFINED INPUT-
#SHOULD BE IN FORMAT SIMILAR TO "M:\GIS\FinalExam"
filename = arcpy.env.workspace = sys.argv[1]

#create list of CSV files in folder
csvs = arcpy.ListFiles("*.csv")

#Create for loop of CSV files
try:
   for csv in csvs:
     #set parameters
        input = csv
        desc = arcpy.Describe(csv)
        output = (desc.basename + ".shp")
        #call in module
        Script1.scriptOne(input, filename, output)

except:
   print "There has been an error processing the scriptOne module."#
#retrieve all messages
print(arcpy.GetMessages())

This is the actual module:

def scriptOne(input, filename, output):
    #import modules
    import csv
    import arcpy
    from arcpy import env

    #environment settings
    env.overwriteOutput = True

    #Open File
    vehicle = open(input, "r")

    #create csvReader object to process the header
    csvReader = csv.reader(vehicle)
    header = csvReader.next()
    latIndex = header.index("LAT")
    lonIndex = header.index("LON")
    out = arcpy.CreateFeatureclass_management(filename +"\\", output, "POLYLINE", "", "", "", 4326)
    spatialRef =arcpy.Describe(out).spatialReference

    # Write the array to the feature class as a polyline feature
    vertexArray = arcpy.Array()

    # Loop through the lines in the file and get each coordinate
    for row in csvReader:
        # Get the lat/lon values of the current GPS reading
        lat = row[latIndex]
        lon = row[lonIndex]
        # Make a point from the coordinate and add it to the array
        vertex = arcpy.Point(lon,lat)
        vertexArray.add(vertex)

        #create insert cursor to make polyline shapefile
    with arcpy.da.InsertCursor(output, ("SHAPE@",)) as cursor:
        polyline = arcpy.Polyline(vertexArray, spatialRef)
        cursor.insertRow((polyline,))


texts: I was wondering if anyone had a good tutorial for Open Trip Planner. I am not an experienced programmer and I am having a lot of trouble understanding what the documentation wants me to do. I got through basic usage, but after that I became very confused.  For example, I think the next relevant section is configuration. The first paragraph is about build directories. What exactly is a build directory? I've looked all over the internet for, but I didn't see any tutorials, and only a handful of stack overflow/exchange questions.


texts: Exporting same attribute values to separate layers and merging them in QGIS

texts:
I am trying to split watersheds of the same attribute table into separate layers.

The Two two repeated values, for eg. SW1 and SW1 both are splits from a same watershed called 'SW1'. I am trying to achieve a single layer with both SW1 and SW2 without getting dissolved. So that I can run calculations for both of the them for my purpose.
By visually I am trying to achieve like this
Original visualization

The repeated labels of SW1 has came together as a single layer

How to do this automatically without selecting each of them?


texts:
I have the following workflow.

generate vrt from imagelist file(with paths to tifs/grids)
check if coordinate transformation is required, if so convert to target SRS
Create tile using gdal_retile based on the warped vrt

Step 1 works fine using the following command:

gdalbuildvrt -vrtnodata 255 -hidenodata -overwrite -resolution user -tr $cell_resolution -te $bBox -input_file_list $imagefile $outvrtfile

In Step 2, I first assign the source SRS to the generated VRT using gdal_translate to generate an intermediate VRT, then I warp this intermediate VRT to the target SRS using the follwoing command:

gdal_translate -co "TILED=YES" -co "BLOCKXSIZE=512" -co "BLOCKYSIZE=512" -a_srs $sourceSRS -a_nodata none -of vrt -projwin $bBox ${source} ${outvrtfile}

gdalwarp -wm 500 -wo "SKIP_NOSOURCE=YES" -wo "NUM_THREADS=ALL_CPUS" -wo "OPTIMIZE_SIZE=TRUE" -of vrt -t_srs $targetSRS ${outvrtfile_warped} ${transformedfile}

The warp operation is meant (from my point of view) to do the actual transformation from source to target SRS.

Question is am, I doing the right thing? I am asking this because sometimes I get weird and incorrectly transformed outputs.


texts: Raster to Point for certain raster values

texts:
I am trying to export points form a raster using the arcpy's RasterToPoint_conversion, but I don't want to export all the cells, only the ones with a value of 1.

This is the code:

import arcpy

rast_tif = ".......\Local\Loc_rast.tif"
out_raster = arcpy.Raster(rast_tif)
arcpy.RasterToPoint_conversion(out_raster, '.......\Local\Max_points.shp', "VALUE" == 1)

But returns the error

  arcgisscripting.ExecuteError: Failed to execute. Parameters are not
  valid. ERROR 001000: Field: Field 0 does not exist


texts: What is the difference between Bing Maps, Here WeGo and Google Earth imagery in this area?

texts:
We are working on a project about the exploration of historical ruins via publicly available satellite imagery. This project is for a high school scientific research competition.

In the photos below you can see that Bing Maps imagery is different from Here WeGo and Google Earth. Burried ruins of a settlement can be clearly seen in Bing Maps image.

In all of these map services the source of the satellite imagery is shown as DigitalGlobe. We want to know what the difference is between the images apart from "season", if there is any.

coordinates for the example site seen in the photos: 39,050433, 44,015129


texts: How can I share with other users a preconfigured algorithm?

texts:
I´m creating a process in QGIS, and I make it work automatic. But, I´ve used a preconfigured algorithm in order to avoid selecting shapes every time you use the process.

So, is there any way I can share with other users, that preconfigured algorithm?

Here´s a image from the "algorithm", because i think the translation of his name is not correct.


texts: Import GPS data (text documents files) to QGis

texts: I want to import GPS data (text document files) in to Qgis, what  file format are you suppose to use ?


texts:
I've already tried How can I add map from website maps-for-free.com as my base map on QGIS 2.10.1-Pisa? with the adjusted filepath for the metadata.ini file but it didn't work. The FreeMap option doesn't show up in the dropdown menu.

Am I missing something or is this service not available anymore for QGS 3.4?


texts: Get length of a line stored in SQL Server Geometry column with srid 4326

texts:
I have some spatial data stored in SQL database geometry column called sp_geography (The field is defined as Geometry). The SQL below:

f.sp_geography.STLength()

returns, for example 0.0895991008150951 when I was expecting 892.9m and 0.00175496077867079 when I was expecting 132.6m

I believe it is returning decimal degrees but I cannot workout how to convert to metres. How can I retun the result in metres?


texts: What is the significance of 591657550.5 as it relates to google maps scaling factor?

texts: According to this gis stackexchange answer and this stackoverflow question the scale ratio for google maps at any given zoom level is 591657550.500000 / 2^(level-1). That's not the circumference of the earth in meters, miles or kilometers. What is that number? Just curious how the calculations came about.


texts: Gtiff2Tiles - How to generate tiles without processing tons of blank tiles

texts:
I am using GTiff2Tiles to create tiles from my TIF. (I am testing on level 11).

However, I am unable to complete the test because it is generating loads of blank tiles, so many that after 20 minutes it is still generating them. See the image below.

My command I am using is :

gtiff2tiles -i "input.tif" -o "C:\outputdir" -t "C:\tempdir" --minz 8 --maxz 11 --a join --tms true --threads 4

What can I do get rid of all the blank tiles so it generates faster? I am not an expert with this, and I am trying to migrate from gdal2tiles as GTiff2Tiles. Here is what the tiles look like in the output directory - tens of thousands of them.

I have used the gdal tag as Gtiff2Tiles is new and does not have a tag yet, but it is essentially a clone of gdal2tiles.


texts: Adding geospatial filters to GeoServer GetFeature or GetMap request?

texts:
Translation from Google Translate:

I am developing a tool to obtain information from the cadastre. I want
to get the number of properties inside a "multypolygon". Can I do it
with parameters by GET URL?

For example:

http://ovc.catastro.meh.es/INSPIRE/wfsBU.aspx?service=WFS&request=GetFeature&typenames=BU.BUILDING&bbox=742438,4046840,742613,4046970&srsname=EPSG::25829&

Can I pass geospatial filters by URL?

Is WMS or WFS better?

Is it better to install GeoServer on my Apache?


texts: How to separate out bands from hyperspectral data using R

texts:
I'm having one hyperspctral data with 426 bands, with dim 250*250*426. I want to seperate each bands.
I tried like this

img <- stack("name.tif")
rr <- array(, dim= c(250,250,426))
for(i in 1:426) {rr[,,i] <- raster(img, layer = i )}

but it's showing error as

number of items to replace is not a multiple of replacement length

I think the problem is, for dim(raster(img, layer = i )), I'm getting as 250*250*1 . How can I make it just as 250*250 ??

Edit : Or, Is it possible to take my hyperspectral image as a D x N matrix where D is the dimension of the data (no of HS bands = 426) and N is the number of HS pixels(250*250)


texts: Opening geodatabase files using Browse Dialog Command (VB.NET/C#) of ArcObjects?

texts:
I am currently using the following method with ArcGIS 10 to open raster files from a form:

Public Function BrowseRaster(ByVal sTitle As String, ByRef sFolder As String, ByRef sName As String) As Boolean

 Dim pGxDialog As IGxDialog = New GxDialog
    Dim pRasterFilter As IGxObjectFilter = New GxFilterRasterDatasets()

    Dim pFilterCol As IGxObjectFilterCollection
    Dim pEnumGx As IEnumGxObject
    Dim pGxObject As IGxObject

    sFolder = ""
    sName = ""

    pFilterCol = pGxDialog
    pFilterCol.AddFilter(pRasterFilter, True)
    'pGxDialog.StartingLocation = ""
    pGxDialog.RememberLocation = True
    pGxDialog.AllowMultiSelect = False
    pGxDialog.Title = sTitle

        If pGxDialog.DoModalOpen(0, pEnumGx) Then
            pGxObject = pEnumGx.Next
            Dim sFile As New FileInfo(pGxObject.FullName)
            'sName = pGxObject.BaseName
            sName = pGxObject.Name
            sFolder = sFile.Directory.FullName
            Return True
        End If

End Function

How can I modify this function to also open files from inside of a file geodatabase?


texts:
I'm working on a Leaflet map, currently adding content to my popups. In .bindPopup ([" "]) is there any way to format a paragraph of writing on the backend into multiple lines of code within those quotes?

In other words, the following poem is one long quote—is there a way to format it into multiple lines so it's easier to organize headers, paragraphs, etc.)?

L.marker([]).addTo(mymap)
        .bindPopup("Three Rings for the Elven Lords under the Sky, Seven for the Dwarf Lords in their halls of Stone, Nine for Mortal Men doomed to die, One for the Dark Lord on his dark throne, In the Land of Mordor where the Shadows lie, One Ring to rule them all, One Ring to find them, One Ring to bring them all and in the darkness bind them")


texts: Styling lines depending on several attributes

texts:
I have a line shape file which I made from DWG file. I have several fields which I want to use for styling. I have colour, width and type. I need to specify colour, width and type independantly.
For example colour=1, width=2 and type="dashed line" or colour=3, width=1 and type="continious line". Colours are not range or random colours, they should be specified manually 1 is red, 2 - is yellow etc. The same with line types. Width can be either range, or specified manually as well as colour and width.


texts: ol3 error with style function -- feature shows half geometry

texts:
I used stylelayerfunction with unique colors for all the different 700 points. While applying it some features are rendered correctly and some features are half render.

While using single style it rendered perfect.

in stylelayerfunction i just set the different color.

The feature geometry in both the cases are same.

I can not understand why it is happen.


texts: Adding free routing services in ArcGIS Web AppBuilder (developer)?

texts: To add free routing services like openrouteservice or ORSM in WAB(developer) can we use the directions widget or do we have to develop a new widget for it?


texts: Trouble with openlayers v4.6.5 uploading geoJSON

texts:
I have been following a few books but they all seem to be about OpenLayers 3.  Things have apparently changed since then.  I am having a lot of trouble viewing a geojson over a basemap layer.  The bing basemap seems to work just fine.... I have gone into about:configure and turned the security to false on firefox.  Someone said that could be interfering with loading something from disk into the web... Is there an easier way?? Can i just upload these into some database like SQL server and have them inserted by command?? I am going crazy.  I can usually follow along with books or websites but everything is outdated.  This is my code:

<!DOCTYPE html>

<meta charset="UTF-8">
<title>First Example</title>
<link rel="stylesheet" href="ol.css">
<link rel="stylesheet" href="features.css">

<meta name="viewport" content="width=device-width">

<select id="layer-select">
<option value="AerialWithLabels">Aerial with labels</option>
<option value="Aerial">Aerial</option>
<option value="RoadOnDemand">Road</option>

</select>

I have also experimented with editing the vector area like this

var vectorLayer = new ol.layer.Vector({
source: new ol.source.Vector({
format: new ol.format.GeoJSON,
url: 'C:\\Users\\------zips.geojson'
})

});

Should I just switch to leaflet.. Documentation seems up to date and I had no problem with that library when I tried doing this...


texts: Is Geoserver mandatory for developing GIS web application in Java using Ol4JSF in NetBeans IDE?

texts:
I am developing a GIS web application which involves querying and creating map layers dynamically (like page getting refreshed for every 5 minutes) and also static querying etc... This application interacts with PostgreSQL for non-spatial data rendering.. I am using NetBeans to develop the application. Is it mandatory to use Geoserver in the project? FYI, I am planning to use OL4JSF with Jsp to render the maps in browser where user interactions will be carried out.

If Geoserver is mandatory, let me know why?


texts: Can you display mxd files on a app for iPads?

texts:
I would like to know if there is any apps that can displaying ArcGIS data?

I have multiple .mxd's that i would like to display on a portable app.


texts: shapely wrong output -- when snapping points to line shapefile

texts:
I have a series of nodes to be connected to a road network, so as to calculate the cost of transportation. I used shapely to snap the nodes to the road network polylines. Most of the connections worked well, however, there are several connecting to very far road segments, instead of a near road. I don't understand what could have gone wrong here.

here is the essential part of the codes:

# points are the nodes
# roadLine is the road network in linestring
ptsPrj = points.copy()   # find the points project to locations on the road network
ptsPrj['geometry'] = ptsPrj.apply(lambda row: roadLine.interpolate(roadLine.project(row.geometry)), axis = 1)

# join the points with its projecting points on the road network

joinPts = points.copy()
for i in range(len(points)):
    joinPts.at[i, 'geometry'] = LineString([points.at[i,'geometry'], ptsPrj.at[i, 'geometry']])


texts: Seeking ArcPy dissolve to multipart workaround for when arcpy.Dissolve fails

texts:
I have an ArcPy script in which I am using arcpy.Dissolve on a single feature class with multiple features to get a single multi-part feature in the output.  The script is successful only part of the time and the problem, I would guess, relates to the function built-in tiling.

Per Esri:

  The availability of physical memory may limit the amount (and
  complexity) of input features that can be processed and dissolved into
  a single output feature. This limitation could cause an error to
  occur, as the dissolve process may require more memory than is
  available. To prevent this, Dissolve may divide and process the input
  features using an adaptive tiling algorithm...

I tried putting the operation in a while loop that exits when the feature count is 1 but as the Esri page states, that rarely results in success.

How do I complete the dissolve in the midst of this problem?

My system has plenty of memory so the problem cannot be corrected by running on a box with more memory.


texts: Sending coordinates via VHF Radios to Computers

texts:
I need to send and receive coordinates in NMEA format from a VHF radio to a computer.
Any suggestions where to look in order to accomplish this task?
The transmission will occur only when the user presses the PTT button.
Thank you all.


texts: Usage of Custom Transparency band for elevation raster in QGIS 3

texts:
I would like to set a raster transparency to represent its elevation, so that the mountain tops are opaque and gradually become transparent toward their foot.
On the Layer Properties > Transparency tab, in Custom Transparency Options, I found Transparency band and would like to use it.

Problem is, I do not know how to control the transparency according to the band I have just assigned (Band1 in the above image).
At default setting it gives me an uncomfortable image (right-hand-side; below).

I also tried simple color ramp such as Reds, but the output is still weird.

Question: How can I use this Transparency band function ?
Alternatively, giving a list of transparent pixel seems to work. But I have to give a color range in stepwise manner, and it is far from gradual change I hope to achieve.

Note: I borrowed the West Timor elevation data from SAGA GIS Tutorials - Downloads.


texts: Changing the Attribute table view in QGIS 3

texts:
How does one change the Attribute table view on QGIS 3?
An initial answer says to use the Switch to Table view button but I can't seem to find this button.


texts: Create GeoDataFrame polygon geometry from coordinates list of lists

texts:
My current dataframe (clean_geoms) contains two columns (field_geom_type, field_coords) that I want to convert to a single geometry (GeoSeries)
This is a representative sample of my data:
print(clean_geoms.field_geom_type.iloc[0])
print(clean_geoms.field_coords.iloc[0])

Polygon
[[[-96.568927, 46.769008], [-96.5689454, 46.7705433], [-96.5726564, 46.7705222], [-96.572638, 46.7689868], [-96.568927, 46.769008]]]

What is the best way to unpack this list so that it can be converted into geometry?


texts:
i'm looking for a way to disable/enable modify OWS Server tab in project properties. Is it possible using just qgis python api ?

Currently i'm building web application, which provides spatial services with qgis server. And i'm looking for a way to manage this services. Is it only way to do this ? Any other ideas ?


texts: Creating Crown Maps from Catalog using lidR

texts:
I have a large catalog and I'd like to create crown maps for the entire catalog or a portion of the catalog.  I've been successful at creating crown maps from single LAZ files by editing/running the script below:
#Script to generate crown map from single file
library(lidR)
set_lidr_threads(6L)
las <- readLAS("c:/myfile.laz", select = "xyzcr", filter = "-drop_z_below 0")
las <- segment_trees(las, li2012(R = 3, speed_up = 5))
crowns <- delineate_crowns(las)
par(mar=rep(0,4))
plot(crowns)
crowns <- delineate_crowns(las, func = .stdtreemetrics)
spplot(crowns, "convhull_area")

How do I modify this code to generate crown maps from a catalog?