texts:
Anyone knows the the cause and how to resolve the "Could not obtain native oracle connection" error?

I'm using GeoTools 2.7.4, a oracle 11g DB, ojdbc14.jar, and a simple java code on a Tomcat 6 like this:

Map map = new HashMap();
map.put( "dbtype", "oracle");
map.put( "jndiReferenceName", "java:comp/env/jdbc/GisDS");

OracleNGJNDIDataStoreFactory factory = new OracleNGJNDIDataStoreFactory();
log.info("Can process? "+ factory.canProcess(map));
store = factory.createDataStore(map);

SimpleFeatureSource featureSource = store.getFeatureSource(viewName);
log.info(featureSource.getFeatures().getBounds().toString());

And my context.xml:

<Resource name="jdbc/GisDS" auth="Container" type="javax.sql.DataSource"
    driverClassName="oracle.jdbc.driver.OracleDriver"
    accessToUnderlyingConnectionAllowed="true"
    username="***" password="***"
    url="***"
    maxActive="25" maxIdle="8" maxWait="10000"
    removeAbandoned="true" removeAbandonedTimeout="60" logAbandoned="true"
    poolPreparedStatements="true" maxOpenPreparedStatements="100"
    validationQuery="SELECT SYSDATE FROM DUAL"/>

And I always get the "Could not obtain native oracle connection" error:

java.lang.RuntimeException: java.io.IOException: Error occured calculating bounds
        at org.geotools.jdbc.JDBCFeatureSource.getBoundsInternal(JDBCFeatureSource.java:485)
        at org.geotools.jdbc.JDBCFeatureStore.getBoundsInternal(JDBCFeatureStore.java:179)
        at org.geotools.data.store.ContentFeatureSource.getBounds(ContentFeatureSource.java:370)
        at org.geotools.data.store.ContentFeatureCollection.getBounds(ContentFeatureCollection.java:274)
        at ***.GisServlet.doGet(GisServlet.java:82)
        at javax.servlet.http.HttpServlet.service(HttpServlet.java:617)
        at javax.servlet.http.HttpServlet.service(HttpServlet.java:717)
        at org.apache.catalina.core.ApplicationFilterChain.internalDoFilter(ApplicationFilterChain.java:290)
        at org.apache.catalina.core.ApplicationFilterChain.doFilter(ApplicationFilterChain.java:206)
        at org.apache.catalina.core.StandardWrapperValve.invoke(StandardWrapperValve.java:233)
        at org.apache.catalina.core.StandardContextValve.invoke(StandardContextValve.java:191)
        at org.apache.catalina.core.StandardHostValve.invoke(StandardHostValve.java:127)
        at org.apache.catalina.valves.ErrorReportValve.invoke(ErrorReportValve.java:102)
        at org.apache.catalina.core.StandardEngineValve.invoke(StandardEngineValve.java:109)
        at org.apache.catalina.connector.CoyoteAdapter.service(CoyoteAdapter.java:293)
        at org.apache.coyote.http11.Http11Processor.process(Http11Processor.java:859)
        at org.apache.coyote.http11.Http11Protocol$Http11ConnectionHandler.process(Http11Protocol.java:602)
        at org.apache.tomcat.util.net.JIoEndpoint$Worker.run(JIoEndpoint.java:489)
        at java.lang.Thread.run(Thread.java:662)
Caused by: java.io.IOException: Error occured calculating bounds
        at org.geotools.jdbc.JDBCDataStore.getBounds(JDBCDataStore.java:1137)
        at org.geotools.jdbc.JDBCFeatureSource.getBoundsInternal(JDBCFeatureSource.java:478)
        ... 19 more
Caused by: java.sql.SQLException: Could not obtain native oracle connection for class org.apache.tomcat.dbcp.dbcp.PoolingDataSource$PoolGuardConnectionWrapper
        at org.geotools.data.oracle.OracleDialect.unwrapConnection(OracleDialect.java:490)
        at org.geotools.data.oracle.OracleDialect.readGeometry(OracleDialect.java:434)
        at org.geotools.data.oracle.OracleDialect.readGeometry(OracleDialect.java:424)
        at org.geotools.data.oracle.OracleDialect.decodeGeometryEnvelope(OracleDialect.java:371)
        at org.geotools.jdbc.JDBCDataStore.getBounds(JDBCDataStore.java:1123)
        ... 20 more

I've tried to connect to a 10g DB, I've tried using the ojdbc6, and I always end with the same result.

The strange is, If I replace the getBounds() with size() it works.


texts: Zebra grid is not displayed in map composer

texts:
in Qgis 2.6, in Map composer, while selecting grid frame style as Zebra, it does not works. Also while drawing coordinates points as exterior in map convas, it is not displayed.
Actually I am facing this problem in Mac OS. while composing a Map, the Grid and coordinates are not displayed, but when I export/ save as the map canvas in as JPEG/PDF, the Grid and coordinates is displayed in output.


texts:
I try to parse a KML file to retrieve the polygon names and -coordinates. I try to do this with fastkml (also tried pykml and element.tree, but no result). Reading the KML seems to be good. However, then creating the string to read (and parse) it, seems to create a problem. I have this code:
from fastkml import kml

with open('Test.kml', 'rt', encoding="utf-8") as myfile:
    doc = myfile.read()

k = kml.KML()
k.from_string(doc)

placemarks = list(innerFeature[0].features())

However, I get this feedback from the terminal:
Traceback (most recent call last):
  File "parse_kml.py", line 17, in <module>
    k.from_string(doc)
  File "/Users/*****/opt/anaconda3/lib/python3.8/site-packages/fastkml/kml.py", line 89, in from_string
    element = etree.fromstring(
  File "src/lxml/etree.pyx", line 3237, in lxml.etree.fromstring
  File "src/lxml/parser.pxi", line 1891, in lxml.etree._parseMemoryDocument
ValueError: Unicode strings with encoding declaration are not supported. Please use bytes input or XML fragments without declaration.

I installed the fastkml library with pip. In the past I also installed some libraries with Anaconda (e.g. geopandas) which created some errors during installation. I am now doubting if this error message is caused by a corrupted library system, or if the problem is in my coding.


texts: Selecting only mixed attributes (house numbers blended with other characters and symbols) in ArcGIS Desktop?

texts:
In ArcGIS, I have an attribute table with house numbers, i.e. numeric and mixed.

I want to select only mixed house numbers.

Mixed means house numbers blended with other characters and symbols, e.g. 65g, 0665, 75-a, a-b, 44-65 etc.


texts: PostGis's ST_MapAlgebra rasters do not have the same alignment

texts:
I'm having a problem with ST_MapAlgebra.

An error stating The set of rasters provided (custom extent included, if appropriate) do not have the same alignment is shown.

I've been trying to know where the error comes from and nothing seems to fix it.

My SELECT is this one:

SELECT ST_MapAlgebra(rast,
        ST_AddBand(
            st_setsrid(
            ST_MakeEmptyRaster(
                 (430276.740772 - 429856.518926)::numeric::integer,
                 (4581915.18987 - 4581587.0391)::numeric::integer,
                  429856.518926, 4581915.18987, 2),
            3857),
            '8BUI'::text, 1, 0),
        '[rast1]', '8BUI', 'SECOND')
FROM bigraster;

SELECT ST_SetBandNoDataValue(rast, 1) FROM bigraster; is '8BUI'
select st_pixelwidth(rast) FROM bigraster is 2

What else could be wrong?


texts: Custom QGIS maptool to return current event

texts:
This is what I am trying to do.
I designed an interface which can switch between selection mode and placement mode.

Selection mode activates the selection tool. Placement mode is what I am struggling atm. when the user selects "placement mode" It activates Custom Maptool which I designed. The maptool will ask the user to click on map canvas. And will wait until something is clicked or Esc button clicked.(I dont know how to force the code to wait for the result from click).
The maptool should be disabled after one click. and return the value of click point to the method which activates maptool.

here is my code for maptool which needs to return clicked point(QgsMapMouseEvent object). and deactivate itself after one click.

class EventHandler(QgsMapTool):

def __init__(self):
    canvas = iface.mapCanvas()
    self.canvas = canvas
    QgsMapTool.__init__(self, canvas)
    self.setCursor(Qt.CrossCursor)
    #self.current_event=""#this  holdes an instance of current event
    self.button_clicked=""

def canvasReleaseEvent(self, event):
    """This method passes the instance of current event"""
    global current_event
    current_event=event
    self.canvas.unsetMapTool(EventHandler())
    Return event #here I want to return value

def keyPressEvent(self,event):
   print event

Here is the method which initiates CustomMapTool

def placement_mode(self):
    canvas = iface.mapCanvas()
    result=canvas.setMapTool(EventHandler())
    #should wait here till user click on map canvas
    print result

my problems are how to return values to the placement_mode method. and how to wait until the result is available


texts:
Is there a way to cut a polyon in equal parts?

Basically, I have several rectangles that I need to divide into equal parts, either 3, 4 or 5 equal parts.

I am using ArcGIS 10.3.1


texts:
When I create a map in print composer the map looks exactly how I want it, but when I save it as a PDF some of the sympbols have changed. Note the roads are still visible in the PDF but some have changed to be a light grey, from the black that they are suppose to be? Why is it doing that and how do I fix it? If I unselect print as a raster, then the colours of the map change all together.
Map as seen in print composer.

PDF as raster.


texts:
I have shapefile with next attribute table :

(name of features, unique number for each feature of layer)

I must split shapefile by NAME field, but in result in attribute tables of layers should be only unique number column (see picture):

Split layer by Name I create in QGIS. Is there any solution to delete column "Name" in many shapefiles (about 800)?

I use ArcGIS, QGIS, PostGIS, OGR...


texts: How to perform a non-linear regression pixel by pixel?

texts:
Using nlsLM from Package ‘minpack.lm’ is straightforward for a simple example like this:

MODEL:y=(exp(a*x+b*z+c)+d)^f
x=c(0.5,0.3,0.2,0.4)
z=c(0.1,0.6,1,0.9)
y=c(0.2,0.3,0.9,0.9)
fit=nlsLM(y~(exp(a*x+b*z+c)+d)^f,start = list(a = -.03, b = 0.5, c = 1,d=0.02,f=0.003))

And no problem with that. My real data are multiple rasters for one year.So I have one year of data for x,y,z but as rasters that contain pixels,I want to do a non-linear regression pixel by pixel(each regression for one pixel is not related to the next pixel). example data:

   r <- raster(nrows=10, ncols=10); r <- setValues(r, 1:ncell(r))
   r1 <- raster(nrows=10, ncols=10);r1 <- setValues(r1, 1:ncell(r))
   r2 <- raster(nrows=10, ncols=10);r2 <- setValues(r1, 1:ncell(r))
st1=stack(r,r1,r2)

   rl <- raster(nrows=10, ncols=10); r <- setValues(r, 1:ncell(r))
   rS <- raster(nrows=10, ncols=10);r1 <- setValues(r1, 1:ncell(r))
   rT <- raster(nrows=10, ncols=10);r2 <- setValues(r1, 1:ncell(r))
st2=stack(rl,rS,rT)

  re <- raster(nrows=10, ncols=10); r <- setValues(r, 1:ncell(r))
  ru <- raster(nrows=10, ncols=10);r1 <- setValues(r1, 1:ncell(r))
  rg <- raster(nrows=10, ncols=10);r2 <- setValues(r1, 1:ncell(r))
st3=stack(re,ru,rg)

so first x would be the first pixel of st1(for the whole temporal period)

so first y would be the first pixel of st2(for the whole temporal period)

so first z would be the first pixel of st3(for the whole temporal period)

do regression for this pixel then do the same for all other pixels

As a result I get a map (raster) 10*10 for each best parameter(a,b,c,d,f)

I hope it is clear and you have any help!


texts: How to add a progress bar to a python plugin?

texts: I would like to add some type of progress bar to my Qgis python plugin.  Can someone point me to an example?


texts: Correlation of yield data within a buffer zone

texts:
Working with QGIS 2.12.1-Lyon, I am trying to compare one point (single shapefile) with multiple others that fall within a buffer zone of 10m i.e. comparing one shapefile with another.

The end product I am trying to calculate is a correlation between the one point and; the nearest point, the average within a 10m radius and the average within a 5 meter row.

The data I am working with is yield data from a wheat field. The correlation is with 1 biomass sample and multiple yield data points.


texts: Mask raster with shapefile in Python using GDAL/OGR

texts:
I'm building a script in Python to mask out land and keep the ocean. As seen in the following I have a shapefile for the coastlines (one feature only).

I want to mask out the the land to black or null in Python using GDAL/OGR.

I'm fairly new in using the gdal librairies and I've tried the following without luck:

    inraster = 'raster.tif'
    inshape = 'shape.shp'
    subprocess.call(['gdalwarp', inraster, outraster, '-cutline', inshape,
                 '-crop_to_cutline', '-cwhere'])


texts: Is it possible to integrate the mapping capabilities of QGIS with Stata?

texts:
I was wondering  of the possibility of integrating the mapping capabilities of Q-GIS with Stata, by using the programming aspect of Stata as opposed to using R or python.

Is this possible?


texts: Identifying FRONT LAND using ArcGIS ModelBuilder

texts:
I work at a town hall.
I need to identify the front of the land. Here in my country we call it of "TESTADA" or "FRENTE DE LOTE". I did not find the exact term in English.

THE FRONT OF A LAND is a line near the entrance to the ground. That is, is the entry of a house. It is where the entrance gate, with the port number.
To exemplify following image below.

Map before:

Map later with the front lines of terrain - As should be:

A LOT has only one FRONT OF LAND. In corners it's main street.
There are old buildings that will stay as they are because it is historic.
Currently we do these lines of FRONT LAND manually.
We would like to automate part of it. That could be using ModelBuilder.

In short:

In: LAYER OF LAND, LAYER STREET, BUILDING LAYER.
Output: LAYER with FRONTS OF LAND LINES.

Another example of front lot:
The term FRONT OF LAND refers to front of property.The input of the house is his front.


texts: Restricting what data gets uploaded to QGISCloud

texts:
I have set up a free account on QGISCloud and installed the Cloud plugin but I am having problems restricting the data that gets uploaded.

I have four (plus a background) layers open in QGIS but when I try and upload data using the plugin it lists a fifth layer which I had deleted from the project.  Unfortunately that layer is full of "non public" information.  I first tried saving the project and restarting QGIS in case it was some how lingering in memory but I had the same result.

Any clues as to what is going on here?

The selected layer (wp_master) is not referenced by this project and, so far as I can tell should not be uploaded.

Some of the data in wp_master is potentially sensitive and I definitely do not have authority to upload this data anywhere remotely public.

QGIS 2.18.9 running on a Mac ( High Sierra )


texts: How to use joinFieldNamesSubset from QgsVectorJoinInfo class

texts:
I have written a small python plugin to refresh a join on a point layer in QGIS.

point_layer=None
for lyr1 in QgsMapLayerRegistry.instance().mapLayers().values():
  if lyr1.name() == "point_layer":
    point_layer = lyr1

text_table=None
for lyr2 in QgsMapLayerRegistry.instance().mapLayers().values():
if lyr2.name() == "text_table":
  text_table = lyr2

POINT_ID='ID'
TEXT_ID='ID'
joinObject = QgsVectorJoinInfo()
joinObject.joinLayerId = text_table.id()
QgsVectorLayer.removeJoin(point_layer, text_table.id())
joinObject.joinFieldName = TEXT_ID
joinObject.targetFieldName = POINT_ID
joinObject.memoryCache = True
point_layer.addJoin(joinObject)

I tried to get a subset of fields to be used from joined layer text_layer.

joinObject.joinFieldNamesSubset = 'TYPE'

There is no error, but it joins all the fields, not only the field 'TYPE' (attribute of layer text_layer).

How can I set the attributes I want to join? And where do I have to call joinFieldNamesSubset?


texts: Missing 0 values of a Rasterlayer prevent correct calculation with Zonal Statistics tool in Qgis

texts:
For my master thesis I want to calculate the mean percentage of sealed surface in some parks of Berlin.

For this purpose I have a rasterlayer containing this data. In a previous step I already created the polygon shapefiles of the parks.

For the calculation of the mean value I used the Zonal Statistics Tool in Qgis.
After checking the results though, it became clear to me that something must have gone wrong. Even small park with very little sealed surface showed up with values of 70 %.

I'm pretty sure the reason for this is that somehow the rasterlayer does not correctly show 0 value pixels. Since only pixels different from 0 was used for the calculation, which could also be seen in the "count" value. I added a screenshot to make this point clear. In the picture you can see that the tool only was calculating with the 8 values of the pixels with values bigger than 0 and completely neglecting the big area with 0 value, resulting in a mean value of 72.25 % which is much higher than the true value.

One step I took to solve this problem was to remove the tick for the "no data value: 0" in the properties section of the raster. In the metadata it still shows "no data value" as 0, though and the results of the calculations did not change. According to the colleague from whom I got the layer, it should contain 0 values, though.


texts:
How can I display frequently traveled routes?

I think this has been asked before in a question I am unable to relocate which was about how to display multiple overlapping routes on a map, so that frequently traveled routes (with many overlaps) are bright, and single routes are dim. I believe they were trying to re-create an effect from a map they had found elsewhere, and there was an image of that map included in the question.

The answer used transparency and blending modes to create a beautiful effect, where overlapping lines appeared to glow.

The solution looked something like this, except that the lines followed roads in the basemap. They may also have been categorized with a red-yellow-green color scheme.

The basemap in the solution (and possibly also in the question) was a dark grayscale map, similar to ESRI Gray (dark) or CartoDB Dark Matter.


texts:
Possible Duplicate:
  How to delete fields in QGIS?

How to delete a column of a shapefile in QGIS?
The "Delete Column" tool is not highlighted even in editable mode.


texts:
I'm trying to automatically measure the distance between 500 points in order of attribute value. So the points have the value 1 to 500 and in that order I want to measure the individual distance of all points from its predecessor.

Is there a tool or plugin to execute this analysis in QGIS?


texts: What is the projection of an Android Location class instance?

texts:
What is the projection of a Location object returned from the LocationClient by a class implementing the LocationListener interface with a fused provider in the Google Play Services SDK? For example, if I performed the following,

1)

Location mLocation = mLocationClient.getLastLocation();
double mLongitude  = mLocation.getLongitude();
double mLatitude   = mLocation.getLatitude();

2) converted those values to String objects,

3) performed I/O to a file/network/database,

4) converted or outputted them to a .csv file of the form,

Longitude,Latitude
-85.26418009999999,35.1050902

5) and then finally imported this file into QGIS, what projection would I choose?


texts: How to set the style of a layer, in OpenLayers 3, without auto refreshing the layer?

texts:
I want to set the style on run time, but I don't want the layer to auto refresh.

I want to be able to perform some changes before I manually refresh the layer.

I looked into the setStyle method and I wrote my own customSetStyle method:

// Exactly like the 'setStyle' method from OL3
customSetStyle = function(style) {
    this.style_ = style !== undefined ? style : ol.style.defaultStyleFunction;
    this.styleFunction_ = style === null ?
      undefined : ol.style.createStyleFunction(this.style_);
    //this.changed();    // I commented this line out
};

By commenting out the this.changed(); line, the layer is not auto refreshed anymore.

I want to know if there is a better way to prevent the layer from auto refreshing AND if the way I found has nasty consequences.


texts: Setting ArcGIS Desktop license to lower level for testing whether Python script still runs?

texts:
I have a Python script that will need to be given to a client.

We are currently evaluating what license level the client will need in order to run the script for their software acquisition purposes.

We all have ArcGIS 10.1 Advanced license in our company.

I checked all the tools and extensions used in the script against the web help and found that they should work on the Standard license level.

Is there a way I can test the script against the Standard license level, when we have no Standard license level in or company (only the higher Advanced level)?


texts: Conditional Statements in QGIS raster calculator?

texts:
What is the equivalent statement of the following expression in the QGIS raster calculator?

    Con(((Abs("DEM" - "FSMin") > 120) | (Abs("DEM" - "FSMax") > 120)), 1, 0)

I have no problem running this command in ArcMap but as soon as I paste it in QGIS the raster calculator there rejects it.


texts: Extract data to file after performing kriging interpolation in ArcGIS

texts:
For 3 meteorological stations, I created maps by Kriging Interpolation tool to extract data for other places. How is it possible to extract a data file (Excel or text) from this interpolation?

(this work is done by ArcGIS9 and Arcmap version 9.3)


texts: Is it possible to read/write macros using PyQGIS?

texts:
Macros are saved in the .qgs project file. However, looking at the QgsProject Class, there is no mention of a method to read or write a macro.

Other than dangerously editing the .qgs project file, is there some other method of accessing a macro?


texts: Can you add the "measure" tool to model builder?

texts: I am wondering if there is a way to add the measure tool, usually residing in the Tools toolbar, into model builder?


texts: How do you create a text-buffer in QGIS' print composer?

texts: When creating a label you can create a buffer of another color around your text. Is there any way to do so in the print composer?


texts: What is the definition of pt_spacing on Point File Information tool?

texts:
I can't find any definition for the outputs of Point File Information tool (3D analyst), specifically the pt_spacing field.

Is this the average distance between points, or the average number of points per square foot?


texts: Where to find Vietnam administrative level 4 shapefile?

texts:
I am new to mapping, and would like to make a custom map using Vietnam's administrative boundaries shape file.  I found a great source at gadm.org, but it is licensed for non-commercial use only.

I am not currently planning to do anything commercial, but didn't know if there might be another source that is not as restrictive just so that I don't have to worry down the line?  I have looked around and it seems really hard to find.

I looked at OpenStreetMap (OSM), but I think they don't get as specific as the gadm.org files, OSM seems to only get down to rural districts and not the commune level.

I might be wrong about OSM though, as I haven't really delved into it.


texts: Plot the longitude and latitudes on map using geotools

texts:
I am new to geotools and i find bit difficult to plot a point using the latitude and longitude values. can any one help me out to solve my problem.
I tried the following code but not helpful.By quickstart i can able to pass a shapefile and getting the map but unable to plot the lat long values on it.

   SimpleFeatureTypeBuilder b = new SimpleFeatureTypeBuilder();

         b.setName( "MyFeatureType" );
          b.setCRS( DefaultGeographicCRS.WGS84 );
          b.add( "location", Point.class );
          //building  the type
          final SimpleFeatureType TYPE = b.buildFeatureType();

        SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE);
        GeometryFactory geometryFactory = JTSFactoryFinder.getGeometryFactory();
        Point point = geometryFactory.createPoint(new Coordinate(longitude, latitude));
          featureBuilder.add(point);
          SimpleFeature feature = featureBuilder.buildFeature( "fid.1" );
          DefaultFeatureCollection featureCollection = new DefaultFeatureCollection("internal",TYPE);
          featureCollection.add(feature);

Adding to map:

List a=new ArrayList();
a.add(-155.2792);
a.add(67.3623);
double longitude=(Double) a.get(0);
double latitude=(Double)a.get(0);
File file = JFileDataStoreChooser.showOpenFile("shp", null);
if (file == null) { return; }
FileDataStore store = FileDataStoreFinder.getDataStore(file);
SimpleFeatureSource featureSource = store.getFeatureSource();
MapContext map = new DefaultMapContext();
map.setTitle("Show Map");

Layer layer = new FeatureLayer(featureSource, null);
        map.addLayer(layer);
    map.addLayer(featureSource, null);
    JMapFrame.showMap(map);


texts: US Census ZCTA shape attribute table: shapeids with suffixes

texts:
What does it mean (i.e. what situation is being addressed) when the attribute table for the ZCTA shapes contains two shapeids for a given GEOID10 (e.g. 10307, Staten Island) each having a numerical suffix, but there is no shapeid for the GEOID10 without a suffix:

10696 -- does not exist
10696.0
10696.1

The data in the rows are identical.

shapeid ZCTA5CE10   GEOID10 CLASSFP10   MTFCC10 FUNCSTAT10  ALAND10 AWATER10    INTPTLAT10  INTPTLON10
10696.0 10307   10307   B5  G6350   S   5085103 17221   40.5091830000   -74.2377853000
10696.1 10307   10307   B5  G6350   S   5085103 17221   40.5091830000   -74.2377853000


texts: Intersect two polygons but preserving the shape of one

texts:
I have two polygons layers; a selection of land titles and a buffer layer. I want to intersect the two layers but preserve the shape of the land title polygon to show all the land plots that touch the buffer layer.
The buffer is the green hatched area and the land titles is the layer below.


texts: How do I create a map view with an insert in QGIS?

texts: I would like to create a map view with an insert in QGIS to show a zoomed view of my project area and its location on the National map of Zimbabwe in one view. I can not seem to find a tutorial on how to accomplish this. May you kindly assist on how I can create such a map. [![The kind of map view that I would like to produce looks like the view on the attached image.


texts:
I work for a small non profit conservation organization  and we have a reasonably high turnover of people... so it stands to reason that some people don't want to change from a preferred program.
this in mind...

Can i install ArcGIS 10.1 and QGIS on the same PC without problems?


texts: Assigning attribute based on selection points on a graph?

texts: I have a bunch of points on a map on a river that branches. I want an easy way to select one branch from the other and was hoping I could just somehow "grab" those points on the map and add a column to the attribute table with those points that would delineate the two branches.


texts: QGIS Modeler: Layer Rename does not output

texts:
I am currently in the creating a process in my QGIS model builder where I am trying to rename a finished process layer based on some best practice rules. Here is the model:

Now here is the expression I am using to rename my layer:
"algplanid" || '_' ||
CASE
WHEN to_string(@benoemingcontouren) ='0'
THEN
'CT'
WHEN to_string(@benoemingcontouren) ='1'
THEN
'ONT_CT'
WHEN to_string(@benoemingcontouren) ='2'
THEN
'ROO_CT'
ELSE
null
End

"algplanid" only has one Value to speak off in its column so it should not matter much.
so in the End it should look like this:
RUP_xxxxx_xxx_xxxxx_xxxxx_CT (xxx are numbers but are censored for privacy)
Now the problem:
Why does my rename layer algorithm do not seem to output, there is not even a temporary output function, what is going on here? And is there even another way to get this done?


texts: Collecting training data for land cover classification

texts:
I want to begin mapping various land cover such as coconut and pine. My issue is I do not have much experience in terms of collecting land use information or analyzing satellite imagery and analysis.
What is the preferred method to collect training data suitable for land cover classification of satellite imagery?


texts: Speed Up WMS of Mapserver in Openlayers 3

texts:
I have a problem with slow rendering of WMS in my map server. It takes a minute when zooming and panning of map.

here is my code to add WMS layer to openlayers. Is there a way to speed up the process ?

   layerSources[i] = new ol.source.TileWMS({
                            url: mapAgenURL,
                            params: {
                                'LAYERS': wms,
                                TRANSPARENT: true,
                                'TILED': true,
                            },
                            serverType: 'mapserver',
                            projection: 'EPSG:32651',
                            crossOrigin: 'anonymous'
                        });

                        layers[i] = new ol.layer.Tile({
                            source: layerSources[i],
                            visible: init,
                            preload: Infinity
                        });


texts: Are folders always part of ArcGIS Server REST API service name?

texts:
I'm using an ArcGIS REST map service, and my question is how do I construct the URL to the service?

Specifically, what relationship is there between the Service Name, the Test Server, and the folder directory structure?

This is the relevant information to the service I am trying to use:

  Service Name: Louisville/LOJIC_LandRecords_Louisville
  Folder: Louisville Test Server url:
  http://sampleserver1.arcgisonline.com/ArcGIS/rest/services

The server send me for the folders:

{"currentVersion" : 100,
  "folders" : [
    ...
    "Fire",
    ...
  ],
  "services" : [
    {"name" : "Geometry", "type" : "GeometryServer"}
  ]
}

and for the service:

{"currentVersion" : 100,
  "folders" : [],
  "services" : [
    {"name" : "Fire/Sheep", "type" : "MapServer"},
    {"name" : "Fire/Sheep", "type" : "FeatureServer"}
  ]
}

You can see, the folder is part of the service name, is this always so? What is, if there is more then one level of folders?


texts: WFS returns 2 records, QGIS shows only one

texts:
When I show a WFS layer (published by GeoServer) in QGIS, I'm missing one of the two records.

What I do:

Select the WFS layer
Setup a query (properties -> query-builder) (eg postalcode = '1234')
Then hit 'test' in the query builder -> "The WHERE clause returned 2 rows"
Then I look at the layer and i'm missing one of two polygons
Also, the attribute table contains only one row.

So, there is a mechanism with which QGIS merges records. In my case the records are, besides 2 columns (geo and one other), completely similar.
There is no primary key (it is a db-view in PostgreSQL) so, maybe, QGIS just takes the first column?

Is there a setting in QGIS to make sure every record is returned?


texts: Vertex Order/Initial Vertex in LineString

texts:
I would like to change the initial vertex in a LineString, I would also like to change the direction it is "drawn". Are methods available to edit this information?
The line is closed.


texts: QGIS 2.14.x missing Globe-PlugIn on Fedora 23/24

texts: QGIS 2.14.3 on Fedora 23 does not provide the Globe-Plugin. The package qgis-plugin-globe like used in *buntu is also not available. Also postinstalling osgearth did not help. Anyone a solution to this ?


texts: Create folders to group KML Placemark in a big file

texts:
We have a huge KML file having Placemarks which contain SimpleData elements.

We need to filter the elements for a value inside SimpleData and create a folder to group all the related Placemarks inside.

Is there a tool or editor which could perform it?


texts: How to center continuous raster layer by subtracting mean and dividing by standard deviation using QGIS raster calculator?

texts:
I'm trying to extrapolate a spatial model to a map in QGIS using the raster calculator. The model predicts the value of a pixel given the environmental conditions at this location (e.g. elevation, slope, ect.). Hence, each explanatory variable is a raster layer. For the model I centered my covariates by subtracting the mean and dividing by the standard deviation. Now I want to do the same to my raster layers. As trivial as this seems, I didn't manage to do so. I tried (among similar attempts - like changing upper case, etc.):

(my.raster-mean(my.raster))/sd(my.raster)

The calculator tells me the expression is valid, however, mean() and sd() don't work.

And I could not find a site with the built-in functions for the raster calculator.


texts: Connecting PostGIS raster data to GeoServer

texts:
The data I have used is DEM. It creates a table and has 2 columns in PostGIS.

Is it necessary to create a new table again and follow the configuration process as mentioned in the Image Mosaic JDBC documentation?

The error -

  INFO   | jvm 1    | 2016/01/20 19:28:24 | java.lang.RuntimeException:
  Could not list layers for this store, an error occurred retrieving
  them: Failed to create reader from file:coverages/oek.pgraster.xml.xml
  and hints null INFO   | jvm 1    | 2016/01/20 19:28:24 |  at
  org.geoserver.web.data.layer.NewLayerPageProvider.getItemsInternal(NewLayerPageProvider.java:151)
  INFO   | jvm 1    | 2016/01/20 19:28:24 |     at
  org.geoserver.web.data.layer.NewLayerPageProvider.getItems(NewLayerPageProvider.java:59)
  INFO   | jvm 1    | 2016/01/20 19:28:24 |     at
  org.geoserver.web.wicket.GeoServerDataProvider.fullSize(GeoServerDataProvider.java:242)
  INFO   | jvm 1    | 2016/01/20 19:28:24 |     at
  ......

I removed the rest of the redundant log entries.


texts: GxFilterWorkspace allows a folder that is not a IWorkspace or IGxDatabase

texts:
I am trying to open a ArcGIS file browser. And For this, I used IGxDialog interface's DoModalOpen method. To specify what the user can select using a GxDialog browser is based on the filters (GxObjectFilter) held by the object. To browse workspace, I set ObjectFilter to GxFilterWorkspaces. But it allows folders that contains layers/feature classes/raster datasets that is not IGxDatabase or IWorkspace.

        IGxObjectFilter gxObjectFilter = new GxFilterWorkspaces();
        IGxDialog gxDialog = new GxDialog();
        gxDialog.AllowMultiSelect = false;
        gxDialog.Title = title;
        gxDialog.ObjectFilter = gxObjectFilter;
        gxDialog.ButtonCaption = "Open";
        gxDialog.RememberLocation = true;

        IEnumGxObject enumGxObject;
        gxDialog.DoModalOpen(0, out enumGxObject);
        while ((gxObject = enumGxObject.Next()) != null)
        {
            if ((gxObject as IGxDatabase2) != null)
            {
                return gxObject.FullName;
            }
        }

Question:

What should I do to select only workspace through this browser? Or is there any other option to browse GIS dataset?


texts: Slow SDE Performance After Append - Geometric Network

texts:
I had an odd thing happen today, I do almost daily appends to our fiber network (fiber, splices, etc) and usually, it goes smoothly and with no complications. For some reason, I had an Append stall for hours on the cable layer (it usually takes no more than ~20 minutes, and this was no bigger than usual). I had no choice but to cancel it as it was tying up the resources during work hours.
Ever since that append has been canceled the database and layers have been so slow that it is almost unusable for people.

I was wondering if this has happened to anyone else, Nothing else has changed in the database and it was loading/selecting/editing with no issues. I had a few people look into it and noticed nothing wrong, and I did multiple "rebuild indexes" and database compressions.

I have an Enterprise ArcSDE Geodatabase on SQL Server 2016 SP2 Enterprise, we are using ArcMap 10.7.1


texts: Publishing PostGIS Data in Portal for ArcGIS Enterprise

texts:
I have a PostGIS table with millions of points (x,y,z,time).

I have to develop a web app with the ArcGIS JavaScript API, to make some statistics like this application:
https://ycabon.github.io/presentations/2019-devsummit/Working-With-Your-Data-Using-the-ArcGIS-API-for-JavaScript-Feature-Layer/demos/4_client-side_statistics/index.html

I cannot find any sharing method in ArcGIS Pro, except to copy as layerPackage.
And, as a layerPackage in Portal fort AE, I can not find any option to create a FeatureLayer to use in my app.

Which is the best method to publish PostGIS data in Portal for ArcGIS Enterprise for using with ArcGIS JavaScript API 4?


texts: I'm having an issue getting the OSM2PGSQL cgywin windows package to stay open. When I click the OSM2PGSQL.exe file after unzipping the command prompt window opens then immediately closes. I have tried opening it from command prompt to see if OSM2PGSQL would stay open but it still closes. Does anyone know how to keep it from closing.


texts:
I'm trying to publish a map service using simple xy data from a table in SDE; however, whenever I try to create the MSD, I'm told "no can do" by ArcMap.

I would like to be able to simply add xy coordinates in SQL and have them update in real-time in the map service.

I've experimented with model builder, but nothing's worked yet.

Any advice?


texts: Showing size of pie chart in legend of QGIS?

texts:
I'm using a pie chart map in which the size of the pie chart convey an attribute value (scaled size).

Is it possible to show this classification in the legend of QGIS 3?


texts: Understanding min and max values with Mosaic to New Raster?

texts:
I have to combine several DEM rasters into one from this source: http://srtm.csi.cgiar.org/SELECTION/inputCoord.asp, so I use Mosaic to New Raster tool.

I set up the tool:

I input the rasters (all from the same source, same size and no projection);
Set up the Pixel Type (32 bit floating point as in the original rasters);
Set up the Cell size as the original rasters;
Number of Bands = 1 as the original rasters; Mosaic Operator, I have done it with BLEND and MEAN (I get the same result).

The problem I have is that the resultant mosaic shows a different range of maximum and minimum values than the maximum and minimum values of the individual raster e.g., raster 1 (-5123.8, 23.25), raster 2 (-5974.6, 40.09), raster 3 (-57770.2, 38), raster 4 (-2534.3, 23.55), and final mosaic raster (-5975.8, 81.1).

I guess this solution is not right, at least I didn't expect to get that.
Anybody has an idea if this is ok, and if this is not, how to solve it and get a proper mosaic raster with right maximum and minimum values?

I am using ArcGIS 10.2.2 for Desktop.


texts: ArcGIS JS API 3 Add protected WMS services from Geoserver

texts:
How to add WMS service protected with basic authentication from Geoserver using ArcGIS Api 3? I have tried authentication with ArcGIS Proxy, something like this:
<serverUrl url="Geoserverurl/wms"
            Username="username"
            Password="password"
            matchAll="true"/>

And with
urlUtils.addProxyRule({
                urlPrefix: "Geoserverurl/wms",
                proxyUrl: "routeToAshxDocument"
            });

It doesn't work. ESRI does not provide username and password parameters for new WMSLayer class. I have also tired to write my own proxy, but it seems wms needs several requests, one for GetCapabilities and then one for reading image. Or something like that. If anyone have experience to let me know what would be the best way to do it.


texts: OpenLayers 2.13.1 WMSGetFeatureInfo not working on mobile

texts:
I'm using OpenLayers in version 2.13.1. On desktops everything is OK, alert is showing when I click on route. But on mobiles it is not working, probably touch event is not triggered.

var infoControl = new OpenLayers.Control.WMSGetFeatureInfo({
    url: '/wms',
    title: 'Identify features by clicking',
    layers: [ncn],
    infoFormat: 'application/json',
    formatOptions: new OpenLayers.Format.JSON,
    queryVisible: true
});

infoControl.events.register("getfeatureinfo", this, function (data) {
    alert('abcd');
});

var map = new OpenLayers.Map(
    {
        div: div,
        projection: "EPSG:900913", // WGS84
        featureEvents: true,
        layers: [mapboxLayer[1], mapboxLayer[0], ncn],
        controls: [infoControl],
        resolutions: [896, 448, 224, 112, 56, 28, 14, 7, 3.5, 1.75, 0.875],
        restrictedExtent: mapExtents,
        fractionalZoom: g_environment.isMobile(),
        isValidZoomLevel: function (zoomLevel) {
            return !isNaN(zoomLevel) && !(zoomLevel < MapViewOSGreenSpaceMap.prototype.MIN_ZOOM_LEVEL);
        }
    }
);

I have no idea what I should do... I never used OpenLayers and I have task to fix this.

As I see getfeatureinfo event should be triggered after getting response from WMS. But on mobiles no request to WMS is sent.


texts: I have a polygon layer which is listed at the top of the layers panel, however it is showing underneath other layers on the map canvas.  It is in the right location but despite being at the top of the layers list is underneath everything on the canvas and I cant figure out why.  What should I do?


texts: Joined attributes not recognized in mean coordinates

texts:
I am trying to find mean coordinates of a shapefile to which I have joined a CSV file for population attributes, but these will not show.

I'm working with county data of each US state based on the year 1860. I have the shapefile for the states which has polygons as counties. It has its own set of attributes but not that of population in the year 1860.

I imported the CSV file and joined the data to the original shapefile of the state. This all works fine.

The problem comes when I try to use the mean coordinate function to make a weighted centroid it does not recognize the joined attributes.

What am I doing wrong?


texts: Split Up Big Polygon Into Smaller Ones By Condition

texts:
I have some big Polygon. Within it collection of leaflet layers (points). Every point have some numeric property. What I want is split up big Polygon into smaller ones.
Every smaller polygon should contain points with ~equal (+-200 ok) sum of points properties. On the left side of my example page I added image of desirable result.
Here is my simplified example with sufficient code and comments.

So my first step is find some start point within big Polygon. It's should be point near Polygon edge - the northernmost point for example.

var nothernmostPoint= 0;
var nothernmostLayer= 0;
L.geoJSON(features, {
    pointToLayer: function (feature) {
        return L.circleMarker(feature.geometry.coordinates.reverse(), defaultPointStyle);
    },
    onEachFeature: function (feature, layer) {
        if (feature.geometry.coordinates[0] > nothernmostPoint) {
            nothernmostPoint = feature.geometry.coordinates[0];
            nothernmostLayer = feature;
        }
    }
}).addTo(map);

Second step is to find next nearest point(s) to my start point.

var geoJ = L.GeometryUtil.nClosestLayers(map, features,
nothernmostLayer.geometry.coordinates, 5);

Then sum up their properties. If sum is smaller than needed I go to step 2 and repeat, if sum is satisfy condition then I draw polygon with chosen points inside it and find nearest point to my last point and repeat search for points for the next smaller polygon.

My current difficulty is to find nearest point(s) to my start point. For that purpose I use GeometryUtil Leaflet plugin. Red dots are points that GeometryUtil found as nearest to my start point (the green one). It's certainly not what I was expected.

What am I doing wrong?

Maybe I should use different algorithm and/or tool for that task.

If it may help - all that data I store in PostgreSQL with PostGIS extension. Maybe this can be done on database level.


texts:
I am seeking advice on working with world data sets where the GCS has a positive range (0 degrees to 360 degrees). The data I am working with are NetCDF oceanographic data and have positive coordinate values as mentioned. When displayed in a normal GCS WGS84 in ESRI's ArcGIS 10, it is offset from other data that exist in the range of -180 to 180. If I reproject it, gaps occur around the prime meridian (mostly due to its origins, sometimes as far off as 25 degrees W). My thought is to create a custom projection/coordinate system that has the positive values and hope that other world data sets reproject on the fly just fine. Does anyone see any issues with this or have other solutions? Hope this makes sense. I appreciate any help in advance.

** Update **
This is a screenshot of a countries data set with graticule on top to show the normal CRS of GCS. Just importing the NetCDF data as-is brings it in but the values don't start until 20.5 degrees east of the Prime Meridian.

Normal GCS http://grafa.co/rnd/img/GCS_Normal.png

If I choose to use the NetCDF's CRS (which is actually the same) it will reproject the world data on the fly if it is defined. Note the graticule is undefined so it does not reproject.
Sort of reprojected GCS http://grafa.co/rnd/img/GCS_Other.png

But then if I try to reproject everything to a normal GCS with negative values, it's like wrapping a flat map of the data around the world and then it disappears when it hits the Prime Meridian.
Sort of reprojected GCS http://grafa.co/rnd/img/reprojected.png

Now, I know there are no values from 0-20.5 as is stated in the metadata. But why can't the values from rest of the data display in the gap? I even tried a Shift in the Raster Tools to no avail.


texts: create route network by connecting points together

texts:
Below is a screendump of my .csv file. It contains 2 pairs of X,Y coordinate columns.

The dataset is about transport of chemicals or dirt by trucks from one destination to the other (one XY coordinate to the other).

I also have a network dataset containing all roads of the Netherlands, including maximum speed per road segment.

I want to map which roads were used during transportation.

I have to create a line-layer of shortest route over my network for every record.

So I want ArcMap to calculate the shortest route from start position to end position for every record in my database and make this a line shapefile.

How can I make such a line dataset based on the .csv file I have and the network dataset.

Do I do this with the network analyst?


texts: How to Add Name of Different Category in Attribute Table

texts: I'm stuck in a problem when I started calculating area in arcmap 10.1. In the picture, I categorize my polygon into four unique labels - i.e. vegetation, sand fill, water body, built-up area. Now I want to calculate the area of these four categories. But when I opened the attribute table and try to calculate the area of each category, there are many polygons having the same gridcode. Also I tried to add another field called Name but can't. Can anyone please tell me how I can do this ?


texts:
I am trying to line up CAD data with a shapefile I have. The shapefile is in NAD 83 PA South feet, and the CAD data is in red below and is in a custom projection NAD 83 Longwood Feet, but the projection is: Rectified_Skew_Orthomorphic_Natural_Origin
As you can see the CAD data comes in much larger and off center. Any ideas how to fix this?


texts: FME: No attributes exported when writing to SHAPE file

texts:
Setup:
I have some input data (SHAPE) which I use a Clipper transformer on. Depending on the result of the Clipper, I assign a value to an attribute. In the end I assign the two results/streams to a SHAPE writer.

Problem:
When I run the project my resulting SHAPE file has no attributes.

Observations:

If I open my output (by clicking on the arrow/triangle next to the OUTPUT) on the transformers before my writer, the attributes are available and the same for both the streams.
If I select Fanout By Attribute in Properties in my writer, all my attributes are available. If I select the User Attributes tab, no attributes are listed.

Through these observations, my attributes seems valid (supported by the transformers and the Fanout) and unique thus I see no reason why the writer shouldn't export them.

Solution:
Like GIS-Jonathan suggest, I have added my writer and selected the static schema. Then I right click my writer and select Copy Attributes From ... and select the desired source. Now I have all the attributes I want.


texts:
I have rasters with discrete values (Landuse categories, boolean values...). I would like to be able to style them in a discrete way:

0 -> Red
1 -> Blue
2 -> Green
...

I wish something like the "Unique Categories" in ArcGIS: http://resources.arcgis.com/en/help/main/10.1/index.html#/representing_unique_categories_such_as_land_use/009t00000074000000/

However, all the options I see involve color ramps... Have I missed something?

I've opened a feature request, since there seems to be only workarounds available, not solutions: http://hub.qgis.org/issues/14845


texts: Saving and re-applying symbology properties for classified raster in ArcMap?

texts:
I have a classified raster which has the desired symbology configured in ArcMap 10.2, and I want to apply this symbology to other similar rasters.

So I use the Save option on the Symbology tab, intending to save a *.lyr file - but the only option is to create a Raster Function Template XML file.

When I attempt to import this symbology back onto the same layer, the only option is to import a *.lyr file:

This is preventing me from applying this saved symbology onto multiple layers.

How can the symbology Save and Import functions show different options for the same layer? This sounds like a bug, but I'm hoping for some workarounds/suggestions.


texts:
does anyone know how DEMs that can be downloaded from LandViewer are generated?
They have higher resolution than SRTM and other free sources.


texts: How to find USA Water Hardness shapefiles?

texts:
I want to find USA water hardness shapefiles. I have searched for it from here: https://www.epa.gov/waterdata/waters-geospatial-data-downloads but without success.
Is there any good data source from where I can download these boundaries?


texts: Averaging polylines from a number of lines

texts:
I have a number of polylines with different length, direction, and curvature. I want to create an average or ensemble of the polylines.

What I previously did in ArcGIS 10.1 is to convert the polylines to points in a cell size (25 km) that is quite larger than the study area (10). This way, I could "average" the points. Then I converted the raster cells back into lines. In a way they were "averaged" but the lines are very straight and formed small box-like in the intersection following the cell size.

I have read some solutions using R. I am interested to be provided with suggestions to do this in R.


texts:
I have a PostGIS table of linestrings representing a road network. I want to simplify the geometries such that only lines that connect road junctions are stored.

At the moment the table has a large number of road segments that make up a single road. I want to join as many of these as possible into a single linestring. Such that:-

x----x----x----x becomes x------------x

and:-

     x                        x
     |                        |
x----x----x----x becomes x----x--------x

It is important that I retain the points in each linestring so I am not looking to simplify to a topology.


texts:
I have some geometry in cgal that I would like to export to postgis.  The kernel I'm using is the same as what sfcgal uses (exact predicate exact construction), but I'm encountering some errors in postgis when I use the standard float output for the cgal points with ST_GeomFromText.  For instance when I try to use ST_3DIntersection with a polygon I get an error that the polygon points aren't coplanar, but I have confirmed that that the points are coplanar in my cgal program.

What I'm hoping for is a way to construct the postgis geometry using the exact fraction output from point().exact() in cgal, e.g. ST_ExactGeomFromText('PointZ(0, 0, 8676011023327629/4611686018427387904)').  I haven't been able to find a function like this in the docs though.  If there isn't such a function does anyone have some suggestions on how I can get valid geometry from cgal to postgis?

EDIT:

I've now tried building postgis and sfcgal from source with SFCGAL_CHECK_VALIDITY=OFF as well as a couple other things, but nothing seems to be working. Here's a reproduction:

You'll notice IsPlanar is true, but ST_3DArea raises an error saying the points aren't coplanar.

CREATE TABLE polygon_test();
SELECT AddGeometryColumn('polygon_test', 'geom', 4326, 'POLYGONZ', 3);
INSERT INTO polygon_test (geom) VALUES(ST_GeomFromText('POLYGON Z ((-2456.65849662024 -494.251393883276 412,-2486.7786660429 -413.696999618586 412,-2486.77945264572 -413.696623334753 354,-2456.65889152129 -494.252065179804 354,-2456.65849662024 -494.251393883276 412))', 4326));
SELECT ST_IsPlanar(geom) FROM polygon_test;

SELECT ST_3DArea(geom) FROM polygon_test;


texts: Seeking QGIS tutorials and web resources?

texts:
I wish to familiarize myself with the capabilities of QGIS up to and beyond QGIS version 3.0.

Can you recommend some good web resources and tutorials that would show me QGIS' capabilities?

I am interested in both using QGIS through its GUI and in developing applications using QGIS.


texts: Splitting Quad to 1:7,500 in ArcGIS Desktop?

texts:
Here I have 1:24,000k quads. I know how to use the word split toolbox. What I would like to split 24,000K to 1:7,500 or 6,000k. Is there a way around for me to do that ?

For example,, if you have one quad that is 1:24,000 and you want to break it down to either 1:6,000 or 1:7,500 ?

I have about 269 quads and need them to break down to where I want 1:6,000 or 1:7,500.

I am using ArcGIS Desktop 9.3.1.


texts: Use spei() function with user specified parameters from a rasterstack in R?

texts:
This is a follow up to Using spei() function on time series from rasterstack in R?

I am trying to use the spei() function from the SPEI package to calculate drought metrics. I have a rasterstack of 30 years of monthly data, so 360 layers. I want to run spei() once on historic climate data to calculate the coefficients for fitting the SPEI index. Then I want to apply those coefficients to new monthly data for recent and future climate to actually calculate SPEI (fitted) values on those data. So, @Kamo helped me in the previous Q&A to run spei() and return a rasterstack.

library(SPEI)
library(raster)

# Generate a sample raster stack time series with 360 layers
r <- raster(nrows=10,ncols=10,vals=rnorm(100))
rstack <- stack(r)
for(i in 1:359){
  rstack <- stack(rstack, raster(nrows=10,ncols=10,vals=rnorm(100)))
  cat(paste("..",round(((i+1)/360)*100,1),"%")) # check progress in %
}

# Change the function SPEI so it outputs a numeric vector of coefficients
# Modified from @Kamo
funSPEI <- function(x, scale=2, na.rm=TRUE) {
  as.numeric((spei(x, scale=scale, na.rm=na.rm))$coefficients)
}

rst.coeff <- calc(rstack, fun = funSPEI)

Now I want to use these coefficients as the param argument in spei(). According to spei documentation (https://cran.r-project.org/web/packages/SPEI/SPEI.pdf), user specified params "should be a three-dimensional array with dimensions (3,number of series in data,12), containing twelve parameter triads (xi, alpha, kappa) for each data series, one for each month".

Instead of 12 parameter triads, what I have is a raster stack of 36 layers (3 parameters times 12 months). How can I use that as input to spei()?

I have tried:

funSPEI.2 <- function (x, scale=2, na.rm=TRUE, p) {
  as.numeric((spei(x, scale=scale, na.rm=na.rm, p))$fitted)
}

rstFit <- calc(rstack, fun = funSPEI.2, params  = rst.coeff))

But I get the error:

#Error in .calcTest(x[1:5], fun, na.rm, forcefun, forceapply) : cannot use this function


texts:
I have numerous population layers for the Baltic Sea. In the attribute table of each layer, I have population densities, name of towns etc, all the same geometry (points).
I'd like to merge the layers together so that I could have all the population densities of the Baltic inside one layer. I tried to used the "merge vector layers" tool, but when I do so, the information gets all mixed up in the output attribute table. For example, the first 11 columns are linked to Germany (picture 1), but then, the next columns show only zeros and no data (picture 2)

Is there a way to make all this tidier?


texts:
I have installed the Cartegraph 8, where I would like to know if it is possible to add the MAPDirector for ArcGIS Extension to ArcMap.

As per the following documentation, it seems like the extension should be available from the Customize > Extensions > Map Director for ArcGIS. Unfortunately I was not able to get the extension displayed in ArcMap.

Cartegraph Documentation:
http://ua.cartegraph.com/ua/Content/MapSolutions/MapArc/ArcMapExtension_Task.htm

I could also find the following ESRI link, however I was not able to get the how to steps to add this extension to ArcMap.

ESRI Documentation:
http://www.esri.com/news/arcnews/spring03articles/business-partners.html

I am suspecting if this is something to do with the registration of the DLL. However when I went to the C:\ProgramFiles (x86)\Cartegraph\, I could only see a system folder, where I can see many folders such as bin, Utils, etc. However I am not sure which DLL that I would need to register to make the extension available in ArcGIS.

ArcGIS Version: 10.1
Cartegraph Version: 8
OS: Win 7 64-bit

Please suggest me on this issue.


texts: Summing attribute of feature in network along route using ArcGIS Network Analyst?

texts:
In the walking network that a friend and I are creating, we want to calculate an overall safety rating for each route. Each segment of our walking network has a certain safety value. We have calculated a "total safety" rating of each segment by multiplying the length of the segment by the safety value of the segment (1 to 5). When creating walking routes using Network Analyst, we want to be given the sum of all the "total safety" ratings for all the segments along the route. Using this, we can normalize this value by the total distance, giving us an overall safety value for the entire route.

Does anyone know how to do this?


texts: Extracting attributes of points within overlapping polygons in ArcPy?

texts:
I have a set of overlapping polygons and a line of points. A single points can fall within any number of polygons, and I need to determine the attribute values of points that fall within each polygon.

My outlined workflow was:

Select each polygon individually using SelectLayerByAttribute_management
Select the points which fall within that selected polygon using SelectLayerByLocation_management
Iterate through the selected points and extract their attribute values using SearchCursor

The issue I'm running into is that according to this post the polygons for step 1 have to be a "layer" rather than a feature class or shapefile. However for step 2, the polygons cannot be "layers" and have to be a feature class or shapefile.

I was relying on the selection of the polygon features in step 1 and 2 to to be the same. This works inside of ArcMap using the Python console, but I am unable to produce more than the first line.


texts: Loading a CSV file into QGIS without longitude or latitude data

texts:
I am trying to load a CSV table into QGIS. There are no longitude or latitude values and the options for X and Y axis look like this:

These are fields in the table I am trying to load. and this is the file:
https://data.london.gov.uk/download/land-use-ward/1f9be3a5-7ea7-4996-8875-fa060d8335e4/Land%20use%20borough.csv


texts: Downloading Met Office Datapoint tiles for specific region

texts:
Does anyone know how to use the Met Office Datapoint API to download tiled imagery?

There are some "detailed instructions" here, but it is not clear how to relate tile coordinates to real world coordinates.

For example, to download a tile of rainfall radar imagery, I can use a request like

http://datapoint.metoffice.gov.uk/public/data/inspire/view/wmts?REQUEST=gettile
&LAYER=Precipitation_Rate&FORMAT=image/png&TILEMATRIXSET=EPSG:4258
&TILEMATRIX=EPSG:4258:0&TILEROW=0&TILECOL=0&DIM_RUN=2013-11-20T03:00:00Z
&DIM_FORECAST=%2B9&key=<API key>

This is fine for zoom level 0 (as there is only one tile at this level), but what I can't figure out is which row and column numbers to use to get a specific tile.

More information about the service can be found using the getcapabilities command (you need to register to get an API key first).

http://datapoint.metoffice.gov.uk/public/data/inspire/view/wmts?REQUEST=getcapabilities
&key=<APIKEY>

It gives information for each available spatial reference and zoom level, for example this is EPSG:4258 level 9:

<TileMatrix>
    <ows:Identifier>EPSG:4258:9</ows:Identifier>
    <ScaleDenominator>128062.35453175165</ScaleDenominator>
    <TopLeftCorner>-10.67 71.05</TopLeftCorner>
    <TileWidth>256</TileWidth>
    <TileHeight>256</TileHeight>
    <MatrixWidth>512</MatrixWidth>
    <MatrixHeight>512</MatrixHeight>
</TileMatrix>

I can turn each tile into a GDAL dataset then use the heights and widths to place it in the grid. All I need to do is add a geotransform - I've got the top left coordinates, but can't work out the pixel resolution.

Presumably I can do something with the scale denominator, but I'm a but stumped.

Maybe I'm barking up the wrong tree and there is a better way to do it. Can anyone help?


texts: Why is Maplex making all my labels disappear?

texts:
I've got a shapefile with 300+ points in it, and am trying to organise the labelling.

I installed a trial version of the Maplex licence, as via some research I've come to the conclusion that this will give me many more options in terms of placement and conflict resolution, as I can't seem to get the desired effect with the standard labelling engine.

However, when I enable the maplex label engine, all my labels disappear, and I can't figure out how to get them back on? I get all the maplex options in terms of placement, but find it odd that I can't see a single label no matter what I click? I'm hoping there's a checkbox somewhere or extra step I have to perform before this will work, can anyone shed any light on this?

My labels are all showing fine using the standard labelling editor (though they all overlap, and I want to be able to fix them without converting to annotation and placing each individually, as there are so many!) But as soon as I enable Maplex, they disappear.

[EDIT: Screengrabs added:
Prior to enabling Maplex, I can see every point labelled:

Then when I hit Enable Maplex, they all disappear?

Hitting View Unplaced Labels doesn't change anything?

And even when I make sure Never Remove Labels is ticked, I still can't see any labels?

Is there really a scale limit to labelling? This is not the actual file I was trying to label initially, but the same thing is happening to every file I have?

[EDIT: On further research, I've noted a similar problem on another Stackexchange post here, with a variety of suggestions, all viable but also not working for me. Not sure if it's reassuring that it's maybe not just me, or disconcerting!]

[EDIT: On contacting Esri Tech Support about this: "In general, we support installation and evaluation issues with trial ArcGIS features, but ArcGIS 9.2 has been retired for over 18 months and we no longer support this version of our product."

...I would guess that this means it's probably a compatibility issue with 9.2 and maplex rather than a standalone maplex issue, as I've just tried a single point at 1:250 scale with 72pt size labelling and again, works prior to enabling maplex but disappears after. Thanks for all comments and suggestions.]


texts: Convert data type of percent fields to number in Carto Builder

texts:
I have a dataset in Carto Builder with a column that lists percentage values with a percent symbol (e.g. "10%", "50%", "35%"). By default Carto has categorized the fields in this column as strings. I want to use this column to apply shading to a choropleth map based on percent. Is there any way to get Builder to process this column numerically?

If I try to convert the data type of the column from String to Number all the fields become null. I can do the conversion in my source program but I am looking for a way to do this all in Carto Builder.


texts: Fiona and Geopandas not applying CRS to shapefile during creation

texts:
I'm trying to create a shapefile of points from a coordinate list. But, when I use fiona to create the shapefile Arc Catalog is showing that there is no CRS. I tried to add the CRS with Geopandas and got the same result. I have my code below with file paths removed but it works with the files I use.

# import
import os
from PIL import Image, ExifTags
import geopandas as gpd
import fiona
from fiona.crs import from_epsg
from shapely.geometry import Point
from shapely.geometry import mapping
import shapely
from collections import defaultdict
# define functions
def dms_to_dd(d, m, s):
    dd = d + float(m)/60 + float(s)/3600
    return dd
pathforresult = r"Result.shp"
NCDict = {}
WCDict = {}
# opening files and extracting exif
Path = r"final project pics"
filelist = os.listdir(Path)
for f in filelist:
    if f.endswith(".jpg"):
        with open(Path + "/" + f, 'r') as I:
            print(I)
            img = Image.open(Path + "/" + f)
            exif = {ExifTags.TAGS[k]: v for k, v in img._getexif().items() if k in ExifTags.TAGS}
            print(exif)
            meta = exif['GPSInfo'][2]
            meta = [x[0] for x in meta]
            d = meta[0]
            m = meta[1]
            s = meta[2]
            NCDict[os.path.basename(I.name)] = dms_to_dd(d=d,m=m,s=s)
            # Ncoords = {I:dms_to_dd(d=d, m=m, s=s)}
            # print(Ncoords)
            meta2 = exif['GPSInfo'][4]
            meta2 = [b[0] for b in meta2]
            d = meta2[0]
            m = meta2[1]
            s = meta2[2]
            WCDict[os.path.basename(I.name)] = dms_to_dd(d=d,m=m,s=s)
print("North Coordinate List " + str(NCDict))
print("West Coordinate List: " + str(WCDict))

# writing xy coords to new file
corddict = [NCDict, WCDict]
dd = defaultdict(list)

for d in (NCDict, WCDict):  # you can list as many input dicts as you want here
    for key, value in d.items():
        dd[key].append(value)

print(dd)

schema = {'geometry': 'Point', 'properties': {'name': 'str:10'}}
# create and save the resulting shapefile
with fiona.open(pathforresult, 'w', crs=from_epsg(4326), driver='ESRI Shapefile', schema=schema) as output:
    for key, value in dd.items():
        point = Point(value[0],value[1])
        prop = prop = {'name': key}
        output.write({'geometry': mapping(point), 'properties': prop})
        print(output)
result = gpd.read_file(pathforresult)
print(result.crs)
result.crs = {'init': 'epsg:4326'}
print(result.crs)
result.to_file(r"Result2.shp", driver="ESRI Shapefile", schema=schema)


texts: 'Fields' variable does not show specified shapefile fields in Arcpy custom script

texts:
I'm creating a script in ArcMap (ArcGIS 10.7.1).

I have created tool parameters to specify the workspace where the data is stored (Workspace), a shapefile displaying a boundary (Field boundary), a shapefile of soil points that will be interpolated (Soil points), and a few others (crop type, yield goal). Also to note, the soil points contain fields such as "Nitrogen", "Phosphorus", "Potassium", and a few others. I have verified their presence.

Below that, there is a field called "Nitrogen field" that should contain a list of all the fields that are found in the soil points shapefile, but they are not there when I open up the drop-down box. When setting up the tool parameters, I selected this parameter as a 'fields' data type, but nothing is there. How can I make these fields visible in the drop-down box?

Here's a picture of the tool parameter box to provide a better visual:


texts:
Is there a way to draw a perfectly vertical line in QGIS's Composer. I wanted to created a visual separator between two elements for aesthetic reasons.

I've tried the arrow tool, removing the arrow head and holding the shift / ctrl / alt buttons to see if it will lock it to horizontal or vertical angle but no luck.

The only other alternative I could think of was making and importing an SVG, but I wanted to check here first, as that seemed a bit ridiculous.

Thanks for the help


texts: Exporting inset map as shapefile from ArcGIS for Desktop?

texts:
I learned how to set the inset map in the data frame properties.

I would like to export a permanent inset map layer or shapefile when I finish my map.

What is best way to do it or is there any alternative method?


texts: Convert features from EPSG:3857 to EPSG:4326 projection not working OpenLayers 4?

texts:
I have a simple task to do, supposedly OL4 and format.GEOJson.readfeatures should do that automatically. But it doesn't. Remembering that is a FeatureCollection and has CRS set.

I have a map with EPSG:4326 projection which works fine with same projection. If I change the map projection to EPSG:3857 and use format.GEOJson.readfeatures it all goes well. But when I try to use the other way, it doesn't.

I tried all properties and methods, searched all over Google, but no luck.

<script>
$(document).ready(function () {
    $('.tap-target').tapTarget('open');

    var typeSelect = 'None';

    <!-- Creating Basic Map -->
      var styles = [
        'Road',
        'Aerial',
        'AerialWithLabels',
        'collinsBart',
        'ordnanceSurvey'
      ];
      var layers = [];
      var i, ii;
      for (i = 0, ii = styles.length; i < ii; ++i) {
        layers.push(new ol.layer.Tile({
          visible: false,
          preload: Infinity,
          source: new ol.source.BingMaps({
            key: '#mykey',
            imagerySet: styles[i]
            // use maxZoom 19 to see stretched tiles instead of the BingMaps
            // "no photos at this zoom level" tiles
            // maxZoom: 19
          })
        }));
      }

      var map = new ol.Map({
        interactions: ol.interaction.defaults().extend([new app.Drag()]),
        layers: layers,
        // Improve user experience by loading tiles while dragging/zooming. Will make
        // zooming choppy on mobile or slow devices.
        loadTilesWhileInteracting: true,
        target: 'map',
        view: new ol.View({
          projection: 'EPSG:4326',
          //center: [0, 0],
          //projection: 'EPSG:4326',
          center: [-48.669159, -26.904534],

          zoom: {{ geo|length > 0 ? 17 : 13 }}
        })
      });

      var select = document.getElementById('layer-select');

      function onChange()
      {
        var style = select.value;

        for (var i = 0, ii = layers.length; i < ii; ++i)
        {
            layers[i].setVisible(styles[i] === style);
        }
      }

      select.addEventListener('change', onChange);
      onChange();
      <!-- Creating Basic Map -->

      <!-- Draw Feature -->
    var source = new ol.source.Vector({
        wrapX: true
    });

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

    var stroke = new ol.style.Stroke({
        color: '#3399CC',
        width: 3
    });

      var sourceStyle = new ol.style.Style({
          image: new ol.style.Icon(/** @type {olx.style.IconOptions} */ ({
              anchor: [0.5, 36],
              anchorXUnits: 'fraction',
              anchorYUnits: 'pixels',
              src: '{{ asset('/images/material-icons/maps/1x_web/ic_place_black_36dp.png') }}'
          })),
          stroke: stroke,
          fill: fill
      });

      var vector = new ol.layer.Vector({
          //projection: 'EPSG:25830',
          style: sourceStyle,
          source: source
      });

      map.addLayer(vector);
      vector.setVisible(true);

      $.ajax({
          type: 'get',
          url: 'http://datosabiertos.malaga.eu/recursos/urbanismoEInfraestructura/sedesWifi/sedesWifi-25830.geojson',
          type: 'json',
          success: function (r) {
              // eval pode incomodar, verificar isso para evitar problemas, mas feito assim para testes

              //EPSG::3857

              //console.log(r);

              //console.log(new ol.format.GeoJSON().readProjection(r));

              var feats = (new ol.format.GeoJSON()).readFeatures(r);
//I've tried to use featureProjection and dataProjection with no luck

              source.reportError = true;
              source.addFeatures(feats);
          }
      });

});
</script>


texts: The Vertex Editor in QGIS is displaying x,y coordinates of the nodes in longitude and latitude. Is there any way to make it display the x,y coordinates in meters?


texts:
I want to rasterize a shapefile (I don't mind whether to use QGIS or python) but I want to rotate my shapefile first and than have my raster aligned to this rotation.
So lets say I rotate my shapefile by 30° (I can do this easily with QGIS), I then want to rasterize that projection of the data into a rectangular grid. Is that possible? I tried a quick search on the internet but I always find it hard to google gdal/qgis functions.


texts: Can you set the PostGIS SRID when creating a Feature class through ArcCatalog?

texts:
I have ArcSDE installed in a Postgresql database and I am using PostGIS geometry. I'm not a DB minded person so feel free to let me know if my process is not right.

I create my feature class in ArcCatalog (no probs) I then go into PgAdminlll to check out the new table. There are two things I notice when I check out my new table.

1) The record for my new table in the public.geometry_columns table has an SRID of -1. Is there any way to set the SRID during the feature class creation? If not how do I do this in PgAdminlll? Is it something like updateGeometry or addGeometry?

2) Also I notice that in PgAdminlll under the newly created tables properties there is no primary key assigned. So what I have been doing is assign the sde created objectid field as the primary key. Is this acceptable or should I be creating a new PK?

Thanks
Ando


texts: Export Overpass-turbo Query Result to OSM File

texts:
I have a tool which has input is an osm file, I ussually get the osm file by exporting map data from http://openstreetmap.org.

But now I have to get only data of roads so I use this query http://overpass-turbo.eu/s/2vM and it show what I need.

My problem is when I select "Export" tab and chose jOSM option for exporting to osm file, the site show me pop up :

  Remote control not found. :( Make sure JOSM is already running and properly configured.

I don't know how to configure JOSM properly.
( Very sorry because my english skill is not good enough)


texts: How to reproject .las file from known EPSG (6499) to a custom coordinate system with WKT?

texts:
I am trying to use PDAL (filters.reprojection) to convert a .las file from a known EPSG:6499 to a custom projection.  Can I just put my WKT string into the pipline?
"type":"filters.reprojection",
"in_srs":"EPSG:6499",
"out_srs":"my_WKT_string_goes _here?"


texts:
QGIS is very powerful when it comes to styling Vector Data. When it comes to Raster data, ArcGIS still has some advantages including Styling based on Quantiles. See my earlier post:  For this I want to create a plugin. I got stuck when I want to determine quantiles for my raster (.tif) in a QGIS python Plugin based on this tutorial.

I got my raster stored as variable "layer". Now I want to determine quantiles. For this I could use a grass command called quantile. (r.quantile) however I don't know how to call this command in my plugin and how to store the output.

Bonus:
Next I want to set my raster classes equal to the percentile score or store the output as a style for QGIS (.qml)


texts:
I quite often get stuck trying to solve particular problems, or getting syntaxes right.

Apart from the details provided in the links on GDAL Utilities, does anyone know of any 3rd party resources/manuals explaining it in a bit more user friendly way?


texts: How to access data from a registered SQL Server database in ArcGIS Server

texts:
Little bit of background on the problem: I have been working on a PHP application with an SQL Server backend. The application works through a web browser, but, I've recently been tasked with allowing data access through ArcGIS Server in order to be able to use Survey123 to update data in the database.

The SQL Server database has been registered with ArcGIS Server and when I go to the Site tab -> Data Stores the database name shows up in the list of data stores and has been successfully validated.

However, my question now is, how do I access the database through ArcGIS Online? More specifically, how can I connect to the database through an online map?

I've looked on Google and this site for information on this but I haven't been able to find any. Does anyone know of a good place to find information on accessing a registered database through ArcGIS Server and displaying the data in a map?


texts: QGIS table join does not work with a join field starting with '0' (zero)

texts:
I'm trying to perform a table join in QGIS (2.14.5-Essen) between a shapefile and a attribute table.
The join field and the target field is the ID of a spatial enumeration unit of 8 numerical characters length. (Examples: 04041133, 12214422, ...).
The join does only work successfully for IDs that do NOT start with '0' (zero). For the IDs starting with '0' (zero) the other attributes are set NULL.

I have found posts of people who seem to have the same problem, but just workarounds where offered as solution. I have read that the plugins "MMQGIS" and "ftools" offer table joins, too, but both are not working in my QGIS-installation.
As mentioned, the join works successfully for a part of the rows; so it is probably not an issue of data types.
The solution i see in the moment is, creating new ID fields without zero as starting character in both tables, but is there a simpler solution?


texts:
Is there any free OSM Routing REST API without limits which supports multiple waypoints between the start point and destination AND is able to optimize these waypoints instead of calculating a route in the given order?
I've only found MapQuest which is limited by 500 Requests/day.

I know that a similar question has been asked 2 years ago but maybe there is an API now.


texts: Save/loading shapefiles with SF in R: no simple features geometry column present

texts:
I have this code:

library('sf')
#init data
xes = c(-122.133533701866, -122.141590665524, -122.382782560179, -122.302227764504, -122.347524548303)
yes = c(47.617342005802, 47.6102707417054, 47.6704037709971, 47.6118279934801, 47.7237286721321)
a = data.frame(id = 1:5, X = xes, Y = yes)

#init sf object
a = st_as_sf(a, coords = c('X','Y'))
st_crs(a) = 4326
st_write(a, "C:/Users/USER/Downloads/aa.shp", delete_dsn=TRUE)
test1 <- st_read("C:/Users/USER/Downloads/aa.shp", stringAsFactors = FALSE)

The final line (test1) returns:

  Error in st_sf(x, ..., agr = agr, sf_column_name = sf_column_name) :
    no simple features geometry column present

Which is odd, because using rgdal and converting to an sf object works just fine:

#Assumes code block above has already been run
library('rgdal')
a1= as(a, 'Spatial')

writeOGR(a1, 'C:/Users/USER/Downloads/aaa.shp', driver = 'ESRI Shapefile', layer = 'aaa', overwrite_layer = T)

#works
test2 = readOGR('C:/Users/USER/Downloads/aaa.shp')
test3 = readOGR('C:/Users/USER/Downloads/aa.shp') #from the earlier block
identical(test2,test3) #returns TRUE

b = st_as_sf(test2)

identical(a,b)

  [1] TRUE

SessionInfo:

  R version 3.5.1 (2018-07-02) Platform: x86_64-w64-mingw32/x64 (64-bit)
  Running under: Windows >= 8 x64 (build 9200)

  Matrix products: default

  locale: [1] LC_COLLATE=English_United States.1252
  LC_CTYPE=English_United States.1252    LC_MONETARY=English_United
  States.1252 [4] LC_NUMERIC=C
  LC_TIME=English_United States.1252

  attached base packages: [1] stats     graphics  grDevices utils
  datasets  methods   base

  other attached packages: [1] rgdal_1.3-6       sp_1.3-1
  data.table_1.11.8 openxlsx_4.1.0    sf_0.7-1

  loaded via a namespace (and not attached):  [1] Rcpp_0.12.19
  lattice_0.20-35 class_7.3-14    crayon_1.3.4    grid_3.5.1
  spData_0.2.9.4  DBI_1.0.0       magrittr_1.5    e1071_1.7-0     [10]
  units_0.6-1     pillar_1.3.0    zip_1.0.0       rlang_0.3.0.1
  rstudioapi_0.8  tools_3.5.1     yaml_2.2.0      compiler_3.5.1
  classInt_0.2-3  [19] tibble_1.4.2


texts:
Leaflet has a markerclustring plugin that clusters markers.

I am asking is there a solution to cluster your (lat/lon) with Geohash on Leaflet map in the a similar way to that on Kibana?


texts: Recreate a Local Custom Coordinate System from Survey Data in GIS?

texts:
I have a .DWG file of a series of survey points that has no associated projection file with it. The surveying group I received the data from states that the system they use is a local coordinate system in a "modified UTM 20N". The data are located on St. Kitts in the Caribbean. Here are some of the points in WGS84 DMS, WGS84 DD, and in the local coordinate location Easting and Northing (meters):

    Point_Id    Latitude    Longitude   Ellip_Height    Lat_DD  Long_DD CustomX CustomY CustomZ GRS
1   17° 16' 53.25342" N 62° 40' 43.25561" W 63.148  17.281459   -62.678682  2796138.580 -5412477.224    1882617.923 WGS84
2   17° 17' 58.85760" N 62° 41' 43.83679" W 82.806  17.299683   -62.695510  2794282.417 -5412782.257    1884549.473 WGS84
3   17° 20' 08.90031" N 62° 42' 26.74262" W -17.577 17.335806   -62.707429  2792566.808 -5412220.661    1888336.213 WGS84
4   17° 23' 29.53690" N 62° 46' 26.37868" W 6.196   17.391538   -62.773994  2785445.293 -5413844.723    1894230.305 WGS84
5   17° 24' 22.92018" N 62° 47' 54.72740" W -21.295 17.406367   -62.798535  2782889.824 -5414577.508    1895788.138 WGS84

Here is the info from the Leica File that the surveyor gave me:

 Ellipsoid A: WGS 1984
    Ellipsoid B - GRS 1980
    Type: Classical 3D
    Transformation type: Molodensky Badekas

    dx: 18.2922m
    dy: -35.5306m
    dz: 12.4868m
    Rx: -1.69221"
    Ry: -1.69378"
    Rz: -1.92504"
    SF: 5.7362 ppm
    x0: 2787863.0804
    y0: -5414455.175
    z0: 1888935.4451

Here is a screenshot I received from the surveyor

How do I define the reference system in ArcGIS Desktop with those given parameters?

I have seen the custom projection tool in ArcGIS, but I don't see any entry fields for the rotations he gave me (Rx, Ry, Rz):

Here is the screenshot from ArcMap


texts: How to make Geoserver dynamically calculate the bounding box of a postgis backed Layer

texts:
I set up a store for a postgis db.

I set up a layer that queries a table.

Geoserver requires that you set the Bounding Box of the layer when you create the layer.  The problem is that I don't know beforehand what the bounding box will be, since someone else will be writing on the table in the future.

Of course I could set a huge bounding box to cover all the expected points, but isn't it possible for Geoserver to compute the bounding box at runtime based on the actual data in the table?


texts: Calculating Z-scores using ModelBuilder in ArcGIS Desktop?

texts:
I am trying to take the natural log and then calculate z-scores on cell values in several composite band GRID rasters in ArcGIS Desktop.

I am not fluent in Python and do not see a way to automate this within ArcGIS Desktop.

I am in need of some coding help or a suggestion on how to automate this within ModelBuilder.


texts: Where to find Documentation for developing QGIS Plugins?

texts:
I thought I would try to write a plugin for QGIS (some Python code that I hoped I could easily port) I started following http://www.qgisworkshop.org/html/workshop/plugins_tutorial.html but that turned out to be hopelessly out of date. I turned to http://www.qgis.org/en/docs/pyqgis_developer_cookbook/intro.html instead and started going through the intro for the Python console and stumbled at the second query.

layer = qgis.utils.iface.activeLayer()

worked fine but then

layer.getLayerID()

threw up

Traceback (most recent call last):
File "<input>", line 1, in <module>
AttributeError: 'QgsVectorLayer' object has no attribute 'getLayerID'

and sure enough code completion in the console would seem to suggest that this doesn't exist. So my question is "where do I find up to date, usable instructions for developing plugins for QGIS 2.2?" It would seem that unless you are already "in the know" there is no way in. Or am I just being paranoid?


texts:
I'm a newbie to QGIS. Initially, I only want to display a map on the screen, nothing more, but I couldn't find a solution for my problem on the web. Also, I want to use that map as QWidget; is it possible or not? If not then I'll quit trying to implement a map with QGIS. (I used QML map but it doesn't provide my requirement.)
I using Qt 5.12
MinGW 73_64
QGis 3.14
My .pro file:
QT       += core gui
QT       += xml

greaterThan(QT_MAJOR_VERSION, 4): QT += widgets

TARGET = ss
TEMPLATE = app

DEFINES += QT_DEPRECATED_WARNINGS

CONFIG += c++11

SOURCES += \
        main.cpp \

HEADERS += \
        mainwindow.h

FORMS += \
        mainwindow.ui

INCLUDEPATH +=  "C:/OSGeo4W/include"
INCLUDEPATH += C:\OSGeo4W64\apps\qgis-ltr-dev\include

LIBS += -L"C:\OSGeo4W64\apps\qgis-ltr-dev\lib" -lqgis_core -lqgis_gui

DEFINES+=CORE_EXPORT=
DEFINES+=GUI_EXPORT=

qnx: target.path = /tmp/$${TARGET}/bin
else: unix:!android: target.path = /opt/$${TARGET}/bin
!isEmpty(target.path): INSTALLS += target

The main.cpp is :
#include <qgsapplication.h>

int main(int argc, char ** argv)
{
  // Start the Application
  QgsApplication app(argc, argv, true);
  QgsApplication::setPrefixPath("C:/OsGeo4W64/apps/qgis-ltr-dev");

  // Start the Application Event Loop
  return app.exec();
}

And the compiler says:
 error: undefined reference to `QgsApplication::QgsApplication(int&, char**, bool, QString const&, QString const&)'
 error: undefined reference to `QgsApplication::~QgsApplication()'


texts: Why hybrid GIS application gets GPS position only from Google & not from the device

texts:
I have built a GIS application using Leaflet on the web, after I put this web application to run on Android as a hybrid GIS application I made the application work with offline basemap (offline usage only) but the GPS finder can't get the position it trys to retreive it fom the Internet and not from the device.
is it because my GIS application isn't Native android ?


texts: Should I set the projection if the data has already been projected?

texts:
I need to plot some data which was produced on a native grid  with a lambert conformal conic projection. I need to create a geotiff (using python and gdal) and was wondering how to go about and do this. I have been told the following:

  Is your data on a native grid? A native grid is a model that was
  designed from the beginning with a particular map projection. In order
  to plot the data exactly as the designers planned, we do not want to
  transform the data to a projection, but simply plot it. This is the
  essence of the examples on this page.

  Just because your data has 2D lat/lon arrays, does not make it a
  native grid. Native grids are plotted differently than other data with
  2D coordinates.

from: https://www.ncl.ucar.edu/Applications/lcnative.shtml

The dataset comprises of a set of gridded data which were produced on a lcc grid. The guide states:

  In order
  to plot the data exactly as the designers planned, we do not want to
  transform the data to a projection, but simply plot it.`

How can I make a geotiff and what spatial and projection information should I use? The image xs and ys I have are in degrees but the grid spacing for xs and ys were difined in metres.


texts: Rectangle Buffer (with dimensions given in metres) around Points on map (lat long). Using QGIS 3.0

texts:
I have points on a map of the globe given with lat|long coordinates. I want to create a rectangular buffer around each one, but I want to be able to set the dimension of this buffer in metres. The boxes will have dimensions of around 100m x 200m so I need my solution to be pretty accurate. The tools I've found so far assume (fairly) that the parameters I pass to it (width, height) are of the same type as the point coordinates. Instead I want to pass parameters in metres whereas my point coordinates are lat long values.

Any ideas?

Converting metres to degrees somehow is an option but I'm worried about accuracy.


texts: Displaying unit as 'text' within label in QGIS?

texts:
When labelling areas in QGIS I can't figure out how to display the unit after the number value from a field.

For example the area of my polygon might be 1.23 hectares, which has been calculated using field calculator $area/10000.

How do I then display, for example, 'Ha' after the value, or do I find the tool to do so in the labelling options on the layer properties?


texts:
I installed on Oracle VM VirtualBox two Windows Server 2008 R2.

•   One with ArcGIS Server and ArcGIS Desktop 10.2.2 10.2.2

•   The second SQL Server 2012 in Mixed Mode

Because I'm not connected to a network, I use the Host-Only Ethernet Adapter Mode to let the host communicating with guests.
My questions are:
On what virtual machine should I have to install ArcSDE and Microsoft SQL Server 2012 SP1 Native Client.

If you could point me to a good tutorial, that would be great.


texts: Deleting definition queries of multiple layers in ArcMap?

texts: I'm searching for a way to freeze/delete definition queries of several layers at once (instead of doing it manually by opening the properties of each layer). I'm working in ArcMap 9.3.1 on an ArcSDE 9.3.1 Geodatabase.


texts: Moving substring from end to beginning of string

texts:
Disclaimer: I'm very new to QGIS and using expressions.

I have a large data set that contains a time range (00:00-00:00) and 1-6 three-day abbreviations of days of the week (MON, TUE, WED, etc.). I managed to get all the time ranges and day abbreviations to be the same format, but in some strings, the time is before the days of the week and in others it comes after.

eg.
"08:00-17:00 MON TUE WED THU FRI"
"MON THU 11:30-13:00"
"TUE 08:30-10:00"

I want to move the time range to the beginning of all the strings. I've gotten as far as filtering the features that end in numbers...

regexp_match("SIGNDESC1",'(\\d)$')

... and selecting the last 11 characters of the string:

right("SIGNDESC1", 11)

... but am stumped about moving them to the beginning of the string.


texts: How to change a marker on click with OpenLayers 2

texts:
I've created a map and tried, but it doesn't change until another event is caught or zooming map or by switching the browser tabs. What is the reason for this?

Please help me, I've fiddled the code here http://jsfiddle.net/HMFLf/.

This is what I've done to change marker:

vectorLayer.events.on({
                featureselected: function (event) {
                    var fixer = confirm("Are you want to change clolor this marker");
                    if (fixer == true) {
                        size = new OpenLayers.Size(21, 25);
                        event.feature.style = { externalGraphic: "http://www.openlayers.org/dev/img/marker.png", graphicHeight: 25, graphicWidth: 21, graphicXOffset: -((size.w / 2) - 2), graphicYOffset: -(size.h - 2) };
                    }
                }
            });


texts: Create a runtime geodatabase in QGIS for an Android app

texts:
I need to generate a runtime geodatabase for map visualization on my android device. Back in the days we used ArcGIS with local file-geodatabases, where we were just able to export a runtime geodatabase. Nowadays we have to use QGIS with an object-relational database and I was not able yet to find any solution to my problem. I found this Is there a way to create a Personal Geodatabase in QGIS? , but this is already 5 years old and does not really fit my use case.

Is it possible to create a runtime geodatabase out of a QGIS project?

Or asked differently:

Is there a way to accomplish in QGIS what ArcMap can do in the picture below?


texts: How do I add .e00 data into my QGIS project?


texts: Why does GDAL open a DTED file as a vector format?

texts:
GDAL version: 2.1.1

Language: C++

Given the path of a geospatial file, I want to know whether the file is a raster or vector format.  The only way I have found so far to make that distinction is to attempt to open the file as vector format first.  If that succeeds, then treat it as a vector file; otherwise, attempt to open the file as raster format.  The following code snippet is an example.

This worked for GeoTIFF and shape files.

However, for a DTED file, opening it with the GDAL_OF_VECTOR flag gives me a valid dataset handle.  So, my code incorrectly decides that it is a vector format.

Why did GDALOpenEx() open the DTED file when the GDAL_OF_VECTOR flag was used?

void test( string const & path )
{
  GDALDatasetH hDataset;
  bool isRaster = false;

  unsigned int openFlags = GA_ReadOnly | GDAL_OF_VECTOR;
  hDataset = GDALOpenEx( path.c_str(), openFlags, NULL, NULL, NULL );
  if ( hDataset )
  {
    printf( "Opened as vector\n" );
    isRaster = false;
  }
  else
  {
    openFlags = GA_ReadOnly;
    hDataset = GDALOpenEx( path.c_str(), openFlags, NULL, NULL, NULL );
    if ( hDataset )
    {
      printf( "Opened as raster\n" );
      isRaster = true;
    }
    else
    {
      printf( "Not opened, skip it\n" );
      return;
    }
  }

  // Process the file differently, depending on its format.
  if ( isRaster )
  {
    ...
  }
  else
  {
    ...
  }

  GDALClose( hDataset );
}


texts: How to change the path to the mysql socket in Quantum GIS?

texts:
I would like to use Mysql in QGIS. But when I want to create a new connection to the database (Add Vector Layer), I can't connect because QGIS is looking for

/tmp/mysql.sock

But I'm using XAMPP for Mac OS X so the socket is in

/Applications/XAMPP/xamppfiles/var/mysql/

Is there a way to change the path to the mysql socket in QGIS ?

Thank you.


texts:
I have a feature layer and some pins on it, with mouseout, mouseover events callbacks on the feature. When the mouse is over the pin icon, that triggers mouseout event (mouseout from the feature layer). I am trying to stop event propagation but nothing works.

How can the leaflet event propagation be stopped?

None of those worked:

DomEvent.stopPropagation(e);
DomEvent.preventDefault(e);
DomEvent.stop(e);
stopPropagation(e.originalEvent);
preventDefault(e.originalEvent);
stop(e.originalEvent);


texts: I've created a freehand shape, but can't delete it. Screen shot attached. I can select it but the trash bin is grayed out, and the delete button prompts the error code "provider does not support delete". How do I get rid of this thing? Or do I have to delete and recreate the whole layer?


texts: Plotting a point within a bounding box as a percentage of X/Y

texts:
I have two pairs of lat/long coordinates that represent the top left and bottom right corners of a bounding box.

I have a third pair of lat/long coordinates representing a point within this bounding box.

How might I go about converting this point to a percentage of the X/Y position within the bounding box? For example, a point at X=50% and Y=50% would be bang in the centre of the bounding box.

It does of course need to take into account positive and negative coordinates.

Any guidance appreciated!


texts: How do I accurately measure a gps track over a long distance?

texts:
I am measuring the length of the Continental Divide Trail that goes from the Mexican border to the Canadian border. I have a .gpx file of the entire route in NAD83 with tracks about every 50 meters. I've tried two methods:

If I convert to North America Equidistant and measure the length
of the line in QGIS and then divide by 1609 (meters in a mile), I get
2845 miles.
If I convert to WGS84, measure the length of each line segment and then multiply by the cosine of the latitude, I get 2952 miles.

Published lengths have the trail at about 3050 miles. Which is most accurate? Is there a more accurate method? Do I need a 3d length measurement since the trail goes over mountainous terrain?


texts: Creating matrices of coordinates for a given geoTIFF

texts:
I have GeoTIFFs that are in a given projection, I'm trying to figure out what the easiest way would be to produce corresponding matrices for the the GeoTIFFs that hold the lat, lon, and meters locations for each pixel using the projections of the geoTIFFs.
Is there a straight-forward way of doing this that anyone can recommend?
I have tried to do this using GDALwarp to assign a coordinate system, and then used geotiffread in MATLAB to read the data, but I can not ascertain whether or not there is a way to use MATLAB (or another program) to build matrices that hold the coordinates information for each pixel.


texts: Convert a SpatialPolygonsDataFrame to raster using rasterize function

texts:
I'm trying to convert my map (a SpatialPolygonsDataFrame) into a raster, using rasterize function, here is a bug that I don't understand...

> class(parcel) # this is my map
[1] "SpatialPolygonsDataFrame"

> attr(,"package")
[1] "sp"

> head(parcel)
AREAHA RECNO      AREA AREASQKM        OS_06 GRD_CAT_06
0 0.0519     0  519.1941        0 chemin prive        way
1 0.0127     1  127.0963        0 chemin prive        way
2 0.0596     2  596.2204        0 chemin prive        way
3 0.0248     3  248.0452        0    route D90        way
4 0.0375     4  375.0141        0 route privee        way
5 0.1260     5 1260.3941        0 chemin prive        way

> dim(parcel)
[1] 20069     6

> object.size(parcel)
88090144 bytes

> r <- raster(ncol=4001, nrow=4001)
> extent(r) <- extent(parcel)
> parcel_ras<-rasterize(parcel, r, "GRD_CAT_06", CRS("+init=epsg:27572"))
Erreur dans data.frame(ID = as.integer(putvals), value = lvs) :    arguments imply differing number of rows: 20069, 13

I just don't understand this comment : what is "putvals" ? There are 20069 polygons in my map, so what is "13" ?

I can't reproduce a smaller example of my map, so here are my data : http://dl.free.fr/fxdT2oczB


texts: ArcGIS Manager fails login once, then succeeds

texts:
Often, when a user attempts to login to ArcGIS manager through online, although the username and password is correct, it will display,

  The username or password you entered is incorrect.

however clicking Login a second time allows user to login without issue.

What could this be a symptom of?

I have seen a similar question that suggests that one should preface a domain, however I am fairly sure that the ArcGIS Server is configured using a standalone server (therefore I'm unsure what domain to preface the username with).