texts: GRASS GIS: r.in.gdal does not find raster data in .hdf file (modis data)

texts:
I'm having trouble in importing MODIS .hdf data in GRASS GIS. When I use r.in.gdal it returns me this error:

ERROR: No raster bands found in </Volumes/Maxtor/TESI_satellite_imagery/IMMAGINI_DI_PARTENZA/serie_storica_estate_dic2015-feb2020/modis/MOD02QKM.A2019044.1415.061.2019045012150.hdf>.

I've tried to uninstall and reinstall gdal on my computer (S.O. MacOS Catalina), and rebuild all the bindings with python, but I didn't manage to solve the issue. I used Homebrew and followed its instructions to do it. I've also tried to use conda previously, but it has issues while "solving the environment" and cannot complete the procedure. What am I doing wrong?
I've also tried the i.modis.import, but I get this error

ERROR: pymodis library is not installed.

I can't understand where this comes from, since if I run pip install pymodis or pip3 install pymodis on the terminal it says:

Requirement already satisfied: pymodis in ./opt/anaconda3/lib/python3.8/site-packages (2.0.9)
Requirement already satisfied: future in ./opt/anaconda3/lib/python3.8/site-packages (from pymodis) (0.18.2)
Requirement already satisfied: requests in ./opt/anaconda3/lib/python3.8/site-packages (from pymodis) (2.24.0)
Requirement already satisfied: GDAL in ./opt/anaconda3/lib/python3.8/site-packages (from pymodis) (3.1.3)
Requirement already satisfied: numpy in ./opt/anaconda3/lib/python3.8/site-packages (from pymodis) (1.18.5)
Requirement already satisfied: urllib3!=1.25.0,!=1.25.1,<1.26,>=1.21.1 in ./opt/anaconda3/lib/python3.8/site-packages (from requests->pymodis) (1.25.9)
Requirement already satisfied: certifi>=2017.4.17 in ./opt/anaconda3/lib/python3.8/site-packages (from requests->pymodis) (2020.6.20)
Requirement already satisfied: idna<3,>=2.5 in ./opt/anaconda3/lib/python3.8/site-packages (from requests->pymodis) (2.10)
Requirement already satisfied: chardet<4,>=3.0.2 in ./opt/anaconda3/lib/python3.8/site-packages (from requests->pymodis) (3.0.4)

How can I solve this problem?


texts: Clipping multiple raster with multiple shapefiles in Python

texts:
I am trying to adopt existing code for clipping multiple rasters with multiple polygons based on this link:
Subsetting multiple rasters with multiple shapefiles using Python
I added the following codes to create the shp_list, shp_names, raster_list, and raster_names.
It is supposed to clip raster by 128x128 dimension, but it is now 129x129. It seems the code changed the float N, S, E , and W to integers. In addition, the output is saved in the main directory not the results folder.
How can I modify the code to result the correct dimension (128x128) and save the clipping of each shapefile into separate folder?
def subs(shp_list, shp_names, raster_list, raster_names, outpath):
    for i, shp in enumerate(shp_list):
        with fiona.open(shp, "r") as shapefile:
            shapes = [feature["geometry"] for feature in shapefile]

            for j, ras in enumerate(raster_list):
                for scene in raster_list:
                    with rasterio.open(ras, "r") as src:
                        out_image, out_transform = mask(src, shapes, crop=True)
                        out_meta = src.meta

                        out_meta.update({"driver": "GTiff",
                                         "height": out_image.shape[1],
                                         "width": out_image.shape[2],
                                         "transform": out_transform})

                        ras_path = f"{outpath}{raster_names[j]}{shp_names[i]}{ras_extension[1:]}"

                        with rasterio.open(ras_path, "w", **out_meta) as dest:
                            dest.write(out_image)

    # number of created subsets
    subset_count = (i+1) * (j+1)
    if len(shp_list) * len(raster_list) == subset_count:
        print(f"Done. \n{subset_count} subsets created")

# calling the function
subs(shp_list, shp_names, raster_list, raster_names, result)


texts: Get QGIS layer symbology as xml string using PyQGIS

texts: Analogue to saving a layers symbology into a *.qml file I would like to have the xml content of a layers symbology in a python string (or QDomDocument). I think the method QgsVectorLayer().writeSymbology() does what I'm after, but I'm unsure how to use its parameters properly. Any suggestions?


texts: QGIS Python variable not updating after new features are inserted

texts:
I have the following script that sets the a feature id automatically after a new feature is digitized.

I added the function inside the project's Macro section inside openProject() function:

The function is being triggered ok every time I open the project. My problem now is that every new feature I create gets the same id (9707).

It seems like max varibles' value is not being updated after a new feature is added and without "Toggle editing" or "Saving layer edits".

Current code:

def openProject():
    from qgis.core import QgsMapLayerRegistry
    import os
    # Set active layer
    #layer = qgis.utils.iface.activeLayer()
    layer=None
    for lyr in QgsMapLayerRegistry.instance().mapLayers().values():
        if lyr.name() == "Pozos":
            layer = lyr
            break
    # Define function to select added feature and add attribute to field "Name"
    def update(featureAdded):
        idx = layer.fieldNameIndex('id')
        layer.changeAttributeValue(featureAdded, idx, layer.maximumValue(idx) +1)
    # Connect "featureAdded" event to "select" function
    layer.featureAdded.connect(update)

def saveProject():
    pass

def closeProject():
    pass


texts: How to know Band overview resolutions available using gdal

texts: I want to know bands resolution. Is there any way using GDAL that we can see the resolution of bands?


texts: Disabling map drag and zoom in OpenLayers 3?

texts:
I need to disable all map interactions, like map zoom, drag, scroll zoom, any kind of interaction.

Is there some easy way to do it?


texts: Value relation autocompletion and dropdown menu in QGIS?

texts:
I use a value relation widget which is connected to a list of streets.

If I tick "use completer" I can use the autocompletion but the dropdown menu disappear. The dropdown menu appears again only if I type a first letter.

I would like to have both autocompletion and dropdown menu even when no  letter is entered because sometimes the surveyor doesn't know by heart the spelling of a street.

Is it already possible or under development ?

I am using QGIS version 2.18.13

Here is a screenshot with widget properties - the function "use completer" is thicked

Here is the menu without completer

Here is the menu with completer (I have to type at least one letter to let the dropdown menu appear)


texts: Adding QgsLayoutItemTextTable to a print layout

texts:
I am writing a script that automatically generates a print layout. In addition to the map, the print layout should contain a table with results of a series of calculations that I have in the form of a python dictionary. While I have gotten into the workflow of adding layout items to the print layout, QgsLayoutItemTextTable seems to work differently.

This works fine:

project = QgsProject.instance()
layout = QgsPrintLayout(project)
layout.initializeDefaults()

title = QgsLayoutItemLabel(layout)
title.setText('MyTitle')
title.adjustSizeToText()
layout.addLayoutItem(title)

The QgsLayoutItemTextTable class, however, is for some reason not based on the QgsLayoutItem class and therefore

project = QgsProject.instance()
layout = QgsPrintLayout(project)
layout.initializeDefaults()

table = QgsLayoutItemTextTable(layout)
table.addRow(['hello', 'world', '!'])
layout.addLayoutItem(table)

results in the following error:

  TypeError: QgsLayout.addLayoutItem(): argument 1 has unexpected type 'QgsLayoutItemTextTable'

Can anyone explain to me why this is the case and what is the intended workflow?


texts: How to get all the FeatureClasses from IWorkspace?

texts:
I have created a Geodatabase in my SQL Server using ArcCatalog and then I created 5 new FeatureClasses.

Now I want to get all the FeatureClasses in this GDB.

This is my function to get all the FeatureClasses from a workspace:

        public static List<ILayer> getLayersFromWokrspace(IWorkspace workspace)
    {
        List<ILayer> layers = new List<ILayer>();
        IFeatureWorkspace fw = (IFeatureWorkspace)workspace;

        // i tried the both ways to get the datasets in the workspace (i tried the type "esriDatasetType.esriDTAny" too)
        IEnumDataset datasets = workspace.get_Datasets(esriDatasetType.esriDTFeatureClass);
        //IEnumDataset datasets = workspace.Datasets[esriDatasetType.esriDTFeatureClass];

        IDataset dataset = null;
        while ((dataset = datasets.Next()) != null)
        {
            if (dataset.Type == esriDatasetType.esriDTFeatureClass)
            {
                IFeatureClass featureClass = fw.OpenFeatureClass(dataset.Name);
                IFeatureLayer featureLayer = new FeatureLayerClass();
                featureLayer.FeatureClass = featureClass;
                ILayer layer = (ILayer)featureLayer;
                layer.Name = dataset.Name;
                layers.Add(layer);
            }
            dataset = datasets.Next();
        }

        return layers;
    }

My problem is that I can't get all the FeatureClasses (5), sometimes I get 2 and sometimes only 1 FeatureClass.


texts: Styling a raster on relief shading analysis

texts:
Possible Duplicate:
  How to create a style file for QGIS from a plain text file?

I have a text file which contains the following data

eleva     R      G       B
0   153 255 153
100 179 255 179
200 237 250 230
300 255 242 242
400 255 242 237

and I need to convert the text file into a saved color style.How best can I achieve this in QGIS


texts: How to export polyline layer from Yandex Maps

texts:
I am interested in exporting a polyline layer(s) of Moscow bike routes from Yandex maps.  The layer does not belong to me but it is open to public use.  I would prefer a KML, but just getting the coordinates in usable shape would be great.  The coordinates seem to be embedded in the source code but I don't know to reshape that easily into a list of coordinates in .csv format  which can then be mapped in QGIS.

https://maps.yandex.ru/213/moscow/?z=11&ll=37.500488%2C55.732220&l=map&origin=jsapi_2_1_34&from=api-maps&um=constructor%3AZ4zDPV2057F7o6aqu4WkQXKpLoiHO9i2


texts:
I've got several layers/feature classes in an ArcMap 10 project.

Thanks to everyone's help, I'm able to get the areas for the polygons contained within each layer.

By opening the attribute tables for each layer, and right-clicking on the SHAPE_area column, I can obtain the total area for all the polygons in that layer.

Is there a way to summarise the total areas for ALL of my polygon layers? Or do I have to do it manually for each one? Thanks!

P.S. Just to illustrate, suppose I have layers A, B, and C, all of which contain multiple polygons. Right now, I have to right-click on each layer, open its Attribute Table, right click on the "SHAPE_area" field, click on Summarize (or Statistics), see the sum of all the cells, and manually copy that into a waiting Excel spreadsheet.

Is there a "global" Summarize function so that I can get output like this (formatted to be comma delimited)?

Layer, area

A, 12

B, 13

C, 5


texts: Why is my coordinate notation different?

texts:
I was working on a project in WGS84, and I was continually adding layers to it. As I was rendering in print composer and enabled grid-coordinates, I noticed that my coordinates were different. It's supposed to be somewhere 14 north and 121 east (Manila, Philippines).

Can anyone explain how and why this happened, as well as a fix?


texts: What is the best way to migrate an SDE Database

texts:
We have a prod database, with an SDE schema.  We want to install a new SDE instance 9.3.1 and link it to a copy of the prod database.

ArcSDE 9.3.1 on Linux SLES 4

Now, what is the best way to do this.  It isn't an upgrade, so you can't really upgrade the SDE install I think?  If you simply install the SDE on top of a pre-existing SDE schema, there are all sorts of links broken, indexes, data issues, etc.

So, if I am correct, then is the only way to go, is to build an empty database (save from the users table) and then import all of the data via sdeimport (obviously using the date from a sdeexport)?

I am not normally involved with SDE, preferring to leave it with those mystical creatures, the DBA's, but I can't find one here (I believe they Apparate and Disapparate without control) so have no option but to do this myself, and find myself without the 100% self belief it is the right way to go.

Any ideas?

EDIT

We have no way of accessing the original data live, so have to export the Oracle data from one area to another, before we can access it e.g. we cannot simply access the live data 'live' as such.


texts:
To publish my map projects I simply used qgis2web or QGIS Cloud but now I need more advanced features for more maps. Now I'm wondering, could the web space I've got from 'All-Inkl.de' be useful to setup a own Map Server?

  The specifications of my webspace: PHP 5.5 at Linux dd1120
  3.13.0-96-generic #143-Ubuntu SMP x86_64 and Postgres without PostGIS

If not, what do I need to publish my maps online in a more advanced way? i.e. more display, layer and attribute options, the possibility to show the attribute table, editing attributes and downloading layers etc. QgisCloud Pro would be fine, but is still too expensive for me.


texts:
I disabled information from a vector layer with a styledefinition.
Is it possible to rasterize only the styled information in QGIS?


texts:
I have created a Leaflet map on my web page and I would like to show/hide it when I click on a button. When I first load the page, I'd like the map to be hidden.

HTML code:

    <button type="button" class="button button--default expandable__button" aria-controls="gw-example-button-dropdown1" aria-expanded="false" id="gw-example-expandable-button1">Select on the map</button>

CSS:

#gw_map1 {
display: none;
}

jQuery:

$("#gw-example-expandable-button1").click(function() {
    $("#gw_map1").toggle();
});

If I click on the button the map div appears and disappears, but the content of the map is not correctly displayed.

In fact, If I do not set the 'display' into css the map is correctly loaded into the page as in the figure:

but if I set the 'display' to 'none' and then I click on the button the map is displayed incomplete and even if I drag on the map it is always incomplete and polygons are not visible.

Any suggestions?


texts: Select by location and copy features to new shapefile

texts:
I have been asked to select point features that fall into a polygon in its own shapefile then copy the points selected to a new shapefile. There are around 300 shapefiles with an indiviidual polygon. Although only a very beginnner in python i thought it would be possible to automate the process.

I based my script of one found on this site but cannot get it too work. Im not sure if it is fit for purpose or not.

We are using arcmap 9.3 so no arcpy unfortunately.

# Import system modules
import sys, string, os, arcgisscripting

# Create the Geoprocessor object
gp = arcgisscripting.create()

# Load required toolboxes...
gp.AddToolbox("C:/Program Files/ArcGIS/ArcToolbox/Toolboxes/Data Management Tools.tbx")

# Script arguments...
Layer_A = sys.argv[1] # polygon
Layer_B = sys.argv[2] # points
outfolder = sys.argv[3]
basename = sys.argv[4]

# create search Cursor that ref the features in layer B
rows = gp.searchcursor(Layer_B)
row = rows.next()

# Create a numbering system for new output files
filenum = 0

while row:
    # Report a message back to the user about progress
    gp.addmessage("Creating Selection #" + str(filenum) + " from " + Layer_A)

    # Create a feature layer to cycle through selection features
    fidvalue = row.GetValue("FID")
    gp.MakeFeatureLayer(Layer_B,"cursor_lyr","FID =" + str(fidvalue))

    # Process: Select Layer By Location...
    gp.SelectLayerByLocation_management(Layer_A, "INTERSECT", "cursor_lyr", "", "NEW_SELECTION")

    outname = outfolder + "\\" + basename + str(filenum)
    # Process: Copy Features...
    gp.CopyFeatures_management(Layer_A, outname, "", "0", "0", "0")

    #iterate the output file number
    filenum = filenum + 1

    #go to next row
    row = rows.next


texts: How do I use Join Attribute by nearest in this case, to avoid referencing an algorithm output in field calculator expression?

texts:
I have a model that I have built to automate the process of assessing wind hazard on forested edges.
Currently I can run it successfully but:

I have to run it once, create the Test layer, and then run it a second time for Offset_Azimuth Calc to complete
Offset_Azimuth Calc expression like so: (referencing the 'Test' layer)

degrees(
    with_variable(
        'line',
        shortest_line (
            centroid (
                array_first(
                    overlay_nearest(
                        'Test',
                        $geometry
                    )
                )
            ),
            $geometry
        ),
        azimuth (
            start_point (@line),
            end_point (@line)
        )
    )
)

I have tried the following, but I think I must be misunderstanding something here - I think it should easier than I am making it.

Separate the model into two, a "prep" model that explodes the lines and assigns azimuth to each segment, and an "assessment" model that has Exploded_Boundary as an input (from prep) and referencing the "Exploded_Boundary" layer within the offset_azimuth calc expression. This will only work if the layer from prep is saved accordingly, and is less ideal for sharing with coworkers.

I have been unsuccessful in my attempt of making it a single model with 2 inputs, and 1 output (wind_hazard). I have tried using Wind_Hazard_Prep within the second Assessment model, but that did not work either - I could not reference the Wind_Hazard_Prep output from within the Offset_Azimuth Calc expression.
I see that there is a workaround by using Join Attribute by nearest, but I cannot figure out how to add and calculate the offset_azimuth field with it, as the problem here is that I cannot reference the exploded lines layer and the Intersection_OUTPUT (point layer) at the same time


texts: ArcGIS: how to obtain an average elevation for grid squares superimposed on a DEM?

texts:
Given the DEM of a certain Country, and a fishnet that is superimposed on it, I need to compute the average elevation and slope corresponding to each square in the fishnet.

[I am using ArcGIS 10.3.1, licence type: Advanced.]

These are the DEM (shades of grey) and the fishnet (red):

The DEM is not a unique raster, but a handful of smaller maps. Their resolution is 10x10 km, whereas the fishnet squares are 2x2 km.

I tried to create a Mosaic Dataset, and then I used Extract -> Clip from ArcToolbox using the DEM as input features and the fishnet as clip features. I thought this wouldhave produced a unique raster map, but I was wrong, and the result is essentially a grid with the DEM resolution (10x10 km) and a bunch of table attributes which do not refer to elevation.

The intended outcome was, instead, to get the average elevation and slope for each one of the smaller tiles in the fishnet.

What is it the best strategy to achieve this?


texts: Looping output as new input many times (50+) in ModelBuilder?

texts:
I want to run focal statistics on a raster, then run focal statistics (same configuration/setup as the first focal stat.) again on the output, and again on the new output, and again, and again.

Maybe 50 times.

Is there a simple way to do this?


texts: ArcPy Calculate Field Management with Variable

texts:
I'm trying to calculate multiple fields that all begin with "B01". This calculation consists of multiplying the fields by a field that contains a percentage ("TractPercent") and then dividing the results of that multiplication by a population field ("BufferPOP").

The problem I'm having is incorporating the variable associated with my "for loop"  with my Python3 expression.

Here is my code:

import arcpy

### set workspaces
arcpy.env.workspace = r"C:\arcGIS_Shared\Python\Demographics.gdb"
relpath = r'C:\arcGIS_Shared\Python'
p = arcpy.mp.ArcGISProject(relpath + r'\Demographics.aprx')
m = p.listMaps('Map')[0]

layers = m.listLayers("DalyCity_HalfMile_Demographics")
for layer in layers:
    print(layer)

    DemoFields = arcpy.ListFields(layer, "B01*", "ALL")
    for DemoField in DemoFields:
        arcpy.CalculateField_management(layer, DemoField, "(+DemoField+ *
        !TractPercent!) / !BufferPOP!", "PYTHON3")
        print(layer.name + ' ' + DemoField + ' Calculated')

When I run the code this way I get the following error, but can't figure out what I need to change.

ERROR 000622: Failed to execute (Calculate Field). Parameters are not valid.
ERROR 000623: Invalid value type for parameter field


texts: Add building height as obstructions to viewshed analysis in ArcGIS Desktop?

texts:
I have LAS data. I have made a DEM from the LAS points. I have a Z-enabled polygon feature representing building footprints.

Somehow I need to incorporate the height of the buildings into my DEM so that I can do an analysis to see what portion of the ground behind the building is not visible from my observer point.

I have attached an image showing the LAS data as elevation along with the line work representing the buildings and the observer point.

How do I perform a view shed analysis in ArcGIS Desktop and incorporate the proposed buildings?


texts: Table with groupings and sums in print composer?

texts:
Assuming I have a basic shapefile with the following columns: name, use, acres

I would like to generate a more visually appealing dynamic table in layout view:

Use         Name    Acres

Building
                Bld A    0.2
                Bld B    0.3
      total                0.5

Parkinhg
                North    0.3
                East      0.5
                South    0.2
                West     0.1
      total                1.1

I've looked around the print composer and found the attribute table, but I can't seem to add sum to the end of each use.  I've also explored the GroupStats plugin, but am having a hard time.

I am wondering if I need to learn Python, or something similar, to get the desired results?


texts: Address for OpenLayer plugin to Google Maps

texts:
In my office, there is a proxy that allows only a few sites to be accessed. My problem is that I need to access through OpenLayers plugin the Google Maps and I cannot.

What could the address be, so that I can ask my network administrator to insert in the exceptions so that I can view those images?


texts: Convert GIS DB PolygonField to WKT or GeoJson in Django

texts:
I am new to geoDjango. I have a polygon field defined in the model like this:
my_poly = models.PolygonField(null=True, blank=True, spatial_index=True, geography=True, srid=4326)

I am getting the record set in view like this:
rset = TestUnit.objects.extra(select={'o':'my_poly'}).filter(reduce(operator.and_, criteria)).values('o')

I see it is returning polygon data looks like this format:
0103000020E610000001000000090000004AB6BA9C127557C060915F3FC4C64440C7BDF90D137557C0187D0569C6C644408690F3FE3F7557C0E35295B6B8CA44403ECA880B407557C05F402FDCB9CA44400F46EC13407557C07D259012BBCA4440F7031E18407557C0C98E8D40BCCA444021B1DD3D407557C03D0CAD4ECECA4440BFB9BF7ADC7457C0D525E318C9D044404AB6BA9C127557C060915F3FC4C64440
I need it to return WKT or GeoJSON format so that I can draw polygons using Leaflet. What are my options?


texts:
We're looking into developing a GEO app for our city. This would be a real time app that displays what's happening in the city- and not meant for data analysis.

The map should convey only the city borders in the app, and anything else will not be available or displayed. The map is small, so the client have a lot of caching for the base layer.

I'm wondering what solution can we use that will fit our needs. Most of the GIS frameworks and APIs, like google's or mapbox are intrinsically world-wide based solutions. Is there any frame work we could use to build a more refined base layer that will be cached in the client? What tools can we use to build our own map?


texts: Image processing using Python, GDAL and Scikit-Image

texts:
I am struggling with a processing and hopefully I will be able to solve here.

I work with Remote Sensing applied to Forestry, especially working with LiDAR data. The idea is to use Scikit-image for tree top detection. Since I'm new in Python, I considered a great personal triumph to do the following:

Import a CHM (with matplotlib);
Run a gaussian filter (with scikit-image package);
Run a maxima filter (with scikit-image package);
Run the peak_local_max (with scikit-image package);
Show the CHM with the local maxima (with matplotlib);

Now my problem. When I import with matplot, the image loses its geographic coordinates. So the coordinates I have are just basic image coordinates (i.e. 250,312). What I need is to get the value of the pixel under the local maxima dot in the image (red dots in the image). Here in the forum I saw one guy asking the same thing (Getting pixel value of GDAL raster under OGR point without NumPy?), but he already had the points in a shapefile. In my case the points were computed with scikit-image (It is an array with the coordinates of each tree top). So I do not have the shapefile.

In conclusion, what I want in the end is a txt file with the coordinates of each local maxima in geographic coordinates, for example:

525412 62980123 1150
...


texts: Styling different shapes in layers - Openlayers-3

texts:
To add styling to Openlayers there seems to be different rules for different shapes.

I am trying to understand when to use what, but it is rather confusing and I am only able to do it with the help of examples.

For example, I have a few points generated on my map. I generated these points with GeoJSON data with the following format:

"type": "FeatureCollection",
    "crs": {
        "type": "name",
        "properties": {
            "name": "urn:ogc:def:crs:OGC:1.3:CRS84"
        }
    },
    "features": [{
        "type": "Feature",
        "properties": {
            "property_a": 4,
            "property_b": 4234,
            "property_c": asdf,
            "property_d": 342342,
        },
        "geometry": {
            "type": "Point",
            "coordinates": [5.53204822887795, 51.609252234043296, 0.0]
        }
    }, etc...

I am trying to create points, polygon's and lines.
What I don't understand is to why every geometry needs it's own style contructors?

I created styling for a filled point with an example:

var _myStroke = new ol.style.Stroke({
  color : 'rgba(255,0,0,1.0)',
  width : 1,
});

var _myFill = new ol.style.Fill({
   color: 'rgba(255,255,255,1)'
});

var circle = new ol.style.Circle({
    radius: 5,
    fill: _myFill,
    stroke: _myStroke
});
var puntObjectenStyle = new ol.style.Style ({
   image:circle
});

var puntObjectenSource= new ol.source.Vector ({
  url: 'puntobjecten.json',
  format: new ol.format.GeoJSON()
});
var puntObjecten = new ol.layer.Vector ({
  source: puntObjectenSource,
  style: puntObjectenStyle
});

The documentation on the Openlayers website doesnt even talk about the contructors ol.style.Circle, ol.style.Fill or ol.style.Stroke?

I have the following style-up for a polygon: (I don't even know if all these settings can be applied!)

 weight: '1.04',
 fillColor: '#da2122',
 color: '#000000',
 dashArray: '',
 opacity: '1.0',
 fillOpacity: '1.0',

How to create the correct functions to format a polygon?
And above all: how can I know?


texts: Converting a modeler 2.x to a 3.x - Done by hand

texts:
After reading here that's no easy way to convert a modeler in 2.18 to 3.10, (model to model3).
I made my model3 by hand based on the model version 2 i had. After some problems with the OUTPUT files, i finally was able to run all algoritm.
But now, all the layers that came as output (and with no errors on the Log) they come always with 0 features.
Basically i do for 5 layers
v.clean Break -> Delete duplicates -> Count Points
And from all i can see in the log, QGIS actually do the all the steps, but the output come out with ZERO features


texts:
Is it possible to filter a map's layer inside Print Composer, like it is possible to filter a layer in Qgis main window?

I have a layer with several points, and I'd like to create several layouts and in each of them have one map displaying only one of those points (filtered by a specific attribute in the attribute table).


texts: How to Edit FileGDB attribute using QGIS?

texts: I am able to attach a file GDB in QGIS (version 1.9). But I cann't edit it's content available in attribute table and 'Toggle editing mode' button was also disabled in Attribute table. Does anyone knows how to edit FileGDB attributes using attribute Table ?


texts:
I have converted a VB6 project to C#, and there is one problem I cannot find documentation on. I am saving an IPropertySet to file by using IMemoryBlobStream and IPersistStream, as seen here:

 public void WriteToFile(string fullPath)
    {
        IMemoryBlobStream pMemStream;
        IPersistStream pPersistStream;

        try
        {
            // Open Stream, must use memory stream instead of FileStream
            pMemStream = new MemoryBlobStreamClass();

            // Write PropertySet to Stream
            pPersistStream = (IPersistStream)m_pPropertySet;
            pPersistStream.Save(pMemStream, 1);

            // Save Stream to File
            pMemStream.SaveToFile(fullPath);
            pMemStream = null;

        }
        catch (Exception ex)
        {
            return;
        }
    }

I get an "HRESULT E_FAIL" with error code "-2147467259" on pPersistStream.Save(pMemStream, 1);

This doesn't throw an error if the PropertySet is empty, only if it has properties within.

When the PropertySet is populated, it is populated by Long data types (int64). Using PropertySet.SetProperty(Name, Value)

Any ideas why this would work great in ArcGIS 8/9 and not in 10 with C#?


texts: Geodatabase Connection to Oracle Active Data Guard Standby

texts:
ENVIRONMENT

Database Server:

Oracle 11.2.0.2.0 Single-Instance Database on Windows Server 2003 R2
Enterprise 64-bit Oracle Active Data Guard for Oracle Database 11g R2
ArcSDE for Oracle 11g 9.3.1 SP2 using SDELOB

Client Machine:

Oracle Client 11.2.0.1.0 on Windows 7 Enterprise
ArcGIS Desktop 9.3.1 SP2 on Windows 7 Enterprise
ArcGIS Desktop 9.3.1 SP2 QIP Patch
ArcGIS Desktop 9.3.1 SP2 Geodatabase Patch
ArcGIS Desktop 9.3.1 SP2 Geometric Network Patch

BACKGROUND

I have successfully created two standby databases from my primary Oracle database; one is a physical standby (with realtime-SQL enabled) and the other is a logical standby. When attempting to connect to either of the standby databases using direct connect in ArcCatalog, I receive the following error:

Error: Operation Failed (-1)

It then occurs to me that both standby databases are running in read-only mode, which is expected when creating standby databases via Oracle Active Data Guard. After some more research, I find that it is possible to modify a standby database such that it can run in read-write mode. After modifying the standby to run in read-write mode, I am then able to make edits to the standby. Afterward, however, I was unable to resume the application of redo logs for the standby database which essentially broke the connection between my standby and my primary. This was evident after receiving a variety of Oracle errors (e.g., ORA-01665 and ORA-01153) when following the Oracle documentation for how to switch from the standby back to the primary. In essence, I wasn't able to re-link the standby and primary databases ever again such that redo logs were applied correctly from the primary to the standby.

The good news is that I can connect to either the logical or physical standby if I perform a Switchover within OEM 12c R3 such that the primary and standby databases swap roles. Therefore, I think it has something to do with the read-only status of the database whereby Data Guard is enabled.

After some discussion with a contact at Esri, part of the solution might reside in the SDE.SERVER_CONFIG table, specifically with the READONLY parameter.  Setting the CHAR field to TRUE and the NUM field to 1 would set the geodatabase to read-only mode (Esri's method of read only).  However, I can't change this value on the standby because it's in read-only mode according to the Oracle software.  I found out that making a geodatabase connection actually performs edits at login.  There is an Oracle package that executes a function or two which performs inserts and object creations within the database.  Setting the geodatabase to read-only mode per the method above may be the way to change this behavior but I'm still waiting for answers at the moment.

Lastly, I believe there are a list of exclusions that are needed when setting up the Active Data Guard replication for logical standby databases.  I'm yet to find those, though.

Has anyone used ArcSDE in conjunction with Oracle Active Data Guard?

If so, how am I able to make a geodatabase connection to the standby?

I am trying to determine the high-level proper workflow for creating a standby database and subsequently configuring the standby to receive archived redo data from the primary database (using either Redo Apply, Real-Time Apply, or SQL Apply).


texts: How to create mask file (txt or netcdf) from shapefile using ArcGIS or QGIS?

texts:
I have a shapefile of el Congo River Basin dy area. Attributes information is just area. Actually I want create a mask with ones in the river basin and ceros outside to just used later in IDL or matlab a mask of the Congo. How can I built a proper mask file (1 inside the congo, o outside) using this shapefile in Arcgis or QGIS ?

I downloaded here: http://hydrosheds.cr.usgs.gov/datadownload.php?reqdata=30bass a shp from the rivers basin of Africa, I selected the congo river basin and I need to convert the whole Africa or just the selected area to a mask of netcdf like a basemap to use later in IDl normally is a matrix of 1 inside the river basin and ceros outside. Thansk !!!!


texts: How to correct depth values with tide data?

texts:
We are trying to derive bathymetry using sentinel images using QGis 2.18.
The training data is from a boat gps track taken in 2 days. The depth value captured needs to be corrected for tidal influences. The track has attributes on date and time, position and depth (sonar).

What is the best way to do it?

Is there any formula or free software to correct the data to obtain the mean lower low water needed for the chart?


texts: Generate subcatchment starting from a DTM

texts: I would like to know if there is a plug-in in Quantum GIS, which can help me to generate subcatchments starting from a DTM.


texts: Excluding nearest points from the SQL query result by distance in Postgis?

texts:
I have a PostGIS 2.4 database with the table that contains more than 300 000 Point features (with big density on the map on high zooms like 10-12 if you will draw them).

So I need to query that table and return the points that are separated from each other by a distance (say, 10 km) to reduce this density. I was wondering if I can do this with the ST_DWithin function or anyone else, but I can not figure out how to do it. Moreover, it needs to work fast (1-2 seconds is the max time for the query). Maybe somebody can help me?


texts:
I have some custom ArcMap addins that I have built over the last couple of years.  In one specific client environment, the addin installs correctly but when the user clicks the tool to execute it, the tool becomes disabled.  I had this issue before at that site but the user that had it has moved on and I don't recall what he did to fix it.  The users there all use remote workstations and do not have admin rights on them.  I've installed the addin in various ways and made sure that they have the appropriate permissions on the various folders where the addin content sits.  I have also made sure the config is specific to their installed version. The same tools work fine in every other environment I have access to. I have suggested running ArcMap with admin rights to see if this fixes the issue but they do not have those permissions and it is taking weeks to get approval to let them test this option.

Has anyone seen this before and know what is going on?


texts: ArcGIS .Net Automation - Error when creating a feature

texts:
I try to create a (point)-feature in an existing Layer via ArcGIS Automation.

I can create the Feature but when I call feature.Store() I get an

  HRESULT: 0x80010105 (RPC_E_SERVERFAULT)

error.

If I just edit an existing feature that I got from a SearchCursor it works flawlessly.

My application is a WinForms .Net 4.5 programm.

I get a connection to a running ArcMap Instance via IAppRot->IAppRef->IApplication->IDocument.

From the IDocument I get a reference to a specific Layer where I need to add a feature.
After that it is very basic stuff:

IFeatureClass featureClass = ((IFeatureLayer)i_Layer).FeatureClass;
IFeature newObject = featureClass.CreateFeature();
// Assign a geometry to the feature
newObject.Shape = geom;
newObject.Store();

Store allways fails with the RPC_E_SERVERFAULT Exception.

I also tried to create the feature via IFeatureCursor->IFeatureBuffer but then it fails with the same error on

featureCursor.InsertFeature(featureBuffer);

Starting an Edit-Session for this operation has no effect. It just fails later when calling

workspaceEdit.StopEditing(true);

I tried with ArcGIS 10.5 and 10.3 but the results are the same.

Other automation tasks like moving the map, setting layer-visibilities or even editing existing features work perfectly.

I read the introduction on ArcGIS automation in the documentation and learned about application's process space and "IObjectFactory".

Do I have to use that somehow to create the new Feature?

How can I do that?

Or is there something else that I am doing wrong?

EDIT:

As requested: The table definition doesn't seem to make a difference. I tried with featureclasses that have multiple fields and assigned values to the fields via

newObject.Value[index] = somevalue;

But now my tests use a featureclass that has only the fields "objectid" and "shape".
The results are allways the same.

EDIT 2:

As requested:
The geometry is created like this:

IGeometry geom = new PointClass();
((IPoint)geom).X = x;
((IPoint)geom).Y = y;

The x and y values are guaranteed to be within the featureclass bounds.

The featureclass is in an file-geodatabase. Z and M values are not activated for the featureclass.


texts: Performing classification and vectorization for batch of rasters

texts:
I have a batch of rasters (circa 8000 GeoTIFF near-infrared ortophotomaps, 12-15 MB size each). I am looking for way to perform bulk classification (supervised or not) and vectorization (to shapefile) of created classes for this batch. I know that it is possible in Erdas Imagine and ENVI, but I am not 100% sure if in open-source software (QGIS, GRASS, SAGA, ILWIS) also.
Does anyone have any experience in this matter?


texts: Combine cluster and categorical mapping on Cartodb

texts:
I came across this map that shows clusters but when you zoom in, the they split into points reflecting categorical data.

How do I replicate this on cartodb?


texts: Convert GPS coordinate to ArcGIS x,y point?

texts:
I'm trying to plot some AGSSimpleMarkerSymbol's in a AGSGraphicsLayer. The issue I'm running into is that the data I'm starting with are GPS coordinates (lat/lon). I'm using AGSPoint to plot the data, however, AGSPoint wants x and y points - not lat/lon.

My question is, how do I convert my GPS coordinates to x and y points?


texts:
I am building a custom map viewer, part of meteorological production system, using Java language and GeoTools library. I have prepared a custom map renderer, using OpenGL, that allows to display maps in different projections and need a bit of help with user camera navigation over maps in Lambert conic projection.

I would like to allow user to navigate the conic map in such way, that map north is always pointing up.
Moving mouse left and right should rotate the map around a pole, and moving up/down is always north-south adjustment:

I am struggling with understanding the mathematical transformation for screen-map coordinates for such use case, could anyone point me to the right direction?

This is a code I am trying right now; It only accounts for x axis mouse movement for now:

origin vector is north pole location
worldPos and prevWorldPos are mouse positions
all vectors are in projection coordinate system

float pdx = prevWorldPos.x - origin.x;
float pdy = prevWorldPos.y - origin.y;

float radius = (float) Math.sqrt( pdx*pdx + pdy*pdy );

float prevAngle = (float) Math.atan2(pdy, pdx);
float da = (float) Math.atan2(radius, worldPos.x - prevWorldPos.x);

float angle = prevAngle + da;

cameraAngle = angle+Angles.PI_div_2;

cameraX = origin.x + (float) (radius * Math.cos( angle ));
cameraY = origin.y + (float) (radius * Math.sin( angle ));


texts: How can I create a merged buffer in QGIS then count points that fall within it

texts:
I have a QGIS plan which shows ponds (as points in a shape file with an OS tile as a background) within a survey area.  I need to calculate the total number of ponds within 1 km of each of the on-site ponds individually.

I have created one layer with 1 km buffers for each pond as a shapefile.

As there are many ponds, this is messy to look at and I would like to create a layer which has the outline of all of the 1 km buffers combined. I have tried selecting the layer with the buffers and vector-geoprocessing tools - dissolve, but the new layer had outlines for each individual buffer.

Question part 1 - how do you create combined buffers, or convert existing buffers into one?

I will then use this outline to manually check the OS plan and add points representing ponds as a new layer "off site ponds".

Question part 2 - How would I then find the number of ponds (onsite and offsite, not counting the pond being queried) within 1 km of each on-site pond (individually).

Sub-question - if one of the buffers dissects a symbol for a pond will it count the pond, or only if the centre point falls within the buffer? (i.e. could the size of the symbol affect the result?)


texts: Is there a performance difference between SQL Server 2008 and 2008 R2?

texts:
We are currently hosting our ArcSDE geodatabase on SQL Server 2008 SP3.

I was wondering if there'd be any type of performance boost if we jumped up to 2008 R2?


texts: Converting OpenStreetMaps PBF to geojson with geocoding

texts:
I am able to convert PBF to geojson, but I lack the additional attributes like state, country, city, suburb, type of business etc which is offered by geocoding service like Nominatim.

The geojson just has a name and type:

{ "type": "Feature", "properties": { "osm_id": "519326004", "timestamp": "2014-10-29T19:46:31Z", "name": "Borgo Agnese", "type": "restaurant" }, "geometry": { "type": "Point", "coordinates": [ 16.6015512, 49.1915338 ] } }

But I need all these additional attributes:

[{"place_id":"4709238","licence":"Data © OpenStreetMap contributors, ODbL 1.0. http:\/\/www.openstreetmap.org\/copyright","osm_type":"node","osm_id":"519326004","boundingbox":["49.1914838","49.1915838","16.6015012","16.6016012"],"polygonpoints":[[16.6015512,49.1915838],[16.601586555339,49.191569155339],[16.6016012,49.1915338],[16.601586555339,49.191498444661],[16.6015512,49.1914838],[16.601515844661,49.191498444661],[16.6015012,49.1915338],[16.601515844661,49.191569155339]],"lat":"49.1915338","lon":"16.6015512","display_name":"Borgo Agnese, 43, Kopečná, Staré Brno, Brno, okres Brno-město, Jihomoravský kraj, Jihovýchod, 60200, Czech Republic","class":"amenity","type":"restaurant","importance":0.201,"icon":"http:\/\/nominatim.openstreetmap.org\/images\/mapicons\/food_restaurant.p.20.png","address":{"restaurant":"Borgo Agnese","house_number":"43","road":"Kopečná","suburb":"Staré Brno","city":"Brno","county":"okres Brno-město","state":"Jihovýchod","postcode":"60200","country":"Czech Republic","country_code":"cz"}}]

I can reverse geocode every entry I have in my geojson created from PBF file but that seems cumbersome to me.

Just to note: I have tried Pelias openstreetmap utilities to upload PBF to Elasticsearch, but the quality of information is not sufficient for me.

Is there any proper way how to do this?


texts:
I have a dojox bar chart on my app, that is populated from a QueryTask everytime the extent changes.

I am making use of the Tufte theme (the only one that has a transparent background), and overriding the bar colours/outlines when i populate the chart:

chart1.setTheme(dojox.charting.themes.Tufte);
    chart1.addSeries("Languages", [
        { y: pctArabic, fill: "#C23A00", stroke: null, tooltip: ctArabic + " in Arabic },
        { y: pctMalay, fill: "#6BAE11", stroke: null, tooltip: ctMalay + " in Malay" },
        ....
    ]);
    chart1.render();

I have also been using the Esri sample as a reference, and have successfuly added a popup to my map which shows some totals per suburb in a dojox piechart:

My codes is identical to the sample, bar the fields have been swapped out.

template = new esri.dijit.PopupTemplate({
          title: "Tweet Languages in {SA2_NAME}",
          mediaInfos: [
            {
              type: "piechart",
              value: {
                fields: [
                  "TotArabic", "TotMalay",
                 ....
                ],
                theme: "Tufte"
              }
            }
          ]
        });

My problem:

I need some control over styling the chart within the popup.

How do I force the pie segments to be specific colours, as per my bar chart example?

When I look at the API reference for Popup, mediainfos does not look like it has this level of control over the chart?

How do I override the tooltips with an alias like I did with my bar chart example?

I will also need to style the popup to have a similar dark grey background to my chart background, but I believe that won't be related to the chart, and should be achievable via styling the popup with CSS.


texts: Agricultural Land Classification from ArcGIS Online not working?

texts:
When I add the Agricultural land classification (ALC) map from ArcGIS Online, the layer in the table of contents loads, however the actual map doesn't.

I really need to ALC for the Humber Region, UK for my final year project, and I can't find any other downloads which actually work.


texts: QGIS Atlas map composition to show only points inside Polygon?

texts:
I want to create an atlas map where I have some 40 polygons, 100 points and 200 lines as different shapefile. Now I want to automate the atlas map preparation for each polygon and needs to show the points and lines which is falling inside of the polygon.

Can anyone help me?


texts: Topography: How to convert raster into vector layers?

texts:
(nb:This project is part of the #Wikipedia #wikimaps project)

1. Question:

How to create vector topographic maps from raster image ?
How to convert raster topographic image into closed vector polygons for the relevant levels as seen in the illustration below ?

Using commands (likely gdal) and a GIS topographic raster image, I want to produce topographic .shp maps with an elegant relief. This .shp should be made of closed polygons for each of the n requested level.

2. Possible direction would be:

Create vector isolines (gdal_contour -a elev -fl -500 0 500 1000 2000 4000 crop.tif contours.shp)
Create a copy of the whole map area into vector polygon,
Cut this polygon with all isolines of level n (i.e. isolines of level 0).
Keep those of these polygons with a value (altitude) superior to the isoline value.
(merge this polygons into one shp layer of level_n.)
Repeat for each level.

I haven't clues how to do that, but it sound promising! (should work) If you can provide the solution to one of the points 2, 3, 4, It will be really appreciated ! and I will be happy to then open a specific question where I will validate your helpful answer.

3. Current workflow: gdal_contours alone is not suitable. I currently crop my source raster ETOPO using gdal_translate then get levels with gdal_contour. Using the libs curl, unzip, gdal, my current worflow is as follow.

# download it:
curl -o ETOPO1.zip 'http://www.ngdc.noaa.gov/mgg/global/relief/ETOPO1/data/ice_surface/grid_registered/georeferenced_tiff/ETOPO1_Ice_g_geotiff.zip'
# unzip it:
unzip ETOPO1.zip
# crop it:
gdal_translate -projwin -005.48 051.30 10.00 041.00 ETOPO1_Ice_g_geotiff.tif crop.tif
# lines for levels:
gdal_contour -a elev -fl -500 -200 0 100 200 500 1000 2000 4000 crop.tif contours.shp

Unfortunally, gdal_contour generates vertor iso lines, which are not suitable. The main trouble is with cropped areas: in the case of a close up on France, the isolines are not closing into a circle as for Islands or continents but open lines, so there are partial parts of the polygons (as the quarter of a circle) whose' relevant fill --a 1/4 pie-- appear as a filled arc (see bug below).


texts: Mismatch in area calculation using Field Calculator vs Group Stats in QGIS 2.0

texts:
I have two polygon files (A being my area of interest and B being a categorised land cover layer). Both layers were in the same projection (projected and not geographic). Using QGIS 2.0, I used the Intersect tool (Vector > Geoprocessing tools > Intersect) to extract and pair the information from my land cover layer (layer B) to my area of interest (layer A).

I then opened the attribute table for my newly created layer and enabled editing. A new column was made (using a decimal field type) and I then used the Field Calculator to calculate the area for each of my table rows. The "preview output" at the bottom of the Field calculator showed a realistic number. However, when I clicked ok the resulting output was generally nonsensical with illogical values produced. I then toggled editing and saved the values.

Then using the Group Stats plugin, I used my intersected layer but used the plugin to calculate the area of my intersected attributes and it produced the correct answer.

Apologies for the long-winded description but can someone tell me why the area calculation using the Field Calculator failed to produce the correct answers? As I said, both layers were projected and matched the projection projection so I don't see this being a projection issue.


texts:
XMLHttpRequest cannot load http://sampleserver1.arcgisonline.com/ArcGIS/rest/info?f=json. No 'Access-Control-Allow-Origin' header is present on the requested resource. Origin 'null' is therefore not allowed access.

I'm getting this error. What might be the problem? I'm beginner to ARCGIS JSAPI,


texts: Update is not a known function on Bounds in version 3.32 of Javascript API

texts:
The typescript code below compiles, but bounds is getting a runtime error suggesting that update is not a known function.
var bounds = self.esriHelper.requireModule('esri/geometry/Extent');
bounds = bounds.update ( -60836536.56026919, -42677544.62462127, 342437.88671754766, 8169589.583117515, self.spatialReference);

This code was originally on version 3.18 of ESRI's Javascript API. I upgraded to 3.32 and added the two lines above to create a map without a basemap
self.map = new Map('map', {
  extent: bounds,
  sliderPosition: 'top-right',
  logo: false
});

I am new to javascript and working on inherited code. I also tried
import {Extent} from 'esri/geometry/Extent';
without success.


texts: Linked point geometries to linestrings in Spatialite

texts:
I have been parsing Canada's radio system frequency open data. It's a single large table per province, containing point geometries for each station callsign. Where a point-to-point link is defined, the record's linkid contains the remote station call sign. Here are a couple of very truncated and simulated records by way of explanation:

callsign    linkid  locname         geometry
==========  ======  ==============  ===========================
CGX934      XLQ644  APSLEY, ON      POINT(-78.090833 44.755556)
XLQ644      CGX934  NEPHTON, ON     POINT(-78.010833 44.622222)
CIN256      CIJ933  MISSISSAUGA     POINT(-79.593611 43.642778)
CIJ933      CIN256  NAKINA ON       POINT(-86.839167 50.170833)

Is it possible to use Spatialite to query the DB to return linestrings of the links? Something like:

linkname            geometry
==================  ======================================================
CGX934-XLQ644       LINESTRING(-78.090833 44.755556, -78.010833 44.622222)
CIN256-CIJ933       LINESTRING(-79.593611 43.642778, -86.839167 50.170833)

Currently I'm manually linking using a Perl program on the text table, which may produce pretty results, but is slow and not very reliable.


texts:
When I click on the point in mapImageLayer the popup is pointing/displaying at wrong place.
var mapDataLayer = new MapImageLayer({
  url: sArcgisURL,
  sublayers: [
    {
      id: sSelectedMapLayer,
      definitionExpression:  "MAP_UID = '"+mapId+"'",
      popupTemplate:{
        title:"Machine Details",
        content:"PIN: 1RW8245GFTX798456
 Dealer: Nortrax
Engine Hours: 365
Last Call In: JUly 15 2015"
      }
    }
  ]
});

Image :  https://i.stack.imgur.com/wajk1.png


texts: Extending lines to Points without altering original geometry using ArcMap?

texts:
Is there an efficient way to extend lines to points in ArcMap?
I have a utilities dataset, where I do not want to change the original geometry, but instead add to it in order for end vertices of lines tor each to the nearby junction points.

I would not want to do this manually thousands of times. Is there a way to extend lines to points? I have tried extending the lines to themselves then snapping all within a certain tolerance, but this alters the original geometry. I have also tried integrating the lines, but again, it alters the original geometry.

After thinking a bit more, I am going to try and create points at the end vertices, then create lines from points using selected end vertices and the junctions, then I will merge the lines with the original polylines. My only concern is ensuring the attributes are appropriately merged and how to obtain an appropriate selection.  I do not want to extend all lines, only certain ones.  The rule is hard to put logically, but some lines used to intersect a cad feature, and others just ended near them...


texts: Select a single entry using SearchCursor and create variable

texts:
I am trying to select a single entry of a floating point entry in a dbf table that contains 9 entries (sampled from a raster image), some of which is NULL using SearchCursor. I want to select one of the entries that is not NULL and assign it as a new variable for use later in my script. I can do a query that excludes NULL values (includes all values greater than 0) but I don't know how to select a single entry that is not NULL and assign as a new variable.

The field name is:

Rt_b1
<Null>
0.138817
0.136027
0.136027
<Null>
<Null>
0.136027
0.138817
0.136027

Any help would be much appreciated. I'm using Windows XP and ArcGIS 10.1

with arcpy.da.SearchCursor(Sample_TIF2, ("Rt_b1"), '"Rt_b1" > 0') as cursor:
     for row in cursor:
           print str(row[0])


texts:
I have a shapefile containing 8,000,000 records, with a file size of 5 gigabytes.  I would like to upload shapefile onto a server, but my file is too big and it gives an error.

How can I reduce the file size without affecting the records contained?

The server only accepts shp & mdb & gdb formats.


texts: Does OSGeo4W install 32 bit or 64 bit version of QGIS?

texts:
If you are installing QGIS via OSGeo installer on a windows 64 bit OS (XP, vista or 7), does it install a 32 bits or 64 bit version of QGIS?

Usually the standalone installer installs a 32 bits version, so I want to use the full potential of my new computer (win 7 64 bits, 16 Mb of Ram).


texts: How to zoom to a GeoJSON feature on map using OpenLayers 3 following a WFS request

texts:
After doing a WFS request, I have fetched all GeoJSON objects and fill in a list. My goal is to zoom specified GeoJSON in fetched all GeoJSON using OpenLayers 3.

In the below code inside the showFeatureByID function, I have reached specified GeoJSON by clicking, after this, I would like to zoom in on map using OpenLayers.

map.on('singleclick', function(evt) {

    var view = map.getView();
    var viewResolution = view.getResolution();
    var source = untiled.get('visible') ? untiled.getSource() : tiled.getSource();
    var url = source.getGetFeatureInfoUrl(
        evt.coordinate, viewResolution, view.getProjection(),
        {'INFO_FORMAT': 'application/json', 'FEATURE_COUNT': 50});

    if (url) {
        //url includes getFeatureInfo request string and use x,y values on map to find clicked object.
        // do wms getFeatureInfo request.
        loadDoc(url);
        interval1 = setInterval(function () {
            if (featureData.status == true) {
                fillDataView(featureData.data[0]);
                showFeatureByID(featureData.data[0]["ID"]);
                clearInterval(interval1);
            }
        }, 100);
    }
});

var featureData;
var allGeoJSONFeatures={};

function loadDoc(url) {
    featureData={status:false,data:[]};
    var xhttp = new XMLHttpRequest();
    xhttp.onreadystatechange = function() {
        if (this.readyState == 4 && this.status == 200) {

             var jsonObject=JSON.parse(this.response);

             for(i in jsonObject){
                 var val = jsonObject[i];
                 if(i == "features"){
                     for(j in val){
                         var geojsonx = val[j];
                         var properties = geojsonx.properties;
                         featureData.data.push(properties);
                     }
                 }
             }
            if(featureData.data.length>0){
                featureData.status=true;
            }
        }
    };
    xhttp.open("GET", url, true);
    xhttp.send();
}

function getAllGeoJSON() {
    var xhttp = new XMLHttpRequest();
    xhttp.onreadystatechange = function() {
        if (this.readyState == 4 && this.status == 200) {
            var jsonObject=JSON.parse(this.response);
            var json = jsonObject.features;
            for(i in json){
                var id = json[i].properties.ID;
                allGeoJSONFeatures["p"+id] = json[i];
            }

        }
    };
    xhttp.open("GET", "http://ip:host/geoserver/xxx/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=BaskentGaz:M_DIZAYN_VANA&outputFormat=application%2Fjson", true);
    xhttp.send();
}

function showFeatureByID(id) {

    var geojson = allGeoJSONFeatures["p"+id];
    alert(geojson.toString());

}


texts:
I am attempting to count how many dissemination areas are in a group of school catchment areas. I am using two shapefiles from different sources and use GeoPandas' overlay function to generate a dataframe with the overlapping areas.
In the attached image, the coloured areas represent catchments, the solid lines represent the boundary of dissemination areas.

This code generates the dataframe with the overlaps.
import geopandas as gpd

schools = gpd.read_file(path to shapefile)
schools=schools.to_crs(epsg=26917)

da = gpd.read_file(path to shapefile)
da=da.to_crs(epsg=26917)

overlaps = gpd.overlay(da,schools, how='intersection')

The code functions properly, but gave me unexpected results: the counts were WAY high. The issue seems to be that each catchment has points where it overlaps minimally with the adjacent dissemination areas. For example, the second image is zoomed in area and you'll notice the purple catchment just barely squeaks over the line which is in the middle of a river. As it stands, the script captures that as being part of the purple and green catchment areas, which would be fine if it were truly part of both catchment areas.
One idea I had was possibly trimming a few meters from the edge of each polygon to reduce the number of "collisions" between adjacent polygons, but I don't know if such a thing is possible. Other suggestions are welcome too.


texts: Is it possible to modify the ArcObjects IProgressDialog2?

texts:
It turns out that the best way to run long tasks in ArcMap is to use the IProgressDialog2. The BackgroundWorker in combination with ArcObjects (STA Threads) is a mess.

However, it seems that the only implementation from ESRI is the IStepProgressor. In my case I do not need and do not want to show a progress bar. Moreover I do not need the cancel button.
Is there a way to remove both things? I want only the text and the circle as you can see in my screenshot

public static void ShowProgressDialog<T>(string message, Action<T> action, T arg)
    {

        ESRI.ArcGIS.esriSystem.ITrackCancel trackCancel = new ESRI.ArcGIS.Display.CancelTrackerClass();

        ESRI.ArcGIS.Framework.IProgressDialogFactory progressDialogFactory = new ESRI.ArcGIS.Framework.ProgressDialogFactoryClass();

        // Set the properties of the Step Progressor
        System.Int32 int32_hWnd = ArcMap.Application.hWnd;
        ESRI.ArcGIS.esriSystem.IStepProgressor stepProgressor = progressDialogFactory.Create(trackCancel, int32_hWnd);

        ESRI.ArcGIS.Framework.IProgressDialog2 progressDialog2 = (ESRI.ArcGIS.Framework.IProgressDialog2)stepProgressor; // Explict Cast

        progressDialog2.CancelEnabled = false;
        progressDialog2.Description = message;
        progressDialog2.Animation = ESRI.ArcGIS.Framework.esriProgressAnimationTypes.esriProgressSpiral;

        action.Invoke(arg);

        trackCancel = null;
        progressDialog2.HideDialog();
        progressDialog2 = null;

    }


texts:
I am working on Lunar SAR data and when opening the raster (.tif) in QGIS shows me distorted image with the following error:
No transform is available between Unknown CRS: GEOGCRS["Moon Spheroid",DATUM["Moon_Spheroid",ELLI… and Unknown CRS: GEOGCRS["unknown",DATUM["unknown",ELLIPSOID["WGS 8….
proj_create_operations: Source and target ellipsoid do not belong to the same celestial body

Zoomed image:

I also have an Excel sheet with Lat, long information of all pixels, but don,t know how can I make use of it.

How to resolve the error and visualize the data properly?
Google drive link to some data: https://drive.google.com/file/d/1uRZNFh3774BxDAE2rEqWthwx0oxuieLz/view?usp=sharing


texts: How to add rules to the topology from a .rul in ArcPy

texts:
I am trying to create a new topology using and ArcPy script and I have a lot of rules inside a .rul archive. The function to do that is:

arcpy.AddRuleToTopology_management(out_topo, rule_type, in_fc1, subtype1, in_fc2, subtype2)

The main fact is that I don't know if I can simply add the .rul archive as a "rule_type" or maybe I have to code a loop to iterate every topology rule and fc related inside the .rul archive.


texts: Not getting updated zoom level after extent is set? (ArcGIS JavaScript API)

texts:
There is a function in my map that works very nicely for zooming to a polygon (parcel) that has been selected:

_zoomToParcel: function(parcel) {
  this.map.setExtent(parcel.geometry.getExtent().expand(1.5));
}

The end user has requested that the zoom level be restricted to no more than 18. If it's a very small parcel, the zoom level may be 19 or 20, and then they have to zoom back out to find enough surroundings to orient themselves. Fair enough. I thought this would be as simple as adding in a get/set zoom:

_zoomToParcel: function(parcel) {
  this.map.setExtent(parcel.geometry.getExtent().expand(1.5));
  if (this.map.getZoom() > 18) {
    this.map.setZoom(18);
  }
}

However, this appears to be getting the zoom level from before the application resets the extent to the parcel. When I console.log the results of this.map.getZoom(), it tells me the initial zoom level (5), not the new zoom level (20) after the extent is reset. Running that function if the map starts off very zoomed-in (e.g. 20) will zoom it back out to 18, meaning the get/set logic works — I'm just not getting the "right" zoom.

How can I catch whether the map has zoomed in too far?


texts: Trying to find historical building footprints or centroids

texts: I'm currently doing a research project analyzing how cities evolve over time and was wondering if there were any good online databases for finding historical building footprints. I've been looking around for a while now but couldn't find what I was looking for.


texts: Offsetting polygon borders where they share border in QGIS?

texts:
Say you have 2 polygons with different coloured borders, where the borders are shared the colour is problematic.

How do you offset/separate them in QGIS 2.10.1?


texts: Inserting an annotation onto all plotted points?

texts:
I have an attribute table of 73 polling stations, and using the labelling function, I have labelled all of the points in preparation for printing on polling cards for our some 60000 electorate.

After playing with the annotation tool, I've found that it is much clearer than any variation of label I can produce - I can also then go through each annotation and re-position it so that it doesn't cover up any nearby street or building names - (as far as I can tell labels cannot be repositioned on an individual basis).

Is there any way to place an annotation onto all of the plotted points whilst retaining the same settings in respect of text font/size and the size of the annotation itself?

Here is an image of the finalised map:

This is an image of the labels, as you can see it does not look as professional as above:


texts:
I have been unable to find any existing repositories or modules for building fixtures and doing unit testing for ArcGIS. I don't think this is a part of computer science that GIS has fully embraced,
Are there any existing tools or python modules for building fixtures/mocking up spatial data/unit testing for Arc?
I know that mock can be used to create anything, but curious to see if there is anything preexisting.
a few resources I was able to find:
https://notesfromthelifeboat.com/post/testing-arcpy-1/
https://notesfromthelifeboat.com/post/testing-arcpy-2/
And some of the related code here from the author:
https://github.com/lobsteropteryx/testing-arcpy


texts: How to avoid duplicate labels on each tile?

texts:
When I display my map using tiles a label gets drawn on each tile - is there any way to avoid this?

I have attached an example image to show the problem:


texts:
I am updating a plugin from QGIS 2 to 3.  I need to iterate through the groups (the names are dynamic) and remove them with their respective child layers.  I need to get a list of groups and then iterate through the groups and delete them.  The original code is as follows:

legend = self.iface.legendInterface()
group_names = legend.groups()
rev = self.groups[:]
rev.reverse()
...
    for i, group in enumerate(rev):
        legend.removeGroup(group)

I replaced self.iface.legendInterface() with QgsProject.instance().layerTreeRoot() but group_names = legend.groups() fails with the message 'QgsLayerTree' object has no attribute 'groups'

What is the proper container of the groups and how do I access the collection; and, once found, how do I delete those groups?  A search for information on how to remove groups from a layer tree in QGIS 3 has proven fruitless.


texts: Using DTED file as raster elevation source at ArcGIS SDK for Qt?

texts:
I am using ArcGIS SDK for Qt 100.7 for development and I was getting elevation from the code

ArcGISTiledElevationSource* elevationSource = new ArcGISTiledElevationSource(QUrl("https://elevation3d.arcgis.com/arcgis/rest/services/WorldElevation3D/Terrain3D/ImageServer"), this);

However I want to gather elevation data when I am offline.

How can I pull elevation data to my GIS when I am offline for .dt2 file format?


texts:
I'm using a big table imported by a *.shp file of the city of San Diego. When I run pgr_dijkstra it takes more or less 5000ms, even for small distances. I was thinking, in order to improve performance, how could I select a subset of the table San Diego? And more important, do you think it would improve the performance?

My idea is to, given a starting and ending point, select only an area (a circle, square,...) that contains both starting and ending point (plus an arbitrary margin). Is it possible?


texts: Pulling date from filename and in proper format for time enabled use?

texts:
My feature attribute table has a folderpath that indludes a date time.
Example:/documents/foldername/xyz2015OCT15063333.shp

It seems very simple, but I am not a programmer and have done my best to calculate the field from examples I have been able to find.

I am using ArcGIS Desktop 10.6.1.

I need this date and time pulled out and placed into a new Date field to enable time.


texts:
I have a dataset that includes populations and subpopulations, so for each main population there may be several subpopulations. I have organised this in the table so that a population has the same code, and then in a separate column if it is a subpopulation then it has been given a letter 'a, b, c, d, etc' with the 'main' populations all having a blank cell under the subpopulation column. For example:

So in the above dataset, row no. 1 is the main population and has the accompanying subpop code column is blank/null and it has 8 subpopulations (a-h).

Is there some way, some expression I can use, under symbology to give all 'main populations' i.e. those with a NULL subpopulation code a particular symbol and all those populations where the subpopulation code column is NOT NULL another symbol? I have many different species, so ideally I would like the main population for each species to be a different colour, but a larger point dot. And all the subpopulations for each species to have the same colour as the main population but be a smaller dot.

So essentially showing something that looks like the above map. I just used filters to get this result and created a duplicate layer so that I could set the symbology differently for each.

Really hope that makes sense! I am still learning QGIS and know the basics, but just starting to get my head around using more advanced features.

Would appreciate any advice you can provide and if there's a better way to tabulate the data.

Thanks a bunch!


texts: Reading WFS layer with sf::st_read. Is it possible to apply a filter along with the request?

texts:
I'd like to read a large WFS layer in R. But since I don't need the entire layer, I'd like to filter it by some value in the attribute table. I know how to read the entire WFS layer with sf::st_read. Is there a way to do this filtering upon layer request? I tried using the where option in st_read but could not figure out how that works.
I can work with other R packages if needed.
Here is a simple reprex of what I'm trying to do.
library(sf)
#> Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 7.0.0
library(ggplot2)

wfs_url <- 'WFS:https://geoservicos.ibge.gov.br/geoserver/ows'
lyr <- 'CCAR:BC250_2017_Unidade_Federacao_A'

br_states <- st_read(dsn = wfs_url,
                     layer = lyr)
#> Reading layer `CCAR:BC250_2017_Unidade_Federacao_A' from data source `WFS:https://geoservicos.ibge.gov.br/geoserver/ows' using driver `WFS'
#> Simple feature collection with 27 features and 9 fields
#> geometry type:  MULTISURFACE
#> dimension:      XY
#> bbox:           xmin: -73.98318 ymin: -33.75118 xmax: -28.84777 ymax: 5.269508
#> geographic CRS: SIRGAS 2000

br_states <- st_read(dsn = wfs_url,
                     layer = lyr)
#> Reading layer `CCAR:BC250_2017_Unidade_Federacao_A' from data source `WFS:https://geoservicos.ibge.gov.br/geoserver/ows' using driver `WFS'
#> Simple feature collection with 27 features and 9 fields
#> geometry type:  MULTISURFACE
#> dimension:      XY
#> bbox:           xmin: -73.98318 ymin: -33.75118 xmax: -28.84777 ymax: 5.269508
#> geographic CRS: SIRGAS 2000

head(br_states)
#> Simple feature collection with 6 features and 9 fields
#> geometry type:  MULTISURFACE
#> dimension:      XY
#> bbox:           xmin: -73.98318 ymin: -33.75118 xmax: -49.69163 ymax: 5.269508
#> geographic CRS: SIRGAS 2000
#>                             gml_id id_objeto              nome nomeabrev
#> 1 BC250_2017_Unidade_Federacao_A.1         1 Rio Grande do Sul      <NA>
#> 2 BC250_2017_Unidade_Federacao_A.2         2          Amazonas      <NA>
#> 3 BC250_2017_Unidade_Federacao_A.3         3              Acre      <NA>
#> 4 BC250_2017_Unidade_Federacao_A.4         4          Rondônia      <NA>
#> 5 BC250_2017_Unidade_Federacao_A.5         5           Roraima      <NA>
#> 6 BC250_2017_Unidade_Federacao_A.6         6             Amapá      <NA>
#>   geometriaaproximada sigla geocodigo id_produtor id_elementoprodutor
#> 1                 Não    RS        43      250000                  NA
#> 2                 Não    AM        13      250000                  NA
#> 3                <NA>    AC        12      250000                  NA
#> 4                <NA>    RO        11      250000                  NA
#> 5                <NA>    RR        14      250000                  NA
#> 6                <NA>    AP        16      250000                  NA
#>                             geom
#> 1 MULTISURFACE (POLYGON ((-53...
#> 2 MULTISURFACE (POLYGON ((-67...
#> 3 MULTISURFACE (POLYGON ((-68...
#> 4 MULTISURFACE (POLYGON ((-60...
#> 5 MULTISURFACE (POLYGON ((-60...
#> 6 MULTISURFACE (POLYGON ((-50...

ggplot(br_states) +
  geom_sf()

# I can filter the data once I've downloaded it
ggplot(br_states[br_states$sigla == 'SP',]) +
  geom_sf()

# How can I filter it before the WFS request?
# Is there a way to filter inside st_read

Created on 2021-03-12 by the reprex package (v1.0.0)


texts: I'm wondering if anyone can point me to resources for making a heatmap entirely within PostGIS. I'm able to create a new raster with an empty band, but it's not clear to me how I can assign values based on the density of nearby vector features.


texts: Extracting values of some columns with ArcPy

texts:
I've been trying to run some code to extract the values from a column and random sample some of them.
The problem is when I try to append values from the cursor to a list.
When I append values to the list, it works, but it appends the value inside a tuple and with a comma.
Here is some sample of the code:
import arcpy
import random

arcpy.env.overwriteOutput = True
ws = r'D:\Projeto_VANT\SIG\proc_parc.gdb'
arcpy.env.workspace = ws

teste = r'teste'
with arcpy.da.SearchCursor(teste, "ID_TALHAO") as cursor:
    list = []
    for x in cursor:
        list.append(x)
    print(list)
    sample = random.sample(list, 2)
    print(sample)

And this is the results (list of values and the sample) I got which were printed:
[(5,), (13,), (20,), (30,), (38,)]
[(30,), (13,)]

How can I extract these values and append them to a list as integer, outside the tuple and without comma?


texts: What are the pros and cons of using the ArcGIS Viewer for Flex vs. creating your own GIS application vs creating your own framework?

texts:
Simple question, but perhaps hard to answer. I have broken the title question into three distinct questions in order to clearly identify what we are considering for options each time we build an GIS web application in Flex.

What are the pros and cons of using the ArcGIS Viewer for Flex application as a starting point for GIS applications?
What are the pros and cons of creating your own application from scratch using the ArcGIS API for Flex? Note: This would be a focused application. In this case it is not intended to be a reusable application from which other applications could be configured / built.
What are the pros and cons of creating your own custom framework using the ArcGIS API for Flex from which multiple applications could be built? While the primary goal would be to build an application that meets the customers needs a secondary goal would be to have a reusable framework from which other applications can be configured / built, much like the Flex Viewer.

Unlike many other decisions the choice to use one method over the other is not an easily reversible one. Assume you will have to live with this decision for a very long time (years).


texts:
I'm having a problem with "flooded" islands when converting OSM natural=water data to shapefiles and opening it in qgis. It doesn't happen consistently; as you can see from the attached screenshot, some islands are fine, but others initially show up, then get covered by the parent water polygon.

Here's how I'm generating my data:

osmconvert north-america-latest.osm.pbf -b=-80,40,-70,50 --complete-ways --complex-ways -o=w080n40.o5m
osmfilter w080n40.o5m --keep= --keep-ways="natural=water" --keep-relations="natural=water" -o=water.osm
ogr2ogr -skipfailures -f 'ESRI Shapefile' water water.osm

Then I open the resulting polygon layer in qgis, with the result illustrated.


texts:
I get the error below when i run the example by snowman2 here Rasterizing polygon grid in Python (geopandas & rasterio)
CRSError: CRS is invalid: <Geographic 2D CRS: EPSG:4326>
Name: WGS 84
Axis Info [ellipsoidal]:

Lat[north]: Geodetic latitude (degree)
Lon[east]: Geodetic longitude (degree)
Area of Use:
name: World
bounds: (-180.0, -90.0, 180.0, 90.0)
Datum: World Geodetic System 1984
Ellipsoid: WGS 84
Prime Meridian: Greenwich


texts: QGIS 3x conditional rotation of the label

texts:
I need to rotate only one label from three seen below

I found some solutions here:
How to rotate labels by the same angle as the map in the composer?
QGIS Labels Rotation
and also here:
Using Trigonometry To Place And Orientate Labels
But it looks like they are referred for all labels belonging to the given layer. I just want to set the condition for for rotating one of them.
My expression could look like this:
 CASE
  WHEN "fid" = 18 then rotation 'degrees 45'
 End

but QGIS says, that my expression is invalid.
Is anyone able to help?


texts: How to detect temporal changes using ArcGIS?

texts:
I want to show urban sprawl in a city over ten year period. I took a MOOC on geospatial revolution, therein a link was shared to do such analysis. I just want to show the urban sprawl, nothing complicated. However, I can't seem to remember that link or technique.

Kindly help me with the quickest technique to show this type of spatio-temporal changes as I have no knowledge of the programming. However, I do remember the technique involved the vegetation, and consequently the built-up space over it over some period of time.


texts: Change the basemap background color from white to black in arcgis js api?

texts:
I have a basemap layer that I published and placed on my application and this layer has a number of blank areas that shows up as white areas since the background of the layer itself is white.

Is there any way of changing this to a different color such as black?

I tried changing the data frame color before publishing the layer however this had no impact on the application itself.


texts: ArcGIS API 4.11 JavaScript Feature Table Alternative

texts:
I cannot find a simple/basic example online of adding a feature table to the bottom of a map in ArcGIS API 4.11 from a Feature layer (using JavaScript). So the functionality would allow you to click on a single symbol on the map and it would bring back that features columns.

Example image:

So a basic code starting point would be:

<link
  rel="stylesheet"
  href="https://js.arcgis.com/4.11/esri/themes/light/main.css"
/>

<style>

  #viewDiv {
    position: absolute;
    left: 0;
    right: 0;
    top: 0;
    bottom: 0;
    height: 100%;
  }

</style>

<script src="https://js.arcgis.com/4.11/"></script>

<script>
  require([
    "esri/views/MapView",
    "esri/WebMap"

  ], function(
    MapView,
    WebMap
  ) {

    const webmap = new WebMap({
      portalItem: {
        id: "83bd8ae0349646a4bd98c1a60e7f835b"
      }
    });

    const view = new MapView({
      map: webmap,
      container: "viewDiv",
      constraints: {
        minScale: 300000
      }
    });

  });
</script>

I believe this functionality has not been carried over to 4.11 yet. The only example on line that I can find is below. But extracting the relevant code from the example I'm finding quite difficult. It seems that the code examples on the ArcGIS website are either very basic, or far more complex and there is nothing in the middle.

ArcGIS Highlight features by geometry

Can anyone provide an example or point in the direction of one which is basic/clear/informative?


texts: Selecting small area from larger map to work within using QGIS?

texts:
I am trying to create a map in QGIS of a specific section of land in Oklahoma that contains subdivisions.  I have been able to create and view a map of the entire county using shapefiles for the sections, subdivisions, blocks, lots and road centers.

Is there a way to single out the data contained in that one section and work with only that?


texts: Modify attributes and field calculator in ModelBuilder?

texts:
I have approximately 200 point shapefiles that have been created from ASCII grid files.

Each one has a Grid_code heading in the attribute table, I need to add in a new heading called height or depth which is a double of length 16, precision 2 that equals the grid code (basically the grid_code but to 2 decimal places). Is there a way to do this in ModelBuilder or another way that saves me from doing them all manually?

Also would be good to delete the grid_code afterwards but this is not essential.


texts: Joining Population Data to Census Tract Shapefile using ArcGIS Desktop?

texts:
I'm trying to create a population density map and I'm struggling with getting the population data to join correctly to a census tract shapefile.

Between the census tract shapefile and population data there is only one attribute shared between the two, the "geoid" column.

However, the amount of numbers in each "geoid" column differs and it's preventing me from creating the join that will allow me to make my map.

I'm pretty new to GIS.


texts:
I'm attempting to do a walking distance analysis using OpenStreetMap data (osm2pgrouting) and PostGIS and following underdark's blog posts, however I need to exclude residential and service roads that cross primary and secondary roads, as they are unsafe to cross. To do this I think my first step is to split edges at intersections. I'm guessing I'll want to use the ST_Split function, but I'm not sure how to put it all together in one SELECT statement.

(A later question will be how to remove such crossing edges based on the highway=* tag, but I need to figure this one out first.)


texts: WMTS in Leaflet - Can't make this code work?

texts:
I know this have been touched upon before, but the other posts did not entirely do it for me.
I'm trying to make a minimal Leaflet .html page with a map, containing a WMTS map. The code looks like this:

    <title>TEST leaflet.TileLayer.WMTS</title>
    <meta charset="UTF-8">
    <meta name="viewport" content="width=device-width">
    <link rel="stylesheet" href="http://cdn.leafletjs.com/leaflet-0.7.2/leaflet.css" />

The resulting page can be seen live here: http://snapseatlas.dk/daftest/daf.html

Not that there is much to see, because only a grey background is shown.

I suspect my problem can be either:
1) I missed something fundamental, because I'm an Leaflet amateur, or
2) Something with the projection (EPSG:25832, is common in Denmark).

So please - Can somebody with a bit more Leaflet/WMTS experience have a look, and tell me what I'm doing wrong.

/M


texts: How to fix problem with ModelBuilder batch processing raster output names

texts: I am using ModelBuilder in ArcCatalog 10.3 to automate a large workflow. As part of this process I am running my inputs as a list (batch processing) and using in-line variable substitution so that I can consistently name the output files. I had the entire model built and everything was working as expected until recently. Now for some raster tools the model is changing the output raster names from what I designate them to a default output name. For example if I have three outputs from euclidian distance all named PRH_Distance_%i% they are automatically being changed to EucDist_PRH_1, EucDist_PRH_2, EucDist_PRH_3 So far just euclidean distance, slope, plus, and aspect are affected. I was wondering if anyone else had encountered a situation like this or had any thoughts on how to get ArcCatalog to use the names I am designating.


texts: Creating a buffer based on area (i.e. sqm) instead of distance

texts:
Is there a way in ArcGIS 10.1 or QGIS (or any add-on) to create buffer based on an area (i.e. sqm) rather than using buffer distance? I'd like to create a buffer of let's say 5% of the area within the buffer. The problem is that the feature is not regularly shaped.
Thanks,
Alex


texts: Inserting zoomed-in plot with full scale plot in sp package in R

texts:
I am using the sp package in R to plot sampled points along a river. However, I would like to zoom in to a sampled area of the entire river, and plot this zoomed area on the full length of the river. Using the spplot function is useful, as I have various variables per point. To "zoom" in, I can simple change the ylim and xlim arguments when using the spplot function.

Using the meuse dataset as an example

library(sp)
library(lattice)
data(meuse)
coordinates(meuse)=~x+y
FUll.study<-spplot(meuse)            #Gives the full study area
ZOOM<-spplot(meuse, xlim=c(170000,  181000), ylim=c(329714,333611))   #Random area to zoom in.

The results look like this:

However the desired result is something like this:

I am trying to find a way to add the ZOOM to an area on the Full.study map.I suspect because spplot is dependent on lattice it may not be possible. In that case, I would have to create a map for each variable independently, and add the respective zoomed area, probably through par(new=TRUE), but this is just a thought. I would still have to figure out how to specify the location and adjust the size of the of zoomed area.


texts: Creating a Legend with arcpy and comtypes

texts:
I am relatively new working with python and arcobjects. Under normal limitations I know that you cannot create a Legend or any arcobject with arcpy.mapping. Is there a generic example of creating a legend or an arcobject using comtypes to get me in the right direction?

I am trying to create as script to automate several maps and export them to jpg after. (All the code is figured out outside the codeblocks with comtypes)

I was able to figure out (with some search) how create a MXD but not sure which arcobjects to be replace for creating a legend with CreateObject.

import arcpy, os
from comtypes.client import GetModule, CreateObject
from Snippets import GetStandaloneModules, InitStandalone

arcpy.env.workspace = r"E:\foo"
wksp = arcpy.env.workspace

# Code that is out of the scope of the question

# Get the Carto module
esriCarto = GetModule(r"C:\Program Files (x86)\ArcGIS\Desktop10.2\com\esriCarto.olb")

for boro in boros:
    # Create a map document object
    mxdobj = CreateObject(esriCarto.MapDocument, interface=esriCarto.IMapDocument)
    # Create a new borough mxd with mapping.MapDocument
    print "Creating Wifi_{}.mxd".format(boro)
    path = r"{}\Wifi_{}.mxd".format(wksp,boro)
    print path
    # Creates mxd
    mxdobj.new(path)
    #delete obj
    del mxdobj

# Rest of the code


texts: Iterating between two layers using QGIS python

texts:
I will try and explain my situation as clearly as possible.

I just created my first plugin that is being used to analyze between various layers and give back quantities to the user. Most of it is straight forward until I have to try and involve two layers.

I have a point layer. With that point layer, I want to create a "buffer" of 10' and use that buffer as an intersect to analyze a polyline layer if it is within range. Once that buffer detects said polyline layer, I want to be able to read from multiple fields attributes to see if they meet a certain criteria. If that criteria is met, then that point will be selected.

To further explain, the point layer represents Assets (eg. Power pole, XFMR, Pedestal, Fuse Pad, etc.). If that asset is a Proposed Pole, I want to create a 10'"buffer". If that buffer intersects one of our Proposed Conductors, I want to to be able to analyze that Proposed Conductor and know whether or not it is a Single Phase or Three phase conductor. If it meets my criteria, I want the code to be able to select that feature (Proposed Power Pole).

What I am creating with this plugin is a Payline schedule. So once the designers are done designing, they will click my magic button for all of their results.

I am running QGIS version 2.0.1

My main goal would be to create the buffer at a point feature and return the polylines within that 10' radius.


texts: Separate polygons based on intersection using PostGIS

texts:
I have a PostGIS table of polygons where some intersect with one another.  This is what I'm trying to do:

For a given polygon selected by id, give me all of the polygons that intersect.  Basically, select the_geom from the_table where ST_Intersects(the_geom, (select the_geom from the_table where source_id = '123'))
From these polygons, I need to create a new polygons such that intersection becomes a new polygon.  So if polygon A intersects with polygon B, I will get 3 new polygons: A minus AB, AB, and B minus AB.

Any ideas?


texts: Is it possible to remove select geometries from GeometryCollection elements within Google Earth Engine?

texts:
I am relatively new to Earth Engine and am trying to understand GeometryCollection elements. I have seen this post, but otherwise haven't found much additional information. Coming from desktop GIS software, the notion of a single Feature with multiple geometric types is unusual. A couple of related questions:
(1) Conceptually: For a given Feature whose geometries are defined by a GeometryCollection, is it possible to remove portions of a GeometryCollection, or does that violate internal topology? For instance, if a single Feature is composed of LineString and Polygon types, does removing the LineString components somehow invalidate the Polygon types, since they are inherently bound together?
(2) Practically: The following example returns a FeatureCollection whose Feature geometries are either Polygon, MultiPolygon, or GeometryCollection. GeometryCollection geometries are composed of LineString and Polygon geometries. If it is possible to remove portions of GeometryCollections per (1) above, how would this be accomplished? The goal is to reduce the parent FeatureCollection into a single geometry type (or set of geometry types) that can be exported to shapefile (it does not appear that Earth Engine has built-in functionality to clean/separate incompatible GeometryCollection elements during export).
(3) There are plenty of operations that could be performed on mixed geometries, and for enormous datasets there may be some non-trivial space/cost reduction by using LineStrings instead of Polygons with near-zero width. However, if there is additional rationale for permitting mixed geometries, I would love to hear it.
// Area of interest
var aoi = ee.Geometry.Rectangle([-82.4685, 23.5529, -75.5471, 27.8867])
// Load and filter Earth Engine asset
var polys = ee.FeatureCollection("FAO/GAUL/2015/level2").filterBounds(aoi)
// Store Feature geometry type as a Feature-level property for filtering
polys = polys.map(function(f) {
  return f.set('geo_type', f.geometry().type())
})
// View full collection
print(polys)
// Isolate and view GeometreyCollection elements
var nonPoly = polys.filter(ee.Filter.eq('geo_type', 'GeometryCollection'))
print(nonPoly)

// Visualize as needed
// Map.setCenter(-78.6, 25)
// Map.addLayer(nonPoly)

// Convert GeometryCollection elements into Polygon or MultiPolygon elements
// ????

How to filter through geometry of feature collection instead of properties?


texts: How to disperse clusters of proportional markers or polygons?

texts:
I am mapping a set of points in ArcMap that are symbolized using proportional marker symbols based on a numeric field value. The scale of the map has lead to alot of the point markers overlapping.

I would like to disperse the markers so that they are not overlapping. A random dispersement would be fine, preferably with control over the spacing between the dispersed markers. It seems that proportional symbols can't be converted to representations (else I could have used the Cartographic Refinement - Disperse Markers tool).

I could buffer the points by the numeric field to create polygons instead of using point markers to get around the representation issue, but is there a way to disperse polygons away from each other so they are not overlapping (similar to the functionality of the Disperse Markers tool for points)?

I'm using ArcMap with an Advanced License. Thanks


texts: Extract a line segment out of a polyline on the basis of numeric values (distances) defined as attributes [QGIS]?

texts:
I have a set of points (red in the picture) in a shapefile containing an ID as well as two numeric values representing distances in Meters for each point.

After creating polylines for each point via r.drain-GRASS-tool representing flow paths to the points (green in picture), I now need to extract a segment of each polyline that is between a distance A and a distance B from the corresponding point.

So what I need is a shapefile containing a set of polylines, one polyline for each point that is between a distance A and distance B from the corresponding point. A and B vary between the different points.

I tried creating two sets of buffers around the points, one for distance A and one for distance B and then afterwards use the vector overlay algorithms of QGIS but the buffers for different points may overlap and a line will then be clipped by the buffers of two different points.

Is there a way to solve this problem with the buffers or even an easier way to do the extraction using QGIS (2.18 or 3)/GRASS?

Oh, and it is not possible do it manually by editing the features, because I want to use it on a set of several hundred points.


texts: Building walls to correct channelized flow in ArcHydro

texts:
My current project entails burning in a known river network into a high resolution (10x10m) DEM. In doing so, I'm attempting to build walls to correct for channelized flow in the DEM, in order to eventually execute the flow direction and flow accumulation tools.
I've identified two areas where there are man-made channels that don't follow the DEM. In these places, I created a polygon shapefile for "Inner Walls" and a polyline shapefile for "breakthrough lines" to put in the Build Walls tool.
In the photo below, I depict an initial flow accumulation DEM in blue just to show why I'm building the walls here.

Yet, I keep getting this error below and am not sure what steps to take to resolve it.


texts: GDAL/Python intersection does not include fields from method layer

texts:
As stated in title. I'm intersecting two layers and expecting fields from both input and intersect/method layers in the result. This is the behavior as describe in API docs (https://gdal.org/python/osgeo.ogr.Layer-class.html#Intersection). The resultant geometries/features are correct, but only the attributes from the input layer exist.
The field names are not the same and being overwritten. I've also tried it with Identity but no difference.
Code really boils down to: input_layer.Intersection(intersect_layer, output_layer)
I could write some manual code to do this by iterating through features and spatial joining to the method layer, but it adds quite a bit of script time. If there's a way to get Intersection to work as I need, that would be most efficient.
Python v3.7.7 and GDAL/Py from GIS Internals release-1911 gdal-3.1.3
FIXED: This was my mistake. The output layer needs to be created as passed as ogr.Layer. If there is no field/attribute schema in this dataset, then GDAL puts attributes of both layers. I had used a shortcut helper function that creates an empty copy of the input layer (so as to automatically set same SRS and such). But while devoid of features, it was copying the field schema, thus overriding default behavior of function.


texts: arcpy: read color value of shapefile loaded in ArcGIS Desktop

texts:
I've got various polyline shapefiles loaded in ArcGIS, each using the default line style with different colors. For a python addin I'm working on I would like to use arcpy to read the RGB value of these colors. Would that be possible?

I noticed this question: Read layer symbology class color with ArcPy
But it deals with raster layers and therefore I assume my question is not a duplicate (please let me know if this is still the case)


texts: Converting CALS type rasters for Arcmap?

texts: I have a bunch (think 400+) scanned old portions of maps that are in .cal format, apparently a CALS type raster format. I am able to open them in Autocad Map 3D, but Arcmap does not recognize them. I'd like to (eventually) get them in a GIS, but not sure how to batch convert .cal to something more appropriate to work with in Arcmap.


texts: Finding a pole of inaccessibility relative to multiple points without using rasters

texts:
I want to find a 'point of inaccessibility' from a series of points, representing towns and villages in Scotland. Using FOSS GIS tools (QGIS, SAGA, GRASS, Postgres/PostGIS ... ).

This is the point furthest from any centre of population. This would be the centre of the largest circle which doesn't contain a town, but it will have a town on its edge.

There's a raster approach, which works - the GDAL Proximity filter.

  Is there a VECTOR based way to do this? One that doesn't require rasterisation?

Here's my attempt with rasters...

Start with a point layer in a suitable projection in meters (I'm using OSGB, 27700)
Add a field called "one", and set to integer with constant value 1.
Vector to raster , editing command to use -a one and using GEO rather than pixel. Make sure pixels are square!
Create contours on the proximity raster (in this case, 1 km wide)

Here's what it looks like. I drew some candidate points by hand  ;-)

I thought that dissolving repeated buffers (using the Multi ring buffer plugin) might work, but I suspect I don't have enough memory for this to work.