texts: Summing attribute values of overlapping polyline portions in QGIS 3.0.3

texts:
Note:
This is basically the same question that I asked before in Summing attribute values of overlapping polyline portions in ArcGIS Desktop?, but this time, I'm trying to solve the same problem using QGIS. The reason is, I need to be able to teach/explain the solution to users in a developing country who cannot afford ArcGIS. I'm not an expert in QGIS either. So here is the question:

I have polylines of 80+ bus routes operating in a city. The polylines have attribute values for the number of passengers, number of trips, etc. Many of the polylines overlap at certain areas of the road network.

For example in the image above, the Green Route has 1,000 passengers and travels 20 trips a day while the Orange Route has 500 passengers and 10 trips a day. I'd like to have a final shapefile with data on 1,500 passengers and 30 trips (overlapping areas), but also data on the 1,000 passengers with 20 trips and 500 passengers with 10 trips (non-overlapping areas). This is just a simple example with two routes but I have 80+ in all that overlap in many locations with each other.

What I have done so far:

  Step 1: Run "Merge vector layers" to combine all the individual polylines into one shapefile.

What I want to do next:

  Step 2: I want to cut/split the overlapping features exactly where the overlaps occur.

  Step 3: I want to fuse/merge/join (whatever the correct term is in QGIS) all overlapping line segments, AND add up their numbers of passengers and trips into the resulting feature. 1,500 passengers and 30 trips in the above example.

  Step 4: Lastly, I want to, somehow be able to find out which routes whose values have been added together. In the above example, a column saying "Green, Orange" would be nice, but doesn't matter if I get "GreenOrange".

What I have tried so far for Steps 2, 3 and 4:

  Multipart to singleparts

Tried to split all overlapping segments, but the splits did not occur exactly where the overlaps start and end.

  Line-polygon intersection

Created a polygon that covers entire area to try to get all overlapping segments split. Splits did not occur exactly where the overlaps start and end.

  Aggregate

Seems to be able to sum up the numbers of passengers and trips, and also able to concatenate the Route Name, but the result is only one feature. In the example, I need to have 5 features as the result.

  Join attributes by location

Referring to Merging overlapping lines into one line. QGIS. Requires an Input Layer and a Join Layer. I only have one layer since Step 1, so this can't be it, unless I'm doing it wrong from Step 1.

  Dissolve

This just returns the separated lines into one, back to Step 1. I was hoping that it dissolve the overlapping segments.

  Collect geometries

Same as above.

EDIT:

Sample Data:

Shp file, Shx file, Dbf file

Result of Step 2: v.clean

Shp file, Shx file, Dbf file

Result of Step 3: Aggregate

Shp file, Shx file, Dbf file

Final Result, Step 4: Delete Duplicate Geometries

Shp file, Shx file, Dbf file

In above image I selected the overlapping line segments. This is how the map looks like:

Two line segments North of the loop are supposed to overlap, but they are now missing. Only a tiny chunk is left.


texts: Is there website with list of common QA/QC errors?

texts:
As I am working on reviewing my soils spatial data, and I am doing my QA/QC  I placed a point layer to check for errors.  I add the field column to mark the error and state the error.

I recently read the articles that came from ESRI on "Starting a Quality Assurance Program", "QA/AC for GIS Data: Initial Quality Control Checks" and "QA/AC for GIS Data: Recording and Tracking Errors" transcripts.

Is there a website of where there are potential list of popular errors that I can take a look at or samples such as pictures of the error?


texts:
I have bands 2, 3, 4, and 5 map.

How can I use these to make a color map in ArcGIS 10.1 for Desktop?

Can you give me steps?

Also I want to make attribute table of colormap.

How this can be done?


texts: Deactivating QgsMapTool after single map click?

texts:
I'm trying that to create a class that inherits from QgsMapTool and deactivates itself after a canvasPressEvent is called. I've tried three options (in the code below) but none of them is working.

class getCoordinateTool (QgsMapTool):
def __init__(self,owner, canvas,onclick,projectInstance):
    """constructor"""
    super().__init__(canvas)
    self.owner = owner
    self.canvas = canvas
    self.onClick = onclick
    self.project_instance = projectInstance

def deactivate(self):
    #option 1:
    #super().deactivate()

    #option 2:
    #super(getCoordinateTool, self).deactivate()

    #option 3:
    QgsMapTool.deactivate(self)
    self.deactivated.emit()

def canvasPressEvent(self, e):
    super().canvasPressEvent(e)
    pos=self.toMapCoordinates(e.pos())
    crsSrc = QgsProject.instance().crs()
    (x,y) = coorTransform(pos, crsSrc,self.project_instance)
    self.onClick(x,y)
    #deactivate the select tool
    self.deactivate()


texts: Using Arcpy to manage creating/deleting/updating features

texts:
Let's say we have a "Land" feature class that has several attribute tables attached to it (Owner, documents, etc), pretty much all of the relationships are Many to Many (With relationship class)

Whenever a user creates a new Land Polygon, he manually types in the primary key (I've found no way of automating this), and sometimes one must cut a Land Polygon into 2 (or more), and when that happens, both Polygons end up with the same primary key.

The primary key ideally should always be unique, but then it would break the relationships when the Polygon is cut.

Our users are not very tech savvy, so my idea was to create a simple GUI form using python and arcpy to populate the related tables, but how can I ideally manage the workflow of creating/cutting the Polygons in this scenario?


texts: How to create text label of different scales?

texts:
I would like to create a different type of map in ArcMAp.I would like to present the frequency of some words. So the dimension of the word (maybe even the color) would symbolize its frequency.

I started by created a polygon for each word and in attribute table I created two columns "word" and "frequency". So how can I use the words from the colomn "word" to create labels and give them size according to the colomun "frequency".

The result I want to obtain is something like this:
http://rising.globalvoicesonline.org/chennai/files/2011/06/figure-1-wordle.png

So if anyone has an idea how to create something like this in ArcGIS, I would be very grateful. Thank you very much.


texts: QGIS not creating new profiles correctly

texts:
QGIS 1.8 is not creating new profiles correctly under Documents and Settings.

QGIS support files are also missing in My Documents for newly created profiles.

I just tried to copy the missing files from an existing working profile to a newly created profile and it was working fine.

However, at my workplace, people don't have a permanent workstation to work on and can log in to any system on the network.

Does anyone have any suggestions?

OS: Windows XP SP3


texts: ArcGIS REST server projection with OpenLayers

texts:
I'm adding an ArcGIS tile server to my openlayers tile server options via this javascript API.  I have it displaying the map correctly, but it seems stretched out horizontally (I have a feeling the default projection may be causing this).

The URL I'm pointing it to is:

http://services.arcgisonline.com/ArcGIS/rest/services/World_Street_Map/MapServer/export

And declaring the layer like so:

var layer = new OpenLayers.Layer.ArcGIS93Rest(null,
        'http://services.arcgisonline.com/ArcGIS/rest/services/World_Street_Map
         /MapServer/export',
        { layers:"0,1,2" } );

Does anyone know the projection I should be using?  The API says I can specify a projection by adding an srs Projection ID property alongside my layers option, but I'm not sure even what the possible values are.


texts: QGIS displays NULL values as zero in the layout tables

texts:
I use QGIS 3.18.2-Zürich to display the table in layout mode. My table contains two fields of type Double. If I set the values to Null (None), a zero is displayed instead of an empty cell.
What should I do to make the Null values appear as empty cells in layout mode?
As a workaround, I can convert Double to Text, but this requires creating an additional column, which is inconvenient.


texts: Conditionally change pixel values in Google Earth Engine


texts:
After running several hundred lines of commands (and getting their subsequent error statements...), it is nice to be able to clear out the console window and start clean.

In the Python Console window in QGIS, I would like to clear it out.  If I close it and then open it again, the buffer is still retained.  If I 'select all' and delete the content, I don't get a prompt and can't do anything to get one back.

Any suggestions on how to clear out the scroll back buffer?


texts:
Where can I find a tool like ArcGis: Raster_to_ASCII

Provided that the tool MUST be accessible through command line interface.


texts: pass selected value from drop down list in updateParameters(self) as input to python script

texts:
I managed to create two drop down lists using the validation code where an user  can select the first parameter from a list and get another filtered list for the second parameter.

However, I am struggling with the passing the selected second parameter as an input to the python script.

def updateParameters(self):
"""Modify the values and properties of parameters before internal
validation is performed.  This method is called whenever a parameter
has been changed."""

accounts = "Database Connections\\someDatabase.sde\\TableToFilter"
rows = arcpy.SearchCursor(accounts)
self.params[0].filter.list = sorted(list(set(row.getValue('USERID1') for row in rows)))

if self.params[0].altered:
    rows = arcpy.SearchCursor(accounts, """ "USERID1" = '""" + self.params[0].value + "'")
    accountList = sorted(list(set(row.getValue('ACCGRPNAME') for row in rows)))
    self.params[1].filter.list = accountList
    self.params[1] = self.params[1].value

#del rows
return

I am using self.params[1] = self.params[1].value to pass the selected value to the rest of the script but I'm not sure I am doing the right thing.

How would you retrieve the selected value and pass it to the script?


texts: Road Graph plugin routing on a virtual numeric field does not work

texts:
I'm trying to calculate the shortest/fastest path based on a 'dynamic' virtual field with the road graph plugin.  The speedvalues are based on a combination of other ppvalues: eg field GFBscore.

But when I want to select GFBscore as the speed field in vector > road graph > settings it does not show up in the dropdown list

Is this due to a bug or do I have to make a feature request somewhere?


texts: Consuming ArcIMS image service in Silverlight application

texts:
I have a vendor who has supplied a secured ArcIMS image service (with raster data) and I am trying to add in a Silverlight based application. When I try to add this as layer in ArcMap and trying to convert this into MSD, it says that this layer is not supported.  I can't really publish this map service as an msd based service. I have been trying lot of different things but haven't found the way to do this. I thought I will ask if other users have tried to do this.

My environment: ArcGIS Server 10

Thanks
Jay


texts:
I'm new to QGIS and I wonder whether there is a tool that helps to adjust the edges between two or more polygons; similar to that topology tool that's used in ArcGIS.

http://desktop.arcgis.com/en/arcmap/10.3/manage-data/editing-existing-features/aligning-an-edge-to-match-another-edge.htm

I’m searching for an equivalent to the following tool (see link):

http://desktop.arcgis.com/en/arcmap/10.3/manage-data/editing-existing-features/aligning-an-edge-to-match-another-edge.htm


texts:
I have a large stream vector file (6GB) which contains four files (shp, prj, shx and dbf) all located together. I cannot display in QGIS. I get this error message from QGIS: The shapefile is not valid or it does not contain valid data source.
When we process the same file with GDAL,Ogr2gr and postgis database In all platforms I am getting the same error message as ERROR 4: Invalid offset for entity 11320128


texts: Why does GeoServer's GetFeature round my numbers, but only for outputFormat=CSV?

texts:
When I make a GetFeature request with outputFormat=CSV  to my GeoServer v2.4.5 I'm finding that the floating point numbers in my returned data have been rounded. When I make the same request with other output formats the floats are returned correctly.

Here's the request returning XML:

 http://geothermal.smu.edu:9000/geoserver/gtda-contributions/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=gtda-contributions:UNDRadiogenicHeat&maxFeatures=1&

This excerpt from the returned XML shows the data as it was entered into the backing database. There are more floating point numbers in the returned data, but I exclude them for brevity:

 <?xml version="1.0" encoding="UTF-8"?>
 <wfs:FeatureCollection xmlns="http://www.opengis.net/wfs" xmlns:wfs="http://www.opengis.net/wfs" xmlns:gtda-contributions="http://scr.siemens.com/gtda-contributions" xmlns:gml="http://www.opengis.net/gml" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:schemaLocation="http://scr.siemens.com/gtda-contributions http://geothermal.smu.edu:9000/geoserver/gtda-contributions/wfs?service=WFS&amp;version=1.0.0&amp;request=DescribeFeatureType&amp;typeName=gtda-contributions%3AUNDRadiogenicHeat http://www.opengis.net/wfs http://geothermal.smu.edu:9000/geoserver/schemas/wfs/1.0.0/WFS-basic.xsd">
    <gml:boundedBy>
        <gml:null>unknown</gml:null>
    </gml:boundedBy>
    <gml:featureMember>
        <gtda-contributions:UNDRadiogenicHeat fid="UNDRadiogenicHeat.2540ea6b-f94a-3088-9b32-762c77a8e57f">
            <gtda-contributions:latitude>47.90752513</gtda-contributions:latitude>
            <gtda-contributions:longitude>-92.79726665</gtda-contributions:longitude>
 ...

The same request with outputFormat=application/json&:

 http://geothermal.smu.edu:9000/geoserver/gtda-contributions/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=gtda-contributions:UNDRadiogenicHeat&maxFeatures=1&outputFormat=application/json&

returns similar (correct) data:

 {
"type" : "FeatureCollection",
"totalFeatures" : 343,
"features" : [{
        "type" : "Feature",
        "id" : "UNDRadiogenicHeat.2540ea6b-f94a-3088-9b32-762c77a8e57f",
        "geometry" : {
            "type" : "Point",
            "coordinates" : [-92.79726665, 47.90752513]
        },
        "geometry_name" : "shape",
        "properties" : {
            "latitude" : 47.90752513,
            "longitude" : -92.79726665,
...

However, if the request is for CSV:

 http://geothermal.smu.edu:9000/geoserver/gtda-contributions/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=gtda-contributions:UNDRadiogenicHeat&maxFeatures=1&outputFormat=csv&

the same portion of data has been rounded off (48 instead of 47.90752513, and -93 instead of -92.79726665):

 FID,latitude,longitude,...
 UNDRadiogenicHeat.2540ea6b-f94a-3088-9b32-762c77a8e57f,48,-93,...

I've looked for configurations that might cause this, but I see none.  I know there are filters that could do this, but I get the impression that these need to be explicitly requested, and the above requests do not include filters.

Is there something I'm overlooking? Is there some way that GeoServer is treating the data so glaringly differently just because I'm requesting CSV?


texts:
I have a base map image with micro regions, however in the attribute table there is no information regarding the long & lat. Is there a way I can create a long & lat for each micro region?

Thanks.


texts: Creating graduated legend for rasters using QGIS?

texts:
I have see a few posts about the difficulties of generating nice looking legends for rasters but they all seem to relate to older versions of QGIS (Las Palmas and older). I'm on 3.4.4 and want to know what is the best way to create a legend entry for a graduated raster with many values such as a DEM.

Are there any plug-ins designed for this that work with 3.4.4?


texts: Numpy raster reclassify outputs only 2 out of 6 classes

texts:
I have a raster (.tif) which I want to classify (6 classes) using gdal, python and numpy. The cell values range from: MIN = -30.9847 to MAX = 7.3505.

From the first two answers in this question, I've come to the following code:

from osgeo import gdal
import numpy as np

driver = gdal.GetDriverByName('GTiff')
input = gdal.Open('D:/test.tif')
band = input.GetRasterBand(1)
matrix = band.ReadAsArray()

# reclassification
matrix[np.where(matrix < 0)] = 1
matrix[np.where((0 < matrix) & (matrix < 0.3)) ] = 2
matrix[np.where((0.3 < matrix) & (matrix < 0.6)) ] = 3
matrix[np.where((0.6 < matrix) & (matrix < 1  )) ] = 4
matrix[np.where((1 < matrix) & (matrix < 1.5)) ] = 5
matrix[np.where((matrix > 1.5)) ] = 6

# create new file
output = driver.Create('D:/test_class.tif', input.RasterXSize, input.RasterYSize, 1)
output.GetRasterBand(1).WriteArray(matrix)

# spatial ref system
proj = input.GetProjection()
georef = input.GetGeoTransform()
output.SetProjection(proj)
output.SetGeoTransform(georef)
output.FlushCache()

The output raster is created, but it only has 2 classes, value 1 and value 6. Why could that happen?


texts:
I have a folder structure /foo/bar/baz/static/tiles/{z}/{x}/{y}.pbf containing many tiles according to a z/x/y file structure.
With /foo/bar/baz/ containing an index.html and a style.json.
I have this nginx config to serve the content of the /foo/bar/baz folder:
server {
    listen 80;
    listen [::]:80;

    server_name my.server.org;

    access_log /var/log/nginx/access.log;
    error_log /var/log/nginx/error.log;

    root /foo/bar/baz;
    index index.html;

    location / {
        try_files $uri $uri/ =404;
    }

}

I thought it would have done the job of retrieving files nested in subfolders according to the /{z}/{x}/{y}.pbf scheme, but when I browse to http://my.server.org the browser displays the index.html fine but the *.pbf files are simply not served (I cannot see them in the network manager, and there's no error actually.)
I'm not used to nginx, how could I load the tiles in the browser using this simple URL pattern: http://my.server.org/static/tiles/{z}/{y}/{x}.pbf?


texts: How to use Field Value from arcpy.GetParameterAsText()?

texts:
I want to to is have a user select a feature class, which then populates the fields that can be selected (this works), create a new field in that feature (works), and then populate the new field with values based on the inputField values.  Right now, it just populates the new field with the field name of the inputField.

It seems like this question has been asked a thousand times, but I still can't get this to work.

import os
import arcpy
from arcpy import env
env.overwriteOutput = 1

inputFC = arcpy.GetParameterAsText(0)
inputField = arcpy.GetParameterAsText(1)
inputFieldLabel = "GeoLabel"
output = arcpy.AddFieldDelimiters(inputFC, inputField)
inputFieldName = arcpy.GetParameterAsText(2)

#Add the GeoLabel field to the Selected Feature Class
arcpy.AddField_management(inputFC, inputFieldLabel, "TEXT")

#Calculate the field values for the selected field, using UniueID field and InputFieldName
arcpy.CalculateField_management(inputFC, inputFieldLabel, output, "Python_9.3")


texts:
There are 5 rasters: r0, r1,r2,r3,r4 with same extent, columns and rows. Random nodata cell disperse in each raster. Now I expect to fill nodata values of r0 with mean values of the other rasters in the same position(the cell of same row and column number). Of course, the mean values of cell exclude nodata values. I can not find function in raster package.
Actually, the basic function I expect is like 'con' in arcmap.


texts: How to align 3D contour file, aerial imagery and design drawings from engineers in microstation v8i series 3.

texts:
I am looking to align a 3D contour dgn file with the design drawings from an engineer in microstation.

Then align an orthorectified aerial image (Scale 1:10 000) to the contour file and design drawing file.

I want to do all this in microstation V8i series 3. And avoid aligning the files in a GIS software.

Please advice on steps to follow.


texts: Transformation matrix error when initializing the map

texts:
I use openlayers 6 in my project.
When initializing map I get this error in browser:

Uncaught AssertionError: Assertion failed. See https://openlayers.org/en/v6.0.1/doc/errors/#32 for details.
    at new AssertionError (http://localhost:1234/mapApp.3f0ed956.js:304:20)
    at assert (http://localhost:1234/mapApp.3f0ed956.js:4706:11)
    at makeInverse (http://localhost:1234/mapApp.3f0ed956.js:9565:23)
    at CompositeMapRenderer.MapRenderer.calculateMatrices2D (http://localhost:1234/mapApp.3f0ed956.js:22861:32)
    at CompositeMapRenderer.renderFrame (http://localhost:1234/mapApp.3f0ed956.js:23178:10)
    at Map.PluggableMap.renderFrame_ (http://localhost:1234/mapApp.3f0ed956.js:17137:20)
    at Map.<anonymous> (http://localhost:1234/mapApp.3f0ed956.js:15948:12)

When I click in error on first row I get description of the error on openlayers site:

Transformation matrix cannot be inverted.

Here is my map initialization:

let extent = boundingExtent([West, South, East, North]);
let mapCenter = getCenter([West, South, East, North]);
const israeliTM = getProjection('EPSG:2039');

const map = new Map({
        target: 'map',
        layers: [
            new TileLayer({ source: new OSM() })
        ],
        view: new View({
            extent: extent,//works ok when commited
            projection: israeliTM,
            center: mapCenter,
            zoom: 10

        })
    });

When I put on comment this row  extent: extent in map initialization above I don't get any errors, so it seems that extent function cause the error.

Any idea why when I try to set extent on map init I get the error?


texts: Copying one or more layers (with properties, style, and path of source data) on ArcMap and paste into different ArcMap project?

texts: I am using ArcMap 10.2.1 and i would like to take all the layers from one plan and bring them into another plan with all of the associated properties.


texts:
I am new to Python and currently trying to get my first script to work. Could someone tell me what is wrong with this code:

total = 0
def accumulate(increment):
 global total
 if total:
  if (total >= 1103):
    total = (total - 1103)
   else:
    total += increment
 else:
    total = increment
 return total

This one (the example from ESRI) works:

total = 0
def accumulate(increment):
 global total
 if total:
  total += increment
 else:
  total = increment
 return total

I am trying to run this code in ArcGIS 10.1


texts: Converting polyline to point (midpoint) including all attributes with Arcpy?

texts:
I have a feature with polylines and want to get the midpoint of every line.
Following pything script I found here works great:

import arcpy, math, datetime, numpy

print ("starting")
start = datetime.datetime.now() # for calculating time of process

#setting the containers
midpoint = #put your file to be populated in here, make sure it already exists
polyline = #put your polyline file in here

#housekeeping
arcpy.DeleteFeatures_management(midpoint)
arcpy.env.overwriteOutput = True

#generating the mid point
with arcpy.da.SearchCursor(polyline, "SHAPE@") as in_cursor, \
arcpy.da.InsertCursor(midpoint, "SHAPE@") as out_cursor:
for row in in_cursor:
    midpoint = row[0].positionAlongLine(0.50,True).firstPoint
    out_cursor.insertRow([midpoint])

#tidy up
del rows, row, updateRows_midpoint, outRow, out_cursor, midpoint, polyline
print "Done in ",datetime.datetime.now() - start, " seconds"

But, all the rest of the attributes are lost. Only the Geometry is written into the new feature.
I'm aware why, but I don't know how to do it.

How do I get them into the midpoint-feature?

Can it be done in the for-loop?


texts: Converting local XYZ to geographic coordinates

texts:
I am working with a computer program called PEP for Windows by Airbus, which is able to calculate exact flight paths of an aircraft, based on many inputs like weight, temperature, wind, etc. The program seems to work with geometrical coordinates, as its output is giving displacement to the east (x) and north (y) in meters, including altitude as z in ft:

0.000 0.000 0.000

0.051 -0.140 0.001

0.232 -0.636 0.003

0.616 -1.692 0.008

1.309 -3.595 0.018

2.423 -6.657 0.033

4.057 -11.147 0.054

6.280 -17.254 0.084

9.124 -25.069 0.123

12.597 -34.610 0.169

16.695 -45.869 0.224

The above example represents an acceleration at an airport during the first 10 seconds of the takeoff roll, where you can see the aircraft moves in southeastern direction on an upslope runway (slope data was given as an input parameter, but as soon as the aircraft lifts off, the z parameter shows either the height above the starting point or altitude above sea level, to be chosen by the user).

My objective is to be able to import this trajectory into Google Earth Pro in order to visualize the flight path, in combination with the elevation profile option to see the terrain clearance of the aircraft at each point of the flight (in the absence of knowing better software to see this). I believe I need to use a program like QGIS to convert the coordinates into Lat Long, right? Could anybody explain me which steps to follow in order to achieve my goal?


texts: Does the Google Earth API offer tools for capturing points from handmade points on paper maps?

texts:
My student team has been tasked with creating a software tool to help a fish biologist import the hand-drawn location markers on his paper maps to a digitized format, simply an output file with the lat/long and some basic ID data on the point. As a bonus, we had discussed outputting the mapped points to Google Earth for viewing.

After researching a little about the Google Earth Engine, I am curious if anyone can suggest whether or not it may offer an interesting solution to our problem. The biologist came up with the idea of using an optical mouse pen device to trace the points on the paper map from reference points (like the map corners). We would find a way to take a little bit of data about his map, trace points with the mouse somehow, and output the point data to a file or mapping software. He envisions a solution he can use on a laptop on a folding table in the field.

We would likely be creating a desktop app, and have pretty much unlimited freedom for selecting the software we incorporate. Does anyone with experience using the API know if it offers features which could help us capture points from a paper map (possibly digital photo upload) or a physical input device as originally specified?

We have a few months to research possible solutions and practice with them before we begin the actual project coding.


texts: GeoServer WFS GetFeature with BBOX returns too many elements?

texts:
My data is stored in PostGIS in EPSG:3035, and I query the data using WFS GetFeature as follows (in CoffeeScript)

  params = {
    service: 'WFS',
    version: '2.0.0',
    request: 'GetFeature',
    typeName: visible_layers.join(','),
    bbox: "#{e.bounds.toBBoxString()},EPSG:4326",
    srsName: "EPSG:4326"
    info_format: 'text/javascript',
    outputFormat: 'text/javascript',
    format: 'text/javascript',
  }

  $.ajax(
    url: GEOSERVER_WFS_URL,
    data: params,
    dataType: 'jsonp'
  )

This builds the following URL:

http://192.168.56.1:8080/geoserver/wfs?
  service=WFS&
  version=2.0.0&
  request=GetFeature&
  typeName=namespace:layer,namespace:layer2&
  bbox=5.18532156944275%2C51.04144786489307%2C5.1875638961792%2C51.04270932015568%2CEPSG%3A4326&
  srsName=EPSG%3A4326&
  info_format=text%2Fjavascript&
  outputFormat=text%2Fjavascript&
  format=text%2Fjavascript&
  _=1479898340470

And this works, apparently, until I want to select the pipeline as well.

My pipelines are all in the same table, but I use a separate layer for each (so same table with different SLD), then I see two things:

For each pipeline layer I get the same pipelines (so it does not use the SLD, which determines what is visible?)
I get pipelines which are not inside the BBOX, but close. Which is weird, because for my other objects this is not the case (only get objects inside the BBOX).

I would want the WFS request to also take the SLD of the layer into account, is that possible? And why does the selection of the BBOX seems more forgiving for some layers?

I am guessing I misconfigured something?

One idea I might have is that WFS is matching on BBOX? So the BBOX of my pipeline falls inside the BBOX, not the actual linestring itself? Because that could explain part of the behaviour. It is still not using the SLD :)


texts: cmv-app 2.0 connect.json is not a function

texts:
I am using the Dev branch of the cmv-app and having trouble getting up and running out of the box. It looks like I am unable to start the web server - error says connect.json is not a function.

I installed the node dependencies as described in the 1.3.4 docs. Not sure if anything has changed in 2.0.0 as far as setting up the dev environment. Are there docs specifically for 2.0.0?


texts: .toList return differences from the code panel and JavaScript API

texts:
Is it possible to remove the user's login requirement and use API only (on row 161 https://pastebin.com/VDmtjdJF)?

When I use print in the panel code, the object returns in json, however when I try the same with JavaScript (console.log) it returns me a completely different object.

CODE LINK:
https://code.earthengine.google.com/7376221e98815345b5f903f32f43dd59
Return from code panel - https://i.stack.imgur.com/AXFiz.png
JS CODE (link for complete code https://pastebin.com/VDmtjdJF):
// Create a list of image objects.
var imageList = collection.toList(collection.size());
//  print('imageList', imageList);
console.log('list');
console.log(imageList);

// Extract the ID of each image object.
var id_list = imageList.map(function(item) {
  return ee.Image(item).id();
});
//print('id_list', id_list);
console.log('list ID');
console.log(id_list);


texts:
I have many files within different folders, for example:
/path/to/dted/e001/n40.dt1
/path/to/dted/e001/n41.dt1
/path/to/dted/e000/n40.dt1
/path/to/dted/e000/n41.dt1

I managed to load single .dt1-file (Python API):
gdal.Open('path/to/dted/e001/n40.dt1', gdal.GA_ReadOnly)

Is there an option to load all the files as one dataset or merge many datasets into one?


texts:
I know, there are some libraries to convert .shp to .json even with javascript, but haven't found a way to convert .json back into .shp that could be implemented in web enviroment.
Do you know or have experience with something? Or what are my possibilities?

Thanks


texts:
I'm using 6s to correct sentinel data, but there is a little deviation that tells me that I'm doing something wrong.

This is the result for example of one scenario:

*******************************************************************************
*                        atmospheric correction result                        *
*                        -----------------------------                        *
*       input apparent reflectance            :    0.100                      *
*       measured radiance [w/m2/sr/mic]       :   37.699                      *
*       atmospherically corrected reflectance                                 *
*       Lambertian case :      0.02879                                        *
*       BRDF       case :      0.02879                                        *
*       coefficients xa xb xc                 :  0.00346  0.10170  0.14026    *
*       y=xa*(measured radiance)-xb;  acr=y/(1.+xc*y)                         *
*******************************************************************************

The problem is that if I calculate the atmospherically corrected reflectance (acr) using the formulas provided with a measured radiance of 37.699, I get a value of 0.02862 which does not match with the presented result = 0.02879. I know it's just a 0,6% error, but it's math... there should be no error.

I guess I'm doing something wrong and need guidance in terms of how to use these values to correct sentinel data, as it seems I am missing a step maybe, which i don't find in 6s' documentation. I need to make the calculation pixel by pixel, so running the algorithm n times (one for each pixel) is not an option.


texts: How to export tracks from QGIS to GPX format including elevation stored in elevation layer created from many hgt files?

texts:
I have many, many, many GPX files (without elevation) imported to QGIS.
I did some operations on them.
I imported hgt files with elevation data.
And I want to export those files from QGIS BUT with included elevation data for each lat,lon points.
How do I do that?


texts: How to add a dynamic timeline to the QGIS canvas?

texts:
I am using the plugin Timemanager which works really great but there is one thing I'd like to improve.

Instead of displaying the date and time on each frame, I would like to display a timeline stretching from *date_start* to *date_end* which shows the current frame time within this timespan.

I understand this possibility is not (yet?) offered by the plugin and therefore I am looking for an amazing idea on how to implement this. A core Qgis functionality that could be a workaround? An already existing plugin that I have missed? A small python script?

Obviously, in order for this to work nicely with the timemanager plugin, this timeline should be displayed in the map canvas, not in the composer.

A solution for Qgis 2.0 would be nice but I could also live with something for 1.8.

If someone has an idea or (even better) has already implemented something similar, please raise your hand!


texts:
In this example, kernel density maps are created.
In this
example, it is shown how to create comparable kernel density raster
maps

Now I would like to understand what are the units of kernel density maps in R.
?kde2d says:

Two-Dimensional Kernel Density Estimation: Two-dimensional kernel
density estimation with an axis-aligned bivariate normal kernel,
evaluated on a square grid.
It returns a list of three components. x, y The x and y coordinates of
the grid points, vectors of length n.
z  An n1 by n2 matrix of the estimated density: rows correspond to
the value of x, columns to the value of y.

I have seen some definitions, for instance from ESRI, but at the end, I have never understood what are the correct units for the legend.
library(MASS)
library(raster)
set.seed(123);

pp <- rnorm(100)
k = kde2d(
  x = rnorm(100) - 43,
  y = rnorm(100) - 23,
  h=0.2 ,
  n=1000,
  lims = c(c(-43,-42),
           c(-23, -22)))
rk <- raster(k)

rk2 <- projectRaster(from = rk, crs = "+init=epsg:31983")

res(rk)
#[1] 0.001001001 0.001001001
res(rk/area(rk))
#[1] 0.001001001 0.001001001
res(rk2)
#[1] 103 111

par(mfrow =c(1,3))
plot(
  rk,
  main = "raster lat long",
  horizontal = TRUE,
  legend.args = list(text='units')
)
plot(
  rk/area(rk),
  main = "raster lat long / raster area(km2)",
  horizontal = TRUE,
  legend.args = list(text='units')
)
plot(rk2,
     main = "raster UTM (epsg 31983)",
     horizontal = TRUE,
     legend.args = list(text='units')
)


texts: Using snap point when using measurecontrol in OpenLayers?

texts:
I have a problem about measurecontrol. I can use measurecontrol and I can show length to use measure control . But my problem is that I didn't use in snap because measurecontrol is a control, not a layer.

    measureControl = new OpenLayers.Control.Measure(
        OpenLayers.Handler.Path,
        {
            immediate: true,
            persist: true
        }
    );

    measureControl.getCustomLength = function(evt) {
        var str = '';
        str = 'Toplam: ' + (evt.units === 'km' ? evt.measure.toFixed(3) : evt.measure.toFixed(1)) + ' ' + evt.units + '&nbsp;&nbsp;&nbsp;';
           // uzunluk =evt.units ==='km';
            uzunluk=evt.measure;
        return str;
    }
    map.addControl(measureControl);
    measureControl.events.on({
        'measure': function(evt) {
            $("#mesafe").html(measureControl.getCustomLength(evt));
        }
    });
    measureControl.activate();

This is my measure control. And I want to snap this to point to my another vector. So when I look openlayer's snap:

    var snap = new OpenLayers.Control.Snapping({
        layer: measurecontrollayer,
        targets: [{
            layer:mylayer,
            tolerance: 15
            }]
        });
        snap.activate();

snap is using only vector but my measurecontrol is not a vector. Please show me some way to achieve this. Also is there another way to use measure in layer linestring layer?


texts: Why mapping a GeoTIFF world map in QGIS at low zoom level (0～2) results in a strange behavior?

texts:
I have created a GeoTIFF world map,
but I notice that it cannot be rendered correctly at zoom level which are low (e.g. zoom level 0, 1, 2).

Is it a predicted behavior or how can I fix this?
I also recorded a Youtube video for demonstration
Update：
There is a screenshot the same GeoTIFF file mapping at zoom level 5 (scale 1：18524092).
I think under this scale, it's mapping well.


texts:
I'm a water resources engineer with literally zero GIS experience. I'm working as a volunteer on a Engineers Without Borders project delivering safe drinking water to a remote community, where my team has provided me GPX files for point of interest, tracks, and LiDAR contours.

In order for me to use the QEPANET plugin (to build a hydraulic water distribution model), I must first add a DEM. The problem: I don't know how to take my LiDAR GPX file and make a DEM out of it, and afterwards, superimpose the points and/or tracks on that DEM.


texts: Probably base64 encoded WKB - Reading in QGis/ArcGIS?

texts:
I have a probably base64 encoded WKB string. How would I load that feature represented by the string into QGis or ArcGIS?

I found the QGis QuickWKT plugin, but that does not handle my base64 encoded string.

I have some Python programming skills.

Get your own WKB from http://nibis.lbeg.de/net3/public/ogc.ashx?Service=WMS&PKGID=22&cardo3SessionGuid=C3_2d80a40c-529e-4574-912b-53a3858b173c&layers=L160&query_layers=L160&width=1&request=GetFeatureInfo&height=1&srs=EPSG%3A31467&version=1.1.1&bbox=3585999%2C5743699%2C3586001%2C5743701&info_format=text%2Fplain&y=1&x=1 if you like.


texts: osm2po > PostGIS > QGIS/pgRoutingLayer Extension Issue

texts:
I've just been testing osm2po and everything worked fine and straightforward until it comes to the QGIS pgRoutingLayer Extension. Selecting the Database and selecting source and target node works fine as can be seen in the screenshot below:

but running the process gives me the following error:

RUN:

EXPORT rsp.:

Layer ist ungültig: Der Layer dbname='routing' host=localhost port=5432 user='postgres' password='admin' key='result_seq'
table="(SELECT mfr_2po_4pgr.*, result.seq AS result_seq, result.node
AS result_node, result.cost AS result_cost FROM mfr_2po_4pgr
JOIN (SELECT seq, id1 AS node, id2 AS edge, cost
FROM pgr_dijkstra(' SELECT id AS id, osm_source_id::int4
AS source, osm_target_id::int4 AS target, cost::float8 AS cost, reverse_cost::float8 AS reverse_cost
FROM mfr_2po_4pgr', 17642444, 10372788, false, true))
AS result ON mfr_2po_4pgr.id = result.edge)" (geom_way) sql=
ist ungültig und kann der Karte nicht hinzugefügt werden

Software config:

osm2po 5.1.0

postgresql 9.4.(6)

postgis 2.2.1

pgrouting 2.1.0

qgis 2.14.0

pgroutinglayer extension 2.0.0

What went wrong here?


texts: Identifying my coordinates from external datasource

texts:
I'm working with data received from an agency outside of my office and two columns in the dataset are labeled 'CompX' and 'CompY'. Here is a sample set pulled from the spreadsheet:

+-----------+------------+
|   CompX   |   CompY    |
+-----------+------------+
| 584490.53 | 4513000.14 |
| 584565.52 | 4513137.77 |
| 583596.45 | 4509846.2  |
+-----------+------------+

This data should all be in New York City, however, using these as decimal lat longs (like 58.449053 & 45.1300014) shows up in eastern Europe. Anyone familiar with this type of coordinate system? I'm trying to create a map based on this data.

Thanks!


texts: Trouble importing qgis.core python module

texts:
I am having issues loading the qgis Python modules in the osgeo4w shell. Ultimately I want to work with the modules in pyscripter. I tried the solution mentioned in Problem with import qgis.core when writing a stand-alone PyQGIS script, by loading the sample my_pyqgis_script.py that was given but got an

  import error: no module named qgis.core.

I also tried importing qgis.core in pyscripter after setting up my python path and path variables but in this case I got another

  import error: no module named sip.

Can anyone help?


texts: How to input Shapefile point data into a raster and have the raster adjust to the new data?

texts:
I have a shapefile with elevation data at specific points. I also have a raster that has incorrect elevation values, but the transition between pixels is what I need. Is there a way to put the values from the shapefile into the raster and have the values around the points adjust to the new values.

Some of the elevations on the raster are correct, so I can't increase all of the values the same amount.


texts: Get measurement values from point along the route Arcmap

texts: I have two shape files, one is a route shapefile with polyline and have measurement values. The other one is point with coordinates but no measurement values. I am trying to snap the points to closest route and get measurement of the route that snapped in(please see picture).


texts: ArcMap Python script for attribute split, replace and concatenate

texts:
I have a problem with attribute table which I'd like to learn to solve without MS Excel.

The Shapefile I have, has an attribute field, which consists of 3 parts, separated by ";" -sign. Here is an example of what the original looks like, how it is split and how the final attribute should look like. Red is the original and Blue is what it looks like at the end.

So what I would like to do, is to edit "RATAOSOITE" with a single python or VBScript in Arcmap. The script should do the following:

a) remove the first number/text before the first ";" sign.
b) change all empty spaces to "Underscore" and change the second ";" -sign to "Underscore".

I have next to no experience with scripting, so this maybe rather basic task is rather hard for me. I have done the same thing with Excel. Copy Attribute table to Excel --> text to columns by ";" --> replace "space" with "_(Underscore)" and "CONCATENATE(a1;Underscore;b1)" or something like that. But it would be great, if I could do this all in Arcmap Field Calculator.


texts: How do you create a local and editable copy of a feature service?

texts:
I am creating an app for a local fire department. I have acquired fire hydrant location data from:

http://gisweb.wsscwater.com/ArcGIS/rest/services/FireHydrant/FireBook/MapServer/0

I have imported it as a feature service in ArcGIS Desktop, however I cannot edit the data. I am trying to make a copy (maybe as a shp file) and start with that data to edit it.

Is this possible? Or is there a workaround since I cannot access it (I have asked multiple times but they have not granted as yet)?


texts: Conditional Statement in QGIS Field Calculator

texts:
I'm using QGIS to calculate total loss in an earthquack attack. I wrote this expression in Field Calculator.

CASE
WHEN  "build_type"  = 1 AND  "MMI"  < 5  THEN  "perc_loss"  = 0
WHEN  "build_type"  = 1 AND  (5 >= "MMI" AND "MMI" < 6) THEN  "perc_loss"  = 20
WHEN  "build_type"  = 1 AND  (6 >=  "MMI" AND "MMI" < 7) THEN  "perc_loss"  = 28
WHEN  "build_type"  = 1 AND  (7 >=  "MMI" AND "MMI" < 8) THEN  "perc_loss"  = 33
WHEN  "build_type"  = 1 AND  (8 >=  "MMI" AND "MMI" < 9) THEN  "perc_loss"  = 35
WHEN  "build_type"  = 1 AND  (9 >=  "MMI" AND "MMI" < 10) THEN  "perc_loss"  = 71
WHEN  "build_type"  = 1 AND   "MMI"  >= 10 THEN  "perc_loss"  = 95 END

There's no error/invalid warning, but no value resulted. Is there any wrong syntax I used? Anyone can help please? :)


texts: Is there a way to edit project variables with modeller in QGIS?

texts: Is there some way to set/edit user defined project variables with a model or a calculator so their values can be read by all the layouts and reports in a project?  I want to calculate a project reference point whose values can appear in all title blocks in various scale maps in the project.  And, the variable names will remain in the project template and recalculated for the next project etc.


texts: Why do the label locations in the Print Composer differ from the pdf?

texts:
When I set up labels in my drawing and view them in the print composer they look the way I want them to. However, when I look at the resulting pdf, it looks as though the labels have shifted on me.

My print composer map is controlled by an atlas. The labels I have denoted in the images below are as follows. Labels 1 and 2 have Not been manually placed and are labeled off the line. Labels 3 and 4 Have been manually placed (controlled by the data defined x/y fields), but are also labeled off the line. As you can see in the PDF Print image, labels 2 and 4 look shifted to the left and label 2 also looks have to shifted slightly down as well. Labels 1 and 3 look to be pretty much where I had them in the drawing, but they too seem to shift ever so slightly.

And the location of where the labels look to move to are not always the same from what I have seen. They seem to shift in all directions.

What I have tried:
Locking layers in the print composer.
Adjusted the DPI.
Exporting as image and straight to the printer (same result).
Print as raster.

Has anyone else experienced this before? And if so, do you have a remedy? I am currently using v. 2.16.2


texts: How to define boundary numbers without decimal points - Geopandas

texts:
Problem
The numbers in a legend has the decimal points (eg. 0.00-4.80 at the top of the legend box), when I want them to be integers without decimal points. All I care is that the boundary values are integers. It doesn't matter whether a certain region in a choropleth is labelled as the category A or B, as long as the A, B and other categories are separated by the integer values so that legend box shows only integer separations.

It's because I am making these maps for publishing articles after the Korean election is over. It would be much easier for the newspaper readers to have the integer numbers in legend!

My current graph

My code

import geopandas as gpd
import matplotlib.pyplot as plt

plt.xkcd()
seoul=gpd.read_file('../json/seoul_municipalities_geo.json')
# You can download the same file from the Github below
# https://github.com/southkorea/seoul-maps/tree/master/juso/2015/json

data_result.plot(figsize=(14,10),linewidth=0.25, edgecolor='black', column='ESRI_PK',cmap='Blues',scheme='quantiles',legend=True)

for index,row in seoul.iterrows():
    xy=row['geometry'].centroid.coords[:]
    xytext=row['geometry'].centroid.coords[:]
    plt.annotate(row['SIG_ENG_NM'],xy=xy[0], xytext=xytext[0],  horizontalalignment='center',verticalalignment='center')
    plt.axis('off')

plt.show()


texts: PostGIS multiple spatial joins in the same query

texts:
While I understand how I would do something like this in non-spatial SQL, I have confused about how I would join data from more than two joins in one query in spatial SQL.
I have a table of zipcode polygons 'Z', a table of restaurant polygons 'F' and restaurant points 'FP'. What I would like to produce as a result is a table with a row for each zip code and columns for the count of the restaurant polygons, the count of restaurant points, and the sum of both restaurant columns within each zip code.
As a test, I am only using 1 zipcode. First just getting the restaurant points:
SELECT z.modzcta,COUNT(f.name)
FROM osm.food_point f
INNER JOIN osm.modzcta z
ON ST_Intersects(ST_Transform(z.geom,3857),f.way)
WHERE z.modzcta = '11238'
GROUP BY z.modzcta;

returns [('11238', 192)]
Then, using the polygons:
SELECT z.modzcta,COUNT(fp.name)
FROM osm.food_polygon fp
INNER JOIN osm.modzcta z
ON ST_Intersects(ST_Transform(z.geom,3857),fp.way)
WHERE z.modzcta = '11238'
GROUP BY z.modzcta;

I get [('11238', 10)]
When I try to combine these two together:
SELECT z.modzcta,COUNT(fp.name) res_points, COUNT(f.name) res_polygons, COUNT(fp.name) + COUNT(f.name) total_rest
FROM osm.food_polygon f
INNER JOIN osm.modzcta z
ON ST_Intersects(ST_Transform(z.geom,3857),f.way)
INNER JOIN osm.food_point fp
ON ST_Intersects(ST_Transform(z.geom,3857),fp.way)
WHERE z.modzcta = '11238'
GROUP BY z.modzcta;

I get [('11238', 1920, 1930, 3850)]
I'm not sure what is going on here.
Also as an aside, is there any way to return the column names with the result of a query in the psycopg2 python module?


texts:
I want to check if schools from my hometown are accessible for students who use public transportation, so I need to perform this analysis in order to see if the schools are located in the range of ~400m, but I want to perform this as a network analysis using the roads, not just the buffer.

I have the following files: bus stations (points)- snapped to the lines, bus lines (vector line), and other roads and park alleys which I got from OSM (vector lines).

Is there a way I can perform this in QGIS? I tried various plugins but they are not working.


texts: How to add information about other tables to a buffer view in Postgresql?

texts:
So basically I'm creating a view of a buffer that fulfills some conditions, and this buffer needs to have data from point attributes in another table

This is my query that doesn't work:

CREATE OR REPLACE VIEW schema1.buffer_50_view AS
SELECT a.siteid,
a.name,
sum(select count from schema2.point_table b where count >= 50 and
st_dwithin(a.geom, b.geom, 100))  as sum_points,
st_buffer(p.geom, 100::double precision) AS buffer_50
FROM (select * from schema1.site_table p where site_progress = 100) a

it doesnt accept the sum line, and I have to add more lines like that, how to make it work?


texts: Assigning values from one polygon to second polygon if located within first using ArcGIS for Desktop?

texts:
I have two layers of polygons. I would like to know if there is a way to assign a value from
a polygon in the first layer to any polygons from the second layer that fall within that polygon?


texts: Setting a theme in the QGIS Print Composer via python?

texts:
I'm using python to command the print composer and export images/pdf (https://opensourceoptions.com/blog/pyqgis-create-and-print-a-map-layout-with-python/)

I'm looping through layers by name to export each layer. It works quite well for all but it's messy if layers are next to each other as I want to only show one. I need and want to use the themes that I've set up but I can't find a way to bring them to Print Composer.

I can change the theme using iface.mapCanvas().setTheme('green'), but this doesn't affect the layout in the Print Composer.

How can I use a different theme for each layouts via Python?
Can it be done through QgsMapThemeCollection and applyTheme()?

from qgis.PyQt import QtGui

layernames = [layer.name() for layer in QgsProject.instance().mapLayers().values()]

for each in layernames:

    layers = QgsProject.instance().mapLayersByName(each)
    print(layers)
    layer = layers[0]
    project = QgsProject.instance()
    manager = project.layoutManager()
    layoutName = each + "_layout"
    iface.mapCanvas().setTheme(each)

    layouts_list = manager.printLayouts()
    #remove any duplicate layouts
    for layout in layouts_list:
        if layout.name() == layoutName:
            manager.removeLayout(layout)

    layout = QgsPrintLayout(project)
    layout.initializeDefaults()
    layout.setName(layoutName)
    manager.addLayout(layout)

    #create map item in the layout
    map = QgsLayoutItemMap(layout)
    map.setRect(20, 20, 20, 20)

    #set the map extent
    ms = QgsMapSettings()
    ms.setLayers([layer])  #set layers to be mapped
    rect = QgsRectangle(ms.fullExtent())
    rect.scale(1.6)
    ms.setExtent(rect)

    map.setExtent(rect)
    map.setBackgroundColor(QColor(255, 255, 255, 0))

    layout.addLayoutItem(map)

    map.attemptMove(QgsLayoutPoint(5, 20, QgsUnitTypes.LayoutMillimeters))
    map.attemptResize(QgsLayoutSize(200, 180, QgsUnitTypes.LayoutMillimeters))

    legend = QgsLayoutItemLegend(layout)
    legend.setTitle("Legend")
    layerTree = QgsLayerTree()
    layerTree.addLayer(layer)
    legend.model().setRootGroup(layerTree)
    layout.addLayoutItem(legend)
    legend.attemptMove(QgsLayoutPoint(245, 5, QgsUnitTypes.LayoutMillimeters))

    scalebar = QgsLayoutItemScaleBar(layout)
    scalebar.setStyle('Line Ticks Up')
    scalebar.setUnits(QgsUnitTypes.DistanceKilometers)
    scalebar.setNumberOfSegments(4)
    scalebar.setNumberOfSegmentsLeft(0)
    scalebar.setUnitsPerSegment(0.5)
    scalebar.setLinkedMap(map)
    scalebar.setUnitLabel('km')
    scalebar.setFont(QFont('Arial', 14))
    scalebar.update()
    layout.addLayoutItem(scalebar)
    scalebar.attemptMove(QgsLayoutPoint(200, 190, QgsUnitTypes.LayoutMillimeters))

    title = QgsLayoutItemLabel(layout)
    title.setText(each)
    title.setFont(QFont('Arial', 24))
    title.adjustSizeToText()
    layout.addLayoutItem(title)
    title.attemptMove(QgsLayoutPoint(10, 5, QgsUnitTypes.LayoutMillimeters))

    #checking theme names can be accessed:
    themes = QgsProject.instance().mapThemeCollection().mapThemes()
    print(themes)

    layout = manager.layoutByName(layoutName)
    exporter = QgsLayoutExporter(layout)

    fn = 'C:/export/' + each + '.png'
    exporter.exportToImage(fn, QgsLayoutExporter.ImageExportSettings())
    #exporter.exportToPdf(fn, QgsLayoutExporter.PdfExportSettings())


texts: How to resolve an InvalidCastException for two Objects of the same type in an ArcMap AddIn?

texts:
I'm trying to work through an InvalidCastException that is being thrown when a form is displayed from an ArcMap Addin. This is a fairly recent occurrence and I haven't identified what has changed that would cause this exception.  The objects are System.Collections.Generic.List (although my gut is telling me that isn't relevant).

[A]System.Collections.Generic.List`1[XXX] cannot be cast to [B]System.Collections.Generic.List`1[XXX].
Type A originates from 'mscorlib, Version=4.0.0.0, Culture=neutral, PublicKeyToken=b77a5c561934e089' in the context 'LoadNeither' at location 'C:\Windows\Microsoft.Net\assembly\GAC_32\mscorlib\v4.0_4.0.0.0__b77a5c561934e089\mscorlib.dll'.
Type B originates from 'mscorlib, Version=4.0.0.0, Culture=neutral, PublicKeyToken=b77a5c561934e089' in the context 'LoadNeither' at location 'C:\Windows\Microsoft.Net\assembly\GAC_32\mscorlib\v4.0_4.0.0.0__b77a5c561934e089\mscorlib.dll'.

I believe this issue is closely related to this one:
https://stackoverflow.com/questions/2500280/invalidcastexception-for-two-objects-of-the-same-type but complicated due to the fact that it is ArcMap that is loading the assemblies. I'm attempted to use the Assembly Binding Log Viewer, but so far, unsuccessful (nothing ArcMap related is being shown).  But I'm guessing ArcMap is loading the same assembly twice.

I've tried clearing my assembly cache and the ArcGIS AddIn folders, but had with no luck.


texts: Problem with downloading dependencies gt-epsg-hsql

texts:
Few days ago I started training with maven repository on www.geotools.org. Right now I have a problem with downloading dependencies gt-epsg-hsql. My Eclipse tell me that kind of Error "Offline/Missing artifact org.geotools:gt-epsg-hsql:jar:13.2". I asked uncle google but he didn't know. So now I'm asking You what can I do? My pom file looks like that:

http://maven.apache.org/xsd/maven-4.0.0.xsd">
    4.0.0
    org.geotools
    tutorial
    0.0.1-SNAPSHOT
    jar
    tutorial
    http://repo.opengeo.org/

<properties>
    <project.build.sourceEncoding>UTF-8</project.build.sourceEncoding>
    <geotools.version>13.2</geotools.version>
</properties>

<dependencies>

    <dependency>
        <groupId>junit</groupId>
        <artifactId>junit</artifactId>
        <version>4.11</version>
        <scope>test</scope>
    </dependency>

    <dependency>
        <groupId>org.geotools</groupId>
        <artifactId>gt-epsg-hsql</artifactId>
        <version>${geotools.version}</version>
    </dependency>

    <dependency>
        <groupId>org.geotools</groupId>
        <artifactId>gt-shapefile</artifactId>
        <version>${geotools.version}</version>
    </dependency>

</dependencies>

<repositories>
    <repository>
        <id>maven2-repository.dev.java.net</id>
        <name>Java.net repository</name>
        <url>http://download.java.net/maven/2</url>
    </repository>
    <repository>
        <id>osgeo</id>
        <name>Open Source Geospatial Foundation Repository</name>
        <url>http://download.osgeo.org/webdav/geotools/</url>
    </repository>
    <repository>
        <snapshots>
            <enabled>true</enabled>
        </snapshots>
        <id>boundless</id>
        <name>Boundless Maven Repository</name>
        <url>http://repo.boundlessgeo.com/main</url>
    </repository>
</repositories>

Could somebody help me?:)


texts:
I'm writing a feature layer editor demo and I'm using the esri.dijit.editing.Editor class.

According to documentation (https://developers.arcgis.com/javascript/jsapi/editor-amd.html), both attribute edits and geometry edits are saved automatically.

Attribute edits are saved automatically when you change selection, change the focus field, or press enter.
Geometry edits are saved automatically.

The editor is working fine but I wanna add a save button to applyEdits only after click it.

Follow my source code and my attempts:

var params = {
    settings: {
        map: map,
        geometryService: geometryService,
        templatePicker: templatePicker,
        layerInfos: [
            {
                featureLayer: featureLayer
            }
        ],
        toolbarVisible: true,
        enableUndoRedo: true,
        maxUndoRedoOperations: 15,
        createOptions: {
            polylineDrawTools: [
                Editor.CREATE_TOOL_FREEHAND_POLYLINE,
                Editor.CREATE_TOOL_POLYLINE
            ],
            polygonDrawTools: [
                Editor.CREATE_TOOL_FREEHAND_POLYGON,
                Editor.CREATE_TOOL_POLYGON,
                Editor.CREATE_TOOL_RECTANGLE
            ]
        },
        toolbarOptions: {
            cutVisible: true,
            mergeVisible: true,
            reshapeVisible: true
        },
        layerInfo: {
            disableGeometryUpdate: true
        }
    }
};

var editorWidget = new Editor(params, "<editor-div-id>");
editorWidget.startup();

I've tried to stop event propagation using dojo/_base/event but the argument passed in before-apply-edits (FeatureLayer) is not a "valid event".

featureLayer.on('before-apply-edits', function(e) {
    // ...
    event.stop(e); // will trow error
});

Tried to clear the edits array too but continue saving:

featureLayer.on('before-apply-edits', function(e) {
    // ...
    e.adds = [];
    e.updates = [];
    e.deletes = [];
});

Someone know how to prevent this behavior?


texts: Making borders between polygons disappear in QGIS?

texts:
I want to have no borders between the polygons of a layer WITHOUT dissolving it. If I set the border color to the same as the fill colour there shouldn't be any borders visible, right? I also tried to set the border width to 0, however, the result looks like this:

you can hardly see the lines but there must be a way to get rid of them like in ArcGIS?


texts: Is there a limit of the number of legend labels passed to ArcGIS server print service?

texts: Is there any limit on the ArcGIS print service for the number of legend labels? I am using ArcGIS server 10.2 and ArcGIS Javascript version 3.10 and am able to print with legend labels (item descriptions) when there are about 60 different items even though the they are not complete. However, when there are more, all the labels disappeared like the image above. Anyone know if there were any limit or anyway to set the limit?


texts: How to make a map in Qgis with a background of google maps?

texts:
I am using Qgis 2.18.9 and I have a very simple request from a client. They have provided me with shapefiles of their mining sites and they want me to produce for them a PDF they can print at a really large scale, about 1 inch equals 2000 feet on a sheet of paper at a size 42 inches x 60 inches with the background looking very close to that of google maps.

I can add a google maps basemap to my drawing but the issue is when I make the PDF all of the labels are so small they are unusable. I am adding the google maps basemap through the use of the openlayers plugin.

Can somebody please guide me in this?

One final note, the final product does not have have the gray color that is seen on google maps used as a background, if it does it does if it doesn't it doesn't.


texts: Leaflet Foreach point in GeoJSON check Distanceto

texts:
So I was wondering how I can get the points from the
GeoJSON located here.

and then for each of them separately check if the distance of Markers is within 5 meters range of these. I have tried but not being able to call points from the GeoJSON and I'm wondering how I would get it to be specific for each of the Markers. So that the if statement return if any of the Markers from the GeoJSON is within 5 Meters of the "Hydrat" Markers below.

 var hydrat1 = [56.41083875205333, 8.864636421203615];
 var hydrat2 = [56.41069629842804, 8.868584632873537];
 var hydrat3 = [56.40941419180167, 8.866181373596193];

Here is the complete code I'm trying to use

function requestGeoJson(){
    fetch('https://laerkeholt.dk/admin/get-features')
      .then(function(response){ return response.json() })
      .then(function(json){
         geoJsonLayer.clearLayers().addData(json);
        // var http = new XMLHttpRequest();
        // var url = "https://maker.ifttt.com/trigger/vandingsmaskine/with/key/";
          //http.open("POST", url, true);
          //http.send(JSON.stringify({"value1": "ID 1"}));
          var hydrat1 = [56.41083875205333, 8.864636421203615];
          var hydrat2 = [56.41069629842804, 8.868584632873537];
          var hydrat3 = [56.40941419180167, 8.866181373596193];

         setTimeout(requestGeoJson, 5000);
      });
};
requestGeoJson();


texts: PyQGIS - Retrieve WMS and ArcGIS services capabilities (layers, styles, crs...)

texts:
I would like to know how to bind the methods used by the "Connect" button which retrieve services capabilities through the OWSConnect to correctly build the URL with paramas (or QgsDataSourceURI) from each layer capacity without opening the dialog or request the GetCapabilities URL by myself.

I'm looking for a "pure" PyQGIS solution, not using a dependency like OWSLib because I need parse ArcGIS Server services (Map and Features).

Any idea?

I'm interested in PyQGIS 2.18.x and 3.x.


texts: ArcMap 10.2.2 Creating two random points approximately 50 meters apart

texts:
I am fairly new to ArcMap 10.2.2 so please bear with me. I am trying to create ten sets of two random points approximately 50 meters apart within a polygon restraining feature.

Edit:

Minimum distance is not a factor, I simply need 10 sets of 2 random points. Each set of two points being exactly 50 meters apart


texts: Checking if a polygon contains multiple point geometries using PostGIS?

texts:
This is not a a duplicate of the "Count Points in Polygon with Postgis", since in the ST_Contains(geometry geomA, geometry geomB), this case has two geometries (start_geom and end_geom) in geomB, while the other question is about only one geometry.

Input polygon table us_county:

Column  |            Type             |                        Modifiers
----------+-----------------------------+--------------------------------
 gid      | integer                     | not null default nextval('us_county_gid_seq'::regclass)
 statefp  | character varying(2)        | |
 geoid    | character varying(5)        |
 name     | character varying(100)      |
 geom     | geometry(MultiPolygon,4269) |

Input point table trips:

Column    |         Type         |                      Modifiers
--------------+----------------------+-----------------------------------
 id           | integer              |
 gid          | integer              | not null default nextval('trips_gid_seq'::regclass)
 start_geom  | geometry(Point,4326) |
 end_geom    | geometry(Point,4326) |

I want to use to check if a polygon in us_county could contain both start_geom and end_geom of a trip points.

I've tried the following queries:

1. this query returns an empty result, which is not right

select us_county.name, count(trips.gid)

from trips, us_county

where st_contains(us_county.geom, trips.start_geom) and
      st_contains(us_county.geom, trips.end_geom)

group by us_county.name;

2. a second thought is to use ST_Collect() function:

select us_county.name, count(trips.gid)

from trips, us_county

where st_contains(ST_Transform(us_county.geom, 4326),
ST_Collect(trips.start_geom, trips.end_geom))

Group by us_county.name;

But this turns out to be very slow, for 3k polygon and 50k points, the query run over 40min.


texts: Add metadata schema to Geonetwork from Github

texts:
I have a local installation of GeoNetwork 3.10.4 where I want to add a metadata schema that is available on Github. The schema can be found here: https://github.com/geoadmin/geocat/tree/geocat_3.4.x/schemas - I want to add the schema iso19139.che.
I tried to copy the folder with its contents to the directory web\geonetwork\WEB-INF\data\config\schema_plugins but the new schema does not show up in the metadata and template section after I restart the application.
Also, I read that it should be possible to upload a ZIP-file with a new schema, but I don't see any place where I can upload a schema in the GeoNetwork admin section.
How could this particular (or any schema) be integrated to another GeoNetwork installation?


texts: Turn disjoint polygons into tessellation

texts:
I have a set of points that have labels. These labels represent geographies; however, some of the labels are wrong. Like in the picture below

So after running some analysis, I end up with a set of disjoint polygons

I would love to be able to turn it into a tessellation similar to this

Obviously I do not need the geometries to line up exactly as above, just an approx, but I was wondering if anyone knew of an algorithm that maybe made a triangulation on the centroids of polygons to get nearest neighbor table and then used a curve smoothing alg coupled with polygon buffering (like STBuffer method on sql spatial types) to change first pic to something like the second pic

so basically I want to turn disjoint polygons into a covering of a space but can't seem to find any algorithms for it


texts:
I have a strange problem that i believe its because the zone pass through international date line.
My goal is to find if a point is within a zone and I use the geoPHP library, everything seems to be working properly expect on cases as the attached zone.
In attached image you can see the zone and the point, but when i check if the point is in zone i get true.

Any ideas about how to handle such cases ?


texts: Change CRS of layer using pyQGIS without saving file and iterating every feature?

texts:
So I have one layerA with CRS:4326 and I can change CRS saving it:

# This layer has crs:4346
layerA = QgsVectorLayer("layerA.sqlite", "layerA", "ogr")

#New crs
crs=QgsCoordinateReferenceSystem("epsg:32629")

#saving layerA to layerB with new crs:32629
error = QgsVectorFileWriter.writeAsVectorFormat(layerA, "layerB.sqlite", "UTF-8", crs , "SpatiaLite")

With this code I can get a new LayerB with another CRS

I have read this: http://docs.qgis.org/testing/en/docs/pyqgis_developer_cookbook/crs

But I think this solution need to iterating over every feature to reproject it.

Is there any solution to reproject a layer with only one instruction (like when saving it) but without need to create a new layer ?

Something like layerA.reproject(newCRS) ?


texts: PyQt5 QGIS Plugin - "native Qt signal is not callable" error

texts:
I'm having fun writing plugins for QGIS, but when I execute the code QGIS send me this message error: native Qt signal is not callable.

This is the "incriminated" code (is a simple horizontal slider that modifies the zoom):

import os
import sys
import inspect
from PyQt5.QtGui import QIcon
from PyQt5.QtWidgets import QAction, QFileDialog
from qgis.core import QgsProject, Qgis, QgsMapLayer
from .change_style_dialog import ChangeStyleDialog
from qgis.utils import iface
from qgis.gui import QgsRubberBand
from PyQt5.QtCore import pyqtSlot
from PyQt5.QtCore import Qt
from PyQt5.QtWidgets import (QWidget, QLCDNumber, QSlider,
QVBoxLayout, QApplication)

# other stuff

def run(self):
  self.dlg.show()

  result = self.dlg.exec_()

  self.dlg.horizontalSlider.setMinimum(10000000)
  self.dlg.horizontalSlider.setMaximum(1000000000)

  if result:
     if self.dlg.horizontalSlider.sliderReleased():
        scale = self.dlg.horizontalSlider.sliderPosition()
        iface.mapCanvas().zoomScale(scale)

This is dlg code:

  class ChangeStyleDialog(QtWidgets.QDialog, FORM_CLASS):
    def __init__(self, parent=None):
    super(ChangeStyleDialog, self).__init__(parent)
    self.setupUi(self)


texts: How to sum a hundred of raster maps in GRASS GIS raster calculator

texts:
I have 120 maps with the same name pattern: arc_01, arc_02, arc_03, arc_(...) and want to sum all of them in one.

Does anyone know how to do it automatically in GRASS GIS raster calculator?


texts: Overlay polygons on raster colour map to extract area values for each colour within each polygon

texts:
Does anyone know an answer to the following:
I have a set of polygons created as the MCP of animal movements. I also have a habitat colour map raster, with 13 different habitat types as different colours. I need to find a way of overlaying one on the other and calculating the area of each habitat type within each polygon. Possible?
I think it has something to do with overlaying, but have not used it before, and usually asking those who know is far faster than self-discovery!
Thanks in advance.


texts: Making a Button in mapinfo with conditional function

texts:
I have following code  :

Sub create_ButtonFromVendor
      Dim i As Integer
      i=1
      Create ButtonPad "allVendor"  As
        ToggleButton
        Icon MI_ICON_ZOOM_QUESTION
          Calling get_Vendor
          HelpMsg vendors(i)
        Show
       For i = 2 TO Ubound(vendors)
          Alter ButtonPad "allVendor" Add
          ToggleButton
          Icon MI_ICON_ZOOM_QUESTION
          Calling vendor_ButtonAction(vendors(i))
          HelpMsg vendors(i)
        Show
      Next
    End Sub

It is small sub procedure. Here I wanted to create multiple buttons with different function in a loop. How would I add function in a loop? In the code I added above does not work as handler must be sub procedure without parameter.


texts: Getting error Route without Measure from Make Route Event tool?

texts:
I'm using make a route event tool, i want to generate lines that have specific measurements.

I made sure that the properties of the polylines were Polyline M as opposed to polyline, as well as my table has a foreign field, which links to the shapefile.

When I use make a route event, the output layer does not show on the workspace.

I open the attribute table of the output that make a route event produces, and it says, "route without measure".

I don't understand what this message means. I made sure that the tabular properties are the same of both the shapefile and the table that contains the measurement of the start to end points. I don't know what else to do to make the output layer display.

Can the "route without measure error" message be attributed to the unique identifier fields in both table or a projection issue?


texts: How could I specify raster's origin, while using arcpy.TinRaster_3d?


texts: Using Maputnik Edior to Style MapBox Vector Tile Data that is Hosted on my Machine without the use of a Server

texts: I'm trying to create a style.json for a MapBox Vector Tile data that is hosted on my machine without the use of a server.  I have a .mbtile data file that I unpacked and placed into a directory on machine using mb-util (https://github.com/mapbox/mbutil).  I'm using the Maputnik CLI because I thought that was the way to proceed with local style development.  However, I don't see a way to specify my data as a source unless I host it with the use of a server.  Is there a way to use the Maputnik Editor outside of hosting my data via a server?


texts: Aggregate heating energy demand based on building data and street segments with PostgreSQL and PostGIS

texts:
I have two datasets / tables:

Georeferenced building data with the yearly heating energy demand per building
Street line segments from OpenStreetMap.

My goal is to calculate the total heating energy demand per street line segment (energy demand per meter of street).
Is there a (ideally performant) way to do the following calculations in SQL?
for every street segment calculate the heating density as follows:
  find all buildings that are closest to this segment:
    sum up the heating energy demand of the matching buildings and divide it by the length of the street segment

So far I could not figure out how to translate the above pseudocode into a SQL query. My dataset contains about 40k buildings and 30k street segments.


texts: Error while unpacking ArcGIS Desktop Tutorial Data?

texts:
Hardware
Hp Compaq dc5700 Microtower
64bit Operating system
Intel pentium dual CPU e2180 @ 2.00 GHz
Memory 4GB

Software
Windows 7 home Premuim Service Pack 1
ArcGIS Desktop 10.3.1

I have just begun training on the ArcGIS Desktop software. As part of the training material I have managed to download the software with little or no problem. However when ever I try to unpack the training files I get the following error message.

Error code 2147467259

I have run a scan on the downloaded file using 7zip and it shows 5 errors,please see below;

While looking for similar solutions I came across reports of people looking to install SQL related softare who fixed their issue. The person who had the issue stated

"I reset the .net permission on the 2008 server buy running : C:\Windows\Microsoft.NET\Framework\v2.0.50727\caspol.exe -machine -reset"

I have attempted this step with no success.


texts: Adding Field To feature class using ArcObjects and C#?

texts:
I want to add a new field to my featureclass.My code for that is:

        public IMxDocument mxDoc;
        public IMap map;
        IFieldEdit2 field = new FieldClass() as IFieldEdit2;
        field.Name_2 = "Parcel_Way";
        field.Type_2 = esriFieldType.esriFieldTypeString;
        field.Length_2 = 50;
        field.DefaultValue_2 = "Parcel";
        IFeatureLayer2 flayer=map.Layer[cmboxLayer2.SelectedIndex] as IFeatureLayer2;
        IFeatureClass featureclass = flayer.FeatureClass;

        try
        {
            featureclass.AddField(field);
        }
        catch (Exception ex)
        {
            System.Windows.Forms.MessageBox.Show(ex.Message, ex.Source);
            throw;
        }

Program return this error :

What do i neglecting ?


texts: Is there any way to import templates in ArcMap?

texts: I have tried to find info on how to import a template for layout in ArcMap. I have a template file pagx file. In Pro there is import by the insert tab, but not in ArcMap. Is there a solution for ArcMap?


texts:
I wanted to know whether it is possible to change the coordinates of a Polygon by typing them in.

I don't want to have to manually drag them into the right position.


texts: How to show Puerto Rico in "albersUsa" projection in Vega-Lite or Vega?

texts:
I have a shapefile that has all the US counties along with Puerto Rico and I want to visualize them in "albersUsa" projection using Vega or Vega-Lite.
Problem is that in other projections, Puerto Rico is visible but in "albersUsa" projeciton it vanishes.
The closest I got to is this example by Mike Bostock.
You can download the counties-10m.json file over there.
Is there a way to render it completely(including Puerto Rico) as an "albersUsa" projection?
My biggest concern is that the "albersUsa" projection itself won't allow Puerto Rico to be displayed.
Any leads on how do I go about using a composite projection in Vega. Since there is an example in d3 and Vega is built on top of d3, I suppose there should be a way.
EDIT
It looks like the file here has a projected version of Puerto Rico in albersUsa projection.
However, when I plot it in Vega-Lite using "albersUsa" projection, I get the following -


texts: How to edit a layer in GeoServer, when we change table schema in a database (SQL Server)?

texts: So, is there a way to change/edit SQL Server layer definition in GeoServer, when a table schema has been changed in a database? Or you have to delete that layer and then build a new one?


texts: How to import ERDAS ECW in QGIS 2.6 using GDAL 1.11?

texts: I am using Lubuntu 14.10, QGIS 2.6, and Ubuntu GIS(unstable). QGIS was used to use GDAL 1.10, and can import ERDAS ecw format. But after GDAL upgraded to 1.11, QGIS couldn't import ecw images. How to solve this problem?


texts:
I have the problem, that my labelled roads (name of the roads) changes when I also change my scale by zooming in or out.
But I want a fixed size and position of the street names independent from the scale.
How can I get to this?

I know there is a possibility to set a X and Y coordinate for each label, but I have really a lot of streets and labels and I want to do this with a rule-based labelling.
So in a perfect way, I just do the rule-based labelling and my street names have a fix position and size.


texts: Use column from Raster Attribute Table as Legend in Leaflet Map

texts:
I have a simple RasterLayer in R. I added a column to the attribute table and would like to use these values in the legend of my leaflet-map.
The RasterLayer looks like this:

> layer1
class      : RasterLayer
dimensions : 38, 40, 1520  (nrow, ncol, ncell)
resolution : 1000, 1000  (x, y)
extent     : 4303000, 4343000, 3365000, 3403000  (xmin, xmax, ymin, ymax)
crs        : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
source     : memory
names      : ue65
values     : 1, 5  (min, max)
attributes :
 ID VALUE percentage
  1     1          0
  2     2       0-15
  3     3      15-20
  4     4      20-25
  5     5        >25

And the code for making the leaflet map like this:

leaflet() %>%
  addProviderTiles('Esri.WorldImagery',group='Imagery') %>%
  addProviderTiles('Esri.WorldStreetMap', group='Streets') %>%
  addRasterImage(layer1, colors = pal, opacity = 0.7, group = "layer1") %>%
  addLayersControl(
    baseGroups = c("Imagery", "Streets"),
    overlayGroups = c("layer1"),
    options = layersControlOptions(collapsed = FALSE)
  ) %>%
  addLegend(pal = pal, values = values(layer1$percentage))

I just cant figure out how to use the percentage-column from the attributes table in my legend. If I use:
values = values(layer1)
it shows 1-5 in the legend. Any hint would be superhelpfull!!


texts: Visualize temperature distribution inside boundary on GEE

texts:
When I run the below code, I am getting only shapefile boundary in color, and inside the border does not visualize any temperature as per the color bar.
/////////////////////////////////////////////////////

// Create a FeatureCollection from the list and print it.
var table = ee.FeatureCollection('users/nidhivermaiiita/district');
print(table);
Map.centerObject(table);
Map.addLayer(table);
var dataset = ee.ImageCollection('MODIS/006/MOD11A1')
.filterBounds(table)
.filter(ee.Filter.date('2020-03-30', '2020-04-29'));
var landSurfaceTemperature = dataset.select('LST_Day_1km');
var landSurfaceTemperatureVis = {
  min: 13000.0,
  max: 16500.0,
  palette: [
    '040274', '040281', '0502a3', '0502b8', '0502ce', '0502e6',
    '0602ff', '235cb1', '307ef3', '269db1', '30c8e2', '32d3ef',
    '3be285', '3ff38f', '86e26f', '3ae237', 'b5e22e', 'd6e21f',
    'fff705', 'ffd611', 'ffb613', 'ff8b13', 'ff6e08', 'ff500d',
    'ff0000', 'de0101', 'c21301', 'a71001', '911003'
  ],
};
//Map.setCenter(6.746, 46.529, 2);
//Map.addLayer(landSurfaceTemperature, landSurfaceTemperatureVis, 'Land Surface Temperature');

// map over the image collection and use server side functions
var tempToDegrees = landSurfaceTemperature.map(function(image){
  var LST=image.multiply(0.02).subtract(273.15).rename('LST');
  return image.addBands(LST);
});
// Create a chart.
var chart = ui.Chart.image.series({
  imageCollection: tempToDegrees.select('LST'),
  region: table,
  reducer: ee.Reducer.mean(),
  scale: 30
}).setOptions({title: 'Land Surface Temperature'});

// Display the chart in the console.
print(chart);
// print and add to the map
//print('image collection in temp in degrees', tempToDegrees);

// Create and print the chart.
print(ui.Chart.image.series(landSurfaceTemperature,table, ee.Reducer.mean(), 30));
//Map.addLayer(landSurfaceTemperature, landSurfaceTemperatureVis, 'LST');
// print and add to the map
//print('image collection in temp in degrees', tempToDegrees);
//Map.addLayer(tempToDegrees, {min: -20, max: 60, palette: landSurfaceTemperatureVis.palette}, 'temp in degrees');
var landSurfaceTemperature = dataset.select('LST_Day_1km').mean();
var landSurfaceTemperatureVis = {
  min: 13000.0,
  max: 16500.0,
  palette: [
    '040274', '040281', '0502a3', '0502b8', '0502ce', '0502e6',
    '0602ff', '235cb1', '307ef3', '269db1', '30c8e2', '32d3ef',
    '3be285', '3ff38f', '86e26f', '3ae237', 'b5e22e', 'd6e21f',
    'fff705', 'ffd611', 'ffb613', 'ff8b13', 'ff6e08', 'ff500d',
    'ff0000', 'de0101', 'c21301', 'a71001', '911003'
  ],
};
{
var image = landSurfaceTemperature.clip(table);
print(image, 'image');
Map.addLayer(image, landSurfaceTemperatureVis);
}
var mean = image.reduceRegion({
  reducer: ee.Reducer.mean(),
  geometry: table,
  scale: 30
});

// Print the result (a Dictionary) to the console.
print(mean);
Export.image.toDrive({
  image: landSurfaceTemperature,
  description: 'LST_Haryana',
  scale: 30,
   fileFormat: 'GeoTiFF',
  folder: 'GEE',
  region:table
});


texts: Add symbols to points in georeferenced images and save the tile with the newly added symbol

texts:
I have photoindex, each tile is a georeferenced jpg and I have a list of coordinates that can be represented as point.

If a point is contained in a tile, I must put a symbol in that tile where the point is to save the tile as jpg using the id of the point as the filename.

Is there an easy way to do this with a lot of points (thousands)?

I know some Python.

Edit: So I have found a way to do it. I'm pretty sure I can still improve some steps but it works.

Make a tile index from the rasters with QGIS tool Tile Index
Convert coordinates to points
Add polygon attributes (From the tile index layer we are gonna take the attribute "name", that is the name of the raster) to points (the coordinates)
Now some python. Iterate through all points

With the new database from the result of step 3. From the newly added
column take the name of the raster and open it with rasterio.open()
function
Then use rasterio method ".index()" this will take the coordinates to
convert them to pixel coordinates (the raster projection and the given
coordinates must be on the same projection)
Use PIL to draw a circle on the image in the given pixel coordinates

Later I'll put the code


texts: How to display a table with more than one geometry column in QGIS?

texts:
In PostGIS, it is possible to use addgeometrycolumn more than once on the same table.

Example:
I have one table "scenes" with many attributes, and 2 geometry columns: the_geom1 and the_geom2.

select addgeometrycolumn('scenes', 'the_geom1', 0, 'POLYGON', 2);
select addgeometrycolumn('scenes', 'the_geom2', 0, 'POLYGON', 2);

In QGIS, when I display "scenes", it only displays the_geom1.

To display the features using the second geometry, I found a solution with :

create view v2 as select id,the_geom2 from scenes;

But I'd like to access all other fields; something like:

create view v2 as select * (except the_geom1),the_geom2 from scenes;

Any ideas?


texts: Using proj to convert EPSG:2263 to WGS84

texts:
I have input data in EPSG:2263 with X and Y coordinates such as:

X,Y
987043,206407
982441,212359
990907,209020

(these locations are within NYC).

I need to convert these to WGS84.  Using QGIS if I import these coordinates and tell it the input data is in "EPSG:2263, NAD83 / New York Long Island (ftUS)" successfully imports the coordinates and renders them correctly (e.g. if overlaid over a WGS84 map).

It gives me this as the conversion options:

+proj=lcc +lat_1=41.03333333333333 +lat_2=40.66666666666666 +lat_0=40.16666666666666  +lon_0=-74 +x_0=300000.0000000001 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs

However, if I try to manually convert a point on the command line, for example:

cs2cs +proj=lcc +lat_1=41.03333333333333 +lat_2=40.66666666666666 +lat_0=40.16666666666666 +lon_0=-74 +x_0=300000.0000000001 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs -f '%.8f'
987043,206407
-73.99000712    40.16666623 0.00000000,206407

The resulting point is obviously not correct and (if shown on Google Earth) is off the coast of New Jersey.

I'll admit projections and coordinate systems and converting between them does still quite baffle me...

Suggestions for how I should be adjusting those cs2cs options to get a correct conversion?


texts:
I am using the following query in my db:

SELECT bg_id, ST_AsGeoJSON(the_geom)
FROM bg
WHERE bg.the_geom && ST_MakeEnvelope(xxx, xxx, xxx, xxx);

It seems to execute alright but when I draw out the geos on my map (I'm using Leaflet for this) it includes a bunch of geos that are not within the bounds I provided.

For a little more detail as to what exactly I'm trying to accomplish... I have a map on a webpage.  As the user navigates throughout the map the bounds of the map are sent back to the db and then the db should send back geos for only what is within the current view.  This seems to only be working if the map is zoomed in very close.  If I zoom out at all it begins to include lots of geos that I am not looking for.  I would only like to include geos that are ENTIRELY within the view.

Thanks.


texts: How to add *.ecw files to QGIS 2.0 on Windows?

texts: How do I add and view *.ecw raster images in Qgis 2.0 ?


texts:
I am trying to change the attributes of a table in pyqgis. However, when finishing the changes, I cannot see them in the layer.

I am using QGIS 3.8.1

Here is the simple version of my code

layer.startEditing()
from qgis.PyQt.QtCore import QVariant
if caps & QgsVectorDataProvider.AddAttributes:
    attrs = { 28 : 60.5, 29 : 70.5} #Changing the values in the
                                    #28th columns and 29th columns to 60.5 and 70.5
    if caps & QgsVectorDataProvider.ChangeAttributeValues:
        layer.dataProvider().changeAttributeValues({ 0 : attrs})
        # changing the first rows values
    layer.updateFields()
layer.commitChanges()


texts: Geographically Semi-accurately georeferencing image with limited information

texts:
It is possible to semi-accurately georeference an image given only the following information?

I'd like to do this programmatically.

Center point of terrain, Lat/Lng (WGS84)
Pitch of camera
FOV
Altitude

My intuition says there's not enough reference points, but I hope I'm wrong.

If it's possible, can you provide any pointers to the math?


texts: Combine multiple images in style of OpenLayers

texts:
is there a possibility to add a second image to this openLayers 3 style?

My code looks as follows:

var klass_red
klass_red = new ol.style.Style({
  image : new ol.style.Circle({
    radius : 11,
    stroke : new ol.style.Stroke({
      color : '#980000',
      width : 3
    })
  }),
  stroke : new ol.style.Stroke({
    color : '#980000',
    lineDash : null,
    width : 3
  }),
});


texts: How to query a data table embedded in ArcPad AXF file package using VBscript?

texts:
I would like to ask for some help on writing some VBscript (Supported in ArcPad) in a user form Command Button where the user activated the form with a unique feature ID(hence named COA_ID) already embeded, to query a data table which is embeded in the AXF file package in ArcPad 10.

The function of the command is basically displaying all records in the data table using the feature ID (COA_ID) as the P. Key. I do not know if a DSN should be established in the system to facilitate this.

Here's my alpha codes:

Option Explicit

Dim cnn
Dim rs

Dim sCOAID

Set cnn = CreateObject("ADODB.Connection")
Set rs = CreateObject("ADODB.Recordset")

cnn.Open "Provider=Microsoft.Jet.OLEDB.4.0;Data Source=C:\saFlush_t2.axf"

''onFeatureChanged = Set sCOAID=ThisEvent.Object

'execute SQL Statements as below,

SQL_query = "SELECT * FROM saFlushing_relate,
      ORDER BY "COA_ID"",
Set rs = cnn.Execute(SQL_query)

RS.Show

RS.Close
CN.Close

Arcpad 's scripting windows detects some syntax errors.


texts:
I was used to develop qgis2leaf on a Ubuntu machine and after each change on the PyQT dialog I've used "make" in the terminal to "compile" the GUI.
How does this work on a Windows machine?


texts: OpenLayers, VectorSource (WFS) adds same features again

texts:
I have some problems with OpenLayers and a WFS vector source. It is a clustered source and it works so to speak that it adds features as it should the first time. But when I pan the map it keeps adding the same features again and again. The WFS request seems to return what it's suposed but it is like it's just adding the same features again.
The result is that if a clusterd point says that it contains 5 features, then I pan the map and the WFS correctly returns these previous 5 points and then another 2 that comes inside the extent, then the map shows 10 + 2. I would expect that it should show 5 + 2 since the 5 old features are the same as the last fetched.
Any suggestion what I do wrong?
this.userMapFeatureSource = new olVectorSource({
            format: new olGeoJSON(),
            url: (extent) => {
                return this.geoserverURL + '/wfs?service=WFS&' +
                    'version=1.1.0&request=GetFeature&typename=cite:usermapsGeom&' +
                    'outputFormat=application/json&srsname=' + this.projstring + '&' +
                    'bbox=' + extent.join(',') + ',' + this.projstring;
            },
            strategy: olBboxStrategy
    });

    this.userMapFeatureClusterSource = new olCluster({
        distance: 36,
        source: this.userMapFeatureSource
    });

    this.userMapFeatureLayer = new olVectorLayer({
        visible: false,
        source: this.userMapFeatureClusterSource,
        style: (feature) => {
            const size = feature.get('features').length;
            const style = new olStyle({
                image: new olIcon({
                    src: '/assets/sysimg/btn_map.png',
                    scale: 0.9
                }),
                text: new olText({
                    font: '11px sans-serif',
                    text: size + ' st',
                    textAlign: 'center',
                    fill: new olFill({ color: 'black' }),
                    backgroundFill: new olFill({ color: '#efeee9' }),
                    offsetX: 0,
                    offsetY: 20
                })
            });
            return style;
        }
    });


texts: Line width proportional to value in attribute table in QGIS

texts:
I am visualizing trade flows in QGIS. I am trying to make the line width of approximately 100 polylines be proportional to the associated value listed for each in the attribute table. Several posts show how to classify them into groups using Natural Jenks through the Layer Properties>>Style>>Graduated field. This is not what I'm looking for. I would like a code where I can say for example line width=value of attribute field. I may also want to do (value of attribute field)/100 or log(value of attribute field) so the scale of the line width is appropriate for the map.

I tried selecting Style>>Single Symbol>>Pen Width>>Variable>>Value as well but that did not do anything.


texts: is there a way to have a routing path plugin or JS to be used for a map by qgis2web plugin of qgis?

texts:
How can I use a JS or plugin on a website, to calculate shortest path between 2 points along a road network?

I made the road network in QGIS, then exported with qgis2web plugin for my website.

I'm not sure but seems like google doesn't go over some lawns and the path will end at the road only.

While making my own network of roads in QGIS and using the road graph I can have what I want. But if I export it using qgis2web, the road graph plugin is left behind. And google APIs only give me the one in the image.

As you can see, the left picture only ends at the road using QGIS. So I made my own paths in QGIS but I can't use it in my website.


texts:
I'm using Leaflet to display map tiles. The map tiles are loaded from a different site using basic authentication.
So the url template looks like this: L.tileLayer('https://username:password@example.com/tiles/{z}/{x}/{y}.png')

Unfortunately this results in map tiles not loading, since e.g. Internet Explorer just ignores the username and password in web site adresses.

Is there a way to succesfully authenticate the Leaflet map tile requests without using a proxy but still using basic authentication?


texts: Decreasing polygon size based on percentage using ArcPy?

texts:
I would like to decrease the size of polygons by percentage and not by distance in order to implement this into ModelBuilder to loop over a number of polygons. To my knowledge this tool is available in QGIS but not in ArcGIS. The scale function does not do the trick in my case.

Does someone know how it is done?

Or may be stumbled over the same problem and could tell how he solved it or has a Python code for implementing into the ModelBuilder?


texts: Free or FOSS Elevation API for Heightmaps (image files) instead of elevation points

texts:
I'm coming from a Blender background, in which simple FOSS extensions exist (like BlenderGIS) that download SRTM/Aster Elevation Data and make a terrain based on this. I'd like to make something similar in the Unity game engine (I can easily implement web requests, it is engine-agnostic), in which the user makes a bounding-box over a map, and then the Elevation Data is downloaded in the form of a Heightmap Image file (I am not sure what this format is called, but the lowest point is black and the highest is white). A simple terrain object is then created and the image file is converted to a file Unity's terrain heightmaps can handle (something like .jpg to .png).
I would also like it so you can feed in the 4 points (x/y values of lat/long) to download the data of a box instead of manually generating a file that has all the values.


texts: convert raster to xyz without no data values

texts:
I want to convert a geotif to an xyz-file by gdal_translate. However, when I use the following statement:

gdal_translate infile.tif outfile.xyz -of xyz

the nodata values are always within the xyz-file. I want only the xyz values of pixels where data is available...

Any suggestions?


texts: Geoserver WPS process list enabling /disabling through functionality code provided by geoserver WEBINF/lib

texts: How can i enable or disable WPS processes through functionality code, i haven't found that related functionality code in geoserver WPS extension jar files those are in WEBINF/lib. Please guide me the specific location of WPS process functionality and enabling and disabling process


texts:
I use a river network on which has a lot of nodes to test pgRouting,

and I run the SQL, searching for the shortest path between nodes.

SELECT * FROM shortest_path('
            SELECT gid as id,
                     start_id::int4 AS source,
                     end_id::int4 AS target,
                     shape_leng::double precision as cost
                    FROM network',
            10, 20, false, false)

But if I want to search for the distance between node#10 and node#21, I need to type the code again and run the SQL.

The result below is parts of distance from node#10 to node#20, but I want is the total distance from starting node to ending node.

In my case, I have 159 nodes and I want to know if there are more efficient way using pgRouting？(like type the SQL code once and I can get all total distance between all nodes)


texts:
I am a new QGIS user, have installed QGIS on a Mac running Snow Leopard OSX 10.6.8 following all the instructions. QGIS seems to work well for everything other than when I boot up I get the message
"Unable to load the plugin: Unable to load required package rpy2.
Please ensure that both R, and the corresponding version of Rpy2 are correctly installed".

I have installed Rpy2 several times and reloaded R. I have R version 2.15 loaded (but also RStudio).

Does anyone have any ideas?