texts:
I am trying to connect from Arc Cataloge (version 10.1) to postgresql 9.2. But every time I try to make "Database Connection" I am getting this error:

Failed to connect to the specified server.
Underlying DBMS error[could not connect to server: Connection timed out
   Is the server running on host"x.y.z.w" and accepting TCP/IP connections on port 5432?

Postgresql is installed on Centos server (PORT 5432) and I can access it using pgadmin.

Using psql --h X.Y.Z.W -U user -d database on another server I am getting the same error.

here is the content of pg_hba.conf:

 local   all             all                                md5
 # IPv4 local connections:
 host    all             all            X.Y.Z.W/32          md5
 host    all             all             0.0.0.0/0          md5
 host    all             all             127.0.0.1/32       md5
 # IPv6 local connections:
 host    all             all             ::1/128            md5

Also in postgresql.conf:

Listen_addresses = '*'
port=5432

in etc/sysconfig/iptables:

 *filter
  4 :INPUT ACCEPT [0:0]
  5 :FORWARD ACCEPT [0:0]
  6 :OUTPUT ACCEPT [0:0]
  7 :RH-Firewall-1-INPUT - [0:0]
  8 -A INPUT -j RH-Firewall-1-INPUT
  9 -A FORWARD -j RH-Firewall-1-INPUT
 10 -A RH-Firewall-1-INPUT -i lo -j ACCEPT
 11 -A RH-Firewall-1-INPUT -p icmp --icmp-type any -j ACCEPT
 12 -A RH-Firewall-1-INPUT -p 50 -j ACCEPT
 13 -A RH-Firewall-1-INPUT -p 51 -j ACCEPT
 14 -A RH-Firewall-1-INPUT -p udp --dport 5353 -d 224.0.0.251 -j ACCEPT
 15 -A RH-Firewall-1-INPUT -p udp -m udp --dport 631 -j ACCEPT
 16 -A RH-Firewall-1-INPUT -p tcp -m tcp --dport 631 -j ACCEPT
 17 -A RH-Firewall-1-INPUT -m state --state ESTABLISHED,RELATED -j ACCEPT
 18 -A RH-Firewall-1-INPUT -m state --state NEW -m tcp -p tcp --dport 22 -j ACCEPT
 19 -A INPUT -m state --state NEW -m tcp -p tcp --dport 5432 -j ACCEPT
 20 -A RH-Firewall-1-INPUT -m state --state NEW -m tcp -p tcp --dport 80 -j ACCEPT
 21 -A RH-Firewall-1-INPUT -j REJECT --reject-with icmp-host-prohibited
 22 -A RH-Firewall-1-INPUT -m state --state NEW -m tcp -p tcp --dport 5432 -j ACCEPT
 23 COMMIT

I also put "postgresql-client-libraries" in /bin directory of ArcGIS for Desktop.

I appreciate any help.


texts: I have a shapefile (points) loaded in PostGIS. What I want to do is execute an ECQL query using geotools filtering points within a certain radius (10 kilometers) from a specific point (user input). Then output a list (table) with the resulting filtered points. I want the output table in a way that I can later call it (as a class) in a Java code for some other algorithms.


texts: Multi ring buffers that are rings/donuts even when overlapping other buffers

texts:
I need to create multi-ring buffers for about 400 points that are quite close to each other. This has required me to set the output buffers to not dissolve in my code as below so that they overlap each other and do not combine with each other.

For each of the buffers, I have an associated set of polygons that fall within the buffers with matching ID's. I need to do a union them and remove the non-overlapping areas of this union.

The current MRB function in arc doesn't allow for the rings to be created if you don't select dissolve and when that happens the overlapping buffers are connected.

This is not ideal for my situation where I am trying to union one buffer to one matching set of polygons of similar ID's.

So, I am wondering if it's possible to create multi-ring buffers with rings even when the buffers overlap each other using arcpy


texts:
I have a python code I'm running in ArcGIS 10.6 that replaces the data source path of all layers in the mxd with a new path. I'm trying to use a wildcard for the current workspace path as many layers in our template mxd have different data source paths.

The part of code in question is this line:

mxd.findAndReplaceWorkspacePaths(r"*", r"P:\JOBS\2018\18-9992\Drawings\GIS")

The r"*" doesn't seem to work as a wildcard.

Here is the full code:

# Importing necessary modules
import arcpy, os

# Setting path to mxd and Data Frame
mxd = arcpy.mapping.MapDocument(r"CURRENT")
df = arcpy.mapping.ListDataFrames(mxd, "UTM_12")[0]

#List layers and replace current path with new path
for lyr in arcpy.mapping.ListLayers(mxd, "*", df):
    mxd.findAndReplaceWorkspacePaths(r"*", r"P:\JOBS\2018\18-9992\Drawings\GIS")

#Delete any layers with broken links (no datasource)
for lyr in arcpy.mapping.ListLayers(mxd, "*", df):
    if lyr.isBroken:
        arcpy.mapping.RemoveLayer(df, lyr)

# Saving the mxd
mxd.save()


texts: QGIS: use graph builder to classify multiple images at once

texts: I am classifying subsets (Sentinel 1, Sentinel 2, NDVI, and Landsat 8 data) in order to see which image as the highest accruacy. These classifications are taking several days for each of the images. Is it possible to use QGIS graph builder to classify the images simultaneously?


texts: Optimal distance and power exponent in IDW

texts:
I want to apply IDW interpolation for my region of interest, and through leave-one-out cross validation I want to find the optimal values for the distance parameter and the power exponent. I know that I can find the optimal distance, running IDW in a loop for many distances and each time for each distance to perform the cross validation. In this case I will have a scatterplot of distance (x-axis) and error (y-axis). The optimal distance value is the one that gives the lowest error. The same I can do for finding the optimal power exponent.
What I want to do is, somehow, run IDW for different distances and different power exponents at once and each time the result of the cross-validation to be stored in a 3D scatterplot in which, in x-axis there will be distance, in y-axis there will be the power exponent and in z-axis there will be the error for each combination of distance and power.

How can I perform this in R?


texts:
I have ASCII data from different parts of water bottoms.
I would like to convert all the parts into shapefiles.

But, when I use the polygonize raster to vector tool in QGIS I get invalid geometries. To solve this I tried using the validation tool and also the v.buffer function in GRASS. After using that it keeps giving new invalid objects.

How can I fix this?


texts: point with multiple texts in open layers 3

texts:
I have a map with vectorial features, represented with a text as icon with font awesome:

new ol.style.Style({
        text: new ol.style.Text({
          text: '\uf005',
          font: 'normal 14px FontAwesome',
          textAlign: 'center',
          textBaseline: 'middle',
          fill: new ol.style.Fill({
            color: 'white',
          }),
          stroke: new ol.style.Stroke({
            color: 'black',
            width: 5
          })
        })
      })

Now i need to add a real text as label, of course with other font.

Is there any way to add 2 texts, with different configuration to the same feature?


texts: Calculating exterior and interior orientation of any drone without camera calibration

texts:
I have a drone, how can I find the projection matrix without calibration?

Since the drone has a GPS, I know some angles and the distance between the two photos in meters. I also know the formula of the projection matrix. P = K[R|t]. I have an R Matrix. I calculated it. I know the focal length of the drone, FOV and the size of the images it takes. How do I calculate K with this information? How do I calculate T?

I have information about how many meters between two photos and gsd(ground sample distance). Should the values in the T vector be in pixels?


texts: QGIS - How to find the distance for the nearest polygon from point for a defined angle?

texts: I'm using QGIS 3.1 and I want to find the distance from a point to the nearest polygon for a given direction. So for example, for each point in my point layer I want to find the distance to nearest polygons located to the south west of each point (225 degrees). Is this possible?


texts: Bottom half of javascript web map application missing

texts: Having an issue where the bottom half of a javascript web application is missing.  I have checked the ArcGIS map services and they are functioning properly.  I'm not a JS developer and am not sure what the issue could be.  Any suggestion on what the cause could be?  I've attached a screenshot of the application where you can see the bottom half missing.


texts: Identify Transformation Used by Append Tool

texts:
In ArcMap, I regularly use the Append Geoprocessing tool to append data from a dataset in one projection into another dataset that is in a different projection.  When you append between projections, as expected, there are slight shifts in the data.  I'm fine with this as the shifts have always seemed to be consistent, meaning it appears to use a consistent transformation in this process.
However, I'm now in a situation, due to some additional datasets in other projections getting involved, that I need to know what transformation is being applied by the append tool.
I know I can use project on the data to a temporary dataset, defining the transformation in the project tool, and then append from that temp dataset. However, that doesn't change the fact that if I do that, without knowing what transformation has already been being used by the append tool, then my data will no longer match where boundaries have historically been according to the append tool workflow (which creates edge-matching issues with other data-sets).
So, does anyone know how the append tool transforms data when it is appending between projections?


texts: I'm trying to create the files in HEC-RAS 5.0.3, but nothing is created. Even though the location of the project is shown in the main window, I can't find the file in this location. That causes this problem?


texts: Building graphs between two point layers

texts:
I am trying to run a script to find the shortest routes between two sets of point layers  (Blue star and Green dots) in QGIS using road network. I am able to use GRASS gis "v.net.distance" for 1 star at a time (see picture below). I am trying for multiple stars at a time but can't iterate more than one layer at using the plugin GUI. I was trying to use network analysis library in PYQGIS using some online scripts but no success so far.
How can I do this?


texts: Raster from a set of Points keeping the resolution

texts:
I have a database with some points covering the earth surface and I want to build a raster (GeoTIFF using gdal_rasterize) which mantains the same spatial resolution of the input points (so without losing any information).

I know that the input points should have a spatial resolution of 7km more or less.
 with resolution of points I mean that there is a point every 7km, or in other words that given a point the closest ones are at a distance of 7km (up/down - right/left). the points represents measurements of a concentration of a certain gas , so I would like to build the image without losing information about the concentration in that point

How can I determine the correct resolution of the image generated by gdal_rasterize (-ts option of the command)?


texts: CAD annotation guide lines do not show in Arcmap

texts:
When I open a CAD file in Arcmap the annotation's guid lines and arrows linking to feature do not appear. see snapshots below.
how can I show thses guides in armap?


texts: How to get feature values from multiple polygons of shapefile while using extract function on raster

texts:
I have a shapefile which contains arround 800 circular polygons with diameter of 60 meters each. Raster is 20 meter resolution so I would have multiple cell values extracted for each polygon.  I am using an extract function on multiple raster layers in a for loop, and am putting the results in the dataframe. I know i can get the cellnumbers within the extract function, but I would also like to obtain the feature values (in my case a column in attribute table containing specific polygon code) directly from the shapefile for each pixel value extracted. Thus I would have in one row a pixel cell number, its extracted value and the polygon code.

my code:

dirs=list.dirs(full.names=TRUE)
dirs=dirs[grepl("R20",dirs)]
dirs=as.list(dirs)

DataToWrite=list()
sumirana=data.frame()
for (j in 1:length(dirs)) {
  setwd(paste("E:/",as.character(dirs[j]),sep="/"))
  files=list.files(full.names=FALSE,pattern="*20m*.tif$")
  for (i in 1:length(files)) {
      curRaster=raster(paste(getwd(),files[i],sep="/"))
      rasterOut=(na.omit(extract(curRaster,azo_plohe, cellnumbers=TRUE,
      weights=TRUE, df=TRUE)))
      if((length(rasterOut) > 1)) {

         DataToWrite=c(DataToWrite,rasterOut[3])

         } else {

         lista=list(rep("x",15))
         tablica_folder=as.data.frame(t(unlist(lista)))
         colnames(tablica_folder)=c(1:15)
         }
   }
   if((length(rasterOut) > 1)) {
      mylist=as.data.frame(DataToWrite)
      naziv_fajla=as.vector(rep((substr((names(DataToWrite)
      [1]),50,55)),nrow(mylist)))
      tablica_folder=as.data.frame(c(rasterOut[2],naziv_fajla[1],
      rasterOut[4],mylist),col.names=c(1:15))
      sumirana=rbind(sumirana,tablica_folder)
      DataToWrite=list()

      } else {
      sumirana=rbind(sumirana,tablica_folder)
      }
}


texts:
How to Add Shapefile to map  in ArcGIS Silverlight application runtime ?

I published the shapefile feature as map service on arcgis server and consumed it my SL Application through the service Url but the shapefile feature does not added to the map ..

what's wrong?


texts: Data not showing in QGIS when connecting to Oracle

texts:
I have an Oracle XE 11g spatial database that is connected to a web map. All the layers of the web map are stored in the database and they work just fine. I sometime need to edit some of the layers in the web map and for this task I use QGIS but the problem is that when I try to add the Oracle tables in QGIS nothing shown on the screen, although the tables have spatial indexes!

Has anyone faced a smimilar problem in QGIS?


texts:
I have seen sample apps for windows phone, IOS, and even a beta for android. does anyone know if there is anything started for Palm, or Blackberry?
Interested in ArcGIS Server mobile applications.


texts: Using Mosaic Dataset on Shared Network Drive?

texts:
I'm working with several colleagues at the same academic institution and we all have access to a shared network drive.  Within this network drive, each of us has our own folder where we store our own geospatial data.

Using ArcMap 10.5, I've created several mosaic datasets - mostly for land cover data - within a colleague's folder (on the shared network) who is the lead on a project, but the source rasters are within my own folder on the shared network drive.  When I view these mosaic datasets on my computer, everything works fine: the overviews are correct and if I bring one of the source rasters into ArcMap, there are no problems.  However, if my colleague (the same person who has the folder where the mosaic datasets are located) attempts to bring the mosaic datasets into ArcMap on his computer, all he sees is a gray rectangle, and if he attempts to load the source rasters individually from the dataset, all that appears again is a gray rectangle.  The appearance is not the only problem; when the source raster is loaded from the dataset, its name is changed and the range of pixel values also changes.

Why is this not loading correctly on his computer, only on mine?  I've checked for broken paths, and there are none. Could this be a problem with using a shared network drive with the source rasters located in my folder, and the mosaic dataset in someone else's folder?


texts: What projection (in meters) should I use for most of the Northern Hemisphere above 50 degrees?

texts:
I am using R to do spatial analysis. I have a river shapefile spanning most of the northern Northern Hemisphere (north of about 50 deg). I am hoping to buffer each river segment by its width in meters. However, the data is currently in WGS84 (aka ESPG 4326), which is lat/lon. The buffer function uses an input distance of whatever the projection is in (if it helps, I am using the 'sf' package's st_buffer function), so I need the projection to also be in meters.

Is there a projection that doesn't distort distances too much that could be used fairly well for latitudes over 50 deg North?

I apologize if this has been asked before. I haven't seen quite the same question, but it is possible that I missed it!


texts: How to join two shapefiles and create a new shapefile using ArcObjects?

texts:
Newbie in ArcObjects programming. Using C#. Developing a desktop application.

I have two shapefiles. One is a point shapefile Called Addresses and another Polygon shapefile which is called Area. The Address shapefile has an attribute called Address_Number while the Area shapefile has a attribute called Area_Number.

I need to select only the features that match "Address_Number" (from Addresses)  to "Area_Number" (from Area) and then create a new shapefile.

How can I do this?

Any help with some code example is greatly appreciated.


texts: trouble with code block execution in python script

texts:
I am trying to calculate values for a field in a python script, but keep getting an error that I don't understand:

ExecuteError: ERROR 000539: Error running expression: Fron000000F1
Traceback (most recent call last):
  File "<expression>", line 1, in <module>
NameError: name 'Fron000000F1' is not defined
Failed to execute (CalculateField).

I can't figure out where the name 'Fron000000F1' is coming from to try and track this problem. I have tested the expression in arcmap using the field calculator with no problems. What is wrong? Code is below:

import arcpy, os, csv
arcpy.env.overwriteOutput = True

out_fc_pm = r'C:\temp\temp.gdb\SignDetect__SOL680_DD_with_PM'

def calc_id(pm,pict_fld):
  import re
  pm_patt = re.compile('\d*\.*\d*')
  pm_match = re.findall(pm_patt,pm[5:])[0]
  pm = pm_match.replace('.','')
  pm_final = pm.zfill(6)
  dir_patt = re.compile('[A-Z]\d*')
  dir_match = re.findall(dir_patt,pict_fld)[0]
  id = pict_fld[:5]+pm_final+dir_match
  if id[-1:].isdigit():
    return id
  else:
    return id+'1'

arcpy.AddField_management(out_fc_pm,'OSMI_ID','TEXT')
arcpy.CalculateField_management(out_fc_pm,'OSMI_ID',calc_id('!POSTMILE!','!Front!'),'PYTHON_9.3')


texts: Interaction with two variables as parameters, saving the raster properties output into a table (.dbf)

texts:
The draft of the model below

TASK: I need to calculate the percentage of river channel that is covered with shadow in each timestep

RESOURCS: I have values for azimuth and elevation of the sun (hourly timestep for daylight hours within a year – over 3000 values)

INPUTS:
1.  A table (.dbf) with azimuth and elevation values
2.  A river chanel mask (single polygon)
3.  A DSM

PROBLEMS

How to make an interaction with two lists as parameters

The elevation value and azimuth value are assingned to a certain time and date,  .I don’t want to calculate hillshade for e.g  elevation row 2 and azimuth row 6 –I want  only elevation row 1 + azimuth row 1, elevation row 2 + azimuth row 2, elevation row 3 + azimuth row 3… etc…

How to save an output – all the means -  back to the input table
– I had an idea to use a calculate field tool, but I don’t know how to “import” the values to the expression field.

I've figured out the solution to my first problem. It was easier than I thought. I Do not need to iterate. The simple for of input variable as a list and specifying the table to derive the values.
But there is still the other problem - how to write the output back to the input table as a next column


texts: Making legend items start in new column in ArcGIS Pro?

texts:
In ArcGIS Pro 1.3.1, I have just looked for the ability to start legend items in a new column.

This functionality is present in ArcMap's Legend Properties on the Items tab which appears below:

Has anyone found the equivalent functionality in ArcGIS Pro, or know when it may be expected to appear?

It is not mentioned in the Work with a legend documentation.

Nor does it appear on the Element: Legend Item pane (see below)


texts: installing exifread for photo2shape plugin

texts:
I've already seen https://gis.stackexchange.com/questions/138207/error-loading-plugin-photo2shape-even-after-installing-exifread-module and Error installing Photo2Shape QGIS 2.6 for Mac, but cannot make it work.

First I installed pip, and then went to c:\python27\scripts and typed "pip install exifread", and it got installed. See

you can see that it is "up-to-date".. but when I try to install photo2shape it insist:

  no module named exifread


texts:
I'm using 2.8.1 QGIS version (Wien) and I have a raster layer with FWI values and a vector with district polygons.

What I want is to calculate the mean pixel FWI value per district.

I know the Zonal Statistics funtion, but when I use her appears this error:

"Algorithm Zonal Statistics starting...
'NoneType' object has no attribute 'exportToWkt' See log for more details".

I don´t know what could be the problem. I'm using the same CRS layers.

Any help?


texts: Seeking solution for PostGIS raster 3d intersection?

texts:
Is there a PostGIS solution to a 3D intersection between a DEM raster (with elevation in band 1) and a 3D line?

Here are some facts:

Rasters only support 2D processing so I can't use ST_3DIntersects directly unless I could convert the rater to a 3D TIN or some other 3D geometry.  If I could do this I could use ST_3DIntersects.
However there isn't a raster to TIN method and the only potential option to convert a raster to 3D Geometry (ST_DRAPE https://trac.osgeo.org/postgis/ticket/2327) is not yet implemented in PostGIS.

Here are two potential solutions that would be desirable but do not seem available:

a custom 3D raster intersection OR
a raster to 3d geomtery conversion so I could use ST_3DIntersects.


texts: Why does Clip analysis generate (some) empty output?

texts:
I am running a script that is working (see below) - it pulls a file from a folder, clips from that file based on an intersection with another file, and creates a new shapefile in a new folder. It has been working, except that it generates empty outputs for only some of the inputs.

With some research I've tried the following fixes to no avail:

Check extent and projection of Clipping feature to match that of the Input features.
Shorten output file name (hence the awkward string parsing in the output)
Standard re-boot of computer, programs.
[edit] Also convert Clipping file from Multi-Part to Single-Part Polygon.

Any ideas? Do you see quirks with my script? Are there quirks of Clip tool?  Thanks!

#Script - Clip selected buffer circles from WEA areas.
#Input: VTRX_1_distX_sel.shp
#Clip: WEA file
#Send output to clip folder

import arcpy import sys import os from arcpy import env

#Set overwrite arcpy.env.overwriteOutput = True print arcpy.env.overwriteOutput

#Set workspace workspace = "C:\\BOEM\\circle_areas\\" arcpy.env.workspace = workspace sel_files = workspace+"selectedcirc\\"

fileList = [] count = 0

directorypath = "C:\\BOEM\\circle_areas\\selectedcirc\\"

for dirname, dirnames, filenames in os.walk(directorypath):
    for filename in filenames:
        if filename.endswith("1_sel.shp"):
         #   filename = "VTR10_1_dist1_sel.shp"
            print "___"
            print "Clipping from " + filename
            print "___"
            in_features = sel_files+filename
            clip_features = "C:\\BOEM\\circle_areas\\WEA_proj_sing.shp"
            out_features = workspace+"clip\\"+filename[0:5]+"_"+filename[11:13]+"c.shp"
            print "Input file " + in_features
            print "Output file " + out_features

            #Clip
            arcpy.Clip_analysis(in_features, clip_features, out_features,"")
            print arcpy.GetMessages()
            print "___"
            print "___"
            print "___"
            ##


texts: Make POST request for WMS in Leaflet for big values of cql_filter

texts:
I am passing a spatial query such as CONTAINS(geom, POLYGON((// geometry))) in the cql_filter, but for big polygons GET request the allowed URI length is violated and hence my requests is getting HTTP 414 status code 'URI too long'.

I looked for the solution so someone suggested me to use POST request and it worked, but I don't know how to modify the WMS GET request in Leaflet to POST?


texts: Determining types of geometry in a single table?

texts:
I need to determine whether any tables within my database have more than one geometry type. For example, I want to know if there are Points AND multipoints within the same table.

Currently, I have been using the command "ST_ASEWKT" to convert the "geom" column into EWKT and scrolling through the results, trying to find more than one geometry type.
But I have table with over 50,000 records.

Is there a better way to do this?


texts: Extending lines to their Bounding Box in QGIS

texts:
How to extend lines to Bounding Box in QGIS as in the following image:

Figure 1: Lines

Figure 2: Lines extended to Bounding Box
Link to Shapefile of Bounding Box: https://drive.google.com/drive/folders/1jsJ8m2dG0DNwF1jepxpWzwfUVI5sVJJe?usp=sharing
Link to shapefile of lines: https://drive.google.com/drive/folders/1BQvrMiPillL5ubYsDHfEyoV9st12z0E9?usp=sharing


texts:
I have several DEM files and I want to merge them to create a mosaic.

I found the St_Union(raster) or St_MemUnion(raster) functions but I don't know how to use them.

Are there any other solutions?


texts: Which spatial statistics tools to use in examining autocorrelation in two layers?

texts:
Possible Duplicate:
  How to compare two spatial point patterns?

I am working on a project and have found myself a bit lost as to which types of Spatial Statistics tools I should use to get the results I am looking for.  In short, my project is looking at two different types of crime and whether they are spatially related to one another.  I am interested in seeing if the presence of one type of crime is related to the incidence of the other type of crime. I pulled in two separate point layers into ArcMap for each of the crimes using XY coordinates.  I am planning to study this on a neighborhood level.  I am unsure which tools or statistical functions to use in order to examine this problem.  Any help is appreciated.  Thanks in advance!


texts: What ArcGIS Online feature service provides European states with unique country codes?

texts:
I'm writing a web page to display ECB statistics on an ArcGIS JavaScript API-based web map.

I'm a bit lost when selecting a feature service so that I can use the UniqueValueRenderer to make the states pretty colors according to their statistics.

Is there a feature service available on ArcGIS Online that I can use for this?


texts: Writing shapefile metadata (shp.xml) programmatically

texts:
I want to write the Metadata for a Shapefile, which are saved in the ".shp.xml"-File. I know, I can look at and write those Metadata with ArcCatalog, but I want to do this programmatically. Is this possible with c# and arcobjects? Or do I have to write directly to the XML-File parsing the Content? If that's the case how can I create the "shp.xml"-File? I have a bundle of Shapefiles, which do not have this File. So maybe I could create the shp.xml-File and write to it using XElement or something like this ...
Do you have any suggestions, how to solve this Problem?
I use ArcGIS 10.1.


texts: Why some coordinate systems define x-axis as northings and some as easting?

texts:
The definition of north and east are pretty straight forward to grasp but only becomes difficult when used interchangably with x-y coordinates which have varying definitions for the direction of the axes. In Mathematics y was always vertical and x was always horizontal so logically I would assume that "up" == "northing" == "y" and "along" == "easting" == "x".

Why is this not the case in GIS?


texts: Efficiently reading and reclassifying many rasters in R?

texts:
I've been tasked to create a suitability analysis of wave conditions in the Gulf of Mexico.  I have 2 thousand or so raster files that are about 8 MB each (2438 columns, 1749 rows, 1km cell size).  The parameter in the raster files is wave period and I'd like to reclassify all of the rasters such that if 4<= wave period <=9 then make cell = 1, else cell = 0.  Then sum up all the rasters into a final raster and divide by the total number of rasters to give a total percentage of suitable observations and export result into some ESRI-compatible format...maybe something that can support floats if need be.  I've not worked much at all with either Python or R, but after searching online it seems to make sense to perform this process in one of those languages.  I've come up with some code so far in R, but am confused about how to make this work.

library(rgdal)
raster_data <- list.files(path=getwd())    #promt user for dir containing raster files
num_files <- length(raster_data)
for (i in raster_data) {                   #read in rasters
   my_data <- readGDAL(raster_data[i])

At this point I'm confused as to whether I should also reclassify and start summing the data within this loop or not.  My guess would be yes since otherwise I think I would possibly run out of memory on my computer, but just not sure.  I'm also not sure about how to reclassify the data.

In researching online I think I would use reclass(my_data, c(-Inf,3,0, 4,9,1, 10,Inf,0)), but does that look right?

And for summing would I use sum(stack(my_data)) and somehow output that.  Also...if this might be more efficiently performed or written in Python I'd be open to that as well.  I truly am a beginner when it comes to programming.


texts:
Selecting a vector layer and running the following script, I expect the total area of my feature-geometries to be the same at the end of execution.

from shapely.wkb import loads
from shapely.ops import cascaded_union, unary_union

def getTotalArea(someLayer):
    area = 0
    for feature in someLayer.getFeatures():
        area += feature.geometry().area()
    return area

layer = iface.activeLayer()
print getTotalArea(layer) # I want this to be the same area at the end of this script !!

# merge geometries
oldGeoms =[loads(feature.geometry().asWkb()) for feature in
layer.getFeatures()]
un = cascaded_union(oldGeoms)
newGeom = QgsGeometry.fromWkt(un.wkt)
print newGeom.area() # still the same area we started with, so far so good

#remove all previous features
layer.dataProvider().deleteFeatures([f.id() for f in layer.getFeatures()])
print getTotalArea(layer) # zero now, because all Featrues are removed

#create and add new Feature, with merged geometry
newFeature = QgsFeature(layer.pendingFields())
newFeature.setGeometry(newGeom);
layer.dataProvider().addFeatures([newFeature])
print newFeature.geometry().area() # just to be sure ... still the same geometry we started with
print getTotalArea(layer) # what the ... !??

However the output (with one of my example layers) is:

5272.52221858

5272.52221858

0

5272.52221858

5272.52221858

5288.52221858

The last line of the script prints a different area than I expect (off by 16). Why ? Is it some kind of bug, or am I missing something important ?

It depends on the layer and the features that are processed, I do not always encounter that strange behaviour. But if I have a layer that gives this "error", I can reproduce it consistently.


texts:
I do not have ArcGIS desktop. And I have installed the ArcGIS Python API.

And I tried the following:

import arcgis
from arcgis.gis import GIS
from arcgis.mapping import WebMap
gis = GIS()
map1 = gis.map('Paris') # Passing a place name to the constructor
                        # will initialize the extent of the map.
map1

No map displayed.

Why is that?


texts:
I have imported OpenStreetMap data for the country Switzerland via osmosis in my MySQL database.

Now I want to look for a specific address in that database, but I can't find a way to get the whole address (street, postcode, town ..).

For  example, the address I'm talking about is: Street: Archstrasse, town: Grenchen, postcode: 2540, county: Amtei Solothurn-Lebern, country: Schweiz.

<!-- language: lang-sql -->
select * from current_way_tags where way_id = 162530695;> &nbsp;

+-----------+----------+-------------+&nbsp;

| way_id    | k        | v           |&nbsp;

+-----------+----------+-------------+&nbsp;

| 162530695 | cycleway | lane        |&nbsp;

| 162530695 | highway  | primary     |&nbsp;

| 162530695 | name     | Archstrasse |&nbsp;

| 162530695 | ref      | 5           |&nbsp;

+-----------+----------+-------------+&nbsp;


texts: Create Custom function to update virtual attribute of point from raster in QGIS?

texts:
I am trying to dynamically update the elevation value of a point based on its location when the user tries to move it in edit mode and after saving it in QGIS. I got a point shapefile and DTM (.asc) raster file.
I want to create a custom function which could be added to the expression of the virtual attribute of the point shapefile.


texts: Delete vertices on polygon feature class, based on the intersection with point shp with python

texts:
I need to delete some vertices from some polygon fc that are unnecessary for the geometry. I have a point shapefile that contains all unnecessary vertices as points - but how can I access the vertices and delete them with a python script or in model builder?

A similar problem has been described here:

Update geometry(polygon) of shapefile based on a point layer
and here: https://community.esri.com/thread/95310 - but I can not install GeoWizard or use any other Software or FME solution because I am no administrator on my pc at work.

I need a python solution ArcGis 10.4.1, probably an UpdateCursor...

I can do this by selecting the single vertices and delete them, but I need an automated solution.

I realized that my adjacent polygons share some vertices (in the second picture the one marked in red) - and my point shp does not give me any information which vertex should be deleted from which of the polygons. That is how some vertices are deleted and the geometry changes (third picture). Probably the red marked vertex should be deleted in the light brown polygon (and not from the one that is a bit darker).

But the script from BERA works perfectly for my initially described problem. I need to find a solution how to find the connection from the vertex to the polygon.


texts: Are there any chances to add grid on Folium map? Something like on this picture.


texts: Replacing Z-field in vector file with attribute table column for export as a Z-enabled KML from ArcMap?

texts:
Using ArcGIS Desktop, I have been unable to take a height field contained within the attribute table of a vector file, and make that height field the z-value of that point.

ArcGIS Desktop has a tool called "Adjust 3D Z" which modifies existing Z values according to a field or a fixed value, but doesn't replace existing values, or add a value if the value is 0 already.

Any ideas?


texts:
Using Qgis 2.8
Trying to label a line/road layer, pull up the dialogue box, everything in the dialogue box seems to work except the labeling doesn't appear on the map. All of the labels don't work for any of the layers, which lead me to believe its something to do with QGIS settings. The shapefiles themselves output just fine.

Simply when I click 'apply changes' nothing happens, but otherwise I receive no errors saying that something on qgis is wrong. What could be the issue?


texts: Using gdal_merge to merge a large amount of files (in pairs)

texts:
I'm new to python and how it can be applied to gis.
I have a large amount of MODIS files to go through and I've already extracted, reprojected and converted to GTiff the bands that I need from each HDF file using convertModisGDAL from pymodis.

Now I need to mosaic pairs of this file list. More specifically, I have this kind of files and I need to merge the h24's with the h25's:

MOD11A2.A2000065.h24v04.005.2007176173328_LST_Day_1km.tiff
MOD11A2.A2000065.h25v04.005.2007176173328_LST_Day_1km.tiff

MOD11A2.A2000065.h24v04.005.2007176173328_QC_Day.tiff
MOD11A2.A2000065.h25v04.005.2007176173328_QC_Day.tiff

Do you have any suggestions on how can I accomplish this with gdal_merge?

I was trying to do so something with this but I have not been very successful:

for root, dirs, files in os.walk(inputFolder):

    for hdf24 in fnmatch.filter(files,pattern_24):
        p1 = os.path.realpath(os.path.join(root,hdf24))

    for hdf25 in fnmatch.filter(files,pattern_25):
        p2 = os.path.realpath(os.path.join(root,hdf25))

    sys.argv = ['-o','outputTif12','-of','GTiff','-ot','Float64','-n','0','-a_nodata', '-9999', p1, p2]
    gdal_merge.main()


texts: QGIS: get unique values of chosen attributes when editing layer style

texts:
In order to style a map comprised by OSM data (shapefiles), I want to use rule-based style.
This is very handy, however: the data has so many attributes with lots of values that I somehow loose the overview.

Suppose I want to style different road types. In the properties window, I open the "fields" tab to see my options. I decide that I want to style according to the field "highway". Is there a way to see the possible values for this attribute without opening the attribute table? What I want to see is something like a list

primary
secondary
tertiary
motorway
residential
...

This could enable making good and fast decisions on layer styles.
If there is such an option, I hope someone can indicate it.


texts: Calculating road density in R using kernel density?

texts:
I have a large (~70MB) shapefile of roads and want to convert this to a raster with road density in each cell. Ideally I'd like to do this in R along with GDAL command line tools if necessary.

My initial approach was to directly calculate the lengths of line segments in each cell as per this thread. This produces the desired results, but is quite slow even for shapefiles much smaller than mine. Here's a very simplified example for which the correct cell values are obvious:

require(sp)
require(raster)
require(rgeos)
require(RColorBrewer)

# Create some sample lines
l1 <- Lines(Line(cbind(c(0,1),c(.25,0.25))), ID="a")
l2 <- Lines(Line(cbind(c(0.25,0.25),c(0,1))), ID="b")
sl <- SpatialLines(list(l1,l2))

# Function to calculate lengths of lines in given raster cell
lengthInCell <- function(i, r, l) {
    r[i] <- 1
    rpoly <- rasterToPolygons(r, na.rm=T)
    lc <- crop(l, rpoly)
    if (!is.null(lc)) {
        return(gLength(lc))
    } else {
        return(0)
    }
}

# Make template
rLength <- raster(extent(sl), res=0.5)

# Calculate lengths
lengths <- sapply(1:ncell(rLength), lengthInCell, rLength, sl)
rLength[] <- lengths

# Plot results
spplot(rLength, scales = list(draw=TRUE), xlab="x", ylab="y",
       col.regions=colorRampPalette(brewer.pal(9, "YlOrRd")),
       sp.layout=list("sp.lines", sl),
       par.settings=list(fontsize=list(text=15)))
round(as.matrix(rLength),3)

#### Results
     [,1] [,2]
[1,]  0.5  0.0
[2,]  1.0  0.5

Looks good, but not scaleable! In a couple other questions the spatstat::density.psp() function has been recommended for this task. This function uses a kernel density approach. I am able to implement it and it seem faster than the above approach, but I'm unclear how to choose the parameters or interpret the results. Here's the above example using density.psp():

require(spatstat)
require(maptools)

# Convert SpatialLines to psp object using maptools library
pspSl <- as.psp(sl)
# Kernel density, sigma chosen more or less arbitrarily
d <- density(pspSl, sigma=0.01, eps=0.5)
# Convert to raster
rKernDensity <- raster(d)
# Values:
round(as.matrix(rKernDensity),3)

#### Results
      [,1] [,2]
[1,] 0.100  0.0
[2,] 0.201  0.1

I thought it might be the case that the kernel approach calculates density as opposed to length per cell, so I converted:

# Convert from density to length per cell for comparison
rKernLength <- rKernDensity * res(rKernDensity)[1] * res(rKernDensity)[2]
round(as.matrix(rKernLength),3)

#### Results
      [,1]  [,2]
[1,] 0.025 0.000
[2,] 0.050 0.025

But, in neither case, does the kernel approach come close to aligning with the more direct approach above.

So, my questions are:

How can I interpret the output of the density.psp function? What are the units?
How can I choose the sigma parameter in density.psp so the results align with the more direct, intuitive approach above?
Bonus: what is the kernel line density actually doing? I have some sense for how these approaches work for points, but don't see how that extends to lines.


texts: Using the jQuery UI Slider to control layer opacity in OpenLayers

texts:
Added some bounty so more people take a look at my question. Thank you.

I have been trying to utilize the jQuery UI Slider to control layer opacity in OpenLayers but I have an issue when trying to set the "value" option to anything other than the "max" value.

What I would like to see is the setOpacity method change my layer's opacity to 20% onload. However, the slider is shown in the correct position but the layer displays the "max" value of opacity (100%/opaque). Currently, a slide event triggers a change in opacity. I need to force the opacity of my WMS layers to be that of the manually set "value" options when the page loads.

Removed the any unnecessary code.
One slider shown below.
Your help is much appreciated!

<body onload="init()">
  <script type="text/javascript">
  $(document).ready(function(){
  $(function() {
      $( "#slider1" ).slider({
      range: "min",
      min: 0,
      value: 20,
      slide: function(e, ui) {
          hii_1.setOpacity(ui.value / 100);
          $( "#amount1" ).val( ui.value );
          }
      });
      $("#amount1" ).val($( "#slider1" ).slider( "value" ) );
  });
  });
  </script>
</body>


texts: Style for OSM in QGIS that can be exported to SLD

texts: I need some production ready styles for OSM maps, which I can export into SLD. I am looking for styles for a long time and unsuccessful, so I decide create custom style by some manual. Can you give me any advice for a book or some tutorial where I proceed how to create nice looking map style from scratch?


texts: Calculating Annual Mean NDVI for set of features using Google Earth Engine

texts:
I'm trying to calculate the annual mean NDVI for a set of polygons that are stored in a table in GEE. I want the annual mean for a series of years to be added to the table as new fields. I have the code to use the MODIS dataset to extract by day but I can't figure out how to get the annual mean. Essentially, this output works, I just need to the annual mean instead of daily value.
var modisNDVI = ee.ImageCollection("MODIS/MCD43A4_006_NDVI").filterDate('2015-01-01', '2015-01-09');

// create dummy image.
 var stacked_image = ee.Image();

// create function to combine all images into one image as bands.
 var combine = function(img, stack){
  stack = ee.Image(stack);
  return stack.addBands(img.select([0], [img.get('system:index')]));
};

var rename_band = function(img){
  return img.select([0], [img.id()]);
};

// stack all the images into a single image.
stacked_image = modisNDVI.map(rename_band).toBands(); //iterate(combine, stacked_image);

// determine scale to perform reduceRegions.
var scale = modisNDVI.first().projection().nominalScale();
    var ndviMean = ee.ImageCollection(years).map(function(img){
      var rR = img.reduceRegions({collection: regions,
                                reducer: ee.Reducer.mean(),
                                scale: 500,
                                tileScale: 4})
      return rR.map(function(feat){
        return feat.set('year', img.get('year'))
      })
}).flatten()

Export.table.toDrive(ndviMean,"NDVI","NDVI","NDVI_zones");


texts: QGIS: Grass Tools: Least cost path using r.drain generating flawed results

texts:
I am trying to use a combination of r.walk and r.drain in the QGIS grass plugin to model the most efficient path between two points, where cost is accumulated principally by movement over deeper water. Ultimately I want to measure the A-B distance whilst sticking within the depth range of the species I am modelling.

I have created a reclassified raster using r.walk, which is coded in 20 equal steps according to cumulative cost travelling from the point in the South West Atlantic.

I then use the r.drain tool to calculate the most efficient route from the red starting point to a finishing point just east of the British Isle for arguments sake. The output is a raster file which I am then attempting to convert to a vector using r.to.vect.line, but this process appears to freeze and I am forced to cancel it. There is no explanation given and so I have no idea why the process is freezing or otherwise not working. Can anyone offer a suggestion?


texts: Outlier filter eliminates parts of grid when applying grid_metrics

texts:
I need to calculate metrics for a normalized and outliers filtered LAScatalog. I am using grid_metrics from the lidRpackage in R. The code for filtering outliers is based in this vignette and, although it seems like it is working, the resulting metrics raster produced lacks some zones of the study area (see Image 1). The eliminated parts seem to be crops, I suspect that the problem comes up when heights are close to 0.
I applied grid_metrics directly to the normalized data to discard errors in the raw data and the resulting metrics covered the missing zones (see Image 2).
The script of the outlier filter is the following:

#Basic function detect outliers is to measure the 90th percentile of height in 10 x 10-m pixels (area-based approach) and
#then remove the points that are above the 95th percentile in each pixel plus, for example, 10% (sensitivity = 1.1).
lasfilternoise = function(las, sensitivity)
{
  p90 <- grid_metrics(las, ~quantile(Z, probs = 0.90), 10)
  las <- lasmergespatial(las, p90, "p90")
  las <- lasfilter(las, Z < p90*sensitivity)
  las$p90 <- NULL

  return(las)
}

#Run the previous function for the catalog clusters
lasfilternoise.LAScluster = function(las, sensitivity)
{
  # The function is automatically fed with LAScluster objects
  # Here the input 'las' will a LAScluster

  las <- readLAS(las)                          # Read the LAScluster
  if (is.empty(las)) return(NULL)              # Exit early (see documentation)

  las <- lasfilternoise(las, sensitivity)      # Filter the noise

  las <- lasfilter(las, buffer == 0)           # Don't forget to remove the buffer

  return(las)                                  # Return the filtered point cloud
}

norm <- readLAScatalog("C:/JAIME/dataLAS_test/norm/") #Directory to normalized data

#Options for the catalog engine
opt_filter(norm)       <- "-drop_z_below 0 -keep_first"
opt_chunk_buffer(norm) <- 10
opt_output_files(norm) <- "C:/JAIME/dataLAS_test/norm_p90/{ORIGINALFILENAME}_p90"

options <- list(
  need_output_file = TRUE,    # Throw an error if no output template is provided
  need_buffer = TRUE,         # Throw an error if buffer is 0
  automerge = TRUE)           # Automatically merge the output list (here into a LAScatalog)

#Apply functionto the catalog
output  <- catalog_apply(norm, lasfilternoise.LAScluster, sensitivity = 1.1, .options = options)
ctg_denoised = catalog(unlist(output))


texts: Transformation on insertion into PostGIS

texts:
If I have a specific SRID set on a geometry column in PostGIS, say 4326; if I then insert geometry with a different SRID, say 4283, what SRID will be stored in the database, will PostGIS transform the geometry to the column SRID on entry or will it just not allow insertion?
INSERT INTO sometable(geom) VALUES ST_GeomFromText('POINT(-38.2481 141.7370)',4283)


texts: add multiple rasters, each located a branch too deep in separate folders - ArcGIS

texts:
I have a year of daily raster datasets (365) that I would like to use the Cell Statistics tool on; however, each raster is located in a separate folder so I cannot add all of them at once.

ie:

-raster_0001
       -raster_0001
            dblbnd.adf
            hdr.adf
            prj.adf
            sta.adf
            w001001.adf
            w001001x.adf
       -info
            arc.dir
            arc0000
            arc0000.nit
            arc0001
            arc0001.nit
        raster_0001.aux
-raster_0002
       -raster_0002
            dblbnd.adf
            hdr.adf
            prj.adf
            sta.adf
            w001001.adf
            w001001x.adf
       -info
            arc.dir
            arc0000
            arc0000.nit
            arc0001
            arc0001.nit
        raster_0002.aux

-raster_0003
       -raster_0003
            dblbnd.adf
            hdr.adf
            prj.adf
            sta.adf
            w001001.adf
            w001001x.adf
       -info
            arc.dir
            arc0000
            arc0000.nit
            arc0001
            arc0001.nit
        raster_0003.aux
...

-raster_0365
       -raster_0365
            dblbnd.adf
            hdr.adf
            prj.adf
            sta.adf
            w001001.adf
            w001001x.adf
       -info
            arc.dir
            arc0000
            arc0000.nit
            arc0001
            arc0001.nit
        raster_0365.aux

I have done a test and if I add a few rasters to a file geodatabase, then I can select these files all at once in the Cell Statistics tool.

My question: Is using python to add these rasters to a file geodatabase the best option? If yes, or no, how come and what would the code/workflow look like?

My guess is "Yes", and would be as simple as iterating through each folder using the "arcpy.RasterToGeodatabase_conversion" method described here raster_to_gdb...but I am not sure how that would look exactly? I am also not so sure that would be the best solution and thought it useful to ask the community for opinions.I would think this is a common problem for people: iterating through multiple folders to add content from each into a geodatabase (for example).

Is a geodatabase even the way to go? I haven't tried an approach using a mosaic raster set yet, not sure if that would work. Catalog? Another tool? Perhaps using Model Builder instead? Thoughts?


texts: Technical comparison of Oracle and PostgreSQL

texts:
I have searched on google and the archives for a technical comparison of Oracle and PostgreSQL (also MySQL), unsuccessfully.

Is there any comparison that covers more than just the basic differences? Maybe you can give me some comparison?


texts: Setting Rasters to the exact same display

texts:
I have been tasked to make a series of .jpeg outputs of different sections of a dike displaying elevation data.
My problem is that the data i have been given is GIGANTIC, almost the whole region made of of lots of small squares of raster - none of which have logical names, just random.

The normal way i would deal with these rasters would be to place them in a raster catalog so that they all have the same display so that the seperation between different rasters will not be obvious. The problem is that there are literally 600 tiny squares of raster making up my dataset and my raster catalog has been loading these rasters for almost an hour now with no sight of an end point!

Is there a way i can manually go though each raster (im very fast at clicking!!) or any otherway anyone can think of to make these rasters uniform in their display within arcmap 9.3.1.?


texts:
Like the author of GRASS not working with QGIS 3.6, I am unable to use GRASS plugins, with the error:

This algorithm cannot be run :-( The specified GRASS 7 folder "C:\OSGEO4~1\bin\bin" does not contain a valid set of GRASS 7 modules. Please, go to the Processing settings dialog, and check that the GRASS 7 folder is correctly configured

This applies whether or not I start QGIS using the "QGIS Desktop 3.6.3 with GRASS 7.6.1" option or not.
I do not have an option to specify the GRASS path in the Settings > Option > Processing menu.


texts: Using state boundary to complete polygon in QGIS

texts:
I'm creating polygons in QGIS using lat/long pairs. These polygons are regions within the state of NY. For each polygon, I know all of the interior points (within the state of NY) and all of the points on the state boundary.

How in QGIS do I use the physical state boundary to complete each polygon (sort of like a convex hull using the state boundary)?


texts: Possible to convert v.net.spanningtree to QGIS custom script so don't need GRASS installed?

texts:
Is it possible (and if yes then how or has this been done) to convert v.net.spanningtree in to a QGIS/python (custom) script that can be accessed in the QGIS graphical modeler in a way that GRASS is not needed to be installed?

Use case:
I have built a model in the QGIS graphical modeler that uses the GRASS spanning tree algorithim, however, I would like this model to work on other people's computers who may not be able to easily install GRASS (or access its plugins through QGIS GUI) so they can use the model, as well as hopefully we could install it on android tablets (that may not have or even be able to install GRASS?).

This model (which currently collates about 50 algorithms which I hope I can simplify a bit!) is currently working well on my desktop but ultimately hoping to get it working on android tablets to be used in the field (though haven't tested QGIS on android yet so not sure of any limitations that may be encountered with that ... though hopefully soon!).

Notes:

The QGIS plugin "Reconstruct Line" (https://github.com/nextgis/ReconstructLine), using the "Insert multiple lines" functionality, also essentially achieves the same result as the GRASS v.net.spanningtree if this is any help as perhaps it would be easier to make this availiable in the graphical modeller/as a custom script.

I am using QGIS 2.16.2 on Windows 7 (but would ultimately need a solution that is workable on other systems)

If the only solution is to create a script from scratch it might be a bit beyond my capabilities at the moment. I do not have python/scripting skills but may be able to work my way through using trial and error in some cases if any proposed solution is close to the mark and not too complex to follow.


texts: Comma delimited Text failed to read my imported CSV file columns

texts:
As you can see after saving my Excel file as CSV (comma delimited) and importing it to the Add delimited layer tool in QGIS, it failed to read the 4 columns and instead read it as 1 single column.
Any idea on how to fix this?


texts: Maximum and mean value is different in Python than QGIS 3.6

texts:
I read a Landsat-8 TIFF file using Python Rasterio to compute the max and mean using the following code:

path = "LC08_L1TP_115036_20181025_20181031_01_T1_B7.TIF"
test = rasterio.open(path).read()
print("value of max and mean", test.max(), test.mean())

value of max and mean 65535 4390.264887378426

But when I opened the same image in QGIS 3.6 it gives the different values as mentioned below:

STATISTICS_MAXIMUM=41986
STATISTICS_MEAN=4629.7234448245
STATISTICS_MINIMUM=0
STATISTICS_STDDEV=3599.2731482793
STATISTICS_VALID_PERCENT=100

Can anyone tell me the reason for the same which one I should consider correct?

Moreover, since 65535 is the highest possible value in 16-bit integer, what does the pixels having value of 65535 in Band 1 to 10 of Landsat-8 denote?


texts: Having changes to a SpatiaLite table column to automatically update GIS layer

texts:
For example, if the length of a line is changed (direct SQL), I would like the layer to reflect the new line length and also update the label for the feature.  How can that be accomplished?
It appears that PostgreSQL Triggers can send a notification to QGIS, but I cannot find the equivalent in for SpatiaLite.


texts: Getting Lat/Long for region and district of Mozambique

texts:
I have a database with region and district names in Mozambique Africa.
I was wondering if it is possible to get lat/Long coordinates of those names so that I can map those locations. In a way going backward?


texts: Raster Calculator - ERROR 000539 - Fail to complete "weighted overlay" calculation

texts:
I have been constantly getting this "ERROR 000539" message when trying to run raster calculator using two rasters. I was having issues running raster calculator initially and made sure that the bit depth and raster type were matching across all rasters. I have been able to get most of my raster calculations to work once I confirmed the bit depth and raster type matched, but am running into this error message on ArcMap and ArcGIS Pro when I try to run this calculation.
Expression: ("[Raster name]" * .5)+("[Raster name]" * .5)

For some process context, I am trying to run the raster calculator tool instead of the weighted overlay tool to create combined rasters. So, by multiplying these two rasters by .5 I am saying that they are receiving equal weight in my combination.
In my research it seems like there is an issue with the python script, but I'm not familiar with python or the process of how to fix a python script within/outside of ArcGIS software.


texts: How do I calculate area error probability of a polygon with knowledge of point errors

texts:
I have a polygon, representing a glacier, which is defined by over 3000 points. This polygon is digitised from an accurate orthophoto. I also have a polygon of said glacier digitised from Landsat. In order to work out the error distribution of the points defining the Landsat polygon, I calculated the shortest distance of each Landsat-polygon point to the extremely accurate orthophoto-polygon. I now know the error distribution of points defining the Landsat polygon (normal distribution centred around a mean close to zero).

My question is:
With knowledge of my point error distribution, how can I estimate an area error probability?

A long search has only shown me how to propagate point errors to area errors for simple shapes (such as squares) with knowledge of x and y errors of points.


texts: Time series condition through field calculator

texts:
I am attempting to filter some time stamped data which I have enabled through the Convert Time Field tool, depending on a true/false condition. I am doing this through the field calculator within Model Builder (Arcmap 10.2).

The data is structured as follows:

I am processing the data to remove points that are less than two hours apart from one another in order to standardize the time intervals. The data is also grouped by series ID (cat) and I have a assigned a sequential field to list the chronological order of the points grouped by ID. I have also added a field which should be populated by 0 or 1 depending as to whether the time condition is met.

How I am trying to implement this condition is through a series of IF ELIF statements. There are three options to this:

Firstly the record is returned true by default as it the start point of the series, this is shown by the sequence value = 0. The timestamp when sequence = 0 is also taken as the variable n.

Secondly, the record is returned as true if the timestamp is atleast 2 hours in the future from the previous record's timestamp (n). n Is then altered to be (n + 2 hours) to be used in the next record's condition.

Thirdly if the timestamp is not 2 hours ahead of the previous record's it is labelled as false and n does not change.

I am attempting this within a python codeblock in field calculator on the conditional field (time_frame_2hr) using the code:

  import datetime

  def fnTime(sequence,vms_time,time_con,order):
    n = df.loc[df[sequence]==0,vms_time]
    previous_time_con = df.loc[df[order],(order-1),time_con]
    if sequence == 0:
      return 1
    elif sequence != 0 and vms_time < (time(n) + (datetime.timedelta(hours=2))):
      return 0
    elif sequence != 0 and previous_time_con == 0 and vms_time >= (time(n) + (datetime.timedelta(hours=2))):
          return 1
    elif sequence != 0 and previous_time_con == 1:
          n = (time(n) + (datetime.timedelta(hours=2)))
          if sequence != 0 and vms_time >= (time(n) + (datetime.timedelta(hours=2))):
            return 1
          elif sequence != 0 and vms_time < (time(n) + (datetime.timedelta(hours=2))):
             return 0
          else:
             return 0

I have set the variables within field calculator using the fields:

fnTime( !sequence! . !date! . !time_frame_2hr! . !OBJECTID!)

The tool runs with no errors however it assigns a false value (0) to all records. I am still learning writing python functions so they're maybe a quicker or efficient way of doing this.


texts: Showing wrong information for Sentinel 5p L2 CO data in VISAN using HARP

texts:
I am using VISAN (Anaconda) IDL for opening Sentinel 5P L2 CO dataset in netcdf format. Whenever I have tried to open the file using the code

S5P_L2_CO = harp.import_product(r"F:\sentinel 5p\New folder\S5P_OFFL_L2__CO_____20200101T055716_20200101T073846_11490_01_010302_20200102T194736.nc")

It works Ok. But when I use

wplot

then it shows "Surface Meridional Wind Velocity". No CO is found. I have tried with three dataset. The visualization is same.

My quesiton is am i doing any wrong or the dataset wrongly submitted in the sentinel5p server.

Here it need to mention, I have used SO2 data, it works like as it should be.


texts: I have a file from MATLAB in (.mat) format. I want to import this in GRASS GIS, but I get the following warning.


texts: Is it possible to add a new line in a non-spatial table using Qgis?

texts:
While using QGIS, and having loaded a non spatial table (from a dbf or a Spatialite table), is it possible to add a new row to its table attributes?

Normally in a spatial layer you have to use "add features button", that obviouslly won't work in this case.

My idea was to use Qgis as somekind of frontend to a Spatialite database, with the ability to edit\add features in spatial tables as well in the non-spatial ones.


texts: Persisting layer tree between QGIS instances using PyQGIS

texts:
I have been working on a PyQGIS project for a bit, and I have been using QgsProject.instance().addMapLayer to add layers to the map.

However, now I need to add layers to my canvas and legend in specific groups, so I changed my logic for adding layers to the map to use the underlying QgsLayerTree object to add groups and add layers to those groups. I got it working, but it also introduced a new problem:

Whenever I add a layer to the map within a layer group using the new logic, the layer is not available whenever QGIS is restarted. The QgsLayerTreeLayer instances of the layers are still available, but it looks like they can't resolve their references to their source or other metadata. From what I can tell, it looks like they only have access to the layer ID.

How do I persist the layer tree between QGIS sessions? Or is there another way to implement adding layers within groups that would avoid this problem?

Adding Layers to Map

layer = QgsRasterLayer(...)
root = QgsProject.instance().layerTreeRoot()
group = root.insertGroup('Layers')
group.insertLayer(0, layer)

Reading Layers after QGIS Restart

root = QgsProject.instance().layerTreeRoot()
layer_tree_layers = root.findLayers()  # Returns list of QgsLayerTreeLayer objects

map_layers = []

for lyr in layer_tree_layers:
    map_layer = lyr.layer()  # Once QGIS has restarted, returns None!
    map_layers.append(map_layer)

I saw that in the QGIS API Documentation, it noted that the QgsLayerTreeLayer.layer() may return None/nullptr "...in the case of a layer node representing a layer which has not yet been fully loaded into a project, or a layer node representing a layer with an invalid data source...."

So why does my QGIS workspace persist normal map layers between instances, but not QgsLayerTreeLayers that were added to the map programmatically or is there an error in my logic that is causing this?

Note: The layers that I am loading into the map are QgsRasterLayers that have been read from geopackage tables. Didn't think it had anything to do with this because I am instantiating the layer objects from the source in the same way as I would any other layers. But just thought I'd mention it.


texts: Does FLEXLM_TIMEOUT environment variable work with ArcGIS Desktop

texts:
In versions of ArcGIS Desktop earlier than 10.0, you could create/set the FLEXLM_TIMEOUT to a higher value for slower connections.  I'm currently trying to troubleshoot a system (one where I have minimal control) with slow connections to the license manager.
I was able to create a FLEXLM_TIMEOUT user variable and fired up ArcMap several times with no timeout error but maybe I'm lucky and the network I'm on is having no issues.
Does anyone know if this environment variable is actually used by ArcGIS (FLEXnet Publisher)?


texts: How to link two different images in eCognition

texts: I'd like to link two different images in eCognition to find the same pixels on both images. The images have same projection (geog lat/long). How can I do this?


texts: How do I get current map rotation value in MapBox?

texts:
I have custom drawn icons that I am rendering on top of a MapBox instance at angles that correspond to their direction of travel, but things start to look at little weird when you right-click and rotate the map.

The icons continue to render at the same angle, while the actual direction is now altered by the map rotation. I've got planes flying sideways...

How do I get the map rotation value from the MapBox instance, so I can automatically adjust the rotation of my icons to account for this?


texts: Getting quantity of intersecting lines of polygons in QGIS?

texts:
im looking for a way to get the quantity of lines from a line-layer that intersect the different polygons of a polygon-layer;

I have:

One line-layer with a clutter of lines, defined by ID, which represent the flight-activity of birds in a specific area. The layer is called "lines".
One Polygon-Layer with a specific grid. The layer is called "grid".

I need for every square-polygon from the "grid"-layer the quantity of different lines from the "lines"-layer, which intersect them. The Value is supposed to show up in a new column of the attribute-table of the "grid"-layer called "LineCount" in this example (the layers in the screenshot are just examples, the original ones have much more lines, so it is not practical to count them by hand ;-)).

The aim is to visualize the density of flight-activity in a raster-like map.

I guess there might be some solution within the refFunctions-plugin, but I'm not able to figure it out.


texts: Point clustering only works for part of area using ArcGIS API for JavaScript

texts:
I am using the 3rd party scripts for point clustering as published by Globoserve, http://globoserve.posterous.com/ .  I followed the instructions to the letter - and just replaced the queryTask variable to point to my point features that are projected in the one of the compatible projection systems.

But the point cluster alogrithm will only cluster a section of points rather than all of my points - see image below.

Can anyone explain why this is happening?

I have tried this with a totally different point layer in US - and it shows a 500 point limit - the same as the pic below!


texts: Merging a raster file to polylines according to their attributes

texts: In ArcGIS Pro, I have a raster layer of slope-aspect of Toronto and I want to merge it another feature which contains bus routes. How can I do it? Join tool is useless since the raster layer is not a feature layer


texts: How to convert 4 band GeoTIFF to PNG or JPEG using gdal_translate

texts:
I have an orthomosaic GeoTIFF and it has 4 bands. I want to convert it to a PNG/JPEG image to show it on Cesium. I've tried gdal_translate like this :
gdal_translate -a_nodata 0 -of PNG 1.tif 1.png -scale -colorinterp_1 red  -colorinterp_2 green -colorinterp_3 blue and it's not working. I'm completely new to QGIS. Can you please guide me what options I should use to get data?
I already have 2D generated picture using photoshop and I think it's not losing any data. I don't want to keep projection or elevation data. I just want a normal image.

Point: it's not possible for JPG because I'm getting this error : libjpeg: Maximum supported image dimension is 65500 pixels

This is gdalinfo result :

Files: /Users/tfiroozian/Downloads/2_111318_SKYG_Nautilus_1_3mi_transparent_mosaic_group1.tif
Size is 41705, 77008
Coordinate System is:
PROJCS["NAD_1983_StatePlane_Texas_Central_FIPS_4203_Feet",
    GEOGCS["NAD83",
        DATUM["North_American_Datum_1983",
            SPHEROID["GRS 1980",6378137,298.2572221010042,
                AUTHORITY["EPSG","7019"]],
            AUTHORITY["EPSG","6269"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4269"]],
    PROJECTION["Lambert_Conformal_Conic_2SP"],
    PARAMETER["standard_parallel_1",30.11666666666667],
    PARAMETER["standard_parallel_2",31.88333333333333],
    PARAMETER["latitude_of_origin",29.66666666666667],
    PARAMETER["central_meridian",-100.3333333333333],
    PARAMETER["false_easting",2296583.333333333],
    PARAMETER["false_northing",9842500],
    UNIT["US survey foot",0.3048006096012192,
        AUTHORITY["EPSG","9003"]]]
Origin = (1285949.954140000045300,10626980.503550000488758)
Pixel Size = (0.081470000000000,-0.081470000000000)
Metadata:
  AREA_OR_POINT=Area
  TIFFTAG_SOFTWARE=pix4dmapper
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  ( 1285949.954,10626980.504) (103d35'10.24"W, 31d46'58.51"N)
Lower Left  ( 1285949.954,10620706.662) (103d35' 8.11"W, 31d45'56.45"N)
Upper Right ( 1289347.660,10626980.504) (103d34'30.89"W, 31d46'59.49"N)
Lower Right ( 1289347.660,10620706.662) (103d34'28.77"W, 31d45'57.43"N)
Center      ( 1287648.807,10623843.583) (103d34'49.50"W, 31d46'27.97"N)
Band 1 Block=41705x1 Type=Byte, ColorInterp=Red
  NoData Value=-10000
  Mask Flags: PER_DATASET ALPHA
Band 2 Block=41705x1 Type=Byte, ColorInterp=Green
  NoData Value=-10000
  Mask Flags: PER_DATASET ALPHA
Band 3 Block=41705x1 Type=Byte, ColorInterp=Blue
  NoData Value=-10000
   Mask Flags: PER_DATASET ALPHA
Band 4 Block=41705x1 Type=Byte, ColorInterp=Alpha
  NoData Value=-10000


texts: Reusing tiles between zoom levels of Leaflet?

texts:
I have a Leaflet map with a GeoJSON Tile Layer and a custom backend for serving geoJSON tiles. The geometries returned from the backend are simplified depending on the zoom level, but within certain bins of zoom levels (1-3, 4-6, etc) the geometries are identical. Right now, new tiles will be fetched for every single zoom level, and I was wondering whether there is a way to reuse tiles when going to a zoom level within the same bin. This would decrease the server load and reduce flickering when zooming in or out.

For example: when going from zoom level 4 to zoom level 5, no new tiles should be fetched and the existing tiles should be stretched. Likewise, when going from zoom level 5 to zoom level 4, the existing tiles should be kept and additional tiles should be fetched to complete the map.

Any thoughts?


texts: Importing CSV file with csvt data in QGIS (2.8) will not assign type "date" to fields

texts:
I have a dataset (in CSV) and am trying to add several fields with date information, so I can filter points based on this data and generate maps (in QGIS 2.8) which will form a chronological series. To accomplish this, I've done the following:

I've taken my original CSV and added four fields (see below for final CSV file, date data is in columns 12-15, titled award1-4).
Data was originally in YYYY-DD-MM format, so I've used the handy command-line csvfix tool to set all the data to be ISO format (YYYY-MM-DD).
To accompany the CSV file, I've generated a CSVT file, with the same name and path to accompany my CSV, following the instructions in the regularly linked blog post by Anita Graser and the related GDAL docs so my CSVT uses only the following types, "String","Integer","Real","Date" (the fourth of these is the crucial one)
I have also followed the helpful suggestion in this thread by Underdark / Kristen to use textwrangler to make sure characters are standard, and removing a second line and converting smart quotes resolved an issue in QGIS with importing.

After following these steps my CSV imports just fine in QGIS using 'Add Delimited Text Layer', but the layer only has Real (as "Double"), Int and String (as "Text"), with the date fields rendered as "String" as well. I'd like to have this CSV file import into QGIS with those date fields as "date" data type in QGIS so I can run queries on them. I'd also like to stick with CSV as the format as this data will be placed in an open repository eventually and I'd like other users to be able to view the data without QGIS, but also be able to import relatively painlessly into their own projects.

Has another user accomplished this in QGIS, and more to the point, what step am I missing here to make these fields import properly as "Date" field data type.

Samples of data are as follows:

ECS-GIS-Locations_2.4.csvt:

"Real","Real","Integer","String","String","String","String","String","String","String","Integer","Real","Date","Date","Date","Date","String","String","String","String","String"

ECS-GIS-Locations_2.4.csv (first few lines):

"X","Y","accuracy code","postcode","unique-id","name","network","denomination","notes","website","awards","icon-size2","award1","award2","award3","award4","wildlife-biodiversity","solar","date-built","churches trust URL","OS ref"
"349344.8977","1053019.497","2","KW17 2BU","2203","St. Ann's, Papa Westray","Orkney Islands","C of S","","","4","3","2003-10-09","2006-08-19","2009-07-07","2013-11-01","","x","1841","http://www.scotlandschurchestrust.org.uk/church/st-anns-papa-westray","HY495516"
"345691","1046915","3","KW17 2DB","2202","Westray Church","Orkney Islands","C of S","","","4","3","2003-10-09","2006-08-19","2009-07-08","2013-11-01","","","1846","http://www.scotlandschurchestrust.org.uk/church/westray-parish-kirk","HY457462"
"272400.702","655262.3425","2","ML3 7DT","2801","St John's Church Of Scotland, Hamilton","South Lanarkshire","C of S","","http://www.stjohnshamilton.org.uk","3","2.5","2003-05-01","2006-03-05","2009-07-09","","","","1835","http://www.scotlandschurchestrust.org.uk/church/st-johns-church-hamilton","NS724553"


texts: Problem creating a simple two-branch if-then-else logic (script tool) in ArcToolbox

texts:
I'm new to Arcpy and script writing and i'm trying to create a two-branch if-then-else logic based on the shapefile type of the input data. What I want to do is to write a script tool that executes action A if the input shapefile is a Polygon and executes action B if the input shapefile is a Line.
I have browsed many websites and the code below, is what I wrote but the ArcMap cannot run the script, showing this message: "SyntaxError: invalid syntax"

the code is as followed:

    # -*- coding: utf-8 -*-
#-----------------------------------------------------------------------------

# Import arcpy module
import arcpy

# Get the feature class from the tool.
#
FL = arcpy.GetParameterAsText(0)

# Determine the shape type of the feature class.
#
dscFL = arcpy.Describe(FL)

try:
  # Set tool output parameters based on shape type.

    if dscFL.ShapeType.lower() == "polyline":
        arcpy.AddMessage("Feature Type is polyine")
        arcpy.SetParameterAsText(1, "true")   # is polyline
        arcpy.SetParameterAsText(2, "false")  # is not polyline

    elif dscFL.ShapeType.lower() == "polygon":
       arcpy.SetParameterAsText(1, "false")   # is not polygon
       arcpy.SetParameterAsText(2, "true")    # is Polygon

I would be thankful if anyone can help me.

I have tried to modify the Parameter tab, which has resulted to this:

the model looks like this:


texts: How to enforce float values for NDVI in ArcGIS?

texts:
I am trying to calculate NDVI from Landsat 4-5 TM imagery (of the year 2008) using the formula NDVI = (band4-band3)/(band4+band3) from Raster Calculator of ArcGIS 10.4.1.
The problem that I am facing is that the resultant raster that I am getting has integer values. I want float values that stretch from -1 to +1.
How can I resolve the issue?


texts: Creating Workspace for Satellite Imagery in ArcGIS Pro?

texts:
I am trying create mosaic from Sentinel – 2 satellite data. To realize this task I’ve used “New Ortho Mapping Workspace” and “Satellite” Type. In next step I’ve chosen appropriate sensor type, the same Spatial Reference like images, and I’ve indicated the path to the folder witch Sentinel – 2 images. (Folder Name: S2A_MSIL2A_20181014T093031_N0209_R136_T34UFV_20181014T111722.SAFE).

The other options have been default.

After pressing the finish button the message appears:

  Queued...

  Running...
  Building image collection...
  Start Time: środa, 13 marca 2019 09:25:40
  2019-03-13T09:25:43.774: Loading raster datasets
  2019-03-13T09:25:45.013: Completed crawling 14 data source items.
  Added 0 mosaic dataset items.
  2019-03-13T09:25:45.022: Error: 8004205f: No new mosaic dataset item
  was added.
  Succeeded at Wednesday, 13 March 2019 09:25:45 (Elapsed Time: 5,40
  seconds)

What kind of data should be imported to receive a mosaic?

Blended RGB images in .tif format or raw data downloaded form the “https://scihub.copernicus.eu/dhus/” website?


texts: Create SDE-file returns warning but creates the file

texts:
When I create an SDE-file using arcpy I get the following warning:

  WARNING 000565: Could not connect to server.

This is the appropriate code

temp = r'C:\Users\admin\AppData\Local\Temp'
host = "myServer"
sdeFileName = "mySdeFile.sde"
service = "sde:oracle11g:myHost/myInstance:" + user
arcpy.CreateArcSDEConnectionFile_management(temp, sdeFileName, myHost, service, username = user, password = password, version = version)

However the file itself was created successfully and I can even put it into ArcCatalogs sde-connections-folder and open it in ArcCatalog. Thus I doubt the error occurs because of any missing privileges.

As this message is quite generic and can be caused by many different problems I don´t want to swallow it.

Has anyone a solution or at least an explanation for this?

I am using ArcGIS 10.2 and Python 2.7 for 32bit.

EDIT: When I use

arcpy.CreateDatabaseConnection_management(temp, sdeFileName, "ORACLE", service, username = user, password = password, version = version)

instead it works without any problems.
However as I also have some customers that still use Application-Server-connections I cannot use this function as there´s no information on the actual host given within the service-paramater as opposed to the DirectConnect-syntax, where the paramater includes both the host and the database-instance.

So what I need is a function that works with all the following inputs:

Application-Server-connection (deprecated)

host = "myServer"
service = "5151:mySdeInstance"
Direct-Connect-connection

service = "sde:oracle11g:tnsName:user"
where tnsName corresponds the name for that connection within tnsnames.ora.
Easy-Connect-connection

service = "sde:oracle11g:myServer/user:user"


texts:
People who already used or are familiar with TerraSAR X, can you point out:

Main applications and methods
Open source software tools


texts: symbol of a newly created point moves to a different place at different zoom level

texts:
I connect QGIS to an Oracle database. I add a point object on my map in QGIS and the new point is shown in the map, so that is great. But once I zoom out the map the symbol of the point moves to another place, when i click the symbol there is no attribute of the point can be shown so i think it is just the symbol becomes crazy.

I have closed everything and loaded the point layer from the database again but the problem is just the same.

I have checked the point in another map software and it looks normal.

Can it be a bug or something i have done wrong?

This problem just occurs in certain area of the map. If I create a point far away from the problematic area then it is normal again.


texts: Save and load the classes definition table in the QGIS reclassify by table tool

texts:
I'm reclassifying an aspect raster which cells values ranges from 0 to 359.999 using the Reclassify by table tool:

Setting the reclassification table manually is painful, especially when you have many classes:

The two first lines in the reclassify table among many others...
But I see no way to either save or load the class definitions. Surprisingly.
Is it possible to save/load a custom definition of many classes to/from a file in order not to have to re-write them all each time this tool is re-open?
Info:
QGIS version: 3.16.3-Hannover
QGIS code revision: 94ac9f21b8


texts: Creating 3D interactive model from lidar data using ArcGIS for Desktop?

texts:
I have a couple of las files for a city and I need to make a 3D model of the buildings, which should display info about the building when clicked on. I am pretty new to ArcGIS so I don't have the best idea to go about this. I have already downloaded LAStools and can use it to view the lidar data.

I have a point cloud set up but I am not sure how to make the model from this point. Would I need different software? (Sketchup, Meshlab, etc..)


texts: Changing basemaps with different spatial references (ArcGIS Javascript)

texts: I've got a new requirement to switch basemaps that have different spatial references. The app is using ArcGIS Javascript, and switching to Silverlight is not an option at this point in development. I saw this question from 2 years ago and was hoping there was an alternative. Is there an easier way to change basemaps that have different spatial references in 3.6?


texts: How to hide attributes in one field via Arcade expression language - ArcGIS Online

texts:
I have a layer and I am configuring pop-ups, but in one specific field - Name, I have blank records which I don't want to include them in the pop -ups.
Does anyone knows how to create the expression, when the field is blank to hide the entire field? Otherwise to show the information about the name.


texts: Where does ArcGIS Online Network Service source came from? Is there any source reference I can refer to? I tried to search on Google but I could not find anything.


texts: How to work with accentuated strings in Python scripts in ArcGIS Desktop Field Calculator?

texts:
`

pre-logic:
def myCalc(name,nameVal):
  if (name == nameVal):
    return 1
  else:
    return 0

var_name=
myCalc(!Epulettipus!,"Lakóépület")

`

This doesn't work with accentuated strings. (But works with unaccentuated strings). What can I do?


texts:
I am fairly new to to QGIS and I am trying to create a contoured map from point data. The points are regularly spaced at 3.5 metres apart and whenever I try and run the IDW interpolation in QGIS 3 it creates a layer that with the minimum and maximum value the same and nothing is visible in the screen.

I have tried the same in version 2.18 and it does create a layer with different values, but it is gridded and not contoured.

I have managed to use these tools before with other data sets, but I am properly stumped. I have tried different settings and values but nothing seems to give any different results in either version.

Picture 1 shows the layer from QGIS 3, picture 2 shows the settings from 2.18 and pictures 3 and 4 show the results from 2.18 being gridded with the red points laid over the top

[


texts: Loading layers in ArcGIS without displaying it?

texts: When loading any layer into ArcMap (with drag-and-drop), it is automatically displayed/rendered into the canvas. I didn't find any option for disabling this behavior. Any suggestion?


texts: Why and How does copying and pasting a "problem" model with errors in it into a new empty model fix the model (ArcGIS, ModelBuilder)

texts:
I often work with ModelBuilder and Models that have previously worked seem to become unstable or develop issues after lots of editing. They return error messages for parts of the model that did not previously return errors.

A solution often recommended is to simply copy and paste the contents of the model into a new empty model.

Does anyone have expertise or advice to share that can explain why this fixes models and if this is likely to be a permanent fix?

I am developing a set of models that I will be sharing and I need them to be robust enough that they work for the model users.

One of the reasons I am asking this is that often when contacting ESRI support they will focus entirely on the reported error messages within ModelBuilder even when you tell them the Model has previously worked OK, or that it runs OK on a different computer or from the Edit vs open window. It seems to me that very often the error messages returned by models are of no use at all, especially when the problem for example is fixed by copying and pasting into a new empty model. It could just be that these error messages are a symptom of a model that is becoming unstable?

Any thoughts or advice would be very gratefully received.


texts: I created a legend for a Map with ArcMap 10.2.2. I converted it to graphic and want to export or save it as an image. How can I do? Can't figure it out.


texts: Map Algebra with Zero-inflated regression weights

texts:
I have a non-normal count variable y and 6 rasters x1 to x6, each mapped onto a (n x n) grid. I have applied a zero-inflated negative binomial regression using the pscl package of R. I now need to generate a risk map using a raster calculator. In this paper, the authors weight their maps after a Poisson regression using the expression

Exp(\beta_0+\sum(\beta_i * x_i))/Exp(\beta_0)

In the current case I have two additional parameters: the log of theta, and the intercept of the zero-inflation model. How can I modify the expression to incorporate my weights? I cannot upload the data as they are large. Any solution using QGIS, the raster package of R, or ArcGIS will be helpful.


texts: How do you compute the earth's radius at a given geodetic latitude?

texts:
(I see there is an equation on wikipedia that does exactly what I am asking but there is no references. I have no way of confirming the validity of this equation!)

I already understand the difference between Geocentric Latitude vs Geodetic Latitude.

Assuming known semi major, a, and semi minor, b, radiuses are given. How do you compute the radius at a given geodetic latitude?

I need some sort of expert confirmation (derivation, link to derivation, confirmation from expert, explanation, etc).


texts: Find closest points in routing network for a set of points to pick shortest route

texts:
I have the following WITH SQL Query which does select a set of filling stations along a route within the distance of 1000m. In the SQL statement I fetch 1 Point from the station and search the closest point in my routing network, but I want to do this for all stations from the WITH Query

WITH stations AS (SELECT id, geom FROM de_tt_stations AS s
  WHERE ST_DWithin(
    s.geom::geography,
    (SELECT ST_LineMerge(ST_union(geom_way))::geography FROM route),
    1000)
)
SELECT 1378549, destination
FROM (
  SELECT id::integer AS destination
  FROM de_2po_vertex ORDER BY geom_vertex <-> (
    SELECT geom FROM stations LIMIT 1)
  LIMIT 1)
AS foo

From the resulting list then I want to calculate the length of each route so that I can choose the closest filling station. Optimal output would look like

start   | destination      | station_id  | length
---------------------------------------------------------
1378549 | de_2po_vertex.id | stations.id | pgr_dijkstra()
...


texts: Persistence in QGIS under Oracle workspace

texts:
I'm using QGIS 2.18.9 to edit my spatial data, my data is saved in Oracle's database.

When I'm editing a spatial layer I can create, edit and delete features without problems, but my database needs work using Workspace Manager, my spatial tables are versioned.
When I created a feature and save in my database I received this error message:

Commit errors:
  ERROR: 1 feature(s) not added.

Provider errors:
  Oracle error while adding features: Oracle error: Could not retrieve feature id -3
  SQL:
    Error: SELECT "CD_AREA" FROM "ALUNO"."VE_AREA" WHERE ROWID=:f

Has anyone received the same error?
Can QGIS and Oracle work together with workspace/versioned tables?


texts: Collect Values + Model Builder ArcGIS 10.1

texts:
I have a question relating to the use of Collect Values in Model Builder.

I've just come across this post (Piping output of Iterate Raster tool into Raster Calculator tool using ModelBuilder?) and it's pretty much the same as what I want to do. I'm also adding in a user defined weighting to the iteration and then attempting to collect the values before running the cell statistics. The problem is that I'm not getting a set of unique outputs at the Output Values stage and it only return results for the last original rasters.

Any idea where I'm going wrong with this?


texts: How to create a loop which creates single mosaic images for each day (from a selected time period) from an ImageCollection using Python?

texts:
I have been trying to create mosaic images from the region of japan from the MCD19A2.006 imagecollection (https://developers.google.com/earth-engine/datasets/catalog/MODIS_006_MCD19A2_GRANULES#description).
I wish to create a single mosaic image for each day of the year from the period 2013/03/01 to 2020/03/01 and be able to download it locally to my computer. I am able to create an image for one day but I am unsure how to create a loop to do this for each day and have it named differently using python. I would also like to know how I can mosaic all the bands together as currently I am only able to do this with one band at a time.
Below is test code that I did for one day (2020/01/01) using only the Optical_Depth_047 band.
roi= ee.Geometry.Polygon(
    [[[122.791488, 29.85464],
      [122.791488, 45.462539],
      [146.780501, 45.462539],
      [146.780501, 29.85464]]], None, False)

collection = ee.ImageCollection('MODIS/006/MCD19A2_GRANULES') \
    .filterDate('2020-01-01', '2020-01-02') \
    .filterBounds(roi) \
    .select('Optical_Depth_047')

mosaic = collection.mosaic()

Map.addLayer(mosaic, {'min': 0, 'max': 500, 'palette': ['black', 'blue', 'purple', 'cyan', 'green', 'yellow', 'red']}, 'mosaic')


texts: Changing Ascii grid file inside a GeoJSON polygon

texts:
I have a Digital Terrain Model ASCII grid file that gdalinfo recognizes as AAIGrid. Its resolution is 2m. I need to alter it simulate the construction of an embankment, so as to feed it to a modelling software.

I can get the GeoJSON polygon for this embankment, as well as height, thickness and so on. What can I use to alter the terrain model file?


texts: pyqgis standalone script fails to execute clip algorithm

texts:
A previous question has established how to invoke QGIS's clip algorithm at the python console. I have done so and demonstrated that the method works satisfactorily on my data.

I would instead like to use this within a standalone script. To do so I've followed this to create the script below.

When I do so I get the following error message:

Unable to execute algorithm
Wrong parameter value: ./path/to/highways.shp

Can anyone outline what is going on here?

import sys
import qgis.core
import qgis.core as qgc
import PyQt4.QtGui

# Initialise the application
app = PyQt4.QtGui.QApplication([])
qgc.QgsApplication.setPrefixPath('/usr/bin/qgis')
qgc.QgsApplication.initQgis()

# Prepare processing framework
sys.path.append('/usr/share/qgis/python/plugins')
from processing.core.Processing import Processing
Processing.initialize()
import processing.tools

vlayer = qgc.QgsVectorLayer("./path/to/highways.shp", "highways", "ogr")

clayer = qgc.QgsVectorLayer("./path/to/boundaries.shp", "boundaries", "ogr")

processing.tools.general.runalg("qgis:clip", vlayer, clayer, "clip_from_script.shp")

qgc.QgsApplication.exitQgis()
PyQt4.QtGui.QApplication.exit()


texts: Style Sheets: Setting up a style sheet to change only styles, no other fields

texts: I'm trying to create style sheets for my projects. I only wish to change the Style, without affecting labelling or the Fields. And I can't, for the life of me, get it to work. I'm using .qml files and editing using notepad/textpad.


texts:
I have some polyline feature objects in python. Now I want to get all the points of polylines.

For example, if a polyline has start point [0,0] end point [5,5].
Result: [1,1];[2,2];[3,3];[4,4];[5,5].

I want to find all the integer points on that line including end points. For straight line this is dead simple, but if polyline has Beizer Curve, Circular Arc, Elliptic Arc geometry types, then how can I do it?

EDIT:

I can only use those tools which are available in all license levels of ArcGIS. For example, ArcGIS Basic.


texts: Least cost path from multiple sources to multiple destinations across resistance layer

texts: I am using ArcGIS Pro 2.7. I have little experience with ArcPy so using geoprocessing tools would be preferable. I would like to determine the number of boats that could potentially visit coastal islands.  I am using a warves point layer to represent the boat launch sites (starting points), and a coastal islands polygon shapefile as the destination islands. I have a wind exposure shapefile to be used as the surface cost with exposure ranked 1-4 (4 is highest cost) that covers the open water between the warves and the islands.  I would like to simulate boats departing from warves as a least cost path line but they can only travel across the surface cost 1 or 2 (avoiding cost 3 and 4) to reach the islands. Then I would like to sum the number of boats (lines) that reach each island from warves within a 10 km total travel distance.


texts:
I'm very new to spatial stuff in R and am trying to do something I imagine is very simple but am having trouble. I have two spatialpolygonsdataframes. 1 is of all US counties and another is of core-based statistical areas (CBSA's). Each CBSA aligns with county borders, but can cover more than one county. In other words, for example, a CBSA can cover all of one county or several adjacent counties, but no CBSA covers part of a county and no county is covered by multiple CBSA's. Attached is a picture of the CBSA's overlapping the counties (Produced in Qgis)

What I'd like to do is basically create a dummy variable called CBSA with a 1 for CBSA counties (red in the map) and a 0 for those not.
My first hunch after some searching was to use over from sp. I tried over(cbsa,county), which spit out a data frame of 929 observations (the number of CBSAs) with attributes from the county dataframe.

This seems right, except quickly looking at the table shows that the counties are not correctly identified. Multiple observations are tied to the same counties and looking in QGIS verifies that many counties that are part of CBSA's do not appear.

Am I using the right command and just doing something wrong or is there another way I should do this?

my code appears below:

library(sp)
setwd("~/Desktop/test_folder")
load("county.rda")
load("cbsa.rda")
test<-over(cbsa, county)

Data:http://www.filedropper.com/testfolder


texts: Using color matching on two geotiff rasters datasets and exporting geotiffs

texts:
I have one geotiff (tiff + tfw file) that serves as a basemap for a global zone. Then, I have several geotiffs, one for each zone of this big orthophotograph (50 zones in total).

I would like to match the colors of each of the zones geotiffs with the base geotiff. Each geotiff was imported as a raster dataset (File -> Add data and chose the geotiff data from the arcCatalog linked directory). Unfortunately, I cannot find the color correction tab in the raster dataset properties.

I can create a mosaic, but since the base geotiff is pretty big, I would like to know the better way to do this. I've searched through the whole ArcGIS doc/google and all docs were not applicable to my case (I'm using raster datasets, ont catalogs, mosaics or workspaces, and cannot find how to match the colors).

PS : Ideally, I would like to automatically create printable documents following a defined layout. This layout is made of
 - a title and inset at the bottom, the inset made of the global geotiff with a red rectangle for the zone, with and north arrow and scale-bar,
 - the detailed zoned with the zoom defined such that the clipping polygon fits the central area.
For each of the 50 zones. I've added clipping polygons for each zone raster, and I would like to automatically create the documents with the big geotiff in an inset and the zone zoomed in. Is it easy to do using ArcGIS ? Again, I could not apply the solutions I've found.


texts:
We are currently working on two cities in the South of Kyrgyzstan (Aydarken and Kadamjay).

They are visible on the map screenshot below (purple dot). For reference, the grid is Lat/Long WGS84 (DD).

Now we received a set of coordinates with an uncommon syntax unknown to us, see screenshot below.

And to make things tastier, EPSG and units used are obviously unknown.

We are currently digging for more details/metadata, but for the moment we simply know they come from a governmental statistics office in Kyrgyzstan.

This syntax is unknown to me. I tried a "brute force" conversion to DMS, DDM and DD, simply by splitting the coordinates and setting them as WGS84.

As you can see, they end up in Kazakhstan (being, in this case, 44°N) - while we know they are supposed to be in southern Kyrgyzstan.

Therefore, there are hypothesises:

Error in source data (wrong coordinates)
The wrong method, these are coordinates in a specific format/unit/EPSG, which I need to identify

Questions:

What means can I use to identify unknown units format/syntax and deduce the right EPSG? I have a look at EPSG.io but the dozen of CRS identified for Kyrgyzstan do not match the syntax and I am lost in the thousands of others.
Has anybody seen similar coordinates before? Data looks quite old, so could be an old Kyrgyzstan or russian CRS.


texts: Referencing file paths stored as attributes in ModelBuilder?

texts:
I'm trying to design a model in ModelBuilder that will derive input files for a process from file paths stored in an index shapefile's attribute field.

The two datasets involved are a shapefile selection grid and a shapefile index of .LAS file footprints. The model will iterate through each grid square and select the index LAS footprints that intersect that square.

Within the index shapefile is an attribute field that has the file path of the corresponding LAS file (i.e. B:\las_files\100.las).  Each grid square will intersect multiple LAS file footprints, so I want the model to collect those LAS files and use them as an input for the Create LAS Dataset tool. The output of the model will be a LAS Dataset for each grid square.

How can I relate the file paths stored in the attribute field to using those external files as model inputs?