texts: How do you explain that band width and relative spectral response of MODIS band do not match?

texts:
To define relative spectral response, I referred to EFFECT OF RELATIVE SPECTRAL
RESPONSE ON MULTI-SPECTRAL MEASUREMENTS AND NDVI FROM DIFFERENT REMOTE SENSING SYSTEMS by David James Fleming, 2006.

  One factor that is often overlooked is the effect of a sensor’s
  relative spectral response (RSR), or spectral response function (SRF),
  on broadband spectral measurements. The RSR describes the quantum
  efficiency of a sensor at specific wavelengths over the range of a
  spectral band. Currently, general descriptors, such as bandwidth and
  average bandpass, are often the only spectral characteristics
  considered in analysis of sensor spectral measurements. However,
  cross-sensor wavelength variations in RSR can lead to measurement
  discrepancies between sensor measurements that make them not directly
  comparable (Teillet et al., 1997). In order to provide consistent
  quantitative spectral measurements of vegetation land cover and
  derived metrics, such as spectral vegetation indices, the effect of a
  sensor’s SRF must be considered and understood.

  [...]

  Spectral bands are often generalized (Pagnutti et al., 2003) in terms
  of full width at half maximum bandwidth and central wavelength
  corresponding to the maximum value of the response function (Liang,
  2004) as shown in Figure 11.

In view of this, I was a bit surprised that, for example, the bandwidth of MODIS band 7 (dotted line) as defined here, seems to not match with the relative spectral response in the way defined above.

Is my understanding of these concepts incorrect?


texts:
I've created a Visual Studio 08 Active X - MXCommands dll for a VBA button within ArcMap 9.3.1.  I'm able to load and use the .tlb on my computer; however when I try it on my test computer the .tlb file is not visible.  The MXCommand is not listed in the Component Category Manager under ESRI MXCommands on the test computer.  I used RegAsm (e.g. C:\WINDOWS\Microsoft.NET\Framework\v2.0.50727\RegAsm.exe C:\Temp\Mydll.dll) to try and register it on the test computer, but this did not help.  Both computers are running MS XP professional.  Any suggestions of what I could do or try next would be greatly appreciated.

I also tried this command to register it:

C:\WINDOWS\Microsoft.NET\Framework\v2.0.50727\RegAsm.exe C:\Temp\Mydll.dll/tlb: Mytlb.tlb/codebase

And received a RegAsm Error  - RA0000 The given path's format is not supported.


texts: Adjecency matrix with road intersaction to the power of 3


texts:
I want to use two different data sources which I fetch from a WMS server, but they do not share any CRS. For example, one of my data data sources comes in either EPSG 6706 or EPSG 4258, 25832, 25833, 25834 and Google Maps is not available in neither of these.

Is there any way around it? From everything I can read, it does not seem possible for QGIS (which is the software I'm using) to reproject on the fly data coming from a WMS, but I wish someone here faced this problem and eventually was able to solve it.


texts: ArcGIS Server Security related justification

texts:
We just installed ArcGIS server(JAVA) and the network security has a few policy that requires us to provide justification on our setting.

Does anyone have any reference document for the following questions?

1) Why must "arcgismanager" have access to the network folder? Basically they believe that noone should have access to the folder.

Where can I find the document that states arcgismanager is required to the folder for configuration of arcgis server.

2) We have a batch job, and it requires ArcGISSOM, ArcGISSOC, ArcGISManager to have access to perform the batch job.

Basically is the same as the above, we need to provide document that it is absolutely necessary for those 3 roles to have their ACL.

where can I find the documents?

3) Remove Excessive Permission from CATALINA_HOME - conf, webapps, bin, logs

Basically we set full control for all the above mention folder to all the 3 roles, they aren't happy with that and just want to bear minimum for the arcgis server(JAVA) to work.

May I know if there is a set of mininum access right for those folders to the 3 roles?

4) 2.1.1 Upgrading of tomcat server in ArcGIS (ArcGIS JAVA). They want to know if it is possible to upgrade ArcGIS server internal tomcat servers -.-

please kindly help me answer the above 4 questions with references to ArcGIS official documentation.


texts: Custom Layer not drawing using ArcObjects?

texts:
I have a Customer layer which inherits BaseCustomLayer.

I add it to my MapControl via addLayer but the function public override void Draw(esriDrawPhase drawPhase, IDisplay Display, ITrackCancel trackCancel) is never hit.

The Custom Layer appears checked in the TOC.

I think I must be missing some initialization one-liner.


texts: Why is my PYTHONPATH setting itself to QGIS?

texts:
I'm trying to cleanup Python on my Mac. I keep getting this message in Terminal:

  Your PYTHONPATH points to a site-packages dir for Python 3.x but you
  are running Python 2.x!
       PYTHONPATH is currently: "/Applications/QGIS3.10.app/Contents/Resources/python"
       You should unset PYTHONPATH to fix this.

What is going on? Why does QGIS keep setting my PYTHONPATH to itself? I keep resetting my PYTHONPATH with unset PYTHONPATH but then this keeps happening anyway.

Update 1: I think Python 2 and Python 3 both use the same PYTHONPATH, so I will have to figure out a workaround. The answers on this SE post may help: https://stackoverflow.com/questions/34066261/how-to-set-different-pythonpath-variables-for-python3-and-python2-respectively. Will post an update if I get it working...

Update 2 (and Solution): The problem is that I set my .bash_profile to QGIS. I forgot I had done this trying to solve another problem I had several months ago, which didn't work, but I never reversed the changes I made. So, problem officially solved.


texts: Creating Raster from Array using GDAL Results in Stripes

texts:
I am using a template raster to create another raster with the same resolution/projection/extent from an array. The original raster looks good, but the new raster results in stripes.

ds = gdal.Open(template_raster)
band = ds.GetRasterBand(1)
arr = band.ReadAsArray()
[rows,cols] = arr.shape
driver = gdal.GetDriverByName('GTiff')
out_ds=driver.Create(output_raster,cols,rows,1,gdal.GDT_Float32)
out_band = out_ds.GetRasterBand(1)
out_band.WriteArray(input_array)
out_ds.SetGeoTransform(ds.GetGeoTransform())
out_ds.SetProjection(ds.GetProjection())
srs = osr.SpatialReference()
wkt = ds.GetProjection()
srs.ImportFromWkt(wkt)
out_ds.SetProjection(srs.ExportToWkt())
out_band.FlushCache()
del input_array
del out_band
out_ds = None
del ds


texts: Download WMS tiles with R and use in tmap

texts:
I have access to a password protected WMS service (https://wms.swisstopo.admin.ch/) and I would like to download tiles for a bounding box with R. Then I would like to use them as background for a map created with tmap. The layer name is ch.swisstopo.pixelkarte-farbe-pk25.noscale. I don't know where to give the authentication information.
I saw this question https://stackoverflow.com/questions/49983053/extracting-data-from-a-wms-layer but using sf::st_read() gives the following error: The file doesn't seem to exist.
I would like to have the WMS data of this extent (CRS is CH1903+):
xmin       : 2581809
xmax       : 2583643
ymin       : 1181706
ymax       : 1182825

Then, I would like to plot the map with tmap::tm_rgb() (https://www.rdocumentation.org/packages/tmap/versions/3.2/topics/tm_raster). How can I do this?


texts: Is there simple way to calculate radius of curvature along polyline?

texts:
I am digitizing  a sinuous linear alignment (in this case a water pipe) and would like to get the radius of curvature of the various bends along the alignment. The main reason is that the  minimum curvature allowable for the exercise is of 200m (or approx 28 degrees along a 200m section).

I am using ArcGIS 10.1 for Desktop.


texts: Imaginary line on the Earth is straight or not?

texts: I was watching Google map. Saw a line which is joinging north to south pole. It is going between middle of China and Mexico. But, the line is not straight. Tried to search in net but not able to get anything. Can anyone check and let me know https://www.google.com/maps/@-19.6123649,-178.5132229,3z or http://pichost.org/image/UJKHD


texts: How to measure the distance between multiple points and multiple polygons?

texts:
I'm trying to find a quick and accuarate way to measure the distance between 3 survey sites and all the designated sites within 10km. The survey sites are on a point vector layer and the designated sites are polygons, both on the same CRS.

I am new to GIS and so far I have only been able to measure the distance to the nearest designated site. As I said I need to know the distance to the nearest boundary of all the sites within 10km, not just the nearest.


texts: PostGIS snap line segment endpoint to closest other line segment

texts:
I have a PostGIS table of singlepart line segments that I've run the ST_ChaikinSmoothing function on. Smoothing all of the line segments worked great, but these simplified geometries created small gaps between the other line segments in the table. (It's a network of flowlines) I need to snap the endpoints of each linestring to the closest segment of another line in the table.

I was directed to this answer about snapping lines, but I think it is snapping features in one layer (line) to features in another layer (input). I always struggle to comprehend PostGIS examples that use abstract text geometries in the queries as opposed to referencing tables of existing features.

My query below is returning an empty table with point geometry. I can successfully run each of these three functions individually to locate the points on the lines, but I can't figure out what i'm missing to snap the lines and create line geometries.

"geometry" is the name of my geometry column in the smoothed_lines table

create table snapped_lines as (
SELECT ogc_fid, name, (ST_LineInterpolatePoint(geometry,ST_LineLocatePoint(geometry,ST_EndPoint(geometry)))) as geometry from smoothed_lines);

UPDATE: Sample dataset: https://www.dropbox.com/s/wz2265symhzsenl/smoothed_lines.pgsql?dl=0

The answer from @gabriel-de-luca below does snap the smoothed lines to the original lines based on matching id's, but ST_Snap moves all vertices in a linestring, not just the end nodes and wrecks the smoothness of the input features.

example after running this ST_Snap query:


texts: PostGIS - query maximum distance between two polygons

texts:
In PostGIS, how can you calculate the maximum distance that one polygon is from the nearest point of another polygon?

I have two tables of polygons, 'red' and 'blue'. They can be joined by field district_code. How do you find the max distance from a red polygon to its joined blue polygon. (Distance indicated by dashed black line in image above).

The following query uses ST_MaxDistance, but this gives distance from furthest point to furthest point, not furthest point to nearest point.

select id, district_code, ST_MaxDistance(red.geometry, blue.geometry) AS max_distance
FROM red
INNER JOIN blue ON red.district_code = blue.district_code;


texts: OL 4.2 WMS getGetFeatureInfoUrl producing incorrect URL

texts:
EDIT: URL works fine, I mistakenly thought i,j, width, and height were sending info specific to the WMS image being checked for features.

I could be misunderstanding how the URL works so correct me if I'm wrong. The I, J, height and width attributes are all used to locate the feature on the image correct? When I provide getGetFeatureInfoUrl with different coordinates and viewResolutions, I always get the same location attributes in my URL.

When I create url using getGetFeatureInfoUrl, no matter what the coordinate or viewresolution, the i and j values are always 50. The height and width are also always 101.

Here are some examples:

"URL: http://localhost/geoserver/workspace/wfs?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetFeatureInfo&FORMAT=image%2Fpng&TRANSPARENT=true&QUERY_LAYERS=workspace%3Atable_points&LAYERS=workspace%3Atable_points&time=1500402731190&INFO_FORMAT=application%2Fjson&I=50&J=50&CRS=EPSG%3A4326&STYLES=&WIDTH=101&HEIGHT=101&BBOX=620954.1142031327%2C-34459035.34341003%2C8526377.327569202%2C-26553612.130043965"
Coordinate: -30506323.736727,4573665.720886167
ViewResolution: 78271.51696402048

"URL: http://localhost/geoserver/workspace/wfs?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetFeatureInfo&FORMAT=image%2Fpng&TRANSPARENT=true&QUERY_LAYERS=workspace%3Atable_points&LAYERS=workspace%3Atable_points&time=1500402501193&INFO_FORMAT=application%2Fjson&I=50&J=50&CRS=EPSG%3A4326&STYLES=&WIDTH=101&HEIGHT=101&BBOX=7044436.526761843%2C9118631.726308389%2C10997148.133444877%2C13071343.332991423"
Coordinate: 11094987.529649906,9020792.33010336
ViewResolution: 39135.75848201024

"URL: http://localhost/geoserver/workspace/wfs?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetFeatureInfo&FORMAT=image%2Fpng&TRANSPARENT=true&QUERY_LAYERS=workspace%3Atable_points&LAYERS=workspace%3Atable_points&time=1500402511193&INFO_FORMAT=application%2Fjson&I=50&J=50&CRS=EPSG%3A4326&STYLES=&WIDTH=101&HEIGHT=101&BBOX=-2822992.6322137713%2C-1506726.7015574016%2C5082430.581152298%2C6398696.511808667"
Coordinate: 2445984.905125633,1129718.974469263
ViewResolution: 78271.51696402048

"URL: http://localhost/geoserver/workspace/wfs?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetFeatureInfo&FORMAT=image%2Fpng&TRANSPARENT=true&QUERY_LAYERS=workspace%3Atable_points&LAYERS=workspace%3Atable_points&time=1500402511193&INFO_FORMAT=application%2Fjson&I=50&J=50&CRS=EPSG%3A4326&STYLES=&WIDTH=101&HEIGHT=101&BBOX=-5405952.6920264475%2C-10429679.635455731%2C2499470.5213396214%2C-2524256.4220896624"
Coordinate: -6476968.028772697,-1453241.085343413
ViewResolution: 78271.51696402048

Here's this HTML if it helps:

var displayFeatureInfo = function(coord, pixel)
{
    var features;
    document.getElementById('info').innerHTML = '';
    var viewResolution = /** @type {number} */ (view.getResolution());
    var url = pointLayer.getSource().getGetFeatureInfoUrl(
        coord, viewResolution, 'EPSG:4326',
        {'INFO_FORMAT': 'application/json'});
    if (url) {
      features = url;
      document.getElementById('info').innerHTML = '<iframe seamless src="' + url + '"></iframe>';
    }
};

map.on('singleclick', function(evt) {
    displayFeatureInfo(evt.coordinate, evt.pixel);
});

Edit:
I believe I found the issue. Here is the relevant openlayers code, line 141(https://github.com/openlayers/openlayers/blob/9fc0fb5e74146cd9a0c99926bed1a0f51fc745ad/src/ol/source/imagewms.js).

When creating the 'extent' openlayers always uses a default image size, which is [101,101].


texts:
I need to do a histogram from a big raster, its a full region in Chile. But I don't like the QGIS's Histograms so I decided to make it in R. I installed sp, raster and rgdal packages, but R brings me this message when I do it:

  "Warning message: In .hist1(x, maxpixels = maxpixels, main = main,
  plot = plot, ...) :   0% of the raster cells were used. 100000 values
  used."

Then I compare both histograms (in QGIS and R) and they don't looks as the same frequency. Any idea how can I fix it?


texts:
I have to have a Shapefile of France of the bigger regions: Ile-de-France, Bretagne, etc. and that is already included in my shape. But its only shown the Arrondissements (smaller regions) in France and their borders.
How can I just show the borders of these bigger regions only?


texts: How to define feature modify event in openlayers 3?

texts:
I'm using OpenLayers-3.8.2. I have modify interaction as follows:

selectedInteraction = ...
modifyInteraction = new ol.interaction.Modify({
   features: selectInteraction.getFeatures()
});

I want to define an event that when a point of feature change, it fires. Something like this, but NOT for drawing, for modifying.

How do I do?


texts: addVectorLayer('filename','layername','ogr') loads layer and names it "filename layername featuretype"

texts:
I have a PyQGIS script that loads shapefile to the canvas with this line:
iface.addVectorLayer('pathto/file.shp','layername',"ogr")

It was working fine and was always loading layer named 'layername', but now it is loading layer and naming it 'layername file Polygon', with exact same script. Do you know why? Same thing happens in both QGIS 2.18 and 3.2.


texts: Showing only one feature per Report page in QGIS 3

texts:
The gift of Reports in QGIS is a massive step in the right direction. I have set up a department template for my team to use and it looks great!

Right now, I am attempting to setup a report that will output one route at a time (I am a transport planner) to allow us to visually inspect how the routes were programmed into the network model. To do this, I used the old method that you could use for atlases by setting up a symbology style that would only show the current atlas feature:

$id = @atlas_featureid

This works great for the first section of the report. Unfortunately, the second section not so much. In the second section, the title matches the zoom (the title gives the route number and description) but the line visible is the wrong route. Is there any reason using the above expression wouldn't also work for reports with multiple sections?

I also tried I also tried "Route" = Attribute( @atlas_feature , 'Route') to see if that would fix it. I had the same problem.

Bonus: I also tried intersects( $geometry , @atlas_geometry ) and now I'm beginning to think this is a bug. See update below.

Update: I have noticed that the layer from Group 1 is still showing even though a different layer is selected for group two (as shown in the pictures below). The setup for Section 2 is provided in the screenshot below, however.


texts:
I have a number of standard projects set up for QGIS users. Some of them take more than 5 minutes to load. Things recently got much worse suddenly - and I want to track down why.
Are there any methods to check what stages of loading a project are slowing things down or/and how long QGIS is taking over each stage of loading a project?
In my particular case the issue may be related to internal network overload, and seems to be related to upgrades (at the same time) from Windows 7 to Windows 10 and QGIS 2.18 to QGIS 3.4. However, I think answers will be more useful for others/Stackexchange if they consider the wider question.


texts:
I'm using ArcMap 10.6.1. I have raster data and census tracts for Los Angeles. My goal is to take the mean value per census tract. I then want to aggregate all the cells in each census tract to be a single cell. So for, instance, I have a few thousand census tracts in my dataset. In the end I'll have a few thousand data points rather than the millions I end up with from 30 m^2 pixels.

I know that I can use the zonal statistics tool to get the mean per census tract. This is what I have been using, however, it keeps the cell size the same so that my analysis ultimately ends up with billions of data points (over 20 years). This is unwieldy to work with.

I thought the solution for this problem would be straightforward, but it doesn't seem to me that there's a tool that can do what I want. I tried using zonal statistics in conjunction with zonal fill, but this didn't work for me.


texts: TAB (MapInfo) Attribute Data not Imported to shapefile using QGIS?

texts:
I am trying to convert some TAB files into Shapefiles by using QGIS. After I import the TAB files, and make sure that these have the correct CRS, it is possible to observe that the spatial data can be exported into a shapefile.

However, and regardless of the tab file imported, the data in attribute table is always missing. Therefore, the only way for me to correctly import a tab file into QGIS, both geometrical and attribute data, is by using Universal Translator (UT).

Is it possible to import both geometry and attribute tables from TAB to SHP, using QGIS?

Here is a sample of the data I am working with (tab files). Sample TAB file


texts: Error on QGIS r.reclass for Mac and Windows

texts:
I am currently trying to use r.reclass from QGIS on Mac computer.

However, if I run r.reclass it keeps getting error message:

Traceback (most recent call last):

File "C:/PROGRA~1/QGIS3~1.8/apps/qgis/./python/plugins\processing\algs\grass7\Grass7Algorithm.py", line 428, in processAlgorithm
Grass7Utils.executeGrass(self.commands, feedback, self.outputCommands)
File "C:/PROGRA~1/QGIS3~1.8/apps/qgis/./python/plugins\processing\algs\grass7\Grass7Utils.py", line 389, in executeGrass
for line in iter(proc.stdout.readline, ''):
UnicodeDecodeError: 'cp949' codec can't decode byte 0xbc in position 11: illegal multibyte sequence
Execution failed after 0.48 seconds
Loading resulting layers
The following layers were not correctly generated.

C:/Users/user/AppData/Local/Temp/processing_173310bd80834b78955ec77e270757ff/cdd1a6c9e96240d48b22fea90a31393b/output.tif

You can check the 'Log Messages Panel' in QGIS main window to find more information about the execution of the algorithm.**

I tried using another computer to solve this problem but it seems it's not working on both Mac and Windows computers.

I wonder whether I HAVE DONE something wrongly? or it this a bug error??


texts: Is there a way to create a new, pre-populated field within a view (Feature Layer), in Arcpy?

texts:
In SQL, users can write a SELECT statement that 'invents' a new field with a default value already populated in it, for all records returned. In the example below, the field Type doesn't actually exist in any table but is made-up as a part of the query:

SELECT
  ID
  , 'Reptile' AS Type
  , Species
FROM AnimalsTable
WHERE Species IN ('Snake', 'Turtle', 'Lizard');

The query above returns a table like this:

ID  |  Type    |  Species
---------------------------
4   |  Reptile |  Turtle
6   |  Reptile |  Snake
8   |  Reptile |  Lizard
12  |  Reptile |  Turtle
13  |  Reptile |  Turtle
18  |  Reptile |  Lizard

Is there a similar way to call a feature class in Python that gives this kind of ability?

Background: in ArcGIS 10.7, I have a script that's designed to truncate an SDE feature class & then append 200,000+ records from three other feature classes (from a mix of SDE, .shp, file-geodatabase types), and all of the resulting records in the updated SDE feature class need to have the exact same value in a single field (my STATUS field needs to be staged with Incomplete, for all of them). I can't add the field & attributes to all three source data sets. The only arcpy options I can see, all require using CalculateField_management, which takes forever! Is there a more efficient way to do this?


texts: How to batch import GDAL supported raster files into GRASS?

texts:
I am trying to import some GDAL supported raster files into a binary raster map layer using the r.in.gdal module from QGis Grass tool.

As I have multiple files to import, I would like to batch import them.

Is is possible or should do I have to repeat the same operation for every single file?

I work with QGis 1.8.0 Lisboa.


texts: How to create atlas along line/route in QGIS

texts:
I'd like to create a mapbook along a given route that is designated by a single linear feature. Is it possible to have QGIS (2.2) automaptically compose a series of pdfs that over lap by a given percent and rotate the map to optimize space?

Thanks.


texts: Enumeration Problem with Eclipse and arcgisscripting?

texts:
i'm recently using Python for my Processing work that I do with ArcGis 9.3.
I figured out, that using Eclipse (with PyDev) could be helpfull for debugging purpose.
The geoprocessor is starting as demanded.

Nevertheless i got into trouble using a Example-Script from Esri:
http://webhelp.esri.com/arcgisdesktop/9.3/index.cfm?TopicName=fieldmappings_properties

fcs = gp.ListFeatureClasses("block*", "Polygon")

for fc in fcs:

This causes an error:

File "C:\Python26\lib\site-packages\win32com\client\dynamic.py", line 228, in __getitem__
    raise TypeError, "This object does not support enumeration"

The Script works fine when starting from ArcCatalog!
It seems that the variable "fcs" is of type "instance" and that python ain't able to enumerate this type.

But why doesn't ArcGIS 9.3 have a problem with this?
Does anybody know this?


texts: What is epsilon and why do we need it in GIS operations?


texts: Jagged Polygon Edge Generalization in ArcGIS Desktop?

texts:
I have a polygon which the edge is quite normal when zoomed out:

But when zoomed in this is how the edge looks like (jagged) because it was converted pixel by pixel from a raster, is there a generalization tool that can make it less jagged?

I am using ArcGIS Desktop.


texts: Conducting operation to shift DEM tilted on Z axis using PyQGIS

texts:
I have two DEMs of a stream channel created with photogrammetry from different times. One of the DEMs used ground control points and is georeferenced correctly. The other model did not use ground control points, but was georeferenced with the gps on the drone used to collect the imagery and was corrected using satellite imagery.
To correct elevation on the DEM I used a DEM of difference between the two, found the mode of the elevation values and added it to the DEM. This worked amazingly well, however there appeared to be a slight, north/south tilt to the DEM. I took several sample from a nearby road and found the mode of those, and then plotted them and found that there is a 1 to 2 meter tilt from the north end to the south end:

I am trying to use something of the form:
new_band =

for [i] in length 1:y[i]{
   zdif_south + dist_north * (zdif_south - zdif_north)/z[y[i]]
    }

I have tried to figure it out in the raster calculator and in python console in QGIS. One problem I am having is I don't know how to index the xy in raster datasets.
Has anyone done anything similar to this?


texts: PostGIS layer attribute table empty in Quantum 2.0.1

texts:
I recently installed QGIS 2.0.1 (Dufour) on a Windows 7 PC (64-bit). Today I tried, for the first time since Dufour was installed, to load PostGIS layers into QGIS. I successfully set up a connection and I could load layers into QGIS (using both DB Manager and the AddPostGIS layers tool). However, although in DBManager the tables show that they contain a number of features, I cannot see them in the Preview tab (of DBManager). Moreover, the attribute tables of the PostGIS layers I loaded into QGIS are empty (0 features, and the "show selected only" box is not ticked).

To make sure it wasn't a problem with my database (PostgreSQL 9.1, PostGIS 2.1.0SVN), I uninstalled QGIS 2.0.1 and re-installed 1.8 (Lisboa). In Lisboa, I can see my layers rendered in QGIS and their attribute tables are not empty.

Is this a known problem with version QGIS 2.0.1 or have default settings changed since version 1.8?

Thanks for your help.


texts:
How can I get the max or min value of one image?
I don't know what's the function is in GEE.
I have to use this for vegetation coverage
For example:
[1,2,3,
 4,5,6,
 7,8,9]
max:9;min:1


texts: Maps zoomed around specific locations (magnifying glass effect)

texts:
How do I create a static map where the scale varies based on location?

Example: to drive from a specific location X in Albuquerque to a
specific location Y in Chicago, I want to "zoom" around X and Y as follows:

Near X and Y, the scale is fairly large, so you can see minor
streets, buildings, etc, near X and Y.
Far away from both X and Y, the scale is fairly small, so you might
just see long lengths of highway, but no details.

In other words, it's like a regular map with a "magnifying glass
effect" around X and Y. The effect should "fade off" smoothly, so
there isn't a break in the map.

Reasoning: when traveling between remote cities, the highway travel is
often easy (in the sense it doesn't require details), but getting from
X to the highway and from the highway to Y can require detailed
instructions, knowing minor streets, etc.

EDIT: I've attached a horrible example of what I want to do. Notes:

This only zooms in on one city, not two.
This is just plain old image zoom. What I would want is more detail in the Albuquerque area.
Notice, however, that the roads/highways are continuous: the magnifying does not have a "cutoff" problem.


texts: PostGIS - Geometry or Geography find points within distance in meters

texts:
I have a requirement where I need to store lat/long coordinates for locations all around the world. I then need to query for those locations that are within a specified distance (in meters) from another location. These distances are usually within the 500-1000 meter range so I don't need crazy precision.

Following various examples I've found, I've been using SRID 4326 with my geometry column in a PostGIS database. The column simply holds lat/long coordinates.

The ST_DWithin method takes a distance in meters for geography types, but I don't know what the unit is for 4326. It doesn't seem to be meters with some queries that I tried. I'm assuming it's degrees.

I've read that for performing Spatial queries on small areas, it's best to use geometries since operations are quicker at the price of precision. How can I use these methods to provide meters for either ST_DWithin or ST_Distance?


texts: Plugin not available : point sampling tool

texts:
I installed the QGIS Standalone Installer Version 2.14 (64bit). I created a map in the program called QGIS Desktop 2.14.3 with GRASS 7.0.4.

I need one tool from Plugin : Point sampling tool. The problem is that I cannot find it anywhere:

When I check the Plugin settings, I can see that:

I am using Windows 7 enterprise 64bit.

Why the repository is unavailable? Is that the problem? How can I solve it?

I checked also the proxy:

Does anyone have the solution?


texts: Stopping GP script from adding temporary data to Table of Contents?

texts:
I have a script that does a series of calculations, producing temporary data along the way. However, at each step of the script, data is added to the TOC and then removed when I call the arcpy.Delete_management function on the data. Before my script reaches the delete function, the temporary data that I've calculated is shown in the TOC. Is there a way to explicitly specify these datasets as temporary, or freeze the TOC during script processing. Example below:

class StationTypology:

    def __init__(self, areaName, bufferDistance):
        self.areaName = areaName
        self.bufferDistance = bufferDistance

    def BusinessDensity(self):
        try:
            # Make a temporary table in memory based on query of feature class
            query = ("""\"AreaName\" = '""" + self.areaName + """' AND
                     \"BufferDistance\" = """ + str(self.bufferDistance) + """ AND
                     \"UseType\" = 'Business'""")

            arcpy.MakeTableView_management("LANDUSE",
                                           "business",
                                           query)

            # Calculate the sum of floor space, from above table
            arcpy.Statistics_analysis("business",
                                      "in_memory\\businessStatistics",
                                      [["FloorArea", "SUM"]])

            # Use SearchCursor to access table records and store rows in a variable
            rows = arcpy.SearchCursor("in_memory\\businessStatistics",
                                      "",
                                      "",
                                      "SUM_FloorArea",
                                      "")

            # Table has one row; access with iteration, grab value of SUM_FloorArea
            for row in rows:
                businessDensity = row.SUM_FloorArea

            del rows  # Housekeeping
            arcpy.Delete_management("business")  # Housekeeping
            arcpy.Delete_management("in_memory\\businessStatistics")  # Housekeeping

            return businessDensity

        except:
            try:
                del rows  # Housekeeping
                arcpy.Delete_management("business")  # Housekeeping
                arcpy.Delete_management("in_memory\\businessStatistics")  # Housekeeping
            except:
                pass

I am not sure if it matters, but the script is being called from a wxPython GUI launched inside of ArcGIS as a Python Add-In. The GUI contains a button with an event that triggers the above script.


texts:
I am trying to covert my feature class to feature layer so that i can pass "IQueryByLayer" over my feature class.

Till now i have my feature class ---

FeatureClass className = new FeatureClass(cc);

want to do something like

IFeatureLayer newLayer = className.as---(FeatureLayer)

Can anybody point out a correct way to impliment it.


texts: WFS Layer does not refresh or give new data when map moved in openlayers

texts:
I am  trying to add a WFS Layer with a BBOX strategy. When I change the extent I want GeoServer to return only the features in that extent. Here is the code I am using to add layer.

wfs = new OpenLayers.Layer.Vector("Fields_WFS", {
    strategies : [new OpenLayers.Strategy.BBOX()],
    protocol : new OpenLayers.Protocol.WFS({
    version : "1.1.0",
    url : GEOSERVER_HOST + ":" + GEOSERVER_PORT + "/geoserver/wfs",
    featurePrefix : "rpid",
    featureType : "FIELDS_SUBSET_BNG_OSGB36_1",
    featureNS : "<URI>",
    srsName : "EPSG:27700" //default
    }),
});

The map is set to a maxextent initially. And I have restricted WFS features to max 100 in geoserver settings.

My problem here is it only shows the same features every time I load and when I move or zoom +- it does not refresh WFS features.

Thanks,
JP


texts: ArcGIS Server - Polygon lines dimenstions

texts:
How to measure the lines of polygon . VB.NET + ArcGIS Server .

I have ArcGIS Server 9.3 - WebADF - VB.NET


texts: OpenLayers Enable KML feature based on user form

texts:
I have a KML layer with over 300 features, but I only want to display features based on a user form.

Each feature has a unique ID and Sector and I have added the entire KML to the map but have the visibility set to false.  How would I go about adding only the features to the map that contain the unique ID and Sector values specified after a user hits submit on a form?


texts: How to specify a port when creating a WFS service?

texts: This is a arcgis forum posting simplified.  When setting up a service using ArcCatalog 9.3.1 the URL port number is not included in the WFS GetCapabilities response.  Is there a way to correct this?


texts: What is the source of horizontal and vertical striping in USGS DEMs?

texts:
When processing 30m and 10m DEM data downloaded through the National Map Viewer from the National Elevation Dataset, we noticed horizontal and vertical striping not only in the produced results, but in merely analytical hill shades of the raw DEMs.  Does anyone know the source?  If not the source, maybe how to remove these artifacts?  These artifacts become very pronounced when using the DEMs to calculate Topographic Indices.  These artifacts remain even after depression filling occurs.

Below are images showing the striping in 30m and 10m data from watersheds in both Pennsylvania and Colorado, and a finished Topographic Index calculation showing the artifacts for a watershed in Syracuse, NY.

Colorado - HUC8 - 10190004 - 10m

Colorado - HUC8 - 10190004 - 30m

Pennsylvania - HUC8 - 02040103 - 10m

Pennsylvania - HUC8 - 02040103 - 30m

Finished TI calculation for Onondaga Creek watershed in Syracuse, NY


texts: Split feature (MultiLine) not working as expected - newist QGIS user (version 2.18.13)

texts:
I would like to split a GPX track (saved as a SHP file multiline feature) at a specific location. The track has several incidents of overlapping line segments (double-backing on trail).
When using the "Split feature" tool, QGIS splits the track at every crisscross in addition to the desired split.  I end up with 330 features, not 2 features). This question appears to have been asked before, but not answered.  What am I doing wrong, or is there a better tool?


texts:
I installed geonode following this instructions: http://docs.geonode.org/en/master/tutorials/install_and_admin/geonode_install/ In the end it say that if you log-in in geonode then geoserver is also authenticated, but in my case that did not happend. So I follow this instructions to configure the OAuth2 http://docs.geonode.org/en/master/tutorials/admin/geoserver_geonode_security/index.html and the test conection worked well. The thing is that when I try to update the style of a layer I get the error: "There was an error saving the style back to the server".

Apache logs says:

  [wsgi:error] [pid 6886:tid 140619259516672] Internal Server Error:
  /gs/rest/styles/line.xml

  ....

  [wsgi:error] [pid 6886:tid 140619259516672]   File
  "/home/geonode/geonode/geonode/geoserver/signals.py", line 133, in
  geoserver_pre_save [wsgi:error] [pid 6886:tid 140619259516672]
  'type': None} [wsgi:error] [pid 6886:tid 140619259516672]
  AttributeError: 'NoneType' object has no attribute 'attribution'

Geoserver logs says:

(When sudo service tomcat8 restart)

  INFO [geoserver.global] - Falling back to embedded data directory:
  /var/lib/tomcat8/webapps/geoserver/data

  WARN [geoserver.config] - Error reading resource URL
  [jar:file:/var/lib/tomcat8/webapps/geoserver/WEB-INF/lib/gs-geofence-server-2.9-SNAPSHOT.jar!/geofence-default-override.properties]
  java.io.FileNotFoundException: URL
  [jar:file:/var/lib/tomcat8/webapps/geoserver/WEB-INF/lib/gs-geofence-server-2.9-SNAPSHOT.jar!/geofence-default-override.properties]
  cannot be resolved to absolute file path because it does not reside in
  the file system:
  jar:file:/var/lib/tomcat8/webapps/geoserver/WEB-INF/lib/gs-geofence-server-2.9-SNAPSHOT.jar!/geofence-default-override.properties

  ....

(When update the style)

  ERROR [geoserver.security] - Could not Authorize OAuth2 Resource due
  to the following exception:
  org.springframework.web.client.ResourceAccessException: I/O error on
  POST request for
  "http://localhost:8000/api/o/v4/tokeninfo/?access_token=gFq6vD2g21eNMAQ9HILSsNEwqXo8Gb":
  Connection refused (Connection refused); nested exception is
  java.net.ConnectException: Connection refused (Connection refused)

  ....

Could you give some recommendations, what things should I check?


texts: Tool or process to automatically snap endpoints of polylines within given proximity using ArcMap?

texts:
I created a network of polylines that were all supposed to be connected.

However, I just realized that many of the endpoints of the lines aren't connected, but are very close.

I think there should be a tool to snap the endpoints of the polylines if they are very close to each other, but I can't seem to find it.

Does anyone know an easy way to do this?


texts: How to calculate the area of multiple polygons and show it in another feature in ArcGIS?

texts: I am trying to calculate the area of vegetation (A) per municipality (B) for all state (645 municipalities). How to do it respecting the boundaries of each polygon (B)? The total area per municipality should appear in the B layer. I tried something using the spatial join tool, but there is no option for areas. I can calculate the number of polygons, but no the area.


texts: OpenLayers: Strange Tiling Issue (Checkered, white and grey) when panning around map

texts:
So this strange tiling issue is only occuring when moving around the map.

Some things to note:

1) You can see the markers, airplane icon, and flight paths. These are all overlays, where the actual world map is the base layer (TMS layer).

2) Before adding these vectors/markers, there was no checkered grey/white tiling issue. You could pan around and everything was super smooth.

Edit: to see a larger version of the map, right click the image, copy image url, and paste it into a new browser window.

Any ideas on why this may be? I understand the tiles (pngs) have to be essentially downloaded, but it seems strange when adding just vectors and markers would cause this to happen.


texts:
Has anyone made an anaglyph photo from Nadir geospatial data?  I have found some excellent examples of making anaglyphs in GIMP from digital photography.  I would like to do a similar process with Nadir geospatial data.  I have a 11km x 11km Nadir image of a TIN that I would like to convert into an anaglyph.  I was wondering if I could apply a shift to one of the images (ie: shift one image X metres west and X metres north) to create the "left" image of the anaglyph.  My question is related to the amount or percentage of the shift needed.  Is there a mathematical way to calculate this, if I know my resolution, extent and scale?

My other thought was to simulate overlap by creating two images that have 60% or 80% overlap between the images.  Has anyone successfully tried this?

Any help would be appreciated.

Thanks...


texts: What are your biggest challenges as a GIS developer?

texts:
What are your biggest challenges when developing GIS software?

Is it coding? Is it understanding cartography/geography/etc concepts (like projections)? Or other difficulties?


texts: Converting GeoJSON to Shapefile using ogr2ogr?

texts:
I have a GeoJson file that I want to convert to ESRI Shapefile using ogr2ogr in FWTools.

I tried using

ogr2ogr -F "ESRI Shapefile" filename.shp geojsonfile.json OGRGeoJSON

I got output

  ERROR 6: Can't create fields of type StringList on shapefile layers. ERROR 1: Attempt to write non-polygon (POINT)
  geometry to  type shapefile. ERROR 1: Terminating translation
  prematurely after failed translation of layer OGRGeoJSON

GeoJSON input we are giving is in standard geojson format.

content of geojson are like  {to share I took only a part from features array }

{"id":1223,
"properties":{"name":"3","level":2},
"map_version":11,"type":"FeatureCollection","x_id":1234,"obj_type":"LevelGeometry",
"features":[

{
"my_area":[131.40093253528894,33.58879853190332,6.515867411008615,3.490549810987222,-0.6237271428108215],
"id":1984791,
"location":{"type":"Point","coordinates":[131.4009325,33.5887985]},
"properties":{"display_name":"xyz"},
"type":"Feature",
"geometry":{"type":"Polygon","coordinates":[[[131.4009508,33.5888314],[131.4008895,33.5887946],[131.4009148,33.5887653],[131.4009761,33.5888020],[131.4009508,33.5888314]]]},
"obj_type":"Geometry"
}

]
}


texts: Conversion differences going from EPSG:4326 to EPSG:26710?

texts:
There are two coordinate systems in question:

"CRS 1" = EPSG:4326 (good ol WGS84)

"CRS 2" = EPSG:26710 (UTM NAD27 Zone 10 (specifically, zone 10S but they don't seem to differentiate)

For the coordinate in question, in "CRS 1" = EPSG:4326:

39.3662N, -120.35175E (i.e. 120.35175W) = NW summit of Castle Peak in Nevada County CA USA

Converting to CRS 2 using various programs (rounding to nearest meter):

All three of these have very good agreement (within 10m, and that's mostly rounding error):

QGIS 2.6:                07 28 241 E   43 60 566 N
caltopo.com:             07 28 239 E   43 60 561 N
Terrain Navigator Pro:   07 28 241 E   43 60 565 N

These two are very far away from the others and from each other:

Locus Pro (android app): 07 28 153 E   43 60 552 N
http://twcc.free.fr:     07 28 145 E   43 60 785 N

Note for Locus Pro I added EPSG 26710 into its config file, which is a documented and supported thing to do.  After emailing with the developer he confirms it is (or should be) using the correct proj4 code:

EPSG:26710 : +proj=utm +zone=10 +ellps=clrk66 +datum=NAD27 +units=m +no_defs

Which is also what shows up when you view the EPSG info in twcc.

Any idea what I'm missing?  Hopefully just something obvious?  Bottom line is that I'd like to get accurate conversion from WGS84 to and from UTM NAD27 CONUS zone 10S on the Android.


texts: How do I separate OpenLayers methods from my own, while using inheritance?

texts:
I am working on a GIS API, that uses OpenLayers. I initialize the map lke this:

xgis.map = OpenLayers.Class(OpenLayers.Map, {
    var mapOptions = {};
    mapOptions.allOverlays = false;
    /*
        more mapOptions here
    */
    OpenLayers.Map.prototype.initialize.apply(this, ['map', mapOptions]);
});

Of course, I do this because I want to create my own methods, like this one, for example.

xgis.map.prototype.zoom(action) {
    switch(action) {
        case 'in:
            this.zoomIn();    // Uses the OpenLayers zoomIn method
            break;
        case 'out:
            this.zoomOut();    // Uses the OpenLayers zoomOut method
            break;
};

The problem is this 'zoom' method I created is not recognized as a function (I get: number is not a function), because OpenLayers has a 'zoom' numeric value where it stores the zoom.

Is there a way to use inheritance and access the parent's methods like this?

this.parent.zoomIn()

This would definitely fix the 'zoom' name problem.


texts: Leaflet error: Uncaught TypeError: Cannot read property '_controlCorners'

texts:
I'm trying to implement the Leaflet choropleth example on a map I already built.

I've followed the instructions and seen many implementations of the i.control and control hover option on Github but I can't figure out what my type error is.

Uncaught TypeError: Cannot read property '_controlCorners'

I dug into the leaflet.js code but frankly can't figure out what controlcorners is doing.

It has a problem with the line "info.addTo(map);" (see below)

I've console.logged info and it is defined.  Where is this undefined part then?

Error in Chrome console:

Uncaught TypeError: Cannot read property '_controlCorners' of undefined
o.Control.o.Class.extend.addTo
(anonymous function)

Code:

var map;

var mapquestOSM = L.tileLayer("http://{s}.mqcdn.com/tiles/1.0.0/osm/{z}/{x}/{y}.png", {
    maxZoom: 19,
    subdomains: ["otile1", "otile2", "otile3", "otile4"],
    attribution: 'Tiles courtesy of MapQuest . Map data (c) OpenStreetMap contributors, CC-BY-SA.'
});

// experiment with info.update

            function style(feature) {
                return {
                    fillColor: getColor(feature.properties.POBTOT),
                    weight: 2,
                    opacity: 1,
                    color: '#666',
                    dashArray: '3',
                    fillOpacity: 0.7
                    };
                }

                var geojson;

                function highlightFeature(e) {
                    var layer = e.target;

                    layer.setStyle({
                        weight: 3,
                        color: '#F00',
                        dashArray: '',
                        fillOpacity: 0.7
                    });

                    if (!L.Browser.ie && !L.Browser.opera) {
                        layer.bringToFront();
                    info.update(layer.feature.properties);
                    }
                }

                function resetHighlight(e) {
                    geojson.resetStyle(e.target);
                    info.update();
                }

                function zoomToFeature(e) {
                    map.fitBounds(e.target.getBounds());
                }

                function onEachFeature(feature, layer) {
                    layer.on({
                        mouseover: highlightFeature,
                        mouseout: resetHighlight,
                        click: zoomToFeature
                    });
                }

                var info = L.control();

                info.onAdd = function (map) {
                    this._div = L.DomUtil.create('div', 'info');
                    this.update();
                    return this._div;
                };

                info.update = function (props) {
                    this._div.innerHTML = '
Población en México
' + (props ?
                        '' + props.NOM_ENT + '
' + props.POBTOT + ' personas'
                        : 'Mueva el cursor sobre el mapa');
                };
                console.log(info);

                info.addTo(map);


texts: Using ogr2ogr to convert GML to shapefile in Python?

texts:
I am trying to convert a GML to an ESRI shapefile using ogr2ogr utility in a python script.

I have successfully installed the GDAL/OGR package via osgeo but I am now struggling to find/understand any details on the syntax for using ogr2ogr in python.

After importing ogr, all I have found is ogr2ogr -f "ESRI Shapefile" output.shp input.gml

How do I use this in python? Is it as simple as assigning the input and output files?

Everything I try I just get syntax errors. Any pointers in the right direction to get me started would be great.


texts: How to get current scale of a Leaflet map?

texts:
I've been looking through the Leaflet API to see if there is a getScale method or scale property of a map, and have't seen anything.

Does anyone know if Leaflet provides a way to get the current map scale?


texts:
According to step 4 of the MapCruncher sample,
http://help.arcgis.com/en/sdk/10.0/arcobjects_net/conceptualhelp/index.html#/Sample_Multithreaded_MapCruncher/000100000mmm000000/

Once the cache is finished processing, open Engine or Desktop, navigate to the cache folder, and add the cache as a raster.

Basically I want to use the precooked maps in a dynamic display environment, but I'm not clear what they mean by "add the cache as a raster".  How can I do this and do you think it would even work?


texts:
I have some tables with columns that I need to join/copy to another table, in case an attribute table. I was using join but I realized that it would be a better solution to add new fields to my initial table (attribute table) and then fill those fields with the columns that I need from the other tables than to join that tables into my initial table.
Is there any tool that does that? Copy a column from one table to another?
My problem is to distinguish the columns, I want to have specific names for each column, therefore I was thinking about create the fields and then fill them with the data from other tables.


texts: Creating Sentinel-2 cloud free, cloud-shadow free composite or scene on Google Earth Engine (GEE)?

texts:
Using this code, I am able to get a cloud-free composite, but it still includes shadows in my image.

var collection = ee.ImageCollection("COPERNICUS/S2")
  .filterBounds(ee.Geometry.Rectangle(-85.14404296875, 46.9502622421856,
  -71.87255859375, 41.60722821271717));
var S2_cloud = collection.min();

How do I create a Sentinel-2 cloud free, cloud-shadow free composite or scene on Earth Engine (GEE)?


texts: Avoid multiple labeling in qtiles / QGIS

texts:
After failing to start the slippy map tiles generator .py module (where to put globalmercator.py...?), i tried this plug-in.

How can i set qgis(2.14), to not move labels, as the tile generator moves? So this way only once will be printed everthig on the features?
Everything is zoom dependant on the map, features were created nicely, but labels are cutted, or multiplicated.


texts: How do I get the centroid of line with OGR

texts:
Can somebody explain how to get the centroid of a feature with ogr? I have no idea what the centroid() function itself does.

>>> linelyr = lineshp.GetLayer()
>>> f = linelyr.GetFeature(0)
>>> g = f.geometry()
>>> g.GetPoint()
(5.85136, 49.974, 0.0) # starting point of the line

c =g.Centroid()
>>> c.GetX()
0.0


texts: Sentinel 1 Extra Wide in Google Earth Engine?

texts:
I have managed to display Sentinel-1 IW images in Google Earth Engine following this example but when I change the instrumentMode to EW nothing shows in the map canvas and it says "No bands to visualise". I have tried changing dates etc but I get nothing. I guess I'm doing something wrong, but I have no clue what. Any ideas?

// Load the Sentinel-1 ImageCollection.
var sentinel1 = ee.ImageCollection('COPERNICUS/S1_GRD');

// Filter by metadata properties.
var vh = sentinel1
  // Filter to get images with VV and VH dual polarization.
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VV'))
  .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VH'))
  // Filter to get images collected in interferometric wide swath mode.
  .filter(ee.Filter.eq('instrumentMode', 'EW'));

// Filter to get images from different look angles.
var vhAscending = vh.filter(ee.Filter.eq('orbitProperties_pass', 'ASCENDING'));
var vhDescending = vh.filter(ee.Filter.eq('orbitProperties_pass', 'DESCENDING'));
var vhDesc2017 = vhDescending.filterDate(ee.Date('2017-02-10'), ee.Date('2017-02-16'));
var vhAsc2017 = vhAscending.filterDate(ee.Date('2017-02-10'), ee.Date('2017-02-16'));

// Create a composite from means at different polarizations and look angles.
var composite = ee.Image.cat([
  vhAsc2017.select('VH').mean(),
  ee.ImageCollection(vhAscending.select('VV').merge(vhDescending.select('VV'))).mean(),
  vhDesc2017.select('VH').mean()
]).focal_median();

// Display as a composite of polarization and backscattering characteristics.
Map.addLayer(composite, {min: [-18, -18, -18], max: [-30, -30, -30]}, 'composite');


texts: Choosing diffusion coefficient and other parameters for r.sim.water in GRASS7?

texts:
I'm trying to run r.sim.water in GRASS7 to simulate overland flow during a heavy rain event.
My DEM has a resolution of 1m. I prepared it by adding values for buildings and included ditches and culverts using r.carve. I calculated xy-derivatives using r.slope.aspect.

My study area is about 8,5 sq km. I was not able to run r.sim.water with 1m resolution so I changed the resolution to 2m using r.resamp.rst. Including NULLS, the input raster size is now about 4.000.000 cells. r.sim.water did still not launch so I changed the number of walkers to the exact number of cells (instead of twice the number).
Now I get values for water depth that are far higher than the actual depression underlying the cells, which is physically not possible. Changing the value for the Manning coefficient even increases those "water mountains".

I think, to avoid this I have to change the diffusion coefficient and/or the values for hmax, halpha and hbeta. However, I have no idea of how they actually influence my result or what could be reasonable values for my purpose. Does anyone have experiences with that or can tell me where I could find information about it?

Besides, the water depth after 10min is far higher if I set the duration to 60 or 120 min as if I chose, for example, only 30 min, which also doesn't seem logical to me. Can anyone tell me why this could happen?

I did not find any information about a maximum number of cells for the input DEM. However, I think r.sim.water cannot handle large areas or high numbers of walkers. Does anyone know if there are limitations for the calculation?


texts: How to get originCoordinate value in "Create Fishnet" tool

texts:
I'm trying to develop a Spatial data processing GUI, one of the step uses the "Create Fishnet" tool. From the toolbox it self, there is a drop down let you select template extent, this way user won't need to manually input the X and Y coordinate for fishnet origin and Y-Axis.

Though, from the script, it doesn't seem like there is a way to browse a reference extent. Does anyone know how to do this?


texts: Transforming data from a rotated pole lat-lon grid into regular lat-lon coordinates in QGIS?

texts:
I am working with data from a numerical weather model, the output of which is given in grid cells where I have the lat-lon coordinates of each cell's center point. The lat-lon specifications are, however, in a rotated coordinate system with the south pole at (lon: 10, lat: -40). The problem I need help with is how one transforms a grid (a vector layer) which is defined in the rotated lat-lon coordinates into the regular lat-lon coordinates inside QGIS. The figure below shows my specific problem of a rectangular grid in pole rotated lat-lon coordinates (off the coast of West Africa) and the target points in regular lat-lon coordinates over Denmark. Notice that the square grid doesn’t simply get translated (moved) and rotated, there is also a lattitude-dependent scaling issue.

A previous thread on SE shows how the transformation can be calculated for points with the help of some simple trigonometry (Lon/lat transformation). These calculations work perfectly (for points) when I implement them in R, but I need to do this transformation for a vector layer inside QGIS. It is common for meteorological offices to only specify the coordinates of the South or North Pole in the rotated lat-lon grid. The QGIS-based transformation therefore has to be done the way simondk does it in his reply to the thread I referred to (he implemented a short piece of Matlab code that does this here.

I've tried to do the transformation with the Affine Transformation plugin where I calculate the parameters for the plugin based on this previous thread. This procedure almost gets it right - the grid is rotated to the right general location (above Denmark), but there seems to be something wrong with the scaling. The center points of the transformed grid cells are not exactly at the correct location - this is clearly seen in the corners of the grid. In the image below, my transformed grid is red, while the correct grid center point locations are the black crosses.

I'm not sure what goes wrong here. Can Affine Transformation even do this lat-lon rotation properly?

I've done extensive Googling and there seem to be many people out there with similar problems, but no solutions. Here is a link to a dropbox folder that contains the data I’ve used for the figures, in case you want to test a solution on the actual rotated grid and target points I’m working with.

I prefer QGIS.


texts: 2D Transformation of points in ArcGIS Desktop?

texts:
I have a table of points with x and y coordinates in UTM along with other attributes. I need to display these points in a local system for which I do not have the parameters. I however have the coordinates of two control points (A1 and A2) in the local system and I have used a GPS to obtain their UTM coordinates as well. Is there a way  to carry out a 2D transformation of the points in ArcGIS 10 using the coordinates of the two control points in both systems?

NB: This can be done in surpac from 2d transformation as str file under file tools but it eventually turns up as  line features whereas I need the output as point features.


texts: How to decide how to "stretch" raster imagery?

texts:
I have to symbolize a digital elevation model. ArcMap offers many different setting for stretching the symbology which all give widely different results. Some stretch methods symbolize the data so that there appears to be more elevated areas than other stretches.

Without knowing too much about the area I am symbolizing, how can I be sure that I am using a correct setting? I don´t want to give the impression that the relief is higher than it really is, and I also don´t want to make the area unrealistically flat.

Here is my Histogram:
"The histogram displays the pixel values on the x-axis and the count on the y-axis.When there is no histogram stretch, the input pixel values equal the output stretched pixel values. The greater the stretch, the larger the difference between the input pixel values and the output pixel values."...ESRI

I can see that there are two peaks. The first peak for the flat areas, the second peak for the higher areas. In the default symbology, the flatter areas are very pale.

What should I be aiming for here? To show as much relief detail as possible? or to give a realistic interpretation of the relative heights? If I stretch the lower values, I can show more detail in the flat areas..but to the detriment of the higher areas. How would other specialists visualize these data?


texts: Can MapQuest tiles be used in an Openlayers 3 map with a non-standard projection?

texts: I have a map in OpenLayers 3 with raster and GeoJSON layers that are in an UTM projection and I would like to use MapQuest tiles as a base map while keeping the map and layers projected as UTM.  Is this possible?  I'm using proj4.js to define the UTM projection for the map if that makes any difference.  The raster layers are being pulled from a WMS service.


texts: Layer with correct CRS still displays misaligned with other layers in the same projection

texts:
I understand this type of question has been asked many times, but I'm still struggling to solve my problem.

I'm trying to import the Global Self-consistent, Hierarchical, High-resolution Geography Database (GSHHG) which I downloaded as a shapefile. Additional info (and the dataset is also available through NOAA.

After downloading and extracting the GSHHG archive, I create a new project in QGIS 2.18 and import the "GSHHS_f_L1.shp" shoreline shapefile. This file should be the full resolution (indicated by "f") version of the data containing the boundary between land and ocean, except Antarctica (indicated by "L1") and is in WGS84 (all this according to NOAA).

After importing into a blank project, the project CRS is EPSG:4326, with on the fly reprojection disabled. Which I think is as expected. If I then add the OSM Standard layer using the "QuickMapServices" plugin, OTF reprojection gets enabled but when I look at the shorelines anywhere in the world, they do not align. They are off some small amount consistently everywhere.

If instead of adding an OSM layer using that plugin, I add a layer I have in the same CRS, I still have the same issue. So, it seems to me that the GSHHG dataset is for some reason misaligned. I've tried setting that layers CRS manually to 4326, but it's still incorrect.

Example:

Obviously, I'm still doing something wrong. How to correct this misalignment and why is this happening?


texts:
I am trying to obtain the shapefile output that I am receiving from the gpTool, which is hosted on the ArcServer, and placing it on the Map using,

gp.getResultImageLayer(jobInfo.jobId,"final_shp", imageParams, function(gpLayer) {
            map.addLayer(gpLayer);
          });

where gp is,

gp = new esri.tasks.Geoprocessor("http://domain/ArcGIS/rest/services/webGIS/Shock_Models/GPServer/Income_Log");

This is how it has been explained in the resource but not sure why it is not working.

The firebug gives me,

  NetworkError: 500 Internal Server Error -
  http://webgis.arizona.edu/ArcGIS/rest/services/webGIS/Shock_Models/GPServer/Income_Log/jobs/j10905b1a50f141e0b62f1ab8ef230f51/results/final_shp?dpi=96&transparent=true&format=png8&imageSR=102100&f=image&bbox=%7B%22xmin%22%3A-14519891.438356264%2C%22ymin%22%3A637618.0139790997%2C%22xmax%22%3A-6692739.741956295%2C%22ymax%22%3A6507981.786279075%2C%22spatialReference%22%3A%7B%22wkid%22%3A102100%7D%7D&bboxSR=102100&size=800%2C600

And when I goto that particular link, I get this error,

  Server Error - Object reference not set to an instance of an object.

Any idea how I can correct it?

UPDATE

I was told to use fiddler and see the error details that is occurring over the network and this is the output that I got,

  SESSION STATE: Done. Response Entity Size: 849 bytes.

  == FLAGS ================== BitFlags: [ClientPipeReused, ServerPipeReused] 0x18 X-CLIENTPORT: 2010
  X-RESPONSEBODYTRANSFERLENGTH: 849 X-EGRESSPORT: 2023 X-HOSTIP:
  128.196.53.161 X-PROCESSINFO: firefox:2248 X-CLIENTIP: 127.0.0.1 X-SERVERSOCKET: REUSE ServerPipe#2

  == TIMING INFO ============ ClientConnected:  15:53:51.383 ClientBeginRequest:    15:53:51.494 GotRequestHeaders: 15:53:51.494
  ClientDoneRequest:    15:53:51.494 Determine Gateway: 0ms DNS Lookup:
        0ms TCP/IP Connect: 0ms HTTPS Handshake:    0ms
  ServerConnected:  15:52:45.077 FiddlerBeginRequest:   15:53:51.495
  ServerGotRequest: 15:53:51.495 ServerBeginResponse:   15:53:51.679
  GotResponseHeaders:   15:53:51.679 ServerDoneResponse:    15:53:51.679
  ClientBeginResponse:  15:53:51.679 ClientDoneResponse:    15:53:51.679

  Overall Elapsed:    00:00:00.1850106

  The response was buffered before delivery to the client.

  == WININET CACHE INFO ============ This URL is not present in the WinINET cache. [Code: 2]
  * Note: Data above shows WinINET's current cache state, not the state at the time of the request.
  * Note: Data above shows WinINET's Medium Integrity (non-Protected Mode) cache only.

But I am still confused as to what the error is? This is the application.


texts:
I have a raster in WGS 1984 coordinate system with the pixel values denoting categorical values i.e. the raster is a land use land cover raster, with each value meaning water bodies, ice, forest etc.

I need to convert it to the unamed ellipse sinusoidal projection (with a linear unit of meters), since my remaining rasters are in this projection. I have tried-

projectRaster(raster in WGS 1984 system, crs(raster in sinusoidal projection), method='ngb' since my raster is categorical as mentioned above).

I get this error message- Error in (function (classes, fdef, mtable)  :
  unable to find an inherited method for function ‘res’ for signature ‘"CRS"’.

In context to this question, what does the error message mean and how do I solve it such that the projection works.


texts: Display polygon colours as percent reflectance of sunlight

texts:
I have a set of polygons which I would like to display in different colours.

The only information I can give on the colours is the percent reflectance of sunlight in different channels, e.g.:

B:  0.45 - 0.51µm, 5.1%

G:  0.53 - 0.59µm, 7.2%

R:  0.64 - 0.67µm, 6.1%

NIR:0.85 - 0.88µm, 38.7%

(Landsat 8 OLI channels)

Is there any chance that I can get an RBG or CIR composite, when the values are stored like this in the attribute table of a shapefile?


texts: Deleting last vertex using ArcMap Keyboard Shortcut?

texts:
Is there a keyboard shortcut to delete the last vertex created when making a line or polygon in ArcMap?

Basically the same functionality as right clicking and then picking "Delete Vertex"


texts: Having two different fields symbolized to give four options in ArcGIS Online, Collector

texts:
I have a integer field for "agency", and a string field for Survey Status "SURV_STATU"
Each field has two options. I used this Arcade line to see the desired symbology in ArcGIS Online.
Concatenate([$feature.agency,$feature.SURV_STATU], '')
Since I used Arcade to make that symbology change, the layer is not visible in Collector. Can anyone explain why, or show another way to symbolize the 4 outcomes of one layer?


texts:
I am trying to convert multiple raster (name of raste: 01_sr, 002_sr....012_sr) to vector through loop function and output of vestor name wil same as of input. but loop function is not working. can anyone guide me on this subject, how to fix this error?

Error Massage:
Message File Name   Line    Position
Traceback
    <module>    <module1>   15
AttributeError: 'module' object has no attribute 'RasterToPoint_Conversion'

Code:

import arcpy, os
from arcpy import env
arcpy.env.overwriteOutput = True
env.workspace = "D:\SOLAR-RADIATION-NOAA"
outpoint = r"D:\SOLAR-RADIATION-NOAA"
field = "VALUE"
Rasters = arcpy.ListRasters ()
for inRaster in Rasters:

     outpoint1 = r"D:\SOLAR-RADIATION-NOAA\\"+ "inRaster"

     arcpy.RasterToPoint_Conversion(inRaster, outpoint1, field)

Updated code, still error exist

import arcpy, os
from arcpy.sa import *
from arcpy import env
arcpy.env.overwriteOutput = True
env.workspace = "D:\SOLAR-RADIATION-NOAA"
outpoint = r"D:\SOLAR-RADIATION-NOAA"
field = "VALUE"
Rasters = arcpy.ListRasters ()
for inRaster in Rasters:

     outpoint1 = r"outpoint"+ "inRaster"

     arcpy.RasterToPoint_Conversion(inRaster1, outpoint, field)

Edited Code , Error still not solved

import arcpy, os
from arcpy.sa import *
from arcpy import env
arcpy.env.overwriteOutput = True
env.workspace = "D:\SOLAR-RADIATION-NOAA"
outpoint = r"D:\SOLAR-RADIATION-NOAA"
field = "VALUE"
Rasters = arcpy.ListRasters ()
for inRaster in Rasters:

     outpoint1 = outpoint+ inRaster
     arcpy.AddMessage(inRaster)

     arcpy.RasterToPoint_Conversion(inRaster, outpoint1, field)


texts: Is there any way to use the geoserver kml reflector in cesium viewer? I have millions of features with height value. The geoserver kml reflector works in google earths, but not in cesium.


texts: PostGIS ST_Buffer with attribute preserve

texts:
I'm playing with PostGIS a little bit but I'm stuck on, I think simple thing - how to create buffer from points with attribute inheritence?
Simple st_buffer code:

CREATE TABLE buffer AS
SELECT ST_Buffer(geom,1) AS geom
FROM ne_110m_populated_places;

Original point layer (ne_110m_populated_places) have 100 columns, result has only geom column.
I'm trying to put rest of attribute table from source to target layer but I'm stuck.
I can easly select columns with st_intersect:

SELECT ne_110m_populated_places.*, buffer.*
FROM ne_110m_populated_places
INNER JOIN buffer ON ST_Intersects(ne_110m_populated_places.geom,buffer.geom);

But can't connect both to do the task.


texts: AWS Sentinel 2 incorrect cloudy pixel percentage

texts:
For some tiles with low data coverage percentage metadata and tileInfo contain cloudy pixel percentage value equal 0.

Examples:

http://sentinel-s2-l1c.s3-website.eu-central-1.amazonaws.com/#tiles/37/U/DP/2017/10/14/0/
http://sentinel-s2-l1c.s3-website.eu-central-1.amazonaws.com/#tiles/37/U/DP/2017/4/17/0

Is cloud coverage calculated only for scenes with data coverage larger than some value?


texts: Adding Image Attachments to ArcGIS Pro Report

texts:
I am trying to create a report that includes image attachments for each point in my feature class using ArcGIS Pro 2.4. I found some documentation explaining how to do it in Desktop using the Report Designer.

Is there similar functionality in Pro or a way to program it in Python?


texts: Postgis query to sort geometries by latitude

texts:
I am looking for an efficient way to get geometries sorted by their latitude (northmost in case of polygons) when querying a postGis database

Something that could work in an order by clause which possibly works independent of the geometry type (as long as the CRS is lat/lon based)

I think the following would be the ideal writing for it

order by ST_North(geometry)

but I am not aware of such an operator


texts: Is it possible to access the Go To XY Tool in ArcObjects?

texts: I've looked around I but can't seem to find a coclass/interface that would allow me to access it. I know I could make something similar but I'm curious about interacting with the built-in tool to save some work.


texts: I am a very green when it comes to publishing features and tools in ArcGIS Server. Do I need to know any other programming language apart from Python. My final aim was to come up with a service that can be accessed via any web browser and used by people who have no GIS background to carry out buffering.


texts: Create shapefile line from a list of point coordinates in a stand-alone python script

texts:
I've used pyqgis to calculate the shortest path through a network between two points. I want to turn this shortest path into a polyline.

I've managed to use the QGIS dijkstra function to calculate a shortest path between two points in a network.
If I print the resulting variable point_route as a string, it prints a list of coordinates.

Where do I go from here?

Is there a(nother) well-documented python library that I can use to turn this list of coordinates into a single polyline shapefile? I can't seem to find/understand the right instructions for pyqgis.

Relevant code
(copied from the shortest path chapter of 'pyqgis developer cookbook'):

(tree, cost) = QgsGraphAnalyzer.dijkstra(graph, startId, 0)

if tree[endId] == -1:
    print "Path not found"
else:
    point_route = []
    curPos = endId
    while curPos != startId:
        point_route.append(graph.vertex(graph.arc(tree[curPos]).inVertex()).point())
        curPos = graph.arc(tree[curPos]).outVertex()
        print "appended " + str(curPos)
    point_route.append(tStart)


texts: Loading different image service functions to QGIS

texts:
If someone has created an image service in ArcGIS with function templates in order to process on the fly... How would you bring that data into QGIS and be able to switch between functions?  In ArcMap you open the properties and select the functions tab... But how would one do this in an OpenSourceGIS utilizing the Image Service?

PS it can be accessed through a WMS/ WFS

But is the 'functions template' an ESRI only product when dealing with image services??


texts:
I am trying to find out the counties that don't contain any stores in ArcGIS using Python. Basically, I have a point layer (representing the stores) and a polygon layer (counties). I have managed to write some code to find out the counties that DO contain the stores (it was pretty simple actually). The code is below.

import arcpy
arcpy.env.overwriteOutput = True
path="C:/Users/XARDAS/Documents/ArcGIS/Packages/Romania1000k_9E5B7FEC-6005-4D3A-81EA-E95FAACEF69E/v101/ro1mil.gdb"
arcpy.MakeFeatureLayer_management(path+"/Counties", "Counties_lyr")
arcpy.MakeFeatureLayer_management(path+"/Stores", "Stores_lyr")
arcpy.SelectLayerByAttribute_management("Stores_lyr", "NEW_SELECTION","Type=1")
arcpy.SelectLayerByLocation_management("Counties_lyr","INTERSECT","Stores_lyr",0,"NEW_SELECTION")

So this gives me the counties that have stores but I would like to somehow inverse the intersection for the program to give me the ones that don't have any stores. I have thought about just deleting the selected counties but I don't think that it would be too nice.


texts:
I have about 300,000 CSV files which I want to convert to shapefiles, but ArcGIS and QGIS can only convert one CSV at a time.

Are there any ways to batch convert my CSVfiles?

x y z1 z2 z3 z4
1 2 3  4  5  6

My csv files just have a few columns like above

Can any one help me?


texts: Add feature manually to a vector layer in ol3

texts:
I'm trying to add a layer manually to a vector layer with javascript.  I can't seem to determine why this is failing:

http://jsfiddle.net/Kieveli/f4t6n6v1/4/

I've tried sane coordinates like 16,22, and big ones to match the view's xy values. I get a javascript error from ol3 : TypeError: b.Q is not a function.

HTML:

Javascript:

var vectorSource = new ol.source.Vector({});

var map = new ol.Map({
  layers: [
      new ol.layer.Tile({
          source: new ol.source.MapQuest({layer: 'sat'})
      }),
      new ol.layer.Vector({
          source: vectorSource
      })
  ],
  target: 'map',
  view: new ol.View({
    center: [-11000000, 4600000],
    zoom: 4
  })
});

var thing = new ol.geom.Polygon( [ [16000000,22000000],[44000000,55000000],[88000000,90000000] ] );
vectorSource.addFeature( thing );


texts: Geowebcache pre-seed creates white and black tiles

texts:
I have published Arial photo on GeoServer using ImagePyramid technique. I have created custom gridset for it.
I have created cached tiles using geoWebcache pre-seeding. When I am accessing it in Qgis using WMTS, from zoom level 01 to 09 its seems correct. but after zoom level 09 on somewhere, I found black and white tiles.

    But When I am accessing it on QGIS using WMS (tiled=true). At the same zoom level, it seems correct.

Can anybody suggest why GWC creates white and black tiles and how to fix it?

Debug log :

  2018-07-09 10:35:25,028 INFO [org.geoserver.wms] -  Request: getMap
    Angle = 0.0     BaseUrl = http://localhost:8080/geoserver/  Bbox =
  SRSEnvelope[137392.943896347 : 138088.0333557299, 123511.53727982522 :
  124206.62673920812]   BgColor = java.awt.Color[r=255,g=255,b=255]     Buffer = 0  CQLFilter = null    Crs = PROJCS["TM75 / Irish Grid",
  GEOGCS["TM75",
      DATUM["Geodetic Datum of 1965",
        SPHEROID["Airy Modified 1849", 6377340.189, 299.3249646, AUTHORITY["EPSG","7002"]],
        TOWGS84[482.5, -130.6, 564.6, -1.042, -0.214, -0.631, 8.15],
        AUTHORITY["EPSG","6300"]],
      PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]],
      UNIT["degree", 0.017453292519943295],
      AXIS["Geodetic longitude", EAST],
      AXIS["Geodetic latitude", NORTH],
      AUTHORITY["EPSG","4300"]],    PROJECTION["Transverse_Mercator", AUTHORITY["EPSG","9807"]],    PARAMETER["central_meridian", -8.0],
  PARAMETER["latitude_of_origin", 53.5],    PARAMETER["scale_factor",
  1.000035],    PARAMETER["false_easting", 200000.0],    PARAMETER["false_northing", 250000.0],    UNIT["m", 1.0],
  AXIS["Easting", EAST],    AXIS["Northing", NORTH],
  AUTHORITY["EPSG","29903"]]    Elevation = []  Env = {}    Exceptions =
  SE_XML    FeatureId = null    FeatureVersion = null   Filter = null
    Filters = null  Format = image/jpeg     FormatOptions = {}  Get = true
    Height = 1024   Interpolations = []     Layers =
  [org.geoserver.wms.MapLayerInfo@edba3cf7]     MaxFeatures = null  Palette
  = null    RawKvp = {FORMAT=image/jpeg, STYLES=raster, WIDTH=1024, HEIGHT=1024, LAYERS=New_IrelandLimerick:Ireland_Limerick_New,
  EXCEPTIONS=SE_XML, GWC_SEED_INTERCEPT=true, REQUEST=GetMap,
  SRS=EPSG:29903,
  BBOX=137392.943896347,123511.53727982522,138088.0333557299,124206.62673920812, VERSION=1.1.1, SERVICE=WMS, TRANSPARENT=true}  RemoteOwsType = null
    RemoteOwsURL = null     Request = GetMap    RequestCharset = UTF-8
    ScaleMethod = null  Sld = null  SldBody = null  SldVersion = null
    SortBy = null   SortByArrays = null     SRS = EPSG:29903    StartIndex =
  null  StyleBody = null    StyleFormat = sld   Styles = [StyleImpl[
  name=raster]]     StyleUrl = null     StyleVersion = null     Tiled = false
    TilesOrigin = null  Time = []   Transparent = true  ValidateSchema =
  false     Version = 1.1.1     ViewParams = null   Width = 1024
  2018-07-09 10:35:25,059 DEBUG [org.geoserver.gwc.layer] - --> GWC
  Seeder Thread-11 submitting getMap request for meta grid location
  [196, 72, 9] on org.geoserver.gwc.layer.GeoServerMetaTile@561daeb8
  2018-07-09 10:35:25,059 INFO [org.geoserver.wms] -  Request:
  getServiceInfo 2018-07-09 10:35:25,059 DEBUG [org.geoserver.ows] -
  Getting layers and styles from LAYERS and STYLES 2018-07-09
  10:35:25,059 DEBUG [org.geoserver.ows] - establishing raster style for
  New_layer:RasterData 2018-07-09 10:35:25,059 DEBUG [org.geoserver.wms]
  - setting up map 2018-07-09 10:35:25,059 DEBUG [org.geoserver.wms.map] - setting up 1024x1024 image 2018-07-09 10:35:25,059 TRACE [org.geoserver.jai] - Recycling tile hit on type:1, banks: 786,432,
  arrayLength: {3}


texts: Why is pgrouting driving distance skipping edges?

texts:
I have a road network where I do a driving distance operation through pgrouting. This is my SQL:

SELECT * FROM driving_distance('SELECT gid AS id,source,target, cost FROM elveg',676663,1000,false,false);

From this query I get a strange result:
http://bildr.no/view/1373984

Black is my road network and green is the result from the above query. From the image you see that some edges are not included(not green) while edges which can only be connected through this not included edge in fact is included?

How can this be? Am I doing something totally wrong here?


texts: QGIS- What is the meaning of the Query Builder operators?

texts:
I am learning to filter my data, and I understand how to use the basic operators (e.g. =, OR, AND, !=, less than, more than etc) however there are some operators that are less obvious.

What do the following mean in the query builder operators section?

LIKE
ILIKE
%
IN
NOT
NOT IN

And how might they be used in an expression?


texts: How Mapserver deals with Huge Aerial Images

texts:
I've been trying to evaluate Mapserver , i'm wondering if there's a way to handle the large size of Aerial image with any kind of Tiling or Caching .

Thanx


texts: Good equidistant Canada-wide projection?

texts:
I'm trying to map several thousand points across Canada, and find nearest neighbours for each point. I've figured out how to do all of this in QGIS using the Distance Matrix tool, but my projection was WGS84 so I'm getting distances in degrees. Some searching of StackExchange has led me to realize that I need an (ideally equidistant) projection that calculates distances in metres.

Anyone have any suggestions for a CRS that would work in this scenario? I realize there's going to be some error, but as long as it's not a crazy amount I should be okay…


texts: How to find overlapping polygons within a multipolygon (Oracle spatial)

texts:
I am working with Oracle spatial and have a table with MULTIPOLYGONS.
How can I find the database entries where the polygons of a multipolygon overlap - within the same geometry?

So, a function that checks the polygons of a multipolygon for overlapping.

I found several examples how to compare polygons in order to find overlapping polygons stored in a different row or in a different table:

SELECT *
FROM
  GEOTAB TabA, GEOTAB TabB
WHERE
  TabA.ID = 1234
  and
  SDO_OVERLAPS(TabA.GEOM,  TabA.GEOM) = 'TRUE'

but what I need is something like:

Select * FROM GEOTAB where SDO_OVERLAPSINTERNAL(GEOM) = 'TRUE'


texts: Is it possible to use SELECT ST_GeomFromText( with other SELECT?

texts:
I'm trying to use SELECT ST_GeomFromText() with a value catched from another query:

SELECT ST_GeomFromText(
  ST_AsText(
    ST_GeomFromGeoJSON(
  SELECT
    column-name
  FROM table-name
  ORDER BY id DESC
  LIMIT 1)
  ),
4326);

The value that I'm trying to use is compatible with the SELECT ST_GeomFromText() and I'm using PostgreSQL/PostGIS.

Other solution that would help would be a INSERT INTO with SELECT ST_GeomFromText() as a value.


texts: Selecting point features on the outside of line curves - QGIS 3

texts:
I need to select points that fall on the outside of line curves using QGIS3.
I have evaluated the curvature of the road and used that to select points that are within 20m of a "curved road" (statistically classified using ZScores). However, I am now interested in selecting the points that fall on the outside of the curve. Ideally, I would like to identify these points by flagging "yes" in an attribute field called "OutsideCurve". This process will be repeated on a wider road network.
The image below shows the data I am working with and I have highlighted the points which I deem to be on the outside of a curve.
Could someone please recommend a method/tool/plugin that may help to identify points on the outside of curves?


texts:
Geopandas read_file reads my time field with less precision. How can I make it the same as my file?

I created a file with geopandas to_file(filename, driver='GeoJSON') and it is something like this:

{ "type": "Feature", "properties": { "LineaRuta": "A", "Time": "2018-06-25 11:57:55.537" }, "geometry": { "type": "Point", "coordinates": [ 505529.0, 4788595.0 ] }

I read it with read_file(filename, driver='GeoJSON') and the time (datetime) losses its millisecond precision:

LineaRuta     Time                 geometry
0     A        2018-06-25T11:48:09  POINT (506375 4790366)

I tried to see if it is changing the format with:

schema = gpd.io.file.infer_schema(data2)
print(schema)

But the output says that it is string type. Then why it reads like a datetime? Output:
{'geometry': 'Point', 'properties': OrderedDict([('LineaRuta', 'str'), ('Time', 'str')])}

I need the millisecond precision in the time field, is there any easy way to do it? For example with fiona's read?


texts: How to get the extent of a group of Raster layers and zoom to it with PyQGIS?

texts:
I want to find a way to zoom to a group within the QgsLayerTreeView. The group only contains QgsRasterLayer. The indirect approach would be to create a QgsRectangle and serve it to the QgsMapCanvas.setExtent() method.

Unfortunately I don't see a way to obtain the extent of a QgsRasterLayer, kind of Bounding Box or "geometry", which I need for the creation of such a rectangle.


texts: Representing NULL data for a numeric field in ArcGIS (symbology)

texts:
I have a tank age data set which I would like to present on ArcGIS.

I would like to represent the different tanks using their tank age. I would like the symbol colour to be dependent on the tank age (think red, green, yellow).

The problem is, there are some tanks which I don't have the tank age for.

How can I represent the tanks that don't have tank age?

Should I:

Put in the tank age as 'N' to give it a character input
I can't put a dummy number such as -9999 since I have other graphical representations that rely on that age field and -9999 would skew the charts.

Edit:
- I could create a different column for Null_values and ask the symbology to draw off that. I'm unsure if we can represent one symbology from different fields.
- Grey of blank would be a great visual but I'm not sure how to map off the NULL data.


texts: What are the differences between the photometric options for GDAL_merge?

texts:
I am using GDAL_merge to merge several geotiff tiles into one large tile. One of the creation options I can add is photometric.
There are several options for this.

MINISBLACK
MINISWHITE
RGB
CMYK
YCBCR
CIELAB
ICCLAB
ITULAB

What are the differences between these different outputs? I understand the MINISBLACK and MINISWHITE options, but the others are less clear.

Thanks in advance.


texts: How to convert polygon shapefile to CSV?

texts:
I can convert shapefile (having points) to CSV file, import in MySQL databases and then i can execute my queries for searching. Now I want to do same for polygon shapefile; that is, I want to convert a shapefile having polygons to CSV and then import in MySQL databases for search purposes.

Is that possible?

Actually i have a geojson file having road polylines. i want to convert this geojson to csv then import in databases and then execute my search queries e-g i want to draw on google map only those roads having length less than 10km etc.


texts: ESRI JS API PopupMobile - view attribute page directly

texts:
I am using the ESRI Javascript API to write a mobile web application.  I have used the PopupMobile widget for viewing/editing attributes.  This widget creates a popup on the map when a feature(s) is clicked and then provides a button to access the attributes in a new page.

I would like to bypass the popup altogether and go directly to the attribute page.  Is this possible?

UPDATE:

I created my own onClick handler instead of using the PopupMobile widget, but I'm getting some very interesting results in my attributes page.

Here is a screenshot of the AttributeInspector widget as it displays on the page I set up in my application:

As you can see, the various dijits that make up the AttributeInspector are not displaying properly on the page.

Here is the code I used to set up the AttributeInspector:

mapDeferred.addCallback(function(response) {
    map = response.map;

    currentItem = response.itemInfo;
    $("#webmapTitle").html(currentItem.item.title);

    resizeMap();

    var fl = currentItem.itemData.operationalLayers[3].layerObject

    var layerInfo = [{
        featureLayer : fl
    }];

    var attInspector = new AttributeInspector({
        layerInfos : layerInfo
    }, "attributeInspector");
    attInspector.startup();

    dojo.connect(map, "onClick", function(evt) {

        var toleranceInPixels = 10;
        var pnt = evt.mapPoint;

        var pixelWidth = map.extent.getWidth() / map.width;
        var toleranceInMapCoords = toleranceInPixels * pixelWidth;
        var ext = new esri.geometry.Extent(pnt.x - toleranceInMapCoords,
            pnt.y - toleranceInMapCoords,
            pnt.x + toleranceInMapCoords,
            pnt.y + toleranceInMapCoords,
            map.spatialReference);

        var q = new esri.tasks.Query();
        q.geometry = ext;

        fl.selectFeatures(q)

        $.mobile.changePage($("#attributePage"));

    })
});

And here is the markup:

Attributes

        Map


texts: Determining min and max values in an ASCII raster dataset using Python?

texts:
I have a raster dataset in ASCII format. Using Python, I need to determine the min and max values inside the dataset. I've been told that the header information is key, which contains things like number of rows/columns, cell size, and etc.

Can't you simply skip the header information and read the entire dataset to determine the min and max values?

This is what I'm trying to do. I'm skipping the first couple of lines which contain the header information, and trying to determine the values from there on. The following is sort of what I have, but need some guidance as I am new to Python.

raster_file = open('data.asc', 'r') # Open the file
data = raster_file.readlines()[4:] # Read the lines in the file, and skip the first six lines

for lines in data:
    print max(data) # Find the max value in data
    print min(data) # Find hte min value in data

Any suggestions?


texts: Performing cluster analysis when separate clusters and point layer using ArcGIS Pro

texts:
I have a garbage patch layer (cluster) and a plastic accumulation points (point). The number of points are less compared to the spread of clusters.  there are many extra clusters. I want to identify the specific cluster each point belongs to. finally I want to make this result in to an output.
Is there anyway to do this using ArcGIS Pro?


texts:
I am trying to import ways from OSM maps and I am trying to classify them by type of roads, for example if highway='service' then 'cat1'. I tried to use a case statement but it doesn't work
For info, I am using ogr2ogr in python. Here is my code:

import os
os.system('ogr2ogr -overwrite -skipfailures -f "ESRI Shapefile"  -progress -
sql "SELECT highway,CASE highway WHEN %s THEN %s ELSE %s END FROM lines" -
lco ENCODING=UTF-8 %s %s' %('service','cat1','cat2','C:/Users/mehdi/Desktop/haiti1.shp','C:/Users/mehdi/Desktop/haiti.pbf'))


texts: Gdal_translate and transparent rasters in PDF

texts:
I am trying to make an application that creates geospatial pdf's from my website, but I am encountering problems with Gdal_translate and transparency of rasters.

The setup that I am working up against is based on Mapserver, creating wms services. These services are created with image/png properties, and a predefined transparency.
I would like to be able to combine these wms serices on top of a basemap in a geospatial PDF, but when I combine the layers with Gdal_translate, the transparency on the individual layers are lost.

When I call the layers individually with Gdal I am quite sure the transparency is present, so I assume that the problem occurs in the process of making the overlays in Gdal.

Have anybody experimented with this usage of gdal_translate against Mapserver - and have succeded in creating transparent overlays from rasters?


texts: Calculating the area of one layer inside area of another layer QGIS

texts:
I am currently looking at some Marine Protected Areas (MPAs) in Madagascar, and am trying to work out the percentage reef cover inside an MPA boundary. For example (in the picture below), I have the outline layer for the Ankivonjy MPA, and then another layer for reef cover. So my questions is: would it be possible to intersect the two layers, and then determine the proportion of the MPA that is covered by reef. If so, please could you explain how this is done.

The end goal is to be able to say that "coral reef comprises x% of the area within the Ankivonjy MPA".

NB: I am doing these analyses in QGIS.


texts: Why does setting extent in OpenLayers end up with an adjusted extent?

texts:
If I run the following code in OpenLayers 3...

console.log(map.getView().calculateExtent(map.getSize()));
map.getView().fit(map.getView().calculateExtent(map.getSize()), map.getSize());
console.log(map.getView().calculateExtent(map.getSize()));

I get the following output logged...

453350.0260566147,5335322.150098615,551079.8960333853,5508792.669307384
404485.0910682294,5248586.890494232,599944.8310217706,5595527.928911767

In my project, I would like to be able to store the current extent, so that I can reuse it later and go back to exactly the same extent.  But it doesn't work as I expect it to.

It appears to be zooming out an entire zoom level.

Why is it behaving this way, and how can it get it to go back to the same extent that I had before?


texts: Interpolating sparse bathymetry to correct global grid using ArcGIS Desktop?

texts:
I am using ArcMap 10.5.1

I am trying to update a global bathymetry grid (at 1 min resolution) with local sounding data (resampled to 1 min resolution). The sounding data is sparse and irregular in coverage (criss-crossing lines across the area). I would like to correct the global grid with the ship sounding data. Rather than just replace the global values by the local ones, I would like to correct the global grid with the local soundings interpolated over a larger scale. I have tried spline in tension interpolation over the entire study area but the output was patchy.

I have also tried to follow the methodology from Jeffrey Evans in this post:
Interpolation of multibeam bathymetry but I was unable to figure out how to filter the data using a Gaussian weighted filter.

My question is - How can I smooth a patchy spline interpolated map with a Gaussian filter (sigma = 2) so that the interpolated region is significant within 2 standard deviations of the original sampling point?


texts: Is there a way to bulk download Landsat data through cmd prompt or Python?

texts:
I know I can download landsat data through earth explorer and they have their own Bulk Downloader here: https://lta.cr.usgs.gov/EEHelp/BulkDownloadApplication

But I am wondering if there is somehting similar to the way Sentinel lets you bulkdownload, through their API in Python? Where you can give it an area of interest and a date range and it downloads anything that meets that criteria.


texts: Setting import to get list of values “constraint” for field, as described in INTERLIS Modell, into FileGDB, using Ili2db (ili2fgdb)?

texts:
I’m trying to import an INTERLIS 2 Dataset (Modell (ili) and Data (xtf)) in a FileGDB with the script ili2fgdb. The goal is to make some edit in the Datas in ArcMap and then export them back in INTERLIS. To avoid errors by the export I’d like to have all the “constraints” described in the INTERLIS Modell in my Geodatabse to prevent any non-compliant values from being entered in a field.

I’m struggling with the attribute with a defined list of values / classifications (“Aufzähltypen”):

When there is a constraint describing a limited list of possible values for a field (in DOMAIN) in the INTERLIS MODELL, I want to have a dropdown list with those values when I edit that field in ArcMap.

Does someone know which setting I need in my import to get those values as “Coded Values” in the “Domains” settings of my Layer in the Geodatabase (in order to have dropdown lists for the concerned fields)?

(I’ve managed to import those lists of values as “Lookup-Tables” with options like --createEnumTabs, but it doesn’t set constraints on the fields in the tables that are supposed to be limited by one of those lists of values.)


texts: How can I do edit-and-continue in ModelBuilder?

texts:
I'm constructing a model in ModelBuilder, and I ran into an error with one of my tools. I fixed the error and want to continue running my model from that step onward without having to repeat the entire model. Is there any way to do this?

Thanks!


texts:
I would like to add my own SVG symbol instead of a dot. I tried the following.

layer = iface.activeLayer()
QgsProject.instance().addMapLayer(layer)
svgStyle = {}
svgStyle['fill'] = '#0000ff'
svgStyle['name'] = path
svgStyle['outline'] = '#000000'
svgStyle['outline-width'] = '6.8'
svgStyle['size'] = '6'
symbol = QgsSvgMarkerSymbolLayer.create(svgStyle)
layer.renderer().setSymbol(symbol)

Unfortunately the following error occurs:

TypeError: QgsSingleSymbolRenderer.setSymbol(): argument 1 has unexpected type 'QgsSvgMarkerSymbolLayer'


texts: Understanding URL that IdentityManager.findServerInfo() method expects?

texts:
Does anyone have experience working with the IdentityManager object in the JS API 3.9?

I'm trying to understand what URL the findServerInfo() method expects.

I've tried the base REST url, and individual service URLs, but so far no luck.


texts: How to add alpha chanel to mbtiles in gdal_translate?

texts:
I have lots of 16bit 3-band .tif from one sattelite. In batch processing i've made them 8-bit 4-band (RGBA) with LZW compression (.tif). Now i'm trying using gdal_translate (gdal 2.1 from QGIS 2.16) generate .mbtiles. But some of the output mbtiles have black and white background although, original raster has alpha channel. Format of mbtiles is png.figure1

I've also tried gdal2tiles.py and mb-util utilities. they made good mbtile from point of transparency but all the pngs are turned upside down see  figure2  (gdalinfo shows that origin has wrong Y coordinate with "-") figure2

gdal_translate with option -a_ullr 0 0 256 256 didn't help me to turn pngs. gdalwarp -to SRC_METHOD=NOGEOTRANSFORM helped to turn pngs, but after mb-util the mbtile became the same as figure1 (without alpha)

After that i've tried global mapper trial and it make normal mbtiles as i want to do (see figure 3). figure3

I've noticed 2 things:

Global Mapper makes all the tiles 32bit. I mean if it is not border tile it anyway has 4 band (32bit). after gdal_translate and mb-util non-border tiles (which are tottaly covered by color pixels) tiles have 3 band (24 bit).
the blank tiles for good rasters (for which output mbtiles has trancparency) have size of 872 bytes and for bad rasters size of 871 bytes.

p.s. The blue background - it is QGIS background to clearly see backgrounds of raster. please do not offer nodata because this makes holes in the bright places in raster like roofs etc.


texts:
I've been trying to replace GIS tools (ArcGIS and QGIS) with R and some packages (leaflet, raster, sp and rasterVis) but I can't replicate GIS layout maps.

Has anyone created layout maps using shapefiles, points or raster layers using only R?

Edit.

I discovered prettymapr, a package that works with simple plots. Playing around with it and simple plots I made this:

plot(r,
     breaks = seq(0, 4800 ,400),
     col = c("#a6cee3", "#1f78b4", "#b2df8a", "#33a02c", "#fb9a99", "#e31a1c",
             "#fdbf6f", "#ff7f00", "#cab2d6", "#6a3d9a", "#ffff99", "#b15928"),
     xlab = 'Longitud [O]',
     ylab = 'Latitud [S]',
     main = 'Medias [mm/mes]')
plot(cuenca2, border = 'mediumblue', lwd = '3', add = TRUE)
plot(limites, border = gray(0.1), add = TRUE)
plot(gridlines(r), add = T, col = grey(0.1))
prettymapr::addnortharrow(scale = 0.6, text.col = 'black', cols = c('black', 'black'))
prettymapr::addscalebar(plotunit = 'latlon', widthhint = 0.25, lwd = 1, padin = c(0.15, 0.08), label.cex = 0.9)

r is a tif raster and both cuenca2 (blue) and limites (gray) are shapefiles.

Also, I tried spplot but despite I like the label axis format I can't add a scalebar and a north arrow, perhaps because of the raster's longlat projection.