texts: GeoDjango: Calculate if someone inside circle from given point

texts: I am working in GeoDjango (with a Django 1.11 and PostGIS 2.1.4 backend). I have a model with shops in my database (i have the address of this shops) and I am trying to find all shops within 100m radius of the given address.


texts: Replacing Negative values with "0" for raster layer

texts:
I am trying to change negative values (raster layer) into "0" but I cannot find a way to do that.

Does anyone know how to accomplish it in QGIS 3.2?

I found an answer somewhere on the Internet to use the Python console entering '0 if !field! < 0 else !field!' but I have no experience using Python.


texts:
I'm working on getting some DEMs ready to bring into Blender to render hillshade. I've done this successfully before with low elevations, but on my current project mountaintops seem to get cut off.

When I mosaic the data or rescale it using the map algebra tool, the highest elevations turn white (in a greyscale scheme), and lose all detail. Does anyone know that causes this or have a possible solution?

EDIT: Changing the stretch type from the default (Standard Deviations) to Minimum-Maximum and changing those values to the ones used in my rescale math did the trick!


texts: I was looking for INSPIRE metadata XSD schema definitions on the Web, but nothing found. Only something like THIS. But that resource is "a little" messy... I wanted to use xsd.exe to autogenerate C# classes that could be easily serialized into XML that is INSPIRE compliant. The problem is in finding that schema file. Could you point me where to download it?


texts: What is the equivalent of ArcGIS's "Feature Class to Feature Class" in QGIS?

texts:
I have a layer in QGIS to which I have joined several csv tables.  I want to save this as a permanent and independent feature class.  However, when I right click on the layer in the table of content and select "Save as", and specify a new feature class, it doesn't write out some of the features.  I have discovered that the reason is because some of the fields from he csv contain entries that are over 255 characters long, and by default, the text fields that are created are initialized with a set length of 255.  Any feature that tries to write in a string longer than 255 characters fails and simply does not get written.

I found out that in ArcGIS, I can use the "Feature Class to Feature Class" tool instead of the "Export" dialog, and this tool allowed me to specify the lengths of the newly-created text fields via a field mapping.

Is there an equivalent to "Feature Class to Feature Class" in QGIS? And if so, what is it called?


texts:
I've configured a WMS service on QGIS Server by using https://docs.qgis.org/3.16/en/docs/training_manual/qgis_server/wms.html and QGIS Server on a Windows machine -where do the qgs files go? but was trying to see how I can also serve WFS.
When I try 'http://localhost:8080/cgi-bin/qgis_mapserv.fcgi.exe' I get

WMS works fine using http://localhost:8080/cgi-bin/qgis_mapserv.fcgi.exe?map=C:%5Cpublic_ftp%5Cftpuser%5Cd%5CMGM%5CMapData%5Canimals.qgs

What do I need to setup to allow it to serve WFS?


texts: Saving MXDs to older version of ArcMap automatically?

texts: I work at a company that has three computers running ArcMap, accessing the same files and databases. Two of them run version 10.3, but we just got a new computer with version 10.6. Saving a document in 10.6 means it can no longer be opened in 10.3, which is a problem. We need all three computers to be able to open the maps. What can we do?


texts: Is there a way to create arrow subSymbol on top of geometry generated symbol using PyQGIS?

texts:
Currently, I am trying to figure out how to set a subsymbol for a PyQGIS generated geometrygeneratorsymbol layer. As can be found in the code I tried multiple things resulting in crashing the kernel or nothing at all. The images represent the manual workflow which works fine. I am interested in transferring this 'Simple line' option to 'Arrow' in an automated way.
from qgis.core import (
    QgsVectorLayer,
    QgsSingleSymbolRenderer,
    QgsLineSymbol,
    QgsFillSymbol,
    QgsProject,
    QgsRectangle,
    QgsArrowSymbolLayer,
    QgsMultiRenderChecker,
    QgsProperty,
    QgsSymbolLayer,
    QgsMapSettings,
    QgsSymbol
)

from qgis.core import QgsVectorLayer,  QgsProject
layer = QgsVectorLayer( "?query=SELECT from_country, ubc_processing_plant, delivered_tonnes_2019, make_line(Centroids_countries.geometry, UBC_processing_2030.geometry) FROM ubc_flows JOIN Centroids_countries ON ubc_flows.from_country = Centroids_countries.NUTS_NAME JOIN UBC_processing_2030 ON ubc_flows.ubc_processing_plant = UBC_processing_2030.name", "vlayer", "virtual" )
QgsProject.instance().addMapLayer(layer)

def test_3(self):
        sym = self.renderer().symbol()

        shape_sym = QgsGeometryGeneratorSymbolLayer.create({'geometryModifier': 'difference(difference(make_line(start_point($geometry), centroid(offset_curve($geometry, length($geometry)/-10.0)), end_point($geometry)), buffer(start_point($geometry), 10000)), buffer(end_point( $geometry), 10000))', 'outline_color': 'black', 'symbolType': '1'})
        shape_sym.setSymbolType(QgsSymbol.Line)

        arrow = QgsArrowSymbolLayer.create(
            {
                "arrow_width": "3",
                "head_length": "3",
                "head_thickness": "3",
                "head_type": "0",
                "arrow_type": "0",
                "is_curved": "1",
                "arrow_start_width": "3"
            }
        )
        fill_sym = QgsFillSymbol.createSimple(
            {"color": "#8bcfff", "outline_color": "#000000", "outline_style": "solid", "outline_width": "0.1"}
        )

        #shape_sym.setSubSymbol(arrow)
        #shape_sym.setSubSymbol(sym_layer)
        sym.changeSymbolLayer(0, shape_sym)

test_3(layer)


texts: Substitute new CRS for one specified in the Training Manual, Module 6, Lesson 6.2.5?

texts:
This question came up on the page at the link https://docs.qgis.org/3.16/en/docs/training_manual/vector_analysis/basic_analysis.html.
6.2.5 Try Yourself: Convert Layers’ CRS
Item 6. Change the [CRS] to WGS 84 / UTM zone 34S
In the Save Vector Layer as … dialog box, the fourth text box down shows
 CRS    EPSG:32734 - WGS 84 / UTM zone 34S

This doesn’t match any of the three choices in the dropdown menu. The second choice is
 EPSG:32700 - WGS 84 / UTM (southern hemisphere)

This CRS also appears in the list of Predefined Coordinate Reference Systems displayed in another dialog box obtained by clicking on the little ?globe icon to the right of the CRS menu. EPSG:32734 - WGS 84 / UTM zone 34S does not appear in the list.
Why does the instruction not match the choices? Will the choice of EPSG:32700 lead to error, or is this just a simplification of the reference labels?


texts: how to get Extent or BoundingBox from geoserver in openlayers3

texts:
I use wms layer from geoserver on my map and i use openlayers3 api.
I want to zoom to my selected layer but for this i have too use extent of my selected layer.

mySelectedLayer=map.getLayers().getArray()[map.getLayers().getLength()-1];
mySelectedLayer.getExtent()// result is undefined

why?


texts: I'm working with multi-temporal datasets from modis and had to do various vegetation indices. I ran into a problem now because there is a significant amount of water too, which skews the resulting map. How can I build a mask to just blend out the water bodies?


texts:
I'm using MS4W on my local server.

This are some codes inside my mymap.map file.

MAP
    etc.....

    OUTPUTFORMAT
        NAME pnggd
        DRIVER "GD/PNG"
        MIMETYPE "image/png"
        IMAGEMODE PC256
        EXTENSION "png"
    END
    OUTPUTFORMAT
        NAME pnggd24
        DRIVER "GD/PNG"
        MIMETYPE "image/png"
        IMAGEMODE RGBA
        EXTENSION "png"
    END
    OUTPUTFORMAT
        NAME AGG
        DRIVER "AGG/PNG"
        IMAGEMODE RGB
    END
    OUTPUTFORMAT
        NAME AGGA
        DRIVER "AGG/PNG"
        IMAGEMODE RGBA
    END

    etc....

and this is my code for my PNG file location to be called.

LAYER
    NAME "map1"
    STATUS DEFAULT
    TYPE RASTER
    DATA "/ms4w/Apache/htdocs/my_map/my_map_image/map1.png"
END

Here is my explorer.js file and some codes inside the file:

var view_map1 = new OpenLayers.Layer.WMS(
    "map1 overlay", "http://localhost/cgi-bin/mapserv.exe?map=c://ms4w/Apache/htdocs/my_map/mymap.map&LAYERS=map1&VERSION=1.1.1&SERVICE=WMS&REQUEST=GetMap&FORMAT=image/png",
     {
         layers: 'map1',
         format: 'image/png'
     },{
        'ratio':1,
        singleTile: false,
        maxExtent: my_extent,
        isBaseLayer: false
      })
 map.addLayer(view_map1);

Is my code inside the .map file and .js file correct? or it should be not like that when to declare the OpenLayers and also for the shapefile.

Can anyone spare me some ideas and suggestions?


texts: GeoServer S3 Support for GeoTIFF REST-API

texts:
I am stuck with S3 Support for GeoTIFF while I add a dataset with the REST-API.

I did not find any documentation for the REST-API... so I tried to edit my local API-POST-Call:

<coverageStore>
  <name>NameOrtho</name>
  <workspace>NameWorkspace</workspace>
  <enabled>true</enabled>
  <type>GeoTIFF</type>
  <url>s3://bucket/DEM.tif?useAnon=true&awsRegion=EU_CENTRAL_1</url>
</coverageStore>

But this returns:

---- Debugging information ----
cause-exception     : com.thoughtworks.xstream.io.StreamException
cause-message       :
class               : java.lang.String
required-type       : java.lang.String
converter-type      : com.thoughtworks.xstream.converters.SingleValueConverterWrapper
wrapped-converter   : com.thoughtworks.xstream.converters.basic.StringConverter
line number         : 6
class[1]            : org.geoserver.catalog.impl.CoverageStoreInfoImpl
converter-type[1]   : org.geoserver.config.util.XStreamPersister$StoreInfoConverter
version             : 2.15.0
-------------------------------

What should I do?


texts: Create a upload button in PyQt QGIS plugin

texts:
I want to create 2 upload file button in PyQt.
One for image and one for .csv file.
If the uploaded file is not in image or .csv file it should thrown an message.

I have search a lot on internet and in internet but couldn't find answer.


texts: Merging List of shapefiles gives ERROR 000400?

texts:
I'm trying to create a list of shp that are within my map and get the pathname of each shp so I can then add it as an attribute. Within the process I want to merge them and that's where I'm getting an error line 16:

  Failed to execute. Parameters are not valid. ERROR 000400: Duplicate
  inputs are not allowed Failed to execute (Merge).

I'm not sure if this due because my shp. have almost the same name or if it's something within my code. I must say I'm quite new to python so not really sure what I'm missing.

This is the script I'm using:


texts: Can QGIS make georectified orthomosaics of group of jpegs?

texts:
I have a file full of photo's taken with a 3DR drone over a small coffee farm in Peru. My intention is to create an Orthomosaic of these photos and than use the mosaic to automate a vegetation analysis all with open source software.

I have tested AirPhotoSE, but I am running into an 'IO Error' after importing my photos and the program crashes constantly.

VisualSFM was similarly unsuccessful. Despite having downloaded the supplementary applications ( SiftGPU, PBA, and VLFeat SIFT), whenever I attempt to 'Compute Missing Matches' I get these error messages:

     ERROR: siftgpu failed the detection test!!
     SIFT: 0000, 4000x3000,     ERROR1
     ERROR: unable to locate sift binary.

I am wondering if anyone can suggest another program or approach I can take. QGIS seems like a good next program to test.

Can QGIS make georectified orthomosaics of a group of jpegs?


texts: Merging FBX tiles in a single model file

texts:
I would like to obtain a single FBX file for a 3D photomesh, with texture. As the original scene is quite large, I use Bentley ContextCapture Master to generate tiles (a dozen), each being in FBX format with the .jpg texture on each tile. I cannot do the model in a single file (no tiling), because this would take +100Go RAM.
As I have to import this model in a software (Lumion) that recognizes .fbx + .jpg, I would like to merge all the tiles in a single FBX file.
Do you know if there is any solution for merging this tiles and textures, in a single FBX + JPG pair?


texts: Can the arrows produced by the d.rast.arrow GRASS script be queried or stored?

texts: I have created a map with the flow directions using the d.rast.arrow. Can i query on the map to know where the direction of the arrows are? I am going to use the arrow directions in an algorithm. Or, is there a way to store the arrow directions in a text file so that I can import it in my code? Thanks!


texts: Georeferencing Museum specimens with poor locality descriptions

texts:
I am currently collating a list of specimens from databases of various natural history museums for further research. However, a well known issue associated with majority of historical data is the lack of appropriate latitude and longitude which prevents one from using that data.

There have been ways to overcome that data - such as drawing a buffer around a region and providing a range of uncertainty associated with that location.

For instance, the function - biogeomancer from the package 'spatial' in R, automates the process of georeferencing, provided there are a few textual descriptions such as "2 miles west of XYZ". See documentation here.

However, my main concern is in using such a protocol for regions as big as 200 square km. Is there a way one can overcome that issue? I would love to use this rich trove of museum data, provided I can handle the uncertainty associated with its location.

An example of some specimens in my dataset is shown below. Please note that many of them come with mentions of elevation, but most of the records are very vague.

EDIT

In the comments section, one of you mentioned the purpose of this question and what I tend to achieve from the same.
1. I am interested in how once can reduce the radius of uncertainty from a really wide polygonal region to a smaller radius of uncertainty (if possible).
2. This information will help me carry out future spatial analysis such as species distribution modeling / occupancy modeling for instance.


texts: How to join weighted attribute from polygons to point

texts: I would like to join attribute from polygons to green point. Are there some methods? Something like "weighted mean by distance" of attribute from polygons or else...


texts: Identifying the position of a building footprint on a land parcel

texts:
I am looking to try classify where a building footprint falls on a land parcel.
I have identified the amount of free space that is available on a land parcel and am looking to try quantify this further.
When referring to the building position, I would like to identify if it is situated at the back of a land parcel or in the front of a parcel, or off to the side of the parcel (relative to the road)


texts: Can not SET a Value in table using Python?

texts:
I'am trying to set values in table, the code run good, i don't got errors, but nothing is done !
I'm using Arcgis10.1, arcpy and Python
This is the code :

 import arcpy
 import math
 import arcview
 from arcpy import env

 arcpy.env.workspace = "D:/Users/saadiya/Documents/ArcGIS/Default1.gdb"
 arcpy.env.overwriteOutput = True

 cur = arcpy.UpdateCursor("Join_Dataset")
 for row in cur :
     SW = row.getValue("RefName_SW")
     S  = row.getValue("RefName_S")
     SE  = row.getValue("RefName_SE")
     new_s = S
     if (SW != S):
         new_s = new_s+", "+SW
     elif (SE != S):
         new_s = new_s+", "+SE
         row.setValue("Sud", new_s)
         cur.updateRow(row)
         print new_s

I got this error :

Traceback (most recent call last):
  File "<string>", line 11, in <module>
  File "d:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\arcobjects\arcobjects.py", line 102, in updateRow
    return convertArcObjectToPythonObject(self._arc_object.UpdateRow(*gp_fixargs(args)))
RuntimeError: ERROR 999999: Error executing function.
Objects in this class cannot be updated outside an edit session [Join_Dataset]


texts: Post processing GPS data with open source software

texts:
I am using Magellan Professional Mapper CX working on ArcPad 7.1 platform. It has sub-meter accuracy after post processing. Without post processing the data collected in the field is off by 10 to 15 feet.

Are there any open source software/Web applications that can post process the GPS data?


texts: Exporting lat/long points from shapefile gives crazy results


texts: Sorting GeoJSON features array by property value

texts:
I am looking to sort a GeoJSON file on a property so I can then slice it keeping only the top 5 features. For example, I wish to take this GeoJSON and sort the file in descending order by the incidents property:
{
  "type": "FeatureCollection",
  "features": [
    {
      "type": "Feature",
      "geometry": {
        "type": "MultiPolygon",
        "coordinates": [
          [
            [
              [
                -75.1972382872565,
                39.9294288475177
              ],
              [
    ...
         [
                    -75.1972382872565,
                    39.9294288475177
                  ]
                ]
              ]
            ]
          },
          "properties": {
            "area": "",
            "location": "24th St. & Wolf St.",
            "perimeter": "81903.64182498",
            "dist_numc": "01",
            "sum_area": "",
            "area_sqmi": "216350124.153",
            "district_i": "",
            "objectid": "321",
            "globalid": "c040618c-ac93-4b22-b174-0ea08e0f805d",
            "district_": "1",
            "_feature_id": "1",
            "_feature_id_string": "Police_Districts.1",
            "dist_num": "1",
            "div_code": "SPD",
            "phone": "686-3010",
            "incidents": 11553
          }
        },
    ...

I'm looking to basically rearrange the order of the features array so that it's in descending order and the highest incident features are first.
How would I go about this?
I am used to creating a compare function and then sorting an array with something like this:
//FUNCTION: SORT
 export const compare = ( a, b ) => {
    if ( a.incidents < b.incidents ){
      return -1;
    }
    if ( a.incidents > b.incidents ){
      return 1;
    }
    return 0;
  }

But in this case I need to sort a JSON and only sort an array inside of it so I am unsure how to proceed.


texts: Join non-spatial CSV to GeoJSON data source using ogr2ogr and SQLite dialect

texts:
I am attempting to do a simple join of a non-spatial CSV file to a GeoJSON data source using the built-in SQL capabilities in ogr2ogr. I have successfully managed to do this using the OGR SQL dialect with a command like this:

ogr2ogr -f geojson -sql "select input.*, join_tbl.* from input join 'join_tbl.csv'.join_tbl on input.id = join_tbl.id" output.geojson input.geojson

To run it using the SQLite dialect I added -dialect sqlite to the previous command:

ogr2ogr -f geojson -dialect sqlite -sql "select input.*, join_tbl.* from input join 'join_tbl.csv'.join_tbl on input.id = join_tbl.id" output.geojson input.geojson

This results in the following error:

ERROR 1: In ExecuteSQL(): sqlite3_prepare_v2(select input.*, join_tbl.* from input join "_OGR_1" on input.id = join_tbl.id):
  no such table: join_tbl

I opened up the CSV in Excel and as expected, the table name is the same as the filename. That makes sense. CSVs are flat files and can't contain multiple tables.

So next, I tried removing the table assignment ('join_tbl.csv'.join_tbl -> 'join_tbl.csv') from the previous command:

ogr2ogr -f geojson -dialect sqlite -sql "select input.*, join_tbl.* from input join 'join_tbl.csv' on input.id = join_tbl.id" output.geojson input.geojson

However, this still resulted in an error:

ERROR 1: In ExecuteSQL(): sqlite3_prepare_v2(select input.*, join_tbl.* from input join 'join_tbl.csv' on input.id = pg.id):
  no such table: join_tbl.csv

There is an example for carrying out a join using the SQLite dialect in the ogr documentation, but it is for a .dbf file:

SELECT p.*, NAME FROM poly p JOIN "idlink.dbf"."idlink" il USING (eas_id)

Based on the errors I've been getting, this syntax does not appear to suggest my path forward.

How should I proceed?


texts: Exporting layer metadata information from QGIS?

texts:
Do you know a method to export information properties (with metadata…), like I see in Vector Properties Dialog in QGIS 3.4.  Presentation is great for me.

I can make copy-paste, but it's tedious and not really optimise for paste in Calc.


texts: How accurate is the Scale Line in OpenLayers

texts:
Is there any source describing the accuracy of OpenLayer's scale line?
You can get an scale line by using the constructor:

map.addControl(new OpenLayers.Control.ScaleLine());

However, it seems to be very inaccurate in relation to my geodesic measurements (0.5 km - 1.5 km difference depending on zoom). Is there any info I can cite to explain this different?

I'm using OpenLayers "Release 2.13.1" and OSM as base layer. I'm calculating the distance between point using WGS84 projection with 3 decimals and geopy distance module.

As an example, the distance between the points -2.416,55.703 and -2.434,55.698 is 2066.599 m according to my calculations and 2.5 - 3 km according to the scale line depending on zoom.

Ideally, the info I'm looking for is something like "the scale is accurate in +- 1 km because it calculates the distance using X method". The official doc (http://dev.openlayers.org/docs/files/OpenLayers/Control/ScaleLine-js.html) doesn't mention accuracy.


texts: Saving an image with PyQgis while adding a scalebar on it

texts:
I'm basically making a script to save an image in .png from Qgis with a script and I would like to add a scale bar on it.
I'm using this code :

    size = self.iface.mapCanvas().size()
    image = QImage(size, QImage.Format_RGB32)

    painter = QPainter(image)
    painter.setRenderHint(QPainter.Antialiasing)
    settings = self.iface.mapCanvas().mapSettings()

    layers = settings.layers()
    settings.setLayers(layers)

    job = QgsMapRendererCustomPainterJob(settings, painter)
    job.renderSynchronously()
    painter.end()

I see that it's possible to do it with a Map Composer, but I had other problems while saving from Map Composer (mainly with the quality of the image output) : https://docs.qgis.org/2.14/en/docs/pyqgis_developer_cookbook/composer.html

Is that possible to add a scale bar on my script ?


texts: SHP file of GDAM polygon of my country to my ASSETS not match the map base borders in GEE

texts:
I upload GDAM polygon of my country as shp file to my ASSETS but as layer this borders not match the map base of GEE. Other source Geofabrik OSM shp files are exactly the same as on GEE, but they are not completed for North Macedonia. I try to save GDAM shp with WSG84 EPSG:3857 and no change.
On picture is GDAM shape over OSM and over basic map on GEE.
Any advice?


texts: Defining distance in minutes from motorway in QGIS

texts: First I have downloaded a shapefile with Austria's transport system.Then I created a layer which now only consists of Austria's motorways. I would like to create a sort of buffer around the motorway (drive time isochrone). The buffer is meant to be in a distance of 10 minutes around the motorway. In the screenshot below one can see that the max speed for each road is given. I am not sure if there is still need for other Information besides the max speed. Is there a chance of creating the things i have in my mind with only this Information given?


texts: Openlayers Overview and OpenLayers versions

texts:
I have installed OpenLayers and the Openlayers Overview plug-ins but I am not sure if they are working. I can create new 'open' layers from the Openlayers plug-in drop down, but when I try to add a map from the 'Openlayers Overview' window it gives an error stating that 'The version  of Openlayers Plugin is 0.93, need be version 0.42!'.

Does this 'Overview Add map' do the same thing as the plugin menu options, in which case I guess there is no problem, I just won't use it or should it something else?

Also, the layers don't have a properties dialog. Is that right?


texts: Multiple selecting and deleting in ArcCatalog?

texts:
I have so many files to delete, to do that I am using ArcCatalog but I can do only one by one.

Is there any way to delete many files at the same time. Also, I would like to know how to select more than files in ArcCatalog as well.


texts: Check Geometry Type of a Shapefile Using GeoTools?

texts: I want to get the type (point, line, or polygon) of shapefile using GeoTools but so far I haven't found any code.  Any other Java library also works for me.


texts: Using if statement in R script within QGIS

texts:
I am creating an R script in QGIS and I am facing an issue with an if statement. First, I declare the input from the user

##Which_S2_bands_do_you_want_to_use=selection 2 3 4 5 6 7 8 11 12;2 3 4 8 11 12

Later I use it in the following condition:

if (Which_S2_bands_do_you_want_to_use == "2 3 4 5 6 7 8 11 12") {
        bands<-c("B((0[2348]_10m)|(((0[567])|(1[12])|(8A))_20m)).jp2$")}

And use it in the following code:

  S2<-list.files(S2, recursive = TRUE, full.names = TRUE, pattern=bands)

When I run it, I get the following error:

Error in list.files(S2, recursive = TRUE, full.names = TRUE, pattern = bands) :
object 'bands' not found
Execution halted

It seems the if statment is not been run. Any idea what could be the reason?

Maybe another approach to use the input from the user as a condition to change a variable in the script?


texts: Pie charts in QGIS: attributes generated with Expressions not working when using aggregate function

texts:
I have a problem to add values generated by expressions with aggregate functions to pie charts (see the screenshot below for better understanding). This question is a follow up question to this one here. In QGIS (using 3.16), you can add pie charts from a vector layer / Properties. This works perfectly as long as you select existing attribute values for the pie "slices" (circle sectors): see the two selected attributes (red rectangle) on the screenshot below - they render as expected on the map canvas.
Further, the documentation also mentions that

Generated attributes with Expressions can also be used.

However, when I try to add such generated attributes with expressions, consisting of an aggregate function (sum) on the basis of all selected features (is_selected), no output is generated: see the entries in the blue box on the screenshot. I tried other aggretage expressions like aggregate, with no success. Using expressions without the aggregate and is_selected criteria, it works as expected.
The expression like sum ( "data_protestant", filter:=is_selected () ) is correct and works in other contexts: I can use it to generate pie charts with geometry generator using the expression wedge_buffer - however, this is cumbersome if you have several values the visualize. And the built in pie charts creator clearly offers the possibilities to use expression based values.
Thus the question is: Is the direct rendering of such pie charts (expression based, aggregating over several features under the condition they are selected) a limitation of QGIS or can this be overcome somehow?
This here is a similar question, however the answer did not resolve the problem.


texts:
I downloaded a GPX file from a tourism web site. However, when I tried to import this file to Google Earth or Garmin's Basecamp, they both complained that the file format was wrong.
I opened the GPX file and noticed that the latitude and longitude numbers where quite high : an exemple is for instance : <rtept lat="6923211.8424401991" lon="367940.01590696711">
I looked at the file header and noticed that there was some kind of translation rule :
<extensions>
    ....
    <ogr:XW84>-1.56894</ogr:XW84>
    <ogr:YW84>49.322</ogr:YW84>
    <ogr:XL93>367940.01591</ogr:XL93>
    <ogr:YL93>6923211.8424399998</ogr:YL93>
  </extensions>

I think that this says that a latitude of 6923211.8424399998 should correspond to 49.322
Can somebody confirm that ?
Do you know if there is any tool to correct such format ?
Here's a link to the complete file.


texts: Converting degrees to X,Y coordinates with OpenLayers?

texts:
I have been tasked to center an openlayer map using the coordinates from an input.

My problem is that the coordinates are given in the degrees form.

For instance : 16° 58' 22" N 60° 17' 35" E

For now I achieved to center the map with X,Y coordinates (lat long I guess ?)

const coords = proj.fromLonLat([xCoord, yCoord]);
map.getView().animate({ center: coords, zoom: 10 });

But what I need is to transform the above degrees coordinates to an X,Y format for the code to work.

(I'm using openlayers v4.6.5)

This question has been marked as duplicate of this one but I can't get the answer to work properly.

For the following example coordonates 48 51' 14" N 2 20' 48" E
I called the getDecimalDeg() method (provided in the answer) with the following parameters : 48, 51, 14 and a second time with 2, 20, 48

I got : 480.850.0038888888888888888 and 9793:55 20.33333333333333330.013333333333333334

I guess I'm using it wrong.


texts: How to set layer symbol and layer name side by side in QGIS Composer?

texts:
When I add a legend in QGIS Composer, the layer symbols are set below each layer name. Is there any way to have them side by side?

This is what I get:

"LayerName"   [layer name]
~         [layer symbol]

This is what I need:

~ "LayerName"


texts: Help! I am downloading waypoints collected with our Garmin Montana 680t Software version 280. That works fine. Get a table of waypoints collected. The user also took photos which are geotagged with the waypoints. In order to link the photos with the waypoints, I need to convert images to waypoints in the DNRGPS software using Edit> Add Image Links. When I go through this process, the identifier field is replaced with the name of the photo starting with the default DSC00000. The comment field disappears altogether and this is the information I need for what was at that location (another identifier). So I have to export 2 different files to shapefiles to use in ArcMap. Does anyone out there have any idea why this happens and what I can do to rectify?


texts: GeoPDF export in QGIS requires GDAL 3, how do I get GDAL 3?

texts:
I would like to use the new GeoPDF export feature of QGIS 3.10. The changelog states that it requires GDAL 3. I have QGIS 3.10 installed on OS X from the binary installer. The export GeoPDF feature is unavailable to me because GDAL 2 is enclosed in the program. Under about QGIS I find this:

Compiled against GDAL/OGR 2.4.1
Running against GDAL/OGR 2.4.1

Here are my questions:

Am I right in my assumption that I cannot 'just' install GDAL 3 to solve this issue, but need to compile QGIS from source against GDAL 3?
If assumption 1 is incorrect: how can install GDAL 3 to use the GeoPDF function?
If assumption 1 is correct: are precompiled installers available with GDAL 3?

I would be surprised if an end-user feature like this would require the program to be self-compiled.

Edit: this question is similar to this one for Windows however, for Windows an advanced installer is available. For OS X it is not.

For now I will:

Try QGIS on Linux (don't have Windows).
Inform on the development of the Mac version or add feature request.
Hope I can have this feature on the Mac soon as well.

If all fails I'll borrow a PC with ArcMap.


texts:
Does anyone know if there are any issues for QGIS in accessing MS SQL Server 2016?
I am currently using 2008, but there are a number of new functions within 2016 that are attractive for the upgrade.


texts: Editing attribute columns to switch properties from text to integer?

texts: Instead of adding an entire new column in my attributes table, is there a way I can edit the properties of an existing (unnecessary) column and change the type from integer to text?


texts: How can I copy archive data from one server to another?

texts:
I have an SDE 10.0 server and an SDE 10.1 server and I am trying to migrate the layers from 10.0 to 10.1. (On a side note the 10.0 server is running SQL Server 2005 and the 10.1 server is running SQL Server 2012.)

Some of the layers have archiving turned on and there seems to be no way to transfer the layer as well as all of the archival information. Is there a way to copy the layer and also its archive?

Alternately, can I just turn on archiving on the new server after transferring the layer, then copy the rows from the old server into the newly created archival table with a SQL tool (Import/Export Data)?


texts: Dynamically adding features with cluster strategy

texts:
I am trying to dynamically add features to a vector layer, and use cluster strategy to group those features accordingly.

I am able to do this for the first set of data that gets loaded on the map.

It looks like below code works only the first time -
         vectorLayer.addFeatures(vectordArray);

when i try to call the same addFeatures again with a new set of vectors, I dont see the new vectors. There are no errors either.

Anyone faced a situation like this? I am looking to dynamically add vectors and use CluseterStrategy.


texts: Labelling pie chart legend entries by alias in ArcGIS 10.1

texts:
When generating pie charts in ArcGIS, it appears that the legend entries generated for the different slices are named based on field name, rather than the field alias.

I've done some research, and found a few historical examples where the user could at least change the labels in the table of contents (example), although it would revert as soon as the layer's properties were changed.  In 10.1, it seems that even this ability has been lost.

Is there any way to change the labels for pie charts in the legend, or to force the legend to use the field alias?


texts: Why do I get correct area and intersect area when use wrong projection?

texts:
I need to calculate areas and intersection areas for polygons (some real geographical objects like lake, city, country, e.t.c.). Polygons located in California, New Zealand, Russia.Anadyr, Sweden
All polygons are in WGS84.
What I did using GeoTool java API:

Project all polygons using EPSG:3488, EPSG:NAD83(NSRS2007) /
California Albers and calculated areas and overlap areas.
Did the same using World_Mollweide and World_Eckert_IV
Picked "local specific projections" for polygons from California, New Zealand, e.t.c.

I assume that #3 is the most accurate result, since I pick projection which covers polygon area
Result:
'#2 showed the worst result comparing to #3
'#1 and #3 areas and intersection areas difference is less than 0.1%
Why?
I pick absolutely wrong projection EPSG:3488 (California) for polygons from Sweden and get pretty the same areas and intersection areas?
UPD:
Looks like I didn't explain my confusion correctly. Here is sample output with explanation
#area_from_new_zealand_1
EPSG_27200 area[11733479] CRS[World_Mollweide] area[11736023] diff[2544] [0.0%]
EPSG_27200 area[11733479] CRS[World_Eckert_IV] area[11736033] diff[2554] [0.0%]
EPSG_27200 area[11733479] CRS[EPSG:NAD83(NSRS2007) / California Albers] area[11736034] diff[2555] [0.0%]

#area_from_new_zealand_2
EPSG_27200 area[2952725]  CRS[World_Mollweide] area[2953281] diff[556] [0.0%]
EPSG_27200 area[2952725]  CRS[World_Eckert_IV] area[2953342] diff[617] [0.0%]
EPSG_27200 area[2952725]  CRS[EPSG:NAD83(NSRS2007) / California Albers] area[2953467] diff[743] [0.0%]

#intersection_area_between_two_new_zealand_areas
EPSG_27200 intersection area[1001857] CRS[World_Mollweide]                          area[1002082] diff[225] [0.0%]
EPSG_27200 intersection area[1001857] CRS[World_Eckert_IV]                          area[1002082] diff[225] [0.0%]
EPSG_27200 intersection area[1001857] CRS[EPSG:NAD83(NSRS2007) / California Albers] area[1002096] diff[239] [0.0%]

#area_from_alaska_1
EPSG_3338 area[56278347]    CRS[World_Mollweide] area[56041510] diff[236837] [0.4%]
EPSG_3338 area[56278347]    CRS[World_Eckert_IV] area[56041585] diff[236763] [0.4%]
EPSG_3338 area[56278347]    CRS[EPSG:NAD83(NSRS2007) / California Albers] area[56278426] diff[79] [0.0%]

#area_from_alaska_2
EPSG_3338 area[17564799282] CRS[World_Mollweide] area[17486015889] diff[78783393] [0.4%]
EPSG_3338 area[17564799282] CRS[World_Eckert_IV] area[17486869816] diff[77929466] [0.4%]
EPSG_3338 area[17564799282] CRS[EPSG:NAD83(NSRS2007) / California Albers] area[17566197286] diff[1398004] [0.0%]

 #intersection_area_between_two_alaska_areas
EPSG_3338 intersection area[43808167] CRS[World_Mollweide] area[45066901] diff[1258734] [2.8%]
EPSG_3338 intersection area[43808167] CRS[World_Eckert_IV] area[45163183] diff[1355016] [3.0%]
EPSG_3338 intersection area[43808167] CRS[EPSG:NAD83(NSRS2007) / California Albers] area[43885182] diff[77015] [0.2%]

My confusion is:
EPSG:3488 designed to be used in California
I pick "wrong" projection EPSG:3488 for Alaska, New Zealand areas and see that resulting calculations don't differ "significantly" from correct projections. EPSG:3488 even performs better than Mollweide, Eckert_IV projections designed to be used around the world.


texts: How do you create Revision Table in Arcmap for OSM contribution?

texts:
I have created the following detailed file for upload to OSM in ArcMap but don't know about creating a OSM Revision Table to be able to put the data into OSM.

At 1:10,000

At 1:1000

The tool

What do I do? I searched but found nothing about how to create the revision file.

http://esriosmeditor.codeplex.com/wikipage?title=Tools&referringTitle=Documentation...mentions

  The second parameter is the local revision table storing the changeset
  information that is submitted to the OSM server. This table also
  stores the upload status, i.e. if the data upload was successful or if
  the upload resulted in conflicts. The conflict editor will use the
  same table to guide you through the conflict resolution. The table
  name has the schema of osmrevision where  is the name of
  the feature dataset where the OSM data resides.

but doesn't mention what fields this table needs, the format etc.

When I try the OSM Download Tool I get


texts:
I am trying to find the shortest path between buildings and any electric plants on a network.

The thing is, when I use the following query, I overload my RAM (network_nodes and network_edges each cointains a bit more than 120 000 rows and 5 columns) and therefore do not get any results.

SELECT start_vid as building, end_vid as plant, agg_cost as junction_length from pgr_dijkstraCost
(
'SELECT gid as id, i as source , j as target , length as cost from   network_edges'
, array(SELECT id from "network_nodes" where buildingID IS NOT NULL)
, array(SELECT id from "network_nodes" where plant_name IS NOT NULL)
, false
)

To avoid that, I thought of trying to use pgr_dijkstraCost() One-to-Many instead. In order to do it, I created a table with 2 columns and as many rows as there are buildings nodes in my network_nodes table (around 85 000). The first column contains the ID of each building and the second contains an array of ids which are the closest plants to the corresponding building.

Now, my problem is : how to write a query which repeat pgr_dijkstraCost() One-to-Many using each time as start_vid the first column of this table (integer) and as end_vid the second columns (array of integers) ?

I suspect that I'll have to do some kind of loop but I do not know how to do that in SQL.


texts: In QGIS, how does one find all US counties which lie within a fixed distance from a target county

texts:
Given the following vector input into QGIS...
US Counties
...is it possible, to generate a report detailing, for each county, all neighbors within 50 miles of that county. In short, the final output would have N number of outputs with each output representing the original county as modified by a distance band and detailing the overlapped counties.

I think I see a way to approach the problem with a database spatial query If I convert the input data to points and implement a per point solution. Maybe that is the correct solution but it seems inaccurate as given a large enough gap between adjacent points of our target county, we may miss a neighbor which approaches between those points.


texts: Dot density map "slicing" dots when exported

texts:
I'm running ArcMap 10.3 and every time I export a dot density map it "slices" long, horizontal gashes in sections of the dots (see screenshot below). The cuts are never vertical, always horizontal, and it doesn't matter if I export as a PDF or image file like PNG or JPEG.

In the example screenshot, I have a dot density layer representing population data on top of a choropleth map showing demographic indices. The only other layers in my TOC are a shaded relief raster and a freeway/roads line layer.

Right now I'm simply manually exporting the map as a PNG at a resolution of 400 dpi, but the slices also show up at a lower resolution (such as 300 or 250). It's worth noting that the slices still show up even when I export the map using the arcpy.mapping.ExportToPng tool in my script.

Wondering if anyone else has noticed this problem? Not sure if it's just an exporting bug or something else.


texts: Plotting contour map of modeled pollution data in QGIS?

texts:
I'm a complete beginner to GIS, and I'm trying to create a filled contour map of pollutant concentration data for an area of NSW, Australia.

I've calculated the modeled concentrations of the pollutant as a value for each xy co-ordinate and need some help getting started plotting this as a contour map?


texts: ArcGIS javascript api: Convert a map service with projected coords to web mercator

texts:
here I see a question about tiled map services.  But what about dealing with an arcgis server map service in a projected coordinate system like say 2278, how would I convert this to web mercator so that it will fall correctly in my map?

What are the choices using dynamic layer or other?

I don't have any control over the map service.  Thanks.


texts: Display M-values in a MultiPointM from a PostGIS database as a label in QGIS

texts:
I have MultiPointMs in a PostGIS (v. 2.2.2) database.
I would like to know if it is possible in QGIS (v. 2.14.1) to label each point with its M-Value.
For instance:

MultiPointM(1 1 0, 2 2 1, 3 3 2);
Coordinates obviously made up

And in QGIS I would like to have the point (1 1) to have a label with "0".
The point(2 2) would be labeled "1" and the point (3 3) would be labeled "2".
I would prefer having a script-free solution (as the final user is not myself).

A friend suggested this : How can I label both coordinates in mmQGIS?
Which gives some hints on the direction to go... The $x and $y work
perfectly... But it doesn't seem to like $m (and I wouldn't know where
to search for to know what would be the equivalent).


texts: Applying unique processing extent to each iteration of tool in ArcMAP Model Builder

texts:
I am trying to create an automated ArcMAP (10.1) tool that will read in a file (ASCII most likely) of maximum and minimum lat/long values along with some additional variables. The tool will iterate over each row of this file, and run the processes on the unique variables each time.

A preliminary tool that only runs on one set of variables works perfectly.

However, the 'Extent' Data Type can't accept variables into it (in any way I have found out at least). So I have no idea how to define a unique processing extent for each run of the tool. Below is an example version of what the model might look like (I have tried several different ways, this is just one).

If it matters, I have never tried working with python scripts in ArcMAP, so I'm hesitant to go down that route unless it's fairly straight-forward. Unless there is no other option. The base DEM I will be using will cover the entire globe, so it really is necessary that I only process a small part of it each time I run the tool.


texts: Display symbology categories at certain scales with ArcGis Pro

texts:
We can since a long time display layers at certain scales but
is it possible to display a specific symbology categories at certain scales ?

For example if I have the following feature layer:

Can I display the tertiary road only until the scale 1:20000 ?

ArcGis Pro have an option "Enable scale-based sizing" but I'm only able to increase the size of one specific categorie at lower scale but i'm not able to prevent arcgis to display this categorie at lower scale.


texts: Exporting spatial data linkage from two inputs in ArcGIS?

texts:
I have a csv file containing of points (house locations) and shapefiles consisting of polygons (congressional districts)

How do I export a csv of the points and what district that point falls in?

Because I don't have any objective id for the points so I was unable to use join data feature in Exporting query to identify what polygon certain points fall in using ArcGIS Desktop?.

I am using ArcGIS Desktop 10.3.


texts:
I am trying to import gdal modules from osgeo installed in my win64x machine but I get the following error:

  Traceback (most recent call last):
    File "", line 1, in
      from osgeo import ogr
  ImportError: No module named osgeo

I am using the python shell from arcgis 10.1 that comes with python 2.7.2


texts: Counting and listing all shapefiles that falls inside layer polygon with ArcPy

texts:
I am new in ArcPy.

As you can see in the image below, I have:

One shapefile with 9 polygons (each polygon illustrates a community)
1200 layers. Each shapefile illustrates a GPS track, so they are point layers. I know there are 1200 layer points because I used the "arcpy.mapping.ListLayers" to figure it out.

For each polygon, I want ArcPy to tell me all the layer names that fall "inside", or "are completely contained" by each polygon.


texts: postgis ST_ConcaveHull from MULTIPOLYGON fails


texts: How to get the absolute raster row / column of rastercell when using ArcObject Raster Cursor?

texts:
When using the ArcObject RasterCursor (see snippet below), how do I get the absolute row / column of the current pixel value?

IRasterCursor rasterCursor = raster.CreateCursorEx(null);

do
{
    IPixelBlock3 pixelblock = rasterCursor.PixelBlock as IPixelBlock3;

    Array pixels = (Array)pixelblock.PixelData[rasterBandIndex];

    for (int i = 0; i < pixelblock.Width; i++)
    {
        for (int j = 0; j < pixelblock.Height; j++)
        {
            // How to get absolute row / column of this pixelvalue
            dynamic pixelValue = pixels.GetValue(i, j);
        }
    }
}
while (rasterCursor.Next());


texts: Distance between reference point to many Points along a line

texts:
I have a Line layer and a Point layer in Qgis. The lines snapped to the points. I need the distance between the reference Point (number 12 in the image) to the other points along the line.
I have tried this V.distance, couldn't help me.


texts: Apply and keep advanced labeling settings in ArcMap

texts:
I want to do a map in ArcMap and keep the labeling that I set in the layer properties for the map legend.
Somehow, when i click 'apply' the values of the advanced labeling appear in the legend, but when I click 'OK' afterwards they disappear again and all my settings are lost.

What is wrong?

EDIT:
Elaborating on @Martin's comment:
I tried to edit the legend entries (and my may be the mistake due to your point!?) via the symbology. I thought this worked, because when I set the interval size in the advanced labeling to 100, click OK in the Advanced labeling window and then Apply ('Übernehmen') in the Layer Properties window, the eight labels from 1000 to 1700 appear at the ticks in the legend. However, as soon as I click OK in the layer properties window, they disappear again and all my settings are lost. Maybe that is the strange thing about it.


texts: Overwriting files in ArcSDE with ArcPy gives ERROR 000601

texts:
I'm having a great deal of difficulty overwriting files when using ArcSDE. The usual 'arcpy.env.overwriteOutput' method doesn't seem to work.

When I run this script:

import arcpy
from arcpy import env

env.workspace = "L:\\DRAINAGE\\D_GROUP6\\STAFF\\CJE\\EDRN_3.sde"

arcpy.env.overwriteOutput = True

selected_link = "L:\\DRAINAGE\\D_GROUP6\\STAFF\\CJE\\EDRN_3.sde\\VALLAHC.EDRN_LINK6"
q_fid = 39
edrn_node = "L:\\DRAINAGE\\D_GROUP6\\STAFF\\CJE\\EDRN_3.sde\\VALLAHC.EDRN_NODE6"
join1 = "L:\\DRAINAGE\\D_GROUP6\\STAFF\\CJE\\EDRN_3.sde\\VALLAHC.join39"

arcpy.MakeFeatureLayer_management(selected_link, "feature_layer10")
arcpy.SelectLayerByAttribute_management("feature_layer10", "NEW_SELECTION", '\"OBJECTID_1\" = {0}'.format(q_fid))

arcpy.SpatialJoin_analysis("feature_layer10", edrn_node, join1, "JOIN_ONE_TO_MANY", "#", "#", "CONTAINS")

I get the following error:

Runtime error <class 'arcgisscripting.ExecuteError'>: ERROR 000601: Cannot delete L:\DRAINAGE\D_GROUP6\STAFF\CJE\EDRN_3.sde\VALLAHC.join39.  May be locked by another application. Failed to execute (SpatialJoin).

I have attempted to get around this in a number of ways:

versioning/not versioning the data.
Checking/unchecking the overwrite geoprocessing operations box in geoprocessing options.
Using 'arcpy.Delete_management()' to delete the join at the end of the script. But ArcSDE won't let me and outputs the same 000601 error.

Does anybody have any ideas on how to get around this?

My script works fine when I use shapefiles or a regular geodatabase. I'm thinking there must be something I've missed related to ArcSDE that I'm not aware of as I don't have much experience with databases.


texts: Measuring Sinuosity of Walking Paths using QGIS?

texts:
I'm currently working on a project that examines the movement of consumers through different kinds of retail spaces in Central Asia – namely, bazaars versus supermarkets. I've equipped dozens of shoppers with GPS trackers that record their movements throughout these different retail spaces, and I'm now working on analyzing the resulting data.

I'm primarily interested in comparing overall walking distances and walking speeds – easy enough. But I've also noticed bazaar paths are highly irregular, with shoppers weaving this way and that through the crowded stalls, whereas the supermarket tracks are highly straight, efficient, and consistent between shoppers. So, I'd like to quantify the sinuosity or efficiency of these walking paths – i.e., how frequently the shoppers turn and/or backtrack.

I'm not sure how to approach this analysis – much less accomplish it in QGIS.

Should I define this "spatial efficiency" index as the ratio of overall distance traveled to overall degrees of rotation, or something?

Or should I use the simple river sinuosity calculation (direct distance vs. curvilinear length) – and if so, how should I determine the critical waypoints?


texts: How to compare a theoretical GPS route and a real one?

texts:
In an iPhone app, I'm trying to generate a roadbook automatically based on a real recorded route. I have figured out an algorithm to use a directions service to do it, but my algorithm needs to be able to tell if the real route follows the one calculated by the directions service for any particular segment of it.

My real route is a series of GPS positions that are pretty regularly spaced.
The theoretical route returned by the directions service gives my less points and they are not spaced regularly since it basically puts a point whenever there is a significant heading change like turns and so on.

And what I would like to do is to make sure both routes are somewhat similar, because if there is a local difference, then it probably means the directions service took a shortcut and I need to recalculate a road for this segment.

I found this question that mentions the Hausdorff distance (How to calculate 'off track' error from GPS data?) and it looks like what I'm trying to do, but implementing an algorithm for it in Swift based on pure mathematical notation in the Wikipedia article is way out there for me.

I have tried to implement a naive version to calculate the maximum of the minimum distance between 2 points (see below) but it doesn't work very well because the points on both paths are not aligned so sometimes I have a big distance because a point on the real route is just too far away from the closest point on the theoretical one.

func hausdorffDistance(from:[CLLocationCoordinate2D], to:[CLLocationCoordinate2D]) -> Double {
   var maxDistance = 0.0
   for point in from {
       let pointLocation = CLLocation(latitude: point.latitude, longitude: point.longitude)
       let sortedTo = to.sorted(by: { left, right -> Bool in
           let leftLocation = CLLocation(latitude: left.latitude, longitude: left.longitude)
           let rightLocation = CLLocation(latitude: right.latitude, longitude: right.longitude)
           return leftLocation.distance(from: pointLocation) < rightLocation.distance(from: pointLocation)
       })
       let closestPointLocation = CLLocation(latitude: sortedTo.first!.latitude, longitude: sortedTo.first!.longitude)
       let distanceToClosestPoint = pointLocation.distance(from: closestPointLocation)
       if distanceToClosestPoint > maxDistance {
           maxDistance = distanceToClosestPoint
       }
   }
   return maxDistance
}

Any idea how I could proceed?


texts: Automatic file name for composer exports

texts:
I have to produce pdf maps everyday (corvid-19 stuff) for my organisation and want to create the file name based on the composer name preceeded by the current date.  The name of the composer is the default but I have to type in the date in the form yyyyMMdd but this is a pain.  So when I press "Export as pdf" button my file name will show "yyyyMMdd-LAyout name.pdf".
I've tried project variables, custom functions etc. but can't get anything to work!

Anyone, any ideas? I'm using 3.12.1-Bucuresti on Windows 10.


texts: How to slow the zoom transition speed in Leaflet?

texts:
I'd like to slow the speed of the zoom transition in Leaflet, but I haven't figured out how.

I'm using map.setView() a few times to change the zoom level.

I've tried the solution suggested here:
Setting a slower zoom speed

Which is to increase the transition time in CSS for

.leaflet-zoom-anim .leaflet-zoom-animated {
    transition: 2s;
}

But, it doesn't seem to work. It starts to slow down, but then just jumps to the zoom level, so it doesn't look smooth at all.

The pan animation, and duration works great.

I've tried all these options:

map.setView([lat, long], 14, {
    pan: {
        animate: true,
        duration: 1.5
    },
    zoom: {
        animate: true
    }
});

and

var map = L.map("map", {
    center: [45.2403, -123.8512],
    zoom: 12,
    fadeAnimation: true,
    zoomAnimation: true
});

But nothing seems to slow or smooth the zoom transition.


texts: Assigning the modal average based on points within a polygon

texts:
I have two layers:

Survey Points
Habitat Multi polygon

The Survey points contain the field sensitivity

I want to assign a number to each polygon based on the modal average sensitivity.

This thread has a similar theme

STIntersects - Update field to give me count of points in polygon

There is no way to do this inside QGIS so I've had to resort to using PostgreSQL. I've found an example to experiment with but it do not know how to incorporate the spatial aspect between datasets.

http://www.udel.edu/evelyn/SQL-Class3/SQL3_Stat.html

In summary I'm looking for help on how to script SQL to find the modal point value per polygon and attach it to the polygons attribute table.


texts: Using OpenLayers - seeking example for saving feature to database using Ajax

texts:
I'm trying to save a point to the database using a balloon that has an input field and a save button. The button triggers the code and sends the feature info to the PHP page.
It works, but the save action refreshes the page. I would like it to save in the background without refreshing the page, like Ajax.

Can someone advise how to do this ?

   var mySaveStrategy = new OpenLayers.Strategy.Save({ auto: true });
   var myProtocolHTTP = new OpenLayers.Protocol.HTTP({
            url: "proc/newPoint.php",
            params: {
                    m: 'r',
                    srs: 'EPSG:4326'
            },
            callback: function() { alert("User callback called!"); },
            format: new OpenLayers.Format.GeoJSON({
                    ignoreExtraDims: true,
                    internalProjection: map.baseLayer.projection
            })
    });

I edited the code to reflect my latest code (16.12.2012):

        function handler(request) {
            alert(request.responseXML);
                        }

        var myProtocolHTTP = OpenLayers.Request.POST({
            url: "api.php?action=insertPOI",
            data: "",
            params: {
                m: 'r',
                srs: 'EPSG:4326'
            },
            headers: {
                "Content-Type": "text/plain"
            },
            callback: handler
        });

    function formAction(action) {
        var feature = mytoolbar.selCtrl.openpopup.feature;
        switch (action) {
            case 'delete':
                feature.state = OpenLayers.State.DELETE;
                vectorLayer.strategies[0].save();
                break;
            case 'cancel':
                break;
            case 'save':
                var attrs = feature.attributes;
                attrs.name = $('fname').value;
                attrs.name_lat = $('fname_alt').value;
                attrs.type = $('ftype').value;
                attrs.notes = $('fnotes').value;

                if (feature.mystate) {
                    feature.mystate = null;
                    feature.state = OpenLayers.State.INSERT;
                    var ctrl = mytoolbar.controls[1];
                    ctrl.layer.addFeatures([feature]);
                    ctrl.featureAdded(feature);
                    ctrl.events.triggerEvent("featureadded", {'feature': feature});
                } else {
                    feature.state = OpenLayers.State.UPDATE;
                    vectorLayer.strategies[0].save();
                }
               break;

        }

    function getPopupForm(feature, context) {

        var gid      = feature.attributes.gid      || '';
        var name     = feature.attributes.name     || '';
        var name_alt = feature.attributes.name_alt || '';
        var notes    = feature.attributes.notes    || '';
        var html = "
" +
            "<form action='' name='baloon'>" +
               "<input id='fgid' type='hidden' value='" +
                  gid + "' />" +
               "<button type='submit'  id='baloon-save' value='Save' onClick=\"formAction('save')\" />"+
               "	" +
               "
"+
               "<button type='submit' id='baloon-cancel' value='Cancel' onClick=\"formAction('cancel')\" />&nbsp;" +
               "

"+
               "<button type='submit' id='baloon-delete' value='Delete' onClick=\"formAction('delete')\" />
"+
               "

</form>" +
           "
";

        return html;
    }

'


texts:
the goal of my program is to retrieve data from a field one after another at the end to compare with other values of another table
when I execute the code that allows me to point to the values of the first field I get a string that starts with u '
What should I do to remove the u'

DataSource= Shapefile

cursor=arcpy.da.SearchCursor("Feuil1$","LIAISON_NOM")
for row in cursor:
print(str(row))

i get this

>>>>(u'ROAD 1',)
    (u'ROAD 2',)


texts: GetRasterBand()Method for GDAL in Python

texts:
I'm working on getting the NDVI value of Landsat imagery 2018 using band 4 and 5 but I have encountered a slight problem.

import osgeo
from osgeo import gdal#gdal, data conversion tool
from matplotlib import pyplot as plt
import sys
import numpy as np#for reading numbers in python
b8a= gdal.Open('C:\\Users\\zack\\GIS PROGRAMMINGWORKS\\LC08_L1TP_168061_20171228_20180103_01_T1_B5.tiff')
b4= gdal.Open('C:\\Users\zack\\GIS PROGRAMMINGWORKS\\LC08_L1TP_168061_20171228_20180103_01_T1_B4.tiff')
d8a = b8a.GetRasterBand(1)
d4 = b4.GetRasterBand(1)

Error returned:

AttributeError                            Traceback (most recent call last)
<ipython-input-22-cb1cfc43017d> in <module>
 ----> 1 d8a = b8a.GetRasterBand(1)
  2 d4 = b4.GetRasterBand(1)
AttributeError: 'NoneType' object has no attribute 'GetRasterBand'


texts: How can I get arcpy.TraceGeometricNetwork_management results?

texts:
I have geometric network and try to launch Python script with tracing

gdb = "C:/Users/my_user/Documents/ArcGIS/gn/gn.gdb"
gdb += "/gn/"
geometric_net = gdb + "gn_Net"
flags = gdb + "flags"
result_name = "gn_Net2"

mxd = arcpy.mapping.MapDocument("CURRENT")
df = arcpy.mapping.ListDataFrames(mxd)[0]

result = arcpy.TraceGeometricNetwork_management(geometric_net,
                                       result_name,
                                       flags,
                                       'FIND_CONNECTED'
                                       )

gl = arcpy.mapping.Layer(result_name)
arcpy.mapping.AddLayer(df, gl, "BOTTOM")

This script do nothing (as I see in ArcMap) - no new group layer in ArcMap, no new layers/feature classes/other instances in current geodatabase.

result variable have some value after last row execution, but I can't determine its type - debugger says thatresult.__class__ returns <class 'arcpy.arcobjects.arcobjects.Result'>, and result returns {geoprocessing server result object} gn_Net2

From documentation I presumed that result is group layer that I must to add to map. But I can't find Arcpy method which can add group layer to map - only layer adding to map or existing group layer (arcpy.mapping.AddLayer).

What really arcpy.TraceGeometricNetwork_management returns? Did I understood correctly that it's group layer? And it's ready for add to ArcMap via arcpy.mapping.AddLayer method?
If no, where is tracing result? How can I see it? What is the arcpy.TraceGeometricNetwork_management result?

Edits: I added

ArcMap document (mxd) and data frame (df) getting
Group layer adding

to script.
I launches my script from ArcMap via toolbox-script object.

Result of my changes is:
Group layer adds to ArcMap document, but without found path (no new layer, selectiong or highlighting). Only geometry network and its layers was added.
When I run geometry tracing manually with same conditions (Utility Network Analyst - Solve), found way highlights at map.


texts: Calculate Values from a raster surrounding a list of addresses

texts:
I have a list of addresses and I have geocoded into a shapefile in ArcMap 10.4. I've also made a raster that stores NDVI data (with water masked out as it will skew the values we are hoping to obtain later) within the area that the addresses lie. What I'm trying to do is take the list of addresses and calculate the mean NDVI value within a 250 m radius of each address.

So far, we have used buffer analysis to create a shapefile containing 250 m radius at each address and used the zonal statistics as table tool. The issue, however, is that some of the radii overlap - which the zonal stats tool doesn't know how to handle.

Additionally, we have tried using ModelBuilder to try to iterate over each polygon and calculate zonal stats individually, but when we do this the file will continuously keep overwriting itself with each input.

Finally, we found the spatial analyst supplemental tools that ESRI created, which includes the zonal statistics as table 2. This tool is supposed to have the same functionality as the original tool, but include a fix for overlapping polygons. We ran this tool and the output was exactly the same as when we ran the original tool - did NOT account for overlapping areas.

Last bit of info, python within our ArcMap does not work. We are unsure why at the moment. So, creating a script to do this is not option at the moment.

Does anyone have any suggestions or tips on how we could go about getting our data?


texts:
Is it possible to use the .NET Entity Framework with PostGIS?

I've reviewed a few options, including hookup it up with Npgsql using http://npgsql.com/index.php/2009/08/how-to-set-up-entity-framework-npgsql-part-1/, but it appears that the data types in PostGIS are unsupported in Entity Framework.

The closest I can find to a comprehensive statement of what's supported is with Devart's dotConnect for PostgreSQL EF data type mapping at http://www.devart.com/dotconnect/postgresql/docs/DataTypeMapping.html. It appears to expose only base PostgreSQL types and neither PostGIS or presumably common addons like ltree.


texts:
I am using OpenLayers 2.13.1. I could get my clustering work for the vector layer but when I try adding a styleMap to the vector layer, neither the clusters nor the individual points are visible on the map. Below is my code.

function PlotExposuresOnMap(selectedExposures)
     {
         if (exposureLayer)
             map.removeLayer(exposureLayer);

     var context = {
         symbol: function (feature) {
             if (feature.cluster) {
                 return '../icons/cluster.png';
             }
             else {
                 return '../icons/marker.png';
             }
         }
     }

     var template = {
         externalGraphic: "${symbol}"
     }
     var style = new OpenLayers.Style(template, { context: context });

     var clusterStrategy = new OpenLayers.Strategy.Cluster({ distance: 20 });

     var geojson_format = new OpenLayers.Format.GeoJSON({
         'internalProjection': new OpenLayers.Projection("EPSG:900913"),
         'externalProjection': new OpenLayers.Projection("EPSG:4326")
     });
     exposureLayer = new OpenLayers.Layer.Vector("Exposures", { strategies: [clusterStrategy] /*, styleMap: new OpenLayers.StyleMap(style)*/});

     map.addLayer(exposureLayer);

     exposureLayer.addFeatures(geojson_format.read(selectedExposures));

 }

I wanted to display different icons for the markers based on whether it is a clsuter or an individual marker. I am really not sure what is missing.


texts: How to get screen coordinate from geographic one?

texts:
I need to know, how can I convert some geographic coordinates to the screen one.
I want to remind you about desktop coordinates system, how it look like:

I've developed soft, which is loading map fragments into matrix ( 10x10 ) in the fullscreen mode. The map images are loading perfectly.
But I have a task to draw some real elements on map, which do have longitude, latitude values.
So... I first of all think to convert geographic data to the web mecator values by using such a formulas:
    internal Models.PointGeospatial ConvertFromMercatorToSpherical(double mercatorX, double mercatorY)
    {
        if (Math.Abs(mercatorX) < 180 && Math.Abs(mercatorY) < 90) return null;
        if ((Math.Abs(mercatorX) > 20037508.3427892) || (Math.Abs(mercatorY) > 20037508.3427892)) return null;

        double subitemFirst = mercatorX / 6378137.0;
        double subitemSecond = subitemFirst * 57.295779513082323;
        double subitemThird = Math.Floor((double)((subitemSecond + 180.0) / 360.0));
        double subitemFourth = subitemSecond - (subitemThird * 360.0);
        double subitemFifth = 1.5707963267948966 - (2.0 * Math.Atan(Math.Exp((-1.0 * mercatorY) / 6378137.0)));

        var model = new Models.PointGeospatial()
        {
            Longitude = subitemFourth,
            Latitude = subitemFifth * 57.295779513082323
        };

        return model;
    }

    internal Models.PointMercator ConvertFromSphericalToMercator(double longitude, double latitude)
    {
        if ((Math.Abs(longitude) > 180 || Math.Abs(latitude) > 90)) return null;

        double subitemFirst = longitude * 0.017453292519943295;
        double x = 6378137.0 * subitemFirst;
        double subitemSecond = latitude * 0.017453292519943295;

        var model = new Models.PointMercator()
        {
            X = 6378137.0 * subitemFirst,
            Y = 3189068.5 * Math.Log((1.0 + Math.Sin(subitemSecond)) / (1.0 - Math.Sin(subitemSecond)))
        };

        return model;
    }

Then I'm getting such values after conversion:

Building Id: 302
Geographic data: { Longitude: 37.5303305475 ; Latitude: 55.7481343179 }
Mercator data: { X: 4177857.2858509403; Y: 7508438.839224923 }

Really... I don't even imagine how could I place correctly the item on map just only knowing its longitude, latitude values and its mercator equivalent with having 640x480 opened window for e.g., on which I must place some pushpins...
The one detail I've foud, that first three digits from mercator values are not changing for the one city (and I understand they must not ), so the task to set a correct place for the item in the map with exact zoom level is rather diffcult for me and I don't understand how to solve this problem.
I'm listening for you suggestions, thanks!


texts: Is it possible to symbolise via a related table?

texts:
I´m using ArcGIS for Desktop with an Advanced license level.

I´ve got a feature class with several points in it; each point contains information for location and Name.
E.g. Location [Berlin] , Name [Deponiestrasse 1]
Each Point describes a disposal site.

Now I want to connect the point information with a Table, which contains several information’s for each point. The Table contains information for Name e.g. [Deponiestrasse 1] and several types of disposal.  This types can be [concrete] ; [gravel] ; [brick] ; [toxic waste] and so on, up to 10 for one point.

E.g. point 1: Name [Deponiestrasse 1] and types are [concrete] and [gravel].
My Idea to connect them, was to use the RELATE function.

But now my Problems are starting. I want to draw around each point a simple “ symbol box”, a small box with different colors e.g.  concrete-blue, gravel-red, brick-green……. and so on. But via using the SYMBOLOGY function I can´t symbolise by using a related table.

So my questions are:

Is there any way to do this in this way? Or is there any possibility to write many types like [concrete] ; [gravel] ; [brick] into one attribute field of one point? Maybe you have another idea how to solve my problem.

I hope you understand my problem.
I used google etc. but couldn´t find a solution.


texts: Grass averaging raster not working with: Could not read georeferenced file error

texts:
I am trying to average geotiff using grass files using r.external and r.series. The geotiffs are georeferenced with wgs84. The command is below.

grass74 g.proj -c proj4="+proj=longlat +datum=WGS84 +no_defs" r.external input="/media/geodev/Seagate Expansion Drive/chlor_a/2017/chlor_a-modis-aqua-01_01_2017-18_30_00.tif" band=1 output=tmp153263774753364 --overwrite -o r.external "/media/geodev/Seagate Expansion Drive/chlor_a/2017/chlor_a-modis-aqua-01_03_2017-19_55_00.tif" band=1 output=tmp153263774753365 --overwrite -o g.region n=31.4648437503 s=17.4068080185 e=-80.4330407378 w=-98.0539218184 res=0.0142218572079 r.series input="tmp153263774753364,tmp153263774753365" method=average range="-10000000000,10000000000" output="output90e1839842f5448890b14b9d0c406f01" --overwrite g.region raster=output90e1839842f5448890b14b9d0c406f01 r.out.gdal --overwrite -c createopt="TFW=YES,COMPRESS=LZW" input=output90e1839842f5448890b14b9d0c406f01 output="/home/geodev/Desktop/ava2ee4.tif"

However, I am getting the error below.

Creating new GRASS GIS location/mapset...
ERROR: ERROR 4: g.proj: No such file or directory
ERROR: Could not read georeferenced file g.proj using either OGR nor GDAL

Exiting...

How could I fix this issue?


texts: Creating a Simple, Free Web Mapping Application

texts:
I have recently left my job as a GIS Specialist that I've held for a number of years. This was my first GIS job and I managed to lean a lot. Now that the holidays are over, I've decided to search for the next position.

Obviously, someone in my position should want to garner a few new skills while I'm on the job hunt. One thing that I didn't have a chance to work with, is web mapping applications in my previous position. It's a skill I'd really like to pick up because it's the present and future of GIS.

My idea is, to start an online resume of sorts. Something I can create and send to prospective employers to show I can work in the GIS environment. Although resumes are great, I feel like this is a lot more tangible than spitting out a bunch of facts about what I do in a resume. When I was hiring people, I had a few people do this sort of thing and I find I have a lot more confidence in people I hire because I can see what they are capable of.

I'm wondering, for all of you who work with web mapping, where should I start? I want to start simple, but have the ability to expand when I get comfortable with it. I'm not willing to fork out money for ESRI products because, well.....I'm unemployed. I want to start working on some open source freeware products. I think Quantum GIS with the OSGEO package seems to be the most popular? I've worked with QGIS in my previous position for desktop and development applications, but I haven't created a web application. What other good, free products are out there that you have used with some success. Is there anything really cool in open source like 3D mapping? I was thinking about creating a 3D fly through tour video using some DEM's or LiDAR if I can get my hands on that sort of data. Any thought, suggestions?

Thanks,
Mike


texts: Can QGIS be installed on Ubuntu 17.04 Zesty?

texts: I cannot find a repository for QGIS (2.18 or LTR 2.14) compatible with the latest Ubuntu distro (17.04) in any of the PPA's listed in the installer docs. Do I have to wait for a 17.04 supported release, or is there perhaps a workaround?


texts: SDE Connections in Geoprocessing Service

texts:
I have encountered something kind of odd and I'm not sure what's causing it. I have published a geoprocessing service and it's duplicating the SDE connection file during publishing. The code is very simple, all it's doing is checking that a database connection can be made and then unzipping the uploaded file. The service works, but I'd like to understand why it's duplicating the SDE connection file.

The python code:

import arcpy, zipfile

db_con = r'Database Connections\HbMonitoringTest_nbcidb_HabitatTestWriter.sde'
arcpy.AddMessage(db_con)
if arcpy.Exists(db_con):
  arcpy.AddMessage("your db connection worked")
else:
  arcpy.AddMessage("your db connection sucks")
try:
  myZipFile = arcpy.GetParameterAsText(0)
  arcpy.AddMessage(myZipFile)
  scratch = arcpy.env.scratchFolder
  arcpy.AddMessage(scratch)
  zip_ref = zipfile.ZipFile(myZipFile, 'r')
  zip_ref.extractall(scratch)
  zip_ref.close()
  arcpy.AddMessage("success")
except:
  arcpy.AddMessage("fail")

When I publish this as a geoprocessing service it does this in the /extracted/v101 directory:

Other info:
When I use the service via javascript API, it uses the duplicate database connection (HbMonitoringTest_nbcidb_HabitatTestWriter1.sde). When I run the service in Desktop, it uses the catalog connection which is to be expected (this is the connection without the 1 on the end). This database is registered with the server under the same name (HbMonitoringTest_nbcidb_HabitatTestWriter). The server and publishing machine are two different machines. Server settings do not allow data to be copied to the server. ArcGIS Server 10.4.1, ArcGIS Desktop 10.6.


texts: Source code for interpolation using a known set of random points in an area of land

texts:
In a fixed area of land random points with their values are generated. I need to get a heat map for this area of land showing the variations using those known point values.
I have found that mostly IDW and Kriging are amongst the interpolation algorithms which are used to do this job. As I don't have the relevant background, I find these pretty difficult to understand and implement in code.
The language C# is going to be used to implement the algorithm. So, are there any source code available with appropriate documentation to help understand it, so that it could be easily re-written/modified for this purpose?


texts: How can ASTER elevation data be converted to a CAD model?

texts: I am trying to convert ASTER elevation data to a 3D CAD model using only opensource software.  I have tried using different file formats but I have been unsuccessful.  I use QGIS for my GIS work and would like to open the model using something like Sketchup or Meshlab.


texts: Using arcpy to set a Python expression for Label Expression

texts:
I have a python script that I input into the label expression by hand which assigns different labels depending on the attribute data, see:

def FindLabel ( [a], [b] ):

if [b] != 'ok' and [b] != 'missing' and [b] != 'not transmitting':
  return "<CLR red = '255' green = '170'>" + [a] + " " + [b] + "</CLR>"
elif [b] == 'missing':
  return "<CLR red = '130' green = '130' blue = '130'>" + [a] + " " + "Missing" + "</CLR>"
elif [b] == 'not transmitting':
  return "<CLR red = '255'>" + [a] + " " + "Not Transmitting" + "</CLR>"
 else:
  return [a]

Is it possible to use arcpy to apply this label expression to a layer using labelClasses.expression, or any other method? The help page hints at the possibility by saying of the .expression method:

  Provides the ability to get or set a layer's individual label class expression. This can be as simple as a single field or more advanced using either a VBScript, JScript or Python expression.

but I've never seen a working version.


texts: Can´t add features to an ArcMap project because a GDB lock of my own computer ArcGIS

texts:
I'm trying to add some data from a GDB to an ArcMap project, but i am not abble to do it because of an existing lock.

I have checked the origin of the lock and it says it´s from my own PC. Actually, when I close the project I am able to access the GDB, the problem is that I am making cartography maps, so I need to add the data to that specific project.

I have dealt with locks due to other users accessing the GBD, but never with a lock of my own PC, I dont have idea of how to solve this problem.

PD: The specific project already has some layers from that GDB, the lock problem appeared right after I set the source of the data to that GDB.


texts:
I have ArcObjects .NET code that is running on a machine with ArcGIS for Server 10.1.

The application worked on the dev server.
But when I moved to Test server, I am not able to initialize the licenses.

This is the code I am using

RuntimeManager.Bind(ProductCode.Server);
bool succes = false;
AoInitialize arcObjectInitEnv = new AoInitializeClass();
if (arcObjectInitEnv.IsProductCodeAvailable(esriLicenseProductCode.esriLicenseProductCodeStandard) == esriLicenseStatus.esriLicenseAvailable)
{
    status = arcObjectInitEnv.Initialize(licenseCode);
}
else
{
    throw exception(...);
}

I am sure I authorized the server right.
and I am able to run ArcGIS server manager.
And test the services and they are working fine.

I compile as "Any Cpu" , which is the same I used for Dev.
Is there tool that give me "availablity" for server, like desktop?
I didn't install ArcObjects SDK on the test server, because I assume it is not required there. Is that true?

How can I troubleshoot this?


texts: Debugging RuntimeError: workspace already in transaction mode from arcpy.da.UpdateCursor and ArcSDE feature classes?

texts:
I am making my first attempt at editing an ArcSDE feature class with python through a da.UpdateCursor. I'm essentially taking code I've written for a file geodatabase feature class and applying it to an SDE feature class. Doing so produces an error, and I'm not sure how to rectify the problem.

I'm using ArcGIS 10.1 for Desktop.

Pertinent code:

Struxfeature = r"DatabaseConnections\PipelinePathways.sde\PPGIS.Strux\PPGIS.StruxPts_v04_all"
P6feature = r"DatabaseConnections\PipelinePathways.sde\PPGIS.Strux\PPGIS.Land_Projects_Parcels_P6"

SDE = r"Database Connections\PipelinePathways.sde"
UserID = "E1B8"
Parent = "SDE.DEFAULT"
version = "change_RW_VC_4447_14_to_C"

#Create Version
print "Creating version"
arcpy.CreateVersion_management(SDE, Parent, version, "PUBLIC")
VersionName = UserID.upper() + "." + version

#Layers
arcpy.MakeFeatureLayer_management (Struxfeature, "Struxlyr")
arcpy.MakeFeatureLayer_management (P6feature, "P6lyr")

#Switch to version
print "Switching version"
arcpy.ChangeVersion_management("Struxlyr", "TRANSACTIONAL", VersionName)
arcpy.ChangeVersion_management("P6lyr", "TRANSACTIONAL", VersionName)

#Start editing
print "Initiating editing"
edit = arcpy.da.Editor(SDE)
edit.startEditing()

# Start an edit operation
edit.startOperation()

#Change P6 project numbers
print "Updating P6.\n"
P6Cursor = arcpy.da.UpdateCursor ("P6lyr", ["P6_NBR", "Name"])
for row in P6Cursor:
    codecodecode

The error comes from the line 'for row in P6Cursor:'

Traceback (most recent call last):
  File "C:\E1B8\ScriptTesting\ScriptsIUse\ChangeP6.py", line 81, in <module>
    for row in P6Cursor:
RuntimeError: workspace already in transaction mode


texts: How can I import an image into Google Colab as an Earth engine asset from Google cloud storage

texts:
When I run the following code in Google Colab, the image gets imported from Earth Engine assets folder into Google Colab such that I can visualize it:

import ee
imageROI = ee.Image('users/username/fileName')

However, instead of referring an asset from the Earth Engine, I want to reference an image from google cloud storage. Is there any method to do so?  An example that has not been working:

import ee
imageROI = ee.Image('gs://bucket/fileName')

I am looking for something on these lines so that I can reference an image directly from Google cloud storage.
Any help is deeply appreciated.


texts:
I would like to learn how to do dynamic/interactive maps for users interested in linguistic issues. In short, the idea is that a user could choose for instance between two syntactic constructions (Put on your coat vs Put your coat on) and the map would show the geographic distribution of the option selected by the user according to a database.

I am new to the GIS world. I just worked a bit with QGIS to do some basic maps with my own data, but I do not know how to do interactive web-based maps. So by now I just would like to get basic information in order to understand well the necessary steps I need to follow.

Can you siggest an introductory book, resource, software, etc?


texts: Setting one class in raster symbology using QGIS?

texts:
Is it possible in QGIS to exclude some values from the symbology and define only 1 class with 1 color for all the other values?

I have hard time to find how to do it with the Singleband pseudocolor render type. It seems QGIS only allows 2 classes minimum (why not 1?) The problem is that the 2 classes appear in the legend afterward while in the legend I want only 1 class. I have a raster with values from 0 to 100 and I would like to exclude values from 0 to 29 and draw value from 30 to 100 in one green color.

In ArcGIS there is a convenient option where one can group certain values into one class at the symbology level (without touching anything in the original data).

Is there a way to produce he same effect in QGIS without reclassifying the data?


texts: Google violates same origin policy using KML with Google Maps?

texts:
Does Google Violate the same origin policy when it tries to access the KML files from a globally accessible server?

The condition for Google maps to display the KML is that it should be globally accessible. Now there are 2 problems with this:

Security of the KML files. (This can be avoided by using geoXML3)
Violating the Same-Origin Policy.

Is there a way around to the same-origin policy?

How does google parse the KML file?


texts: Using MMQGIS plugin through Python console?

texts: I know how to make Hub Lines using the GUI, I want to know how to do the same using the Python Console, a tutorial would be great.


texts: Combining duplicate points while preserving their associated attributes

texts:
I'm working with building routes for a lit of technician's stops. In some cases, they have two "stops" at the same customer location -- turning on two types of service at a house, for example, is listed as two jobs even though it's at the same house.

Each point has a start time (when they're supposed to arrive), a duration (how long the work should take to complete), and task type. I need to remove the duplicate geographies, because the resulting route maps aren't as clear as they could be -- the frequent question in meetings is "why isn't there a Stop #1?" and answering "because it's at the same place as Stop #2 and that marker is on top of it" tends to confuse non-GIS colleagues. However, I also:

want to preserve the earlier start time for Task 1, not the later start time of Task 2
need to sum the task durations, so if Task 1 takes 5 minutes and Task 2 takes 10 minutes, the combined single-stop task takes 15 minutes
would like to keep task types around (concatenating the strings could work)

So something like the following example, Joe's first two stops need to be combined into one:

1  8:00  8:24  24  -81.12346  31.54642  DNP  Joe Smith
2  8:24  8:37  13  -81.12346  31.54642  LCT  Joe Smith
3  8:45  8:55  10  -81.65487  31.45121  ABC  Joe Smith

The ideal result is:

1  8:00  8:37  37  -81.12346  31.54642  DNP; LCT  Joe Smith
2  8:45  8:55  10  -81.65487  31.45121  ABC       Joe Smith

I considered dissolving on a lat/long and technician fields, but that wouldn't preserve the other attributes.

I can find the duplicates and build a list of OBJECTIDs that need to be merged with a search cursor, but I am confused how to proceed after that. I've drawn out a sequence of search and update cursor loops which would probably work, but seems inelegant and inefficient (requires looping through the entire feature class many times). Is there a better approach?


texts:
Using v.dissolve in GRASS, I'd like to get a multipolyon as a result. However, I end up with distinct polygons with the same 'dissolve value.'

For example; lets assume I have two polygon features which are not touching each other spatially. Both have a 'DISS" value of 1. If I run v.dissolve on this DISS value, I get the same thing in the output as the input (i.e. two seperate polygons with a DISS value of 1). What I would like to be able to get is an output that contains one single multipolygon. How do I go about this?


texts: I have a lot long/lat data in MS SQL Server, i want to do reverse geocoding by using openstreetmaps or some other free web service, and on the end i want to display all of that in report builder, does anyone tried something similar? How to pass my coordinates to web server, and how to receive answer and put it in some fields?


texts: PostGIS 3D bounding box query that envolves all 3 axes (xyz)

texts:
I have some experience with PostgreSQL, but not with PostGIS, and for the life of me I can't make my bounding box take the Z axis into account.

I have a points table with XYZ columns. I add another one:

ALTER TABLE foo ADD COLUMN the_geom Geometry(PointZ,4326);

Then I populate the_geom column with my XYZ coordinates:

UPDATE foo
SET the_geom =  ST_SetSRID(ST_MakePoint (x,y,z), 4326)

I confirmed my 3 dimensions with:

SELECT f_geometry_column As col_name,
    type,
    srid,
    coord_dimension As ndims
FROM geometry_columns

And finally, I can query my data against the bounding box with:

SELECT *
FROM foo
WHERE  the_geom &&  ST_3DMakeBox( ST_MakePoint(0,0,0), ST_MakePoint(6600, 15500, 0)  )

Which should return no rows since in the bounding box the Z axis is 0 and all my Z data is greater than zero, but the query only filters the X and Y axes, not caring about the Z axis.

What am I'm doing wrong, please?


texts: Create an attribute table with coordinates vertices in QGIS

texts:
I need to create, in a map layout, an Attribute Table item containing the Map Sheet vertices values in different coordinate reference systems, controlled by the atlas (Coverage Layer = Map sheet), as shown in the images.
Any advice to solve this problem?

I tried using vector-> geometry tools-> add geometry attributes but only area and perimeter values were added to the table, so I am considering using the vector-> geometry tools-> extract vertices first and then using the add geometry tool. I do not know if this is the right way to do this and I do not know how to get/add the converted coordinates (UTM and/or geographic) for all map sheets (the red rectangles) according to the other CRS.


texts: Exporting Data Driven pages as JPEG with World File

texts:
I have a script that works to export my pages as jpegs, but I want to add the world file to them. I don't want to have to do all the {df} args and whatnot, I just want to skip to the world file. How can I do this?

Script I currently have:

mxd = arcpy.mapping.MapDocument("CURRENT")
for pageNum in range(1, mxd.dataDrivenPages.pageCount + 1):
    mxd.dataDrivenPages.currentPageID = pageNum
    arcpy.mapping.ExportToJPEG(mxd, r"D:\Matt_work\jpgs" + str(pageNum) + ".jpg")
del mxd

EDIT: If I use this script:

mxd = arcpy.mapping.MapDocument("CURRENT")
for pageNum in range(1, mxd.dataDrivenPages.pageCount + 1):
    mxd.dataDrivenPages.currentPageID = pageNum
    arcpy.mapping.ExportToJPEG(mxd, r"D:\Matt_work" + str(pageNum) + ".jpg", world_file=True)
del mxd

I get this runtime error:

Runtime error
Traceback (most recent call last):
  File "<string>", line 4, in <module>
  File "c:\program files (x86)\arcgis\desktop10.3\arcpy\arcpy\utils.py", line 182, in fn_
    return fn(*args, **kw)
  File "c:\program files (x86)\arcgis\desktop10.3\arcpy\arcpy\mapping.py", line 1026, in ExportToJPEG
    layout.exportToJPEG(*args)
  File "c:\program files (x86)\arcgis\desktop10.3\arcpy\arcpy\utils.py", line 182, in fn_
    return fn(*args, **kw)
  File "c:\program files (x86)\arcgis\desktop10.3\arcpy\arcpy\arcobjects\mixins.py", line 568, in exportToJPEG
    return self._arc_object.exportToJPEG(*args)
TypeError: PageLayoutObject: Error in executing ExportToJPEG


texts: How to use QgsVertexMarker in a Python plugin?

texts:
I am trying to display a vertex marker for a QGIS Python plugin I am writing. I am getting a the following error...

'NameError: global name 'QgsVertexMarker' is not defined'

canv = self.iface.mapCanvas()
canv.setExtent(rect)
marker = QgsVertexMarker(canv)
marker.setCenter(QgsPoint(x,y))
canv.refresh()

What am I doing wrong?


texts: How to filter legend by map content in composer with raster layers?

texts:
Does the "Filter Legend by Map Content" button in composer of QGIs 3.2 version work only for vector layers? It seems it's a bug reported in 2.18.x version and should be fixed for 3.x versions (https://issues.qgis.org/issues/14194). However it seems that in 3.2 version the problem remains. Is there workaround to address this issue?
I've a project with several raster layers and I want to generate a map of each layer, with the specific legend items of each layer. This works very well with vector layers, including automatic maps generation through the Atlas feature. However in the case o rasters, it has no effect, and the legend items of all layers are always displayed in the composer (or layout), as you can see in the figure bellow. I tried also to set the "only show items inside current atlas feature" check box and it didn't work too.


texts: What factors need to be considered to decide the radius of a bubble map or is there a generic formula?

texts:
I recently started with GIS and I am learning as I go about implementing things.
Take the case of JHU Dashboard or the first world map in this NYT article or the maps in this Guardian article, what I do not yet understand is how do they decide the radius of the bubble on the maps?

In some maps I can see a discrete nature of radius like -

But what I do not understand is that we do not know the upper limit of the numbers we would end up mapping. Moreover in most of the maps it is not discrete but continuous.

Especially for me, I made a website too where I display India's infections as an iframe. I arbitrarily multiply the cases with 2000 to make them look good on the screen. But then a few days ago the numbers went up a lot and now everything looks out of proportion even after reducing the multiplying factor from 2000 to 1500.

I just feel very bad about not knowing how to correctly represent this data.

Additional Info -
I use folium to make an html which I embed as iframe. Folium uses Leaflet internally.


texts: Seeking global projection suitable for raster preserving equal area

texts:
I have a global tif file in epsg:4326 and a shapefile of points distributed around the global.
I want to retrieve the sum of the value of pixels around 10km of each point.

One obvious thing to do is to convert tiff to points and intersect in an equal area projection. Unfortunately conversion to points takes ages so i abandoned it.

Then focal statistics with the corresponding buffer and getting value for each point is an option, but what projection do you think is best for this?

I used mollweide and it seems that when I measure the pixels in any area of the world they have the same size.


texts: Leaflet bookmarks with a server-based bookmarks store

texts:
I'm setting up leaflet bookmarks with a server-based bookmarks store as per 'Custom storage(e.g AJAX)' in the Bookmarks documentation (https://github.com/w8r/Leaflet.Bookmarks).

I'm using a perl PUT handler running under Apache and I'm able to save a bookmark just fine into my '/bookmarks/' folder ('204 OK'), the problem is that I can't get a list of boomarks. It always says 'No Bookmarks Yet' - even if I create a new bookmark in the same session.

The bookmark files themselves seem OK with all the required fields. I have verified that I can retrieve the bookmark files with a GET request with no problem.

function getAllItems works OK and gets a json list of filenames from my server like so:-

{
  "bookmarks": [
    { "file": "filename1" },
    { "file": "filename2" },
    { "file": "filename3" }]
}

Javascript is definately NOT my mother tongue and I'm not sure if the above format is correct or even if reading bookmark files has been fully implemented.


texts: rasterio mask aborts kernel. Potential issue with format for geojson feature?

texts:
I have a pandas dataframe where each row is a bounding box. i'd like to crop my raster by each feature.

       xmin       xmax        ymin        ymax  \
1  398296.6  398301.46  3280997.15  3280998.96

                                                 box
1  NEON_D03_OSBS_DP1_398000_3280000_classified_po...

I first created a function to create a geojson polygon from each row.

def data2geojson(df):
    features = []
    insert_features = lambda X: features.append(
            {"type": "polygon",
                 "coordinates":
                 [[(float(X["xmin"]),float(X["ymin"])),
                     (float(X["xmax"]),float(X["ymin"])),
                     (float(X["xmax"]),float(X["ymax"])),
                     (float(X["xmin"]),float(X["ymax"])),
                     (float(X["xmin"]),float(X["ymin"]))]]}
        )

    df.apply(insert_features, axis=1)
    return geojson.FeatureCollection(features)

such that the output of each row looks like

{'type': 'polygon', 'coordinates': [[(398534.16, 3280995.59), (398539.84, 3280995.59), (398539.84, 3280998.99), (398534.16, 3280998.99), (398534.16, 3280995.59)]]}

The idea would be then to iterate through features to create the crops

out_image,out_transform=mask(img.tile,[features[0]],crop=True)

This yields no error message, but just kills the python kernel.

My initial thought is that my format of the geojson is incorrect? I was following rasterio masking with GeoJSON from file issue

System info

3.6.5 | packaged by conda-forge | (default, Mar 30 2018, 00:00:55)
[GCC 4.2.1 Compatible Apple LLVM 6.1.0 (clang-602.0.53)]

rasterio.__version__
'0.36.0'

Edit #1: Looking at other answer, it appears we are looking for a geojson polygon, not feature. Same behavior, simplified the function.

Edit #2: Closed the polygon by providing the initial point as the last coordinate.

Edit #3, in case there was worry that my installation was invalid. I cloned the rasterio repo and ran this example. No issues. From that example, it appears that the features should be a dict, not a geojson polygon directly, which makes sense looking at the source. Revised. No improvement.


texts: create folder with random number for each session using arcpy

texts: I am processing data and it is required to have same name of outputs [standard names] for different geographical regions. Therefore it is necessary to save outputs in different folders each time. Is there a way to create a folder with different name each time when i run script?