texts: What is the difference between spatialreference.org catalogs: EPSG, ESRI, and User-defined?

texts:
I'm trying to understand EPSG codes.  I've been to the spatialreference.org site, and in its "about" link, it states:

"There are currently three distinct catalogs of spatial reference systems on this website -- EPSG, ESRI, and User-defined (called sr-org)."

What is the difference between these catalogs?  Maybe a better question is why are there separate catalogs?

I'm moving from the ESRI desktop world to QGIS, and I want to make sure that I have a firm understanding of how QGIS handles coordinate systems.


texts:
EPSG 25832 = ETRS89/UTM zone 32N…

Wouldn't it be more correct to call it ETRF89/UTM zone 32N and, if so, what's the equivalent EPSG code for ETRF2000/UTM zone 32N?


texts: How would I have an improper DMS value input return a message to the user in a Python 3 function?

texts:
My question is about input functions in python 3 specifically for a latitude conversion function from dms to dd.

I'm new to python3 and having some trouble understanding how I would allow a user to input their latitude DMS coordinate and have them limited by a 0-90 range, so if they put something like 100, or -1, the code would return a statement asking them to input a valid latitude DMS format.

Below is my code - I think my if statement is where the trouble is happening. When I'm indicating if variable(x,x,x) - would that indicate the format the user would need to input their coordinate in?

def dmslat_convert_to_dd(degrees, minutes, seconds):
    dd = degrees + float(minutes) / 60 + float(seconds) / 3600
    return dd

dmslat = input('Please enter your DMS coordinate for Latitude: ')
if dmslat_convert_to_dd(degrees, minutes, seconds) in range(0, 90):
    print(dd)
else:
    print("Please enter correct DMS Latitude coordinates")


texts: Is it possible to render POIs/landmarks from Postgres/postgis DB using google Map like icons? Any readily available icon sets?


texts: Building disconnected web based mapping application?

texts:
I have been tasked with building a disconnected, web-based mapping application.  It must run in a browser, but it will run on a disconnected Windows based laptop.  It needs to be able to do simple map displays as well as display of objects on map and historical displays of the route of an object.  Does anyone have any suggestions about a particular software solution?  I am guessing that I will have to install some sort of mapping server software that has all of the mapping information and then also build a web-based client that accesses the local mapping server software.  Does that make sense?

Here are more details about the rough system requirements.

The finished application will be web-based (runs in a browser)
The entire system will be disconnected from the internet. So it can't rely on any internet-based mapping system like Google Maps.
The main screen will display a map of the earth. It needs to show basic country boundaries but there is no requirement for any other information.
The user should be able to zoom in and out. There will be some zoom limits although I don't know exactly what they are. The application is displaying flight paths for airplanes so hopefully that gives you some idea of the zoom needs.
The default view will keep the icon for the airplane centered and the map will move beneath.
The user should also be able to pan the map in any direction and then return to having the airplane icon centered.
The system will need to be able to dynamically load some overlay maps that display some polygon type shapes.
The main screen will also have text areas that display other types of information such as current airplane location, etc.
Ideally the mapping software would be easily integrated using C# and ASP.NET

These are the basic rough requirements. I have done a LOT of searching on the web for a solution. Most of the mapping software applications that I have found seem to be overkill. They are complex GIS systems that cost lots of money. I just need a simple mapping application like Google Maps that works offline and allows for the integration functions that I have outlined above. Does anyone have any recommendations?

(I mistakenly started a new question that is basically the same thing as this.  Here is a link to the second question: Looking for web-based mapping software recommendation)


texts: Python script to create a point on a line with attribute data

texts:
I have a line feature class with 2 lines (Engineering and Continuous Station series) co-incident on each other. I have created a toolbar(ArcGIS Desktop add-in 10.1) with 3 tools on it, add control point, move control point and delete control point. It's just an empty toolbar with no business logic.

When one selects the add control point tool and clicks anywhere on the line feature class, a point needs to be created on the line and a few attributes from the line needs to be brought in to the point. Also the distance of the point location from the beginning of the line needs to be calculated and attributed as well. The Move control point tool should move the points to a different location by specifying with a click on the line. Delete control point should delete both the control points.

I would like a python script to put them in the individual tools. Any ideas? A pseudocode would also help me work from my end. Thank you!!


texts: Script giving error accessing field name?

texts:
In the code below, the line 'alias = dict[fld.name]' always throws an exception. I don't see why; printing out fld.name works just fine. Also, the exception message is the field name printing correctly.

Printing the dictionary keys and values also works fine. The ultimate goal is to get a string to set as the field alias using AlterField.

The dictionary seems to be built ok. I can print out key/value pairs with iteritems(). The problem is in accessing the values, using field names returned by arcpy.ListFields, which also print on on their own. I've tried setting the key (fld.name) to a variable, enclosing it in quotes, to no avail. I keep getting an error in the form 'KeyError: u'B01001m1'. B01001m1 is a valid key.

with arcpy.da.SearchCursor(metatable, metaflds) as rows:
        for row in rows:
            if not moe in row[1]:
                fn = row[1].replace(' ', '')
                fn = fn.replace(est, '')
                dict[row[0]] = fn

    # Get list of tables to add alias to
    tbls = arcpy.ListTables()
    for tbl in tbls:
        tblflds = arcpy.ListFields(tbl)
        for fld in tblflds:
            try:
                alias = dict[fld.name] # exception here
                print(alias)
            except Exception as e:
                import traceback
                import sys
                tb = sys.exc_info()[2]
                print('Oh no!')
                print("Line {0}".format(tb.tb_lineno))
                print(e.message)


texts: Computing area of given polygon (in vector format) lying between specified lines of latitude?

texts:
I have species distributions (in vector format) and I would like to calculate how much of the area for a given species distribution is within certain lines of latitude, such as temperate and tropic zones.

Using the Wikipedia definition: The north temperate zone extends from the Tropic of Cancer (approximately 23.5° north latitude) to the Arctic Circle (approximately 66.5° north latitude). The south temperate zone extends from the Tropic of Capricorn (approximately 23.5° south latitude) to the Antarctic Circle (at approximately 66.5° south latitude).

Thus, the tropics would be between 23.5° south latitude and 23.5° north latitude.

For example, using this shapefile of the Atlantic Ocean (choose shapefile in the drop-down menue) plotted on a worldmap, one could easily calculate the total area of the Atlantic Ocean;

require(sp)
require(ggplot2)
require(rgeos)
require(rgdal)
require(maps)

setwd("~/test/iho")
ao <- readOGR(getwd(), layer="iho")
aof <- fortify(world, region="name")

# Not necessary for the calculation per se, but still nice. Although not the best looking map
world <- map("worldHires", col="gray90", fill=T)

# Plot
pp <- ggplot(data = world, aes(x = long, y = lat, group = group)) +
geom_polygon(fill = "grey50") +
geom_polygon(data = AO, aes(x = long, y = lat, group = group), fill = alpha("cyan", 0.5)) +
coord_equal()

# Total area
gArea(SpatialPolygons(ao@polygons))
[1] 7512.821

But how could we restrict the calculation of the area to specified lines of latitude, e.g. tropic and temperate zones as defined above?

Plot


texts: Setting "Map View" instead of "Split View" in Leaflet javascript

texts:
I am using geojson-dashboard framework template for Leaflet (http://fulcrumapp.github.io/geojson-dashboard/).

By default, the initial view of the screen is "splitted" or "divided" (half view is the map and the other half is a geojson table). What I'm trying to do is to modify the function in order to set the initial view of the screen only with the map.

The function that acts on that is the following:

function switchView(view) {
  if (view == "split") {
    $("#view").html("Split View");
    location.hash = "#split";
    $("#table-container").show();
    $("#table-container").css("height", "55%");
    $("#map-container").show();
    $("#map-container").css("height", "45%");
    $(window).resize();
    if (map) {
      map.invalidateSize();
    }
  } else if (view == "map") {
    $("#view").html("Map View");
    location.hash = "#map";
    $("#map-container").show();
    $("#map-container").css("height", "100%");
    $("#table-container").hide();
    if (map) {
      map.invalidateSize();
    }
  } else if (view == "table") {
    $("#view").html("Table View");
    location.hash = "#table";
    $("#table-container").show();
    $("#table-container").css("height", "100%");
    $("#map-container").hide();
    $(window).resize();
  }
}

$("[name='view']").click(function() {
  $(".in,.open").removeClass("in open");
  if (this.id === "map-graph") {
    switchView("split");
    return false;
  } else if (this.id === "map-only") {
    switchView("map");
    return false;
  } else if (this.id === "graph-only") {
    switchView("table");
    return false;
  }
});

I've tried in a thousand ways but I don't know how to solve in order to establish "Map View" as the initial view.


texts: Why is there an upper limit to the SRID value in the spatial_ref_sys table in PostGIS 2.1?

texts:
I was trying to insert a reference system in PostGIS from http://spatialreference.org/, and got a constraint violation for spatial_ref_sys_srid_check. I checked that constraint and found that, when inserting in spatial_ref_sys, the srid must never be greater that 998999.

I understand that the root cause of the problem is that spatialreference prepend a 9 to the srid. Probably to avoid collision srid in the srid column, or something of the sort. I don't really care about that, its easy to solve: removing the 9 make it work most of the time.

But why that arbitrary upper limit on the srid number? It is way below the maximum value for a Postgresql integer, so it doesn't seem to be a technological limit, and the number itself look quite arbitrary. Introducing that limit broke the code in http://spatialreference.org/ so there is an obvious downside.

So what was the rational for introducing such an upper limit, at that specific value? Is that a limit imposed by a standard body?


texts: Creating Conical Polygon around multiple circular buffers?

texts:
I have multiple circular buffers in a feature class. My goal is to create a polygon around these buffers, i.e create a conical polygon between 1 and 2 buffer, then the same with 2 & 3  and so on and finally dissolve all the borders to make a final polygon as shown in Fig2.

I tried convex hull by using minimum bounding geometry tool, but it basically covers a larger area if there are more than 2 records.It does not retain the schema and attributes.

I couldnt find any other Geoprocessing tool.I couldnt figure out how to go about it using ArcPy.


texts: Fail to export shapefile created from X,Y table

texts:
I have a dbf file with X,Y coordinates. I chose to display X, Y data in ArcGIS and got a temporary shapefile. However, when I try to export the temporary shapefile, I got an empty one.

There is nothing wrong with the titles of the variables in the table. And I tried to re-save the table in the other format and re-create the shapefile. It does not work.

Did anyone have this issue before? Thanks for your help!


texts:
Trying to add a map service form arcgis server to arcgis desktop, but what surprises me that I can add any layer " service " to the table of content except the Arial photos layer, see the screenshot .

What might be the issue here ?


texts: Need help regarding customizing Dojo/digit tool for ArcGIS Java Script API 2.8

texts:
I am developing this tool (Similar tool for the India like ‘State’ and then          ‘district’ level) and I am trying different methods so please let me know which method is suitable or not suitable to develop similar tool.

Tool Details:

According to the tool the all data and extent is stored in .json file and its reference given in the zoom_to.js file.

Methods that I am following:

1.Instead of .json file I am storing my data in sql server (where column type is geography so that its storing the extent i.e. xmin,ymin etc.) and then assigning those result to that tool.

2.I am modifying the .json file and replacing the values but getting some errors.

3.Created the two combo-box  and assigned  for “State” and for “District” respectively but not getting ‘zoom to‘ effect .

So Please let me know whether I am going to correct direction or not

Thanks in advance : )


texts: I am very much new to InSAR data. At the moment, i'm stuck at this website https://ursa.asfdaac.alaska.edu. Just wanna know, is it possible to get multi-temporal InSAR data for a particular region from this site or from any other sources. Will be OK even if it's a paid service.


texts: Polygon misalignment with OpenLayers & Ordnance Survey data (possible projection issue)

texts:
I'm trying to plot polygons on an OpenLayers map of Scotland which represent the area covered by UK Ordnance Survey maps. These polygons should cover the whole of Scotland and have no gaps between them. However when the polygons (with yellow edges) appear on the map they appear to be misaligned with slight gaps between some of them. The polygons should overlap in places (some locations appear on two maps) but there should be no gaps.

I've checked the co-ordinates (in both their original OS X & Y co-ordinate form and their converted latitudes and longitudes, done using PHPcoord http://www.jstott.me.uk/phpcoord/) and they seem to be accurate and fine. Which makes me think it might be a projection issue, but I can't figure out what I'm doing wrong or how to change it.

I wasn't able to create a JSFiddle using my actual GeoJson file so the JSFiddle below shows my map code but with the GeoJson file missing.

https://jsfiddle.net/MrAndySweet/e4qrs42w/

I've validated the GeoJson file and it seems to be ok.


texts: How do I zoom/Pan to a Leaflet Map, such that the given point is off-center?

texts:
In my Leaflet Application, on some event, I want to zoom to a particular marker, and then show it's popup. I have certain other custom components on top of the map, and the Popup that appears, is hidden behind these components.

This could be resolved, if I could zoom/pan to the map, such that the given point (i.e. location of marker) is not in the center, but  slightly lower.

How do I do this?


texts: Copy and Paste Feature Selection to new project - ArcMap

texts:
The reason for this question is I don't know the terminology to describe what I want to do!

Basically in the attached image I have take a selection of features. I want to, effectively, copy and paste that selection into a map blank map so I can zoom in just one those features without the context of the rest of
the map.

I just want to isolate these features alone for study, however copying the selection and pasting from clipboard gives me the entire map.

Is this doable?


texts: Opening multiple attribute tables using .showAttributeTables()

texts:
I've got a question about .showAttributeTables() - i'm trying to open two attrbute tables at the end of a query (selection of lines within a polygon layer) - I thought I could simply do somthing like this

# Spatial selection
cm_line = QgsProject.instance().mapLayersByName("t_line")[0]
poly = QgsProject.instance().mapLayersByName("t_poly")[0]

for i in cm_line.getFeatures():
    for l in poly.getFeatures():
        if i.geometry().within(l.geometry()):
            # add the selection to the layer
            cm_line.select(i.id())

# open tables
layer = QgsProject().instance().mapLayersByName('t_poly')[0]
iface.showAttributeTable(iface.activeLayer())

layer = QgsProject().instance().mapLayersByName('t_line')[0]
iface.showAttributeTable(iface.activeLayer())

But this opens the same table twice (the first one) - But I thought i'd changed the active layer here :
layer = QgsProject().instance().mapLayersByName("t_line")[0]  ??

Well, at least I didn't get an error message!

I should note that I've already seen this :
Opening attribute table using PyQgis?

and this tutorial here:
pyqgis-101-viewing-vector-layer-attributes


texts: My input join field and output join field are text type, and my problem is that when I use "Add Join" all joins ok, but if I want to use "Join field" there are some records with Null, it seems like it doesn't recognize that the text are the same in fields from both tables.


texts: Script for Sum of Rows in field calculator - VB or Python?

texts: I wrote the TOTAL_FREQUENCY values manually. But I want to automatically write this values with a script. Namely, I want to write sum of column 3 to all of the the column 5 rows. How can I do it with a script?


texts:
I have a working OpenLayers 2 (OL2) application currently. To create the WMS layer in OL2, I use:

layer = new OpenLayers.Layer.WMS('WMS', '../../WMS/floorplan?', { layers: '@Model.LayerName', version: '1.3.0', format: "image/png", shapefile: '@Model.ShapeFile'  }, { isBaseLayer: true, singleTile: false, yx: [] });

My equivalent on OpenLayers 3 (OL3) is:

var layers = [
      new ol.layer.Tile({
          source: new ol.source.TileWMS({
              url: '../../WMS/floorplan?',
              params: {
                  'LAYERS': '@Model.LayerName',
                  'SHAPEFILE': '@Model.ShapeFile'
              }
          })
      })
    ];

var map = new ol.Map({
        layers: layers,
        target: 'map',
        view: new ol.View({
            extent: ext,
            projection: parameters.projection,
            center: [0, 0],
            zoom: 0.5
        })
    });

However, the OL3 one only gets a single tile from the WMS call while the OL2 one gets all the tiles and renders it correctly. What could the issue be?


texts: How to save a vector layer from QGIS to Sql Server 2008 or 2012?

texts:
using QGIS, i've loaded up some administrative boundary data of New Zealand. I wish to now import this vector data from QGIS into my Sql Server 2012.

Is there any way to do this with QGIS 1.8.0-Lisboa?


texts: How can I map transit access without crashing Network Analyst?

texts:
I want to make a nuanced heat map of transit access, but my computer cannot handle the processing.

I'm using a GTFS table with a each trip in a city's bus system. Each trip has a unique trip_id. I've used the lat/long to place a point for every stop on every bus trip.

This looks like a map of bus stops in the city, but in fact there are dozens or hundreds of points on top of each other -- one for each time a bus arrives at that stop.

I want to create the 0.5-mile service area using Network Analyst around each of these points. My computer crashes here. If I could complete this step, I would then merge the resulting polygons by trip_id. This would generate each trip's walkshed as a single polygon. I would then spatial join those polygons and rasterize the count of overlapping polygons to show how many distinct transit trips you can walk to from a given point.

This is superior to an easier method: using the number of times a bus arrives at a stop. If you live between Stop A and Stop B, your location would get "credit" for when Bus X arrives at Stop A and again when it arrives at Stop B, even though in reality that isn't materially increasing your transit access.

How can I change this analysis so that my computer doesn't crash? Or how else can I capture this nuance in a simpler method?


texts: Plot a bubble chart in d3.js by giving an local geojson file as an input

texts:
This is the content in my geojson file

{
"type": "FeatureCollection",
"features": [{
    "type": "Feature",
    "properties": {
        "profit": 326,
        "npa": 174.000000
    }
}, {
    "type": "Feature",
    "properties": {
        "profit": 1762,
        "npa": 1683.000000
    }
}]
}

I am storing the data in a taluks.geojson. I want to give this local geojson file as an input and display the bubble chart based on the value profit.

My bubble chart is based on this : http://bl.ocks.org/mbostock/4063269


texts: Difference in CRS when reading with read_stars und read_ncdf

texts:
I have a NetCDF file that I would like to load into a stars object.
The two functions st::read_stars and st::read_ncdf both load it in, but the data ends up in different locations when I plot it using mapview.
I think I'm not getting something really basic here and just wanted to ask what the best way is to read NetCDF-files into a stars object. (The one with read_stars ends up in the right location).
I'm using the exactly same path.

The description of two objects look like this:
Read with read_stars
stars object with 3 dimensions and 1 attribute
attribute(s):
 DAILYPCP_198001.nc [mm]
 Min.   :0.00
 1st Qu.:0.00
 Median :0.00
 Mean   :0.07
 3rd Qu.:0.00
 Max.   :4.10
 NA's   :194505
dimension(s):
     from  to                                           offset                                       delta  refsys point values
x       1 641                                           605720                                         250 unnamed    NA   NULL
y       1 643                                          5220237                                        -250 unnamed    NA   NULL
DATE    1   1 3652.5 [(days since 1970-01-01 00:00:00 +00:00)] 1 [(days since 1970-01-01 00:00:00 +00:00)]      NA    NA   NULL

x    [x]
y    [y]
DATE

And the one read with read_ncdf
stars object with 3 dimensions and 1 attribute
attribute(s):
 precipitation [mm]
 Min.   :0.00
 1st Qu.:0.00
 Median :0.00
 Mean   :0.07
 3rd Qu.:0.00
 Max.   :4.10
 NA's   :194505
dimension(s):
     from  to         offset  delta                       refsys point values
x       1 641         605720    250 +proj=tmerc +lat_0=0 +lon...    NA   NULL [x]
y       1 643        5059487    250 +proj=tmerc +lat_0=0 +lon...    NA   NULL [y]
DATE    1   1 1980-01-01 UTC 1 days                      POSIXct    NA   NULL


texts: Hardware for Surveying and mapping total station

texts: Is "Microsoft Windows CE 6.0 / Magnet Field" operating system the sole property of TOPCON total station or it can be found in other total stations as well?


texts: pyqgis : Load WMS layer and make him invisible (qgis 3.4.3 )

texts:
I am loading WMS layer in Python console with simple code. Then I try to make the layer invisible in TOC, but not success. I run the code in QGIS 3.4.3 in Python console in Windows 10.

from qgis.core import QgsLayerTreeLayer,QgsProject

from PyQt5.QtCore import Qt

url_cur= "url='https://ags.iplan.gov.il/arcgis/rest/services/PlanningPublic/Xplan_2039/MapServer' layer='0'"

layer = iface.addRasterLayer( url_cur, 'a1', "arcgismapserver")

QgsProject.instance().addMapLayer(layer,False)

treeLay0 = QgsLayerTreeLayer (layer)

treeLay0.setItemVisibilityChecked (Qt.Unchecked )


texts: How can I work with very large shapefiles (~1 GB) in QGIS?

texts:
I'm performing some basic spatial analyses on a large database of traffic injuries in California between 2003 - 2011. Because the dataset is large (nearly 1 GB of points), I'd like to first cut it down to a specific geographic region using a spatial query, but I find that QGIS consistently freezes and hangs if I try to use a spatial query or filter the layer.

What are some ways I can work with this dataset? Are there more efficient data formats I can use?


texts:
I wish to edit and remove qgsfeatures in qgsvectorlayer. Vectorlayer contains postgis table.

QgsVectorDataProvider::capabilitiesString() returns the folowing string:
Add Features, Add Attributes, Delete Attributes, Fast Access to Features at ID, Simplify Geometries, Simplify Geometries with topological validation

So it looks like I can only add new features, not edit/delete already created.
Why?


texts:
I've done this:

while (feature != null)
            {
                if ((feature.get_Value(feature.Fields.FindField("NAME")) != null) && (feature.get_Value(feature.Fields.FindField("NAME")).ToString() != " "))
                {
                    itemList.Add(feature.get_Value(feature.Fields.FindField("NAME")).ToString());
                }
                feature = featureCursor.NextFeature();
            }

Now I have a listBox with fields values. How can I get a feature back by its name?

The idea was:

create a list with objects names,
find feature on map, clicking by its name on list.

Features (objects) represent as points.


texts: Selecting polygons within another polygons

texts:
I have the same city map (BH city) for two different shapefiles with different locations. The first one is the "shared management", this one have bigger polygons, as we can see

shared_management <- read_sf(dsn = 'TGC_POP_DENGUE_00_11_12_TAXA.shp')
plot(st_geometry(shared_management))

And I have the smallers ones ("sectors")

sectors <- read_sf(dsn = 'SETOR_CENSITARIO_2010.shp')
plot(st_geometry(sectors))

Then I need a list of polygons of "sectors" which is inside of each polygon of "shared_management".

I know the number of polygons of "sectors" inside of which polygon of "shared_management". This number is given by

geo.int <- st_intersection(sector, shared_management)
geo.agg <- aggregate( ID_CS_2010~ID_TGC, data = geo.int, FUN = length)
geo.agg1 <- aggregate( ID_CS_2010~IDENT_TERR, data = geo.int, FUN = length)
head(geo.agg1)
IDENT_TERR ID_CS_2010
1          B1        112
2          B2         98
3          B3         85
4          B4        123
5          B5        127

Any idea of how can I obtain which polygons of "sectors" are inside of each polygon of "shared_management"?


texts: MFC Dialog Application using ESRI ArcEngine 10.0

texts:
I'm new to ESRI ArcObjects and I'm trying to create a basic MFC dialog application. I add the License and MapControl to the dialog resource. I compile and run the app, but I get the following runtime error:-

Could not bind to a valid ArcGIS installation.
CoCreateInstance of OLE control {D7EAD71E-B663-4DCE-82C7-AEE2793DF5F3} failed.

      Result code: 0x80040111
      Is the control is properly registered?
      Warning: Resource items and Win32 Z-order lists are out of sync. Tab order may be not defined well.
      MFCDlg2.exe has triggered a breakpoint

I have registered the OCX's and I can see them in the registry. I'm running Win7 (64 bit) and Visual Studio 2008.

Please can anyone help

Thanks.


texts: Is there any way to configure MapFish config.yaml from Javascript?

texts:
I am using the mapfish print module in my webmapping application.In that I need to have various mapfish print layouts with different background pdf options. For that, am I need to create each config.yaml? or can we configure it from client side (from javascript)?

System: Windows XP / Tomcat7 / Geoserver 2.4.0/printing-2.1-SNAPSHOT.jar

For Example:

 #===========================================================================
 A4 portrait:
#===========================================================================
mainPage:
  pageSize: A4
  landscape: false
  rotation: true
  backgroundPdf: '${configDir}/MyCompany.pdf' #Here I need to change as per the client template

  items:
   - !columns
      absoluteX: 50
      absoluteY: 690
      width: 502
      config:
        borderWidth: 0.3
        cells:
          - padding: 2
      items:
        - !map
          width: 500
          height: 600

and my javascript code for layers:

    var printurl = 'http://xxx.xxx.x.xx:8080/geoserver/pdf/print.pdf?spec={' +
'"layout":"A4 portrait",' +
'"srs":"EPSG:404000",' +
'"units":"m",' +
'"dpi":300,' +
'"outputFilename": "mapPrint",' +
'"layers":[{"baseURL":"http://xxx.xxx.x.xx:8080/geoserver/cite/wms",' +
            '"opacity":1,' +
            '"singleTile":false,' +
            '"customParams":{},' +
            '"type":"WMS",' +
            '"layers":["Oxford_Map"],' +
            '"format":"image/png",' +
            '"styles":[],' +
            '"overview":true' +
          '}],' +
  '"pages":[{' +
            '"center":[' + osMap.getCenter().lon + ',' + osMap.getCenter().lat + '],' +
            '"mapTitle":"Sample Title",' +
            '"comment":"This is the map comment",' +
            '"scale":200000,' +
            '"rotation":0}' +
           ']}';

          window.open(printurl);

Kindly suggest me the solution to achieve this. Thanks in advance


texts: Convert bounding box for OL2 to extent for OL3?

texts:
I have a database with cities and values for: bb_left, bb_top, bb_right, bb_bottom, locate_x, and locate_y.

Here's an example:

city    -   bb_left -   bb_top  -   bb_right    -   bb_bottom   -   locate_x    -   locate_y
geneve  -   493485  -   124865  -   506522      -   110870      -   499910      -   117983

The data was used for an older application using OpenLayers 2 which supports a bounding box containing just 4 values, however in OpenLayers 3 - I must create an extent with co-ordinate points but I am unsure how to create the co-ords based on the bounding box data.

Is this even possible to convert ?


texts: QGIS how obtain multicolour lines from polygon cell grid multicolour

texts:
In QGIS I have different lines with same colour that pass through a grid polygon cell (each cell of the grid have a colour).
When a line cross a cell, it must colour with the relative cell colour. Obviously if the same line pass through more than one cell it must colour with multiple colour.

Example:
cell1=yellow, cell2=red, cell3=green
when the line pass through cell1, cell2 and cell3 than the respective part of the line it will colour by yellow, red and green.


texts: postgis How can I obtain directional buffer?

texts: I have a geometry (polygon) which I want to extend with a buffer while testing on it in a where clause. The problem is that I would like the buffer to be only applied according to one or more directions (e.g. North and East) and not in all directions. Is there a way to do it combining operators into a single query? If not, how would it be feasible with a stored procedure?


texts:
For an app I'm developing, I want to use the phone's GPS to locate it. I then want it to list streets that are nearby. The customer might be in his back garden for example, and this may be closer to another street, but not be the street he actually lives on.

I am thinking, I can ask the API several times and skew the coordinates slightly in each request. But that might be expensive in money and time.

Is there an API that can do such a thing in one request?


texts: Positive Topographic Openness - do any FOSS tools (other than SAGA) provide this?

texts:
I'm looking for a FOSS GIS tool to provide Positive Topographic Openness on 32-bit Ubuntu. (which rules out ArcGIS)

Normally I use SAGA GIS for this, but my Ubuntu setup is a bit messed-up after a recent update, and SAGA is crashing on opening (I think it's down to a clash in GDAL versions, but that's a question for another day!)

Are there any other FOSS tools out there which provide this algorithm? I've looked in GRASS and Whitebox GAT but can't find this algorithm.

I understand Relief Visualization Toolbox offers this, but the Linux version only supports 64 bit, and I'm still on 32 bit..


texts: Is it possible to interact with a running instance of QGIS from the command line?

texts:
I just started looking at QGIS on MacOSX. I installed the MacPorts package and found the executable:

/Applications/MacPorts/QGIS.app/Contents/MacOS/QGIS*

Unfortunately it appears that the command line interface can only launch a new QGIS instance and cannot interact with an already running instance. I hope that I'm mistaken. Is there some way to interact (e.g. open files) with a running QGIS instance via a command line interface?


texts: How do I successfully install the LecoS plug-in for QGIS 1.8

texts:
I have been trying to follow the instructions given here, but repeatedly get an error message pop-up titled "Mount" stating Access is Denied.

Does anyone have a solution please?


texts:
Trying to make a simple class that will read shapefiles and eventually be able to pull out the coordinates/attributes for each feature (working with points here). But can't quite figure out geotools.

File file = new File(filepath);
if(!file.exists() || !filepath.endsWith(".shp")) {
    throw new Exception("Invalid shapefile filepath: " + filepath);
}
ShapefileDataStore dataStore = new ShapefileDataStore(file.toURL());
ContentFeatureSource featureSource = dataStore.getFeatureSource();
ContentFeatureCollection featureCollection = featureSource.getFeatures();

SimpleFeatureIterator iterator = featureCollection.features();
while(iterator.hasNext()) {
    SimpleFeature feature = iterator.next();
    List attributes = feature.getAttributes();
    Collection properties = feature.getProperties();
}

Fails at dataStore.getFeatureSource() with the error: java.lang.NoSuchMethodError: org.geotools.data.shapefile.ShapefileDataStore.getFeatureSource() Lorg/geotools/data/store/ContentFeatureSource;

The imports for every geotools class mentioned above are there. As some of you know, geotools comes with a crapton of jars. As I'm unfamiliar with Maven and this is an existing project that I don't know how to convert over, I'm trying to do the library imports the old fashioned way. Right now my imports are: gt-api-11.1.jar, gt-data-11.1.jar, gt-main-11.1.jar, gt-opengis-11.1.jar, and gt-shapefile-11.1.jar. I figure it may be a missing jar that could be the problem but I'm not sure which as every imported class is in one of these jars.


texts: Visual variable renderer doesn’t work on mobile

texts:
I don’t even know where to start with this one. I have developed this simple application: http://matthewiannowlin.com/school-tax-referendum/

If you try this on a desktop, the layers will render correctly, but on a phone the color ramps are not accurate. In the picture the highlighted tract is dark purple, even though it’s median assessed value is $74k and purple is for assessed values of $300k+.

What about the mobile version of this site could change the way the color ramp is rendered? On my end, I have made this responsive through css, so there is nothing different about the script when you view it on a phone.

I’m using ArcGIS JavaScript api 4.9.

EDIT: Here is the code to define and change the renderer.

    var totAssdRenderer = { //Define Renderer: Total Assessed Value
    type: "simple",
    symbol: tractSymbol,
    label: "Neighborhoods (Census Tracts)",
    visualVariables: [{
      type: "color",
      field: "RESULT_1",
      stops: [
      {
        value: 30000000,
        color: "#FFFCD4",
        label: "$30M and below"
      },
      {
        value: 300000000,
        color: "#350242",
        label: "$300M and above"
      }]
    }]
};
...//I establish a few more of these renderers
var tractLayer = new FeatureLayer({//I build a layer using one of my renderers
    url: "https://services.arcgis.com/tKsJAIiLjd90D5q2/arcgis/rest/services/Join_Features_to_SchoolTaxReferendum/FeatureServer/0",
    renderer: medAssdRenderer,
    outFields: ["*"],
    popupTemplate: template,
    opacity:0.5
});
var map = new Map({ //build a map using that layer
    basemap: "gray",
    layers: [tractLayer, districtLayer]
});
changeRenderer = function(chosenRenderer,title,idString) {//on button click, I call this function to apply the renderer defined by the button
    tractLayer.renderer = eval(chosenRenderer);
    document.getElementById("tools").className = "";
    document.getElementById("chooseVariable").innerHTML = title + "Change";
    chosenRendererLabel = title;
    revealMenuToggle = 0;
    var i;
    legendButtonsArray = document.getElementsByClassName("legendButton");
    for (i=0; i < legendButtonsArray.length; i++){
        legendButtonsArray[i].className = "legendButton";
    }
    document.getElementById(idString).className = "legendButton selectedLayer";
};


texts: Get drainage basin of a polygon (with QGIS, GRASS or other FOSS)

texts:
Using a polygon shapefile and a raster DEM, I'd like to find a way to get the drainage basin of a polygon area — that is, the area of all cells which drain into or through the cells contained in the polygon. This is essentially the function provided by r.water.outlet in GRASS except that it only looks at the drainage basin of a single cell.

I figure there must be some way to do this, if only by summing together the drainage basins of all cells within the polygon. However this would require a means of identifying and iterating through all cells contained in a polygon and running r.water.outlet with each cell as a unique input. Ideally I would be able to automate this process with python scripting. Ultimately, I'd like to run it for a large number of polygons in a shapefile layer.

Is this possible?

EDIT:
This work relates to semi-automated assessment of sinkholes, karst features and groundwater susceptibility with GRASS GIS & other FOSS. To illustrate what I am hoping to accomplish here, I've taken some screenshots.

A slopeshade raster overlaid on a map of low-lying areas and depressions (shown in blue, derived with r.neighbors). You can see a number of sinkholes here.

Given a vector layer containing polygon features representing the outlines of sinkholes, I'd like to calculate the geometry of the basin that drains into the sinkhole, which will generally be a very small fragment of the larger watershed basin.

Here is a map showing drainage basins for each sinkhole, calculated manually with r.water.outlet

Here you can see a split basin. The sinkhole contains two cells which are separate drainage endpoints in the DEM. I want to get the 'complete' drainage basin for a sinkhole polygon feature, that is the composite of the drainage basin for all points within the sinkhole polygon.


texts:
I'm using the FeatureTable widget to add a table to map app that I built in the 3.17 JS API. As it's currently set up, the table loads more content as you scroll. My trouble is that my table has tens of thousands of records and the user could effectively be scrolling forever. Is there a way to paginate the table to something like 1000 records per page? Something that looks like the attached image? Can anyone point me to a sample that I can look at that accomplishes this or suggest a methodology?

require([
        "esri/map",
        "esri/layers/FeatureLayer",
        "esri/dijit/FeatureTable",
        "dojo/dom",
        "dojo/parser",
        "dojo/ready",
        "dojo/on"
    ], function (
        Map, FeatureLayer, FeatureTable, dom, parser, ready, on
    ) {
        parser.parse();

    var mapMain = new Map("mapDiv", {
            basemap: "topo",
            center:  [-121, 38],
        zoom: 6,
    showLabels: true
        });

        var myTableLayer = new FeatureLayer("https://data.farwestern.com/arcgis/rest/services/scratch/test_SD_RelatedResources/MapServer/2", {
    mode: FeatureLayer.MODE_ONDEMAND,
    visible: true,
    outFields: ["*"]
    });

        mapMain.addLayer(myTableLayer);

        myFeatureTable = new FeatureTable({
    featureLayer : myTableLayer,
    outFields:  ["*"],
    map : mapMain,
    gridOptions: {
        allowSelectAll: true,
        allowTextSelection: true,
        selectionMode: "extended"
    },
        syncSelection : false,
    }, 'tableDiv');

        myFeatureTable.startup();
    });


texts:
Having problems installing QGIS 2.18 on Fedora 27. Used Copr repo but QGIS simply won't fire up when i try to launch it. Installation seems to go ok without any errors.

Anyone else having issues?


texts: Scale issues with Google Roads XYZ tiles in QGIS layout

texts:
I am trying to include Google Roads as XYZ Tiles in QGIS in a map that I want to export as PDF as a Layout.

Google roads features and the text looks ok on QGIS and in the Layout, however, after I export the map, the text is very small (tiny) in the PDF. I played with the map scale, and when I exported the map with a 1:100 up to 1:1000 scales it is ok, but it is quite small for 1:1000 to 1:3000 scales and almost impossible to read for anything larger. I have tried to change the properties but nothing works for the exported PDF.

I looked for info online, like the following posts

OpenLayers low resolution and/or shifts in QGIS print composer?

Print composer problems Qgis 2.8 Wien with Google Satellite

Problem with print composer when exporting a map with Google Maps Plugin in the background

but most of them are related to issues when the maps are brought through the OpenLayers Plugin or QuickMapServices. I imported Google Roads as XYZ Tiles, as I find it easier, but do not know how to adjust the scales in the output in this case. Does anyone know how to do this?


texts: Is there a 3rd-party Google Maps API for OpenLayers 4?

texts:
Recently, I immigrated to the OpenLayers 4 and I've just wondered if there is any 3rd-party or plugins in order to add Google Maps layers such as Satellite, Hybrid and etc in my project. How should I achieve this?

Why isn't it support anymore?

Edit:

I see the below link in StackExchange. Someone says that we can use leaflet instead of OpenLayers, which isn't my goal and someone else introduced the ol3-google-maps plugin which integrated Google Maps into the OL3. but it's named OL3 and I confused about it. Can I use it with OpenLayers 4 and is it compatible well?

Openlayers v4.0.1 support Google Maps Javascript API?


texts: Use field calculator to calculate geometries with a different SRS than the dataset itself

texts:
I have a polygon dataset containing hundreds of rectangle map 'frames' that represent map extents across the CONUS. The dataset containing the frames is in EPSG:3857, but the frames themselves represent printed map extents in various UTM zones across the country. Each polygon 'frame' feature has an attribute containing the proj4 string for the printed frame, for example EPSG:32614 is:

+proj=utm +zone=14 +datum=WGS84 +units=m +no_defs

I want to calculate the rotation angle of each polygon frame, and I found this excellent post that does just that. But, of course this:

(180/pi())* atan2($y_at(1)-$y_at(0), $x_at(1)-$x_at(0))

calculates the rotation angle of the feature based on the EPSG:3857 SRS. For example, one of the output angles was -55.5 degrees, but when I extract this feature and run the calculation with the dataset defined as 32614, the angle is calculated at -53.9 degrees.

Is it possible to use this proj4 string attribute to calculate the rotation angle of the feature with the EPSG:32614 instead of the EPSG:3857 of the table itself?


texts: How do I access graphic information from an element in an ArcGISDynamicMapServiceLayer? (Flex)

texts:
I've got a Flex application that shows pipeline data on a touchscreen map (for an exhibit). The spec mandates that tapping an element on the map (such as a power plant or substation) will trigger a popup with related data (age of the facility, name of the facility, owner of the facility, etc.)

I'm not sure how to go about this. I've tried multiple approaches including using a MapMouseEvent on my map and trying to catch MouseEvents on the individual layers. I only ever get either the whole map returned or the top layer (not even the layer I'm trying to touch). It's like layers are large fields that act as one single entity and I can't access individual elements of them or other layers that lie beneath them.

Here's the MXML code I use to access the layers from the server:

<esri:Map level="3" wrapAround180="true" id="myMap">
    <esri:center>
        <esri:WebMercatorMapPoint lon="-98.35" lat="39.5" />
    </esri:center>
    <esri:ArcGISTiledMapServiceLayer url="http://server.arcgisonline.com/ArcGIS/rest/services/World_Street_Map/MapServer" id="baseStreet" />
    <esri:ArcGISTiledMapServiceLayer url="http://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer" id="baseImagery" />
    <esri:ArcGISDynamicMapServiceLayer url="http://localhost:6080/arcgis/rest/services/Pennwell_Substations/MapServer" id="Substations" />
    <esri:ArcGISDynamicMapServiceLayer url="http://localhost:6080/arcgis/rest/services/Pennwell_Facilities/MapServer" id="Facilities" />
    <esri:ArcGISDynamicMapServiceLayer url="http://localhost:6080/arcgis/rest/services/Pennwell_Interconnects/MapServer" id="Interconnects" />
    <esri:ArcGISDynamicMapServiceLayer url="http://localhost:6080/arcgis/rest/services/Pennwell_Pipelines/MapServer" id="Pipelines" />
    <esri:ArcGISDynamicMapServiceLayer url="http://localhost:6080/arcgis/rest/services/Pennwell_PowerLines/MapServer" id="PowerLines" />
    <esri:ArcGISDynamicMapServiceLayer url="http://localhost:6080/arcgis/rest/services/Pennwell_PowerPlants/MapServer" id="PowerPlants" />
    <esri:ArcGISDynamicMapServiceLayer url="http://localhost:6080/arcgis/rest/services/Pennwell_HydroelectricPlants/MapServer" id="HydroelectricPlants" />
    <esri:ArcGISDynamicMapServiceLayer url="http://localhost:6080/arcgis/rest/services/USGS_WindTurbines/MapServer" id="WindTurbines" />
</esri:Map>

Here's the Actionscript 3.0 code I'm using to try and access the elements, where log(str:String):void is my equivalent of trace() but prints arguments to a text field in the application.

myMap.addEventListener(MapMouseEvent.MAP_CLICK, function(meEvent:MapMouseEvent):void {
    log(meEvent.toString());
});

I don't want to display the data to an InfoPanel, I just want to be able to access it (for now). I'll be implementing a custom UI for this.


texts: Emulating Identity Tool in ArcGIS Desktop with Basic level license?

texts:
I work at several different companies each with different license levels of ArcGIS.  I got used to using the Identity tool at a company which had access to the ArcGIS Desktop  Advanced license.  However, after switching to another company with only the Standard license, where the Identity tool is not available, I'd like to build a model or script which has the same functionality.  Is there a workflow which will give me the same results as the Identity tool in the Standard license?

I'd like to stay away from third party tools.  I wish to incorporate this into a python script, or ModelBuilder model.


texts: How to add markers to OSM base map in Google Maps API?

texts:
I am using the example from this webpage to display my maps:

http://wiki.openstreetmap.org/wiki/Google_Maps_Example#More_examples

(Example - Using Google Maps API v3 setting OSM as a base map layer)

However, I need to add markers but I get a blank page. See code below:

<!DOCTYPE html>

    <meta charset="utf-8" />
    <meta name="viewport" content="initial-scale=1.0, user-scalable=no" />
    <meta name="apple-mobile-web-app-capable" content="yes" />
    <title>OpenStreetMap with Google Maps v3 API</title>

Please what am I doing wrong..


texts: Publishing error from ArcGIS API for Python?

texts:
I was trying to upload a CSV file on ArcGIS with Python.
I wrote the code :

import pandas as pd
from arcgis import GIS
import os
import datetime
import shutil

URL = "https://www.arcgis.com/"
gis = GIS(url=URL, username="nimish", password="############")
csv1 = '/home/nimish/PycharmProjects/freelancer/project6/csv_files/capitals_1.csv'
cities_df_1 = pd.read_csv(csv1)
print(cities_df_1.head())
print(cities_df_1.shape)

data_pth = "/home/nimish/PycharmProjects/freelancer/project6/csv_files/backup/"
now_ts = str(int(datetime.datetime.now().timestamp()))  # current timestamp

# copy the file, appending the unique string and assign it to a variable
my_csv = shutil.copyfile(os.path.abspath(csv1),
                         os.path.join(data_pth, 'capitals_1_' + now_ts + '.csv'))

item_prop = {'title': 'Sample spreadsheet ' + now_ts}
csv_item = gis.content.add(item_properties=item_prop, data=my_csv)
print(csv_item)

cities_item = csv_item.publish()
print(cities_item)

item_prop = {'title': 'USA Capitals'}
cities_item.update(item_properties=item_prop, thumbnail='data/updating_gis_content/capital_cities.png')

# Uploaded first CSV uptil here

But it displays an error
RuntimeError: Unable to publish item.
User 'nimishbansal' does not have publishing privileges.
(Error Code: 400)

What can I do?


texts: Linear referencing points when the corresponding line ID is known

texts:
I'm trying to use linear referencing on a set of points of which the corresponding line ID is known.

The points are GPS measurements along a railway. Each measurement contains the ID of the railway that is measured. Since the GPS was not very accurate, most measurements are not placed on the railway but a few meters off, or even closer to another parallel railway. I want to reference the measurements to the railway that is predefined per GPS point, and not to the closest railway or line (which is default when using ArcGIS' Locate Features Along Routes tool).

Any suggestions on how to solve this?

I'm working with ArcMap 10.3.1 (Advanced license available). Preferred solution within ArcGIS, but QGIS is also fine.


texts: Adding information from excel into ArcGIS

texts:
I am working in ArcGIS and have a excel sheet with data pertaining to properties.  I want to add this data to a map in ArcGIS of which i have property boundaries as a shapefile?

I have a common folio number between the two datasets.  However I am a novice and all answers have been over my head.  Please point me to a step by step process to inserting databases.


texts: Explaining ArcGIS Raster to Polygon output?

texts:
I want create a new shapefile with selected fields found in a raster which depicts landcover.

I first attempted to run Raster To Polygon.  However, some class names no longer appear, as shown here:

I need help understanding what's going on and the best way to extract ALL the class names into a shapefile.  I'm using ArcMap 10.2.2 and selecting by "CLASS_NAME", Simplify. I have also attempted using "VALUE" and No_Simplify, but get the same results. I added the raster data if that helps


texts: OpenLayers GeoLocation.getPosition() is inconsistent

texts:
In my mobile application I am using OpenLayers and trying to capture current position of the user by using GeoLocation.getPosition() method.

The method sometimes returns position and sometimes return undefined. The behaviour  is strange ! The GeoLocation object has value in it. Below is my code:

Geolocation object creation :

this.GeoLocationObj = new GeoLocation({
        projection: this.map.getView().getProjection(),
        tracking: true,
        trackingOptions: {
            enableHighAccuracy: true
        }
    });

Use of this to get the position:

  var position = this.GeoLocationObj.getPosition();
  this.map.getView().setCenter(position);
  this.map.getView().setZoom(18);
  this.positionFeature.setGeometry(new Point(position) )

What should I change?


texts: Is it possible show the feature attributes dialog when I add a feature in code?

texts:
I want to know if it is possible to make the feature attributes dialog in Qgis show up when I add a feature by code using:

QgsVectorLayer::addFeature( tFeature )

This is available when I add features using Qgis' digitizing tools but not when I add features by code!


texts: Extract raster values gdal looping rasters

texts:
I have this code, most of them getting from this forum. I need to loop over a folder and sub-folders, or just to run a list (putting previously all the raster inside, I don't mind one way or the other).

But I don't find out how to do the loop, I have done it with arcpy in just 30 minutes, but I am fighting with gdal for a couple of days and nothing...

import os, sys

try:
  from osgeo import ogr, gdal
  from osgeo.gdalconst import *
  os.chdir('/home/digd/Desktop/puntos')
except ImportError:
  import ogr, gdal
  from gdalconst import *
  os.chdir('/home/digd/Desktop/puntos')

# open the shapefile and get the layer
driver = ogr.GetDriverByName('ESRI Shapefile')
shp = driver.Open('/home/digd/Desktop/puntos/centroids2.shp')
if shp is None:
  print 'Could not open puntos.shp'
  sys.exit(1)
shpLayer = shp.GetLayer()

# register all of the GDAL drivers
gdal.AllRegister()

# open the image
img = gdal.Open('b5.TIF', GA_ReadOnly)
if img is None:
  print 'Could not open landsat.TIF'
  sys.exit(1)

# get image size
rows = img.RasterYSize
cols = img.RasterXSize
bands = img.RasterCount

# get georeference info
transform = img.GetGeoTransform()
xOrigin = transform[0]
yOrigin = transform[3]
pixelWidth = transform[1]
pixelHeight = transform[5]

# loop through the features in the shapefile
feat = shpLayer.GetNextFeature()
while feat:
  geom = feat.GetGeometryRef()
  x = geom.GetX()
  y = geom.GetY()

  # compute pixel offset
  xOffset = int((x - xOrigin) / pixelWidth)
  yOffset = int((y - yOrigin) / pixelHeight)

  # create a string to print out
  s = feat.GetFieldAsString('Name') + ' '

  # loop through the bands
  for j in range(bands):
    band = img.GetRasterBand(j+1) # 1-based index

    # read data and add the value to the string
    data = band.ReadAsArray(xOffset, yOffset, 1, 1)
    value = data[0,0]
    #print value
    s = s + str(value) + ' '

  # print out the data string
  print s

  # get the next feature
  feat.Destroy()
  feat = shpLayer.GetNextFeature()

# close the shapefile
shp.Destroy()

Uhmm no, I was trying to do a better explanation of my question.

Anyway, I eventually could do it using numpy and GDAL. I needed to get the pixels values within some polygons, so I convert the polygons to rasters and then I could easily perform different masks and calcs.

I find this process really easy with arcpy, but no idea of how to do this with other tools (well, it's easier and faster and better (at least for me) with the rsgislib, but I have to work in a window enviroment and it's not easy to install this library in windows (at least for me).

Again, I don't know of this is an answer or whatever it's (and truly I don't care), you asked me something and I am just trying to give you the best answer that I can.


texts: Difference between "Geographic Transformation" Options Using ArcGIS Project Tool

texts:
In using Arc's Project Tool my Input Feature Class is in the "WGS_1984_Web_Mercator_Auxillary_Sphere" coordinate system.  Once I set the desired "Output Coordinate System" to "GCS_North_American_1983_CSRS" the tool auto-completed the "Geographic Transformation" field to:

NAD_1983_CSRS_To_WGS_1984_2

But I notice there are other options in the dropdown, most notably:

NAD_1983_CSRS_To_WGS_1984_1

What is the difference between these two apparent verions (1 & 2) of transformation?


texts:
I would like to thank's everyone first, i learn a lot here (gis.stackexchange) about GIS. I have one question about the use of GRASS within R session, i'm able to do some basic task, but i have one important issue is to use GRASS graphics through spgrass6::execGRASS.

I'm aware of the capabilities of R in graphics (lattice, ggplot and others), but i would like to be able to use GRASS graphics because of some option like nviz.

Here is a minimal exemple

require(spgrass6)
require(spatstat)
# to use SpatialGridDataFrame in ggplot
fortify.SpatialGridDataFrame <- function(data) as.data.frame(as.im(as.image.SpatialGridDataFrame(data)))

# initialisation and the use of spearfish60 data
initGRASS(gisBase = "/usr/lib/grass64", home = tempdir(), gisDbase = "/home/dickoa/GIS/grass/grassdata/", location = "spearfish60", mapset = "usr1", override = TRUE)

# verif
gmeta6()

execGRASS("g.list", parameters = list(type = "vect"))

parseGRASS("d.mon")

# and this command gave me errors
execGRASS("d.mon", parameters = list(start="x0"))

and the error message is :

ERREUR :Unable to get GIS_LOCK enviroment variable value

And generaly, i really like the R workflow using emacs (ESS) and the idea of keeping everything in script, it's one of the main raison for me to use grass within R,
so i want also to ask for other R and GRASS users if this integration through spgrass6 (Bivand 2011) can allow me to do everything as in interactive session of GRASS.

here is my R session info (i use Kubuntu 10.10 as OS)

R version 2.12.2 (2011-02-25)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
 [1] LC_CTYPE=fr_FR.UTF-8       LC_NUMERIC=C
 [3] LC_TIME=fr_FR.UTF-8        LC_COLLATE=fr_FR.UTF-8
 [5] LC_MONETARY=C              LC_MESSAGES=fr_FR.UTF-8
 [7] LC_PAPER=fr_FR.UTF-8       LC_NAME=C
 [9] LC_ADDRESS=C               LC_TELEPHONE=C
[11] LC_MEASUREMENT=fr_FR.UTF-8 LC_IDENTIFICATION=C

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

other attached packages:
 [1] spgrass6_0.6-26 XML_3.2-0       rgdal_0.6-33
 [4] sp_0.9-79       spatstat_1.21-5 deldir_0.0-13
 [7] mgcv_1.7-4      ggplot2_0.8.9   proto_0.3-9
[10] reshape_0.8.4   plyr_1.4        MASS_7.3-11

loaded via a namespace (and not attached):
[1] lattice_0.19-17    Matrix_0.999375-48 nlme_3.1-98


texts:
I'm writing a script to treat some rasters. One of steps is a reclassification. So I use r.reclass from Grass. But unfortunately, it return a weird error message: IOError: [Errno 13] Permission denied: u'tmp15066743306927' like it's impossible for it to access this temporary file/folder.

Here is my script :

##rasterndvi=raster
##sortieR=output raster

from qgis.utils import iface
import sys, processing, os, glob, numpy
from qgis.core import *
from PyQt4.QtCore import *
from PyQt4.QtGui import *

layer = QgsRasterLayer(rasterndvi)

extent = "%f,%f,%f,%f" %(layer.extent().xMinimum(), layer.extent().xMaximum(), layer.extent().yMinimum(), layer.extent().yMaximum())

print extent

outputs_GRASS7RNEIGHBORS_1=processing.runalg('grass7:r.neighbors', rasterndvi,0,25.0,True,False,'',extent,0.0,None)

layer = QgsRasterLayer(outputs_GRASS7RNEIGHBORS_1['output'])

extent = "%f,%f,%f,%f" %(layer.extent().xMinimum(), layer.extent().xMaximum(), layer.extent().yMinimum(), layer.extent().yMaximum())

print extent

outputs_GRASS7RRECLASS_1=processing.runalg('grass7:r.reclass', layer,'','0 thru 85 = 1\n85 thru 170 = 2\n170 thru 255 = 3',extent,0.0,sortieR)

Somebody can help me ?


texts: Seeking documentation for Netezza Spatial Functions?

texts:
I am looking for some documentation or manual or reference guide for Netezza Spatial Functions.

Any idea where can I get the same?


texts: QGIS 3.6.0 when shutting down hangs with 'Running checks' message

texts: When I have used the Layout Manager to print several different layouts, then go to exit from QGIS (Windows 10), I get a message that the shutdown cannot complete because QGIS is 'Running Checks' and the program just hangs. Eventually, I kill it.


texts: How to calculate the square root of a feature attribute?

texts:
I have a polygon shapefile with over 16k features. Each has an attribute (Double, Precision(10) Scale (2)) that I am trying to calculate the square root of, and store in new field.  Is there a tool, field calculator code block, or script I can leverage?

Example:
Square root of 16 is (4).  (not looking for root mean square).

I see there are functions in Spatial Analyst for square root calculations on rasters, but I can't seem to locate any for shapefiles.


texts: SRTM 30m hillshade output low resolution

texts:
I am trying to hillshade 30m SRTM geotiff. The result is terrible. I know that these rasters normally come with geographical coordinates and i have to reproject it before i apply hillshade because there may be problems with the elevation units. This is the result. The elevation units are the same (meters).
My projected coordinate system is WGS 1984 Mercator. As you can see the resolution is awful.I have done before this step with other rasters and results have been correct.


texts: Leaflet not showing popup depending on zoom

texts:
I've a map with a leaflet popup showing if a click is performed on apolygon.

But for some reason, the click to display popup only works at a certain level of zoom, if it's at a lower level than the one needed nothing happens (well, internally a GET exception is raised about not finding a resource).

camping_vectorgrid.on('click', function(e) {
   displayBase64Image(e.layer.properties.naamnl,e,map);
})
.addTo(map);

With displayBase64Image being:

function displayBase64Image(base64Image, e, map) {
  var img1loaded = false;

  var imagePlaceholder = document.createElement("image-placeholder");
  var image = document.createElement("img");

  image.onload = function() {
    imagePlaceholder.appendChild(image);

    L.popup()
      .bindPopup({ maxWidth: 600 })
      .setContent(imagePlaceholder)
      .setLatLng(e.latlng)
      .openOn(map);
  };

  image.src = base64Image;
  imagePlaceholder.innerHTML = '

Seleccione foto a subir
<input type="file" font-family: monospace; id="inputFileToLoad" align="center" accept=".jpg,.jpeg,.png" onclick="fileClicked(event)" onchange="fileChanged(event);setFile(this);" width="48"/>
<button type="button" onclick="checkandpost();">Subir foto</button>

Fotos ya subidas

';
}

I've checked that the error is due to that:

On displayBase64Image(e.layer.properties.naamnl,e,map);
       })

Depending on the value of zoom e.layer.properties either takes the value it should have, or another one if the image is done with less zoom, one that I'm not even able to find where it comes from.

I honestly don't see anything there that could cause this change of behavior according to the level of zoom.

Any idea of what might be happening?


texts: Howto change OpenLayers BBOX request projection?

texts:
My map uses EPSG:900913 projection on the map.
All coördinates in the database are stored in the 'normal' google format / GPS coordinates (WGS84). For instance, if you save a new location..Before I submit coordinates back to the server they get transformed like

mapPosition = OpenLayers.Projection.transform({
        x: lon,
        y: lat
    }, new OpenLayers.Projection("EPSG:900913"), new OpenLayers.Projection("EPSG:4326"));

The OpenLayers BBOX strategy makes request like:

?bbox=497590.72212427,6707253.6848049,748304.17486475,6814265.5243893,EPSG:900913

How do I change this so the server gets the bounding coordinates in ESPG:4326??
That would be much better then transforming the coordinates serverside.

EDIT
I tried to create a function that appends an extra set of bounds to the request URL:

this.layer = new OpenLayers.Layer.Vector("POI",{
        strategies: [new OpenLayers.Strategy.BBOX({resFactor: 1})],
        protocol: new OpenLayers.Protocol.Script({
            url: 'http://localhost/tests/poi', //this.get_poi_url,
            params: {
                bounds: this.get_poi_url()
            },
            callbackKey: 'jsoncallback',
            srsInBBOX: true
        }),
    });
    tbMap.map.addLayers([this.layer]);

and the function:

get_poi_url: function(bounds) {
    var bounds = tbMap.map.getExtent();
    console.log(bounds);
    bounds = bounds.clone().transform(
        config.projection, config.displayProjection
    );
    console.log(bounds); // OK, it works, I want these coordinates
    var url = ""
        + bounds.top + ','
        + bounds.left + ','
        + bounds.bottom + ','
        + bounds.right;
    return url;
}

The problem is that this set of extra parameters (or if I would build the URL like this) doesn't get an update since it is set with the layer creation. So If you move the map, it would still transmit the original bounds.


texts: I'm attempting to calculate land use mix and in order to do so I need the proportion of each land use within each polygon zone. When i use the tabulate intersection tool I get proportions that are more than 100% not sure what the issue is could anyone help?


texts: Calling SAGA modules through GRASS run_command in a script PYTHON

texts: Does someone have an experience|example of calling modules of SAGA  with the grass.run_command method in a script written in PYTHON?


texts: Legal standpoint of using open source map data in proprietary software?

texts:
I am planning to use an API call to something like http://overpass-api.de/api/interpreter?data=... to gather raw GIS data for an application I am working on.

Can anyone tell me if that affects my ability to keep my code proprietary, or am I compelled to make my code open source (in a similar way to if I was using someone's GPL code, for example)?


texts: I would like to create a bivariate map, and need to create a 2 dimensional legend to accompany it.  How do I create a bivariate map legend (3x3 cells) in QGIS?


texts: Problem of reversed coordinates (XY) when importing WFS data with st_read

texts:
This question was originally posted in StackOverflow.

I am trying to use French spatial data distributed using WFS. I want to use these data in a Shiny app that allows the user to update the data, so I would really prefer to use WFS and not only to embed the data from a shapefile in the app.

I can get the data but there is something off with it in the sense that the coordinates seem to have been reversed (X instead of Y and Y instead of X) resulting in a strange picture of France. I am not very familiar with the details of spatial analyses and data and I am certainly missing something.

I imported the data using the function st_read which I fed with the link to the WFS service. Then I get the upside-down map of French hydrographic network. When I connect QGIS using the same WFS link, however, the map is correct so I guess that sf and QGIS do not exactly handle things in the same way, but investigating why is beyond my knowledge.

When looking at the bounding box we see that the values returned for the x correspond to what we expect for the y and vice versa.

I had the problem with two different data sets distributed using WFS by the same provider (http://www.sandre.eaufrance.fr/). The first one was displaying points and the data contained XY coordinates so I removed the geometry and transformed the data frame using st_as_sf with the names of the X and Y fields in the coords parameter. This is a dirty workaround but it does not work with my second data set that has a MULTILINESTRING geometry. So I would prefer a clean way of importing/transforming the data in the correct format.

Here is a minimal code example to reproduce the unexpected behavior.

library(dplyr)
library(sf)
library(ggplot2)

data <- st_read("http://services.sandre.eaufrance.fr/geo/zon_FXX?Request=GetCapabilities&SERVICE=WFS&VERSION=1.1.0",
                layer = "SegClassContinuiteEco_Liste2",
                stringsAsFactors = FALSE)

st_bbox(data)

ggplot() +
geom_sf(data = data)

The code above gives this map:
wrong map

The correct map should look like that (I rotated and flipped the original):
correct map

My sessionInfo:

R version 3.5.3 (2019-03-11)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows >= 8 x64 (build 9200)

Matrix products: default

locale:
[1] LC_COLLATE=French_France.1252  LC_CTYPE=French_France.1252
[3] LC_MONETARY=French_France.1252 LC_NUMERIC=C
[5] LC_TIME=French_France.1252

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

other attached packages:
 [1] sf_0.7-3        forcats_0.4.0   stringr_1.4.0   dplyr_0.8.0.1   purrr_0.3.2
 [6] readr_1.3.1     tidyr_0.8.3     tibble_2.1.1    ggplot2_3.1.1   tidyverse_1.2.1

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.1       cellranger_1.1.0 pillar_1.3.1     compiler_3.5.3
 [5] plyr_1.8.4       class_7.3-15     tools_3.5.3      jsonlite_1.6
 [9] lubridate_1.7.4  gtable_0.3.0     nlme_3.1-137     lattice_0.20-38
[13] pkgconfig_2.0.2  rlang_0.3.4      DBI_1.0.0        cli_1.1.0
[17] rstudioapi_0.10  yaml_2.2.0       haven_2.1.0      e1071_1.7-1
[21] withr_2.1.2      xml2_1.2.0       httr_1.4.0       fs_1.2.7
[25] generics_0.0.2   hms_0.4.2.9001   classInt_0.3-1   grid_3.5.3
[29] tidyselect_0.2.5 glue_1.3.1       R6_2.4.0         fansi_0.4.0
[33] readxl_1.3.1     modelr_0.1.4     magrittr_1.5     usethis_1.5.0
[37] units_0.6-2      backports_1.1.4  scales_1.0.0     rvest_0.3.2
[41] assertthat_0.2.1 colorspace_1.4-1 utf8_1.1.4       stringi_1.4.3
[45] lazyeval_0.2.2   munsell_0.5.0    broom_0.5.2      crayon_1.3.4


texts: Refreshing combo box in Python add-in of ArcPy?

texts:
I've got an add-in with a combo box and want to refresh the list in the combo box after an event in one of the other tools in the same add-in toolbar.

Can I call the combo box somehow to refresh it outside of the object's onSelChange event?


texts: Calculate values from specific fields for each row in a featureclass

texts:
I have a large script in which I now need to sum the values of a specific list of fields. I have a list of fields:

ave_fields = arcpy.ListFields(dg, "*Ave_STG*")

This list of fields may vary, hence why im using a list, and not explicit field names.

Now, all I need to do is sum the values in these fields for each row and return this value into another field ("weekly_ave").

Ive tried using SearchCursor:

with arcpy.da.SearchCursor(dg, ave_fields) as ucursor:
    total = 0
    for urow in ucursor:
        for i in range(0, len(ave_fields)):
            total += urow[i]
    print total

And ive also tried using NumpyArray to get this, but its not returning the right value:

tbl = [row[0] for row in arcpy.da.TableToNumPyArray(dg,ave_fields)]

This is my first attempt at using Numpy so I may not have got this right.

EDIT:

I have managed to get the numpy array to work properly, however i now need to return this value in the new field

tbl = arcpy.da.FeatureClassToNumPyArray(dg, ave_fields)

print sum(tbl[0]) # this allows me to sum a row within the array

#now i loop through each row and sum the values within the array, and return this value in the weekly_ave field for each row
with arcpy.da.UpdateCursor(dg, "Weekly_Ave") as ucursor:
    for urow in ucursor:
        print sum(tbl[0]) # this successfully sums the row, but i need to do this for all rows

I get the following error:

TypeError: unsupported operand type(s) for +: 'int' and 'numpy.ndarray'


texts: Adding icon in particular places using ArcGIS and show data based on the location

texts:
I am very new to ArcGIS. Currently I am working on a project where I have to show positions of some sensors(placed in fixed locations) which are frequently(every 30 minutes) sending data to a server. I also need to show the data (2/3 values, like: temperature, PH lavel and so on) in map when anybody click to the icon of the sensor or it will automatically show the updated data table.

I have basically almost no idea about ArcGIS. So it would be better if anybody give some step by step instruction or any way so that I can accomplish the task.

Thanks


texts: Clip or extract a raster by row/column number in QGIS or ArcMap?

texts:
I have a georeferenced raster. Instead of extracting or clipping it by a clip feature or co-ordinate, I need to clip it by row and column number.

i.e. I want BX.tif to be columns 2150:2400 and rows 1250:1650 of B2.tif.


texts: Is there automated way of releasing all ArcSDE layers locks?

texts:
We have an issue where users do not close ArcGIS and layers get locked, and sometimes this affects our nightly updates for some of the layers. I know you can individually remove a lock using the ArcSDE command (which we do), but does anyone know of a  way I can create a script that can be run every night to get rid of the locks?

We are running ArcSDE 9.2 and ArcGIS Desktop 9.3.1.


texts: Is it possible to create MBtiles from a TIF?

texts: I have a geo-referenced tif (it's really big around 1GB) and I'd like to convert it to mbtiles so I can serve them and view the map. What is the best tool to do this conversion?


texts:
With QGIS it´s easy to load non-spatial CSV Data (and use it for joins).

This works fine as long the CSV resides on the local file system as this answer suggests. Is there some way to import a CSV from an online resource like https://example.com/myfile.csv without first downloading it? (Maybe by use of some plugin or Python snippet?)


texts: Hack to change data source of a layer file in ArcGIS 9.3.1 with python?

texts: I'm looking for a way to change the source of my layer file using python. Actually I want to do this because I need to re-use a XSL template which is defined under the tap HTML popup of the layer properties dialog. As I can't use the MakeFeatureLayer_management- command to create a layer file that contains the HTML popup information I'm aiming at letting the layer file point to a different data source every time my python code is being executed. I know it is a breeze using arcpy but I'm bound to use the gp-object. This there a workaround anyone knows of?


texts: How can you identify loam/clay settlements with Sentinel-2 Images when using eCognition?

texts:
I am doing my thesis research about analyzing settlements in Mali with Sentinel2 images with the use of the program eCognition (object-based). We have trouble with determining the settlements because they have almost the same spectral values as the surrounding (the houses are from loam, the same as the ground, see picture..).

A few things I already tried / analyzed for inserting conditions for classification:

NDBI
NDVI
Means from all the Sentinel-2 layers
Standard deviations from all the Sentinel-2 layers
Brightness

In the picture beneath, the area and two examples of settlements are visible.

I hope you could give some suggestions.


texts: Project one layer in QGIS using projection CRS of another layer

texts:
I am using QGIS 1.6 and have imported a set of latitude, longitude points from a CSV file as a layer. I also have another layer from a shapefile. When I go to layer->properties for the second layer I see that the "Layer Spatial Reference System" is

+proj=lcc
+lat_1=42.68333333333333
+lat_2=41.71666666666667
+lat_0=41
+lon_0=-71.5
+x_0=200000
+y_0=750000
+ellps=GRS80
+datum=NAD83
+units=m
+no_defs

I would like to reproject the first layer using the second layer's CRS. How can I do this?

When I enable on the fly projection it reprojects the second layer into latitude/longitudes. What I want is the opposite: I want the lat/longs to be projected using the CRS of the second layer. I know that I can choose a CRS when I am saving the first layer but I do not know what CRS I should select.


texts: How to set opacity for GRASS layers using command line?

texts:
In the process of learning GRASS, I am finding that I can learn many of the command-line methods that I need for scripting tasks from using the wxGUI, as it conveniently shows most of the underlying commands that it runs. This has not turned out to be the case for setting the opacity of layers. A little Googling turned up a module called g.pnmcomp, which may be responsible for controlling layer opacity. Unfortunately, that module appears not to be intended for end-users and I have not found documentation that is sufficient for me to make use of it; I'd need to know how to generate the PPM/PGM pairs that it stacks, at least.

There to not appear to be arguments for opacity in d.rast or d.vect, and I suspect that it's the png or tif driver that would actually be handling the opacity of layers. That's the background, but my question is quite simple. How can I set the opacity of layers in maps created via the command line?


texts: Configuring authentication with Google Maps Engine driver from GDAL?

texts:
I'm using last GDAL from trunk, and trying to read from Google Maps Engine (GME) vía GME driver

I already exported

GME_API_KEY: My API Key, for the app I registered
GOA2_CLIENT_ID:  Generated with the "ClientID for installed applications" option in Google OAauth2 API Console
GOA2_CLIENT_SECRET: Idem

Also ran swig/python/scripts/gdal_auth.py script this way:

python swig/python/scripts/gdal_auth.py -s gme

Obtaining a GME_REFRESH_TOKEN, valid for 1 hour.

The problem is I'm getting a HTTP 401 error (Invalid credentials) with requests like this one

ogrinfo "GME:tables={my_table_id}"

Debugging, I've found the Authorization: Bearer HTTP header is incorrect. But I've used the same client id and client secret for QGIS GME plugin and it works.

What should I do to make the driver work?


texts: PostGIS - How do I create a line that directly connects two earth points without travelling over a geodesic?

texts:
In PostGIS 2.1.1, I've noticed that if I create two distant earth points with equal altitude and then create a line connecting them, as I gather interpolated Z coordinate data along the line my Z (height) results remain equal to my original altitude_in_meters.

I've confirmed that it doesn't simply ignore the Z value because I created a line between an sp and ep with different heights, and I get different Z values along the line, indicating slope from sp to ep.  I've also confirmed that these lines are also still curved such that they will always lie above the geoid and never cut through it, which should happen at very great distances with an actual straight line.

This tells me that the line being projected between sp and ep is actually a geodesic.  How can I create an actual straight line connecting these points in 3D space?  My desired result is to see the relative decrease in altitude along the line connecting two distant earth points, which would be most prevalent at its midpoint.

Here's some SQL to show how ST_MakeLine() is creating a geodesic:

WITH
points AS(
 SELECT
   ST_SETSRID(ST_MAKEPOINT(-74.5, 40.5, 50), 4326) AS sp_equal_height
   ,ST_SETSRID(ST_MAKEPOINT(-74.45, 40.45, 50), 4326) AS ep_equal_height
   ,ST_SETSRID(ST_MAKEPOINT(-74.5, 40.5, 100), 4326) AS sp_diff_height
   ,ST_SETSRID(ST_MAKEPOINT(-74.45, 40.45, 50), 4326) AS ep_diff_height
 )
,line AS(
 SELECT
   ST_MAKELINE(sp_equal_height, ep_equal_height) AS line_equal_height
   ,ST_MAKELINE(sp_diff_height, ep_diff_height) AS line_diff_height
 FROM
   points
 )
SELECT
  ST_Z(ST_LineInterpolatePoint(line_equal_height, 0.5)) AS mid_equal_height_z
  ,ST_Z(ST_LineInterpolatePoint(line_diff_height, 0.5)) AS mid_diff_height_z
FROM
  points
  ,line

Results: mid_equal_height_z = 50 and mid_diff_height_z = 75


texts: Exporting PostGIS Geometry POINT as Label text to DXF or DGN

texts:
I need to export, by the ogr2ogr command or similar, text information related to a geometry position from my PostgreSQL(with PostGIS) database to an AutoCAD DXF layer file or to a Microstation DGN layer file, so I can later show that information in the AutoCAD/Microstation GIS.

I have been looking through internet for so many days but all I can get is a layer with geometry points, I always lose the text label information. I know there is the "Text" element type in DGN and the "MTEXT/TEXT" element type in DXF, so I need to generate that kind of features.

In Shapefiles I know that I can export new attribute fields with that information but it have to be obtained in DXF/DGN format.

I need to made it automatic with some kind of server command like GDAL ogr2ogr, but I have tried with the "-select" (getting the text SQL SELECT attribute) or "-lco" or "-dsco" options without any result...

Those are my commands:

ogr2ogr -select "Text" -f DXF text_layer.dxf "PG:host=192.168.1.* user=**** password=**** dbname=****" -sql "SELECT position AS the_geom, name as Text FROM the_geometry_table"

ogr2ogr -select "Text" -f DGN text_layer.dgn "PG:host=192.168.1.* user=**** password=**** dbname=****" -sql "SELECT position AS the_geom, name as Text FROM the_geometry_table"

Any help?


texts: Is there a method to smooth polygons and force all angles to be 90 degrees?

texts:
Is there a way to smooth jagged polygons and get nice 90 degree angles like in the image below? Are any algorithms or methods available?
I've tried GeoPandas' .simplify() function, which uses the  Douglas-Peucker algorithm, but the results aren't as smooth and the quality of results varies depending on the shape of the polygon.
Are there any existing algorithms or approaches to achieve this?

Image Source here.


texts: I'm trying to find a multipoint shapefile to do some testing for an application I'm working on, but I cannot seem to find one (and do not have ArcGIS to create one). Can anyone suggest a place to find a sample multipoint shapefile?


texts: Extracting the neighbouring cluster to Excel

texts: As a beginner in QGIS, I have a question. I have this shapefile which is a map and then divided into clusters. So what I want to do is find the four clusters that share a border with the chosen cluster. So for example, If I pick cluster 91, I would want some kind of automatic process to return cluster 61, 9, 27, and 79. Once, I have these neighbouring clusters then I would want to export the results to excel. Currently, this is a manual process for me but is there a way to automate it?


texts: Determine on which side of a line points are located

texts:
I have a point shapefile (point.shp) and a line shapefile (line.shp). I want to determine on which side of the line each point is located.

The best would be to use a pyqgis command to add a field to the point.shp with value 0 for left and 1 for right.

An alternative that would work is to determine the nearest distance between each point and the line with distance taking negative values for points located on the left and positive value for points located on the right of the line.


texts: Finding distance between point and polyline using ArcGIS Desktop

texts:
I am new to ArcGIS Desktop.

I have two layers in ArcGIS Desktop (i.e. police station (green point as shown in picture) and road networks).

The police station has a shape type point whereas the road network has shape type polyline. The police station has got an ID that consists of two alphabet letters. Road network has also got an ID that starts with two alphabet letters followed by 4 digit number. Road network that comes under the jurisdiction of a police station has the same two alphabet letters as the police station ID has. For example, police station has ID 'AB', so all roads ID that come under the same police station starts with 'AB'.

How do I calculate the distance between the police station and roads that have the same two alphabet letters in their ID?


texts:
I am new to google map api as well as fusion table. I want to make a dynamic query using fusion table. for that i am using the following piece of code:

var layer3
var yr = 2015;
layer3 = new google.maps.FusionTablesLayer({
          query: {

            select: 'Name_Farmer',
            from: '19ZEAwp0-rJNI1alzTRzgjJBqT1SKR3Ij_N_LfAxR',
             where: 'Year > 2015'
// where: 'Year > yr'
          }
        });

    layer3.setMap(map);

If the where clause 'year>2015' works, then why is 'year>yr' not working when I have assigned a global variable for yr=2015?


texts:
I am using 3 feature classes in ArcGIS 10.7

an irregular spaced set of data points.
a set of regular spaced points 100m apart that covers the data area of 1, generated using Fishnet.
a polygon of circles around each of the points generated in 2 (50m radius).

I want to assign the closest data value to each generated point, within the 50m radius. If there is no data point within the radius, I want a null value assigned (not zero - the data can be positive or negative).
I could code this but is there an existing ArcGIS function that does so.
Failures I have tried:

Spatial Join statistics assigned statistics, not closest value.
Spatial Join closest not statistics option claims it assigns the first value it comes across, but appears to have assigned values outside the 50m radius.

The attached pictures show how "Statistics" Join is similar to the data input, but "First" Join has assigned values everywhere in a Nearest Neighbour looking style.


texts:
I am using GPS coordinates and I am trying to create lines (routes) from lon,lat.
I create a geography column from my point and then I create the line as following:

ST_MakeLine(point_geom order by timestamp asc)

Although I have an issue on -180, 180 producing the following:

Does anyone have an idea how I can solve this?


texts: Postgis TopologyException: side location conflict at 226.001 39.5158, what's the problem?

texts:
after st_union() PgSQL returns a note:

NOTE:  TopologyException: side location conflict at 226.001 39.5158

And in the results there are two rows with empty geometry and other rows are ok.
What causes that problem? How to solve it?


texts: How to get the geometry type FIELD value of attribute table?

texts:
The problem is:

I need to copy every cell's value of attribute table of some feature to string array.

I do it by assigning cell's value to string variable StrBuf one by one.

To check the result the MsgBox will show StrBuf value.

Below is a part of code and result MsgBox window:

But when I assigne Shape field's value to StrBuf - the error is: "Run-time error '13'  Type mismath"

My question is which type (instead of String) should I choose for StrBuf variable to avoide this error?

I use arcgis 10.2.2 , scripting in VBA.


texts: Using FME to extract and store data in geodatabase in different feature datasets to create seamless dataset?

texts:
Original data structure as supplied is 124 folders with identical feature data sets:

There is a common attribute in all datasets (metadata) that I need to inspect and extract the data flow for different zones, After manual inspection and extraction, I divided the data in 3 zones (e.g., Z54, Z55 & Z56)

What I want is to extract all features from each folder and place in a feature dataset in a geodatabase as follows:

where Z54_Grids, Z54_Pages etc. etc. are the containers for features extracted as such from individual folder and Z54_UBD_Data_Combined keeps only one dataset for grids from all folders and Pages, Roads, Suburbs & tiles from all folders.


texts: Adding inward facing triangles on a polyline/polygon SLD Geoserver

texts: How do I create inward facing triangles on a polygon/polyline such as below in Geoserver SLD?


texts: Pratical Grid Index for Data Driven Pages

texts:
I'm trying to automatically generate a smart grid index for data driven pages, though frustratingly, the grid index created by the Grid Index Features tool, seems to start the grid at the bottom left corner (see below)

And isn't mart enough to move re-position itself to remove redundant grids, ideally, if it were intuitive it would shift and thus could remove grid one from the equation.

Has anyone come into this issue? Or found a workaround to create a nice grid index? Would love to get your insight to this problem.


texts: 999999 Generic Error when converting CSVs to dBase with TableToTable_conversion()

texts:
I am attempting to convert 5 fairly large CSV files into DBFs using ArcPy.

My code is below.

Currently when I run this script I am met with the following error:

ExecuteError: ERROR 999999: Error executing function.
Failed to execute (TableToTable).

Code below:

# Convert CSVs to DBase tables in preparation for join to FeatureClass

dir_list = os.listdir(some dir)

boro_list = ["Qn", "Mn", "Bk", "Bx", "Si"]

dir_list = [x for x in dir_list if x.endswith("csv") and x[-6:-4] in boro_list]

print(dir_list)

arcpy.env.workspace = "some dir"

for csv in dir_list:
    arcpy.TableToTable_conversion(csv, "some gdb", os.path.splitext(csv)[0] + "_dbf")
    print("Successfully converted" + os.path.splitext(csv)[0] + ".csv to " + os.path.splitext(csv)[0] + ".dbf")

Printing dir_list prior to defining my env workspace yields a list of the 5 desired filenames for the CSVs. When done in ArcMap I am able to convert without problems using the manual version of the TableToTable_conversion tool.


texts: How to change the data type of columns of a vector layer in QGIS?

texts:
Is it possible to change the data type and attributions of data in QGIS?
Almost all of the variables which contains numbers should be changed in all of my vectors.
In ArcGIS the data type is recognised as the right one...

I tried all of the options in the options of the layers, but couldn't find any solution...


texts:
I have a dll and I'd like to call it in Python.  The dll needs some parameters.

Do you know how to pass a radio button parameter and a checkbox parameter?


texts:
I did not find much documentation about GDAL GetAxisOrientation.

I saw an example :

from osgeo import osr, gdal

myproj = osr.SpatialReference()
myproj.ImportFromEPSG(4326)
xaxis = myproj.GetAxisOrientation(None,0)
xaxis = myproj.GetAxisOrientation(None,1)

and the axis orientation for SR=4326 output is 1 3.

Can someone give me an example of the correct use of the GetAxisOrientation and correct meaning?


texts: Add PostGIS database information to XML generated by CartoDB

texts:
I have installed the latest Mapnik renderer following the information on their website and exported a .pbf file to a local PostGIS DB using this command:

osm2pgsql --style ../openstreetmap-carto/openstreetmap-carto.style --create --database gis -U maps -W -H 192.168.220.137 bucharest.osm.pbf

After this I generated a XML style sheet using carto -a "3.0.0" ../openstreetmap-carto/project.mml > mapnik.xml

Now I would like to use this mapnik.xml in a Python script to render a part of the exported map from PostGIS as a test.

#!/usr/bin/env python

from mapnik import *

mapfile = 'mapnik.xml'
map_output = 'mymap33.png'

m = Map(2048, 2048)
load_map(m, mapfile)
bbox=(Box2d( 26.08901,44.420234,26.125488,44.435067))

#m.zoom_to_box(bbox)
render_to_file(m, map_output)

The problem is that the mapnik.xml has no information about my PostGIS DB (I looked in the file and I saw a lot of queries in  element but no host name, user name or password). How do I add more info about my PostGIS DB to this style sheet?


texts: Geoserver output is too small also in WMS

texts:
I encountered a problem in geoserver. My shapefile output is very small even when I output in my WMS using GeoExt tree-legend.

Here is my Layer Preview in Geoserver:

This is my code in GeoExt tree-legend:

Ext.application({
    name: 'Tree Legend',
    launch: function() {
        var mapPanel = Ext.create('GeoExt.MapPanel', {
            region: "center",
            center: [280858.84899, 1594924.75501],
            zoom: 10,
            layers: [
                new OpenLayers.Layer.WMS("TALABA IV",
                    "http://localhost:8080/geoserver/Bacoor/wms",
                    {
                        layers: 'talaba iv',
                        format: 'image/png',
                        transparent: false
                    },
                    {
                        //displayInLayerSwitcher: true,
                        //isBaseLayer: true,
                        //transitionEffect: 'resize'
                        //singleTile: true
                    }
                ),
            ]
        });


texts: How to get extent of raster in QGIS Modeler?

texts:
In my journey to solve a different question, I've decided to create a QGIS processing model to automatically correlate rasters to a projection (UTM) grid.

The inputs are a raster file and the UTM grid (shapefile).

The process is: get raster extent as a polygon, do a join by location with the raster extent on the UTM grid, get the correct zone, and pipe it into GDAL Warp method to output a correctly projected raster.

There is a function in the processing modeler called Raster layer bounds that I can use on the raster input, but cannot use as an input to anything else.

My question is: how do I get from the Raster layer bounds to a polygon in modeler?  I see many tools for making a polygon into other things, but nothing to make a series of points into a polygon.

Perhaps I am going about this wrong and there is something that will help me extract the extent of the raster file more easily, but searching by "extent" yields nothing.  There appears to be a very helpful Polygon from layer extent function, but alas, it is only for vector not raster.


texts:
I have a web application that calls an ArcGIS 10.2.2 secure print service (accesses secured map services). I am calling the print service using a proxy page. Generally this works, but after roughly 30 minutes the print request fails for no apparent reason, and the ArcGIS Server logs contain the following:-

    Error executing tool.: Layer "internal": Unable to connect to map server at
    http://myserver:6080/arcgis/rest/services/Internal/Cadastral/MapServer.
    Layer "jobs": Unable to connect to map server at
    http://myserver:6080/arcgis/rest/services/Internal/General/MapServer.
    Layer "wells": Unable to connect to map server at
    http://myserver:6080/arcgis/rest/services/Internal/Oil_and_Gas/MapServer.
    Failed to execute (ExportWebMap). Failed to execute (Export Web Map).

The error reported back from the http request is:

  Error: Unable to complete operation

The rest of the application works perfectly. If I try printing later (say 20mins), or if I restart the Print Service Tool, it happily prints again for another ~30 mins. It is as if a connection is timing out. Token expiring? Has anyone had this problem before, or can someone suggest where I should direct my attention?


texts:
I would like to save a 3-band raster to file using rasterio.

The official docs have examples for saving a single band image in both the quickstart and in the advanced topics.

What I would like to do is something like..

with rasterio.open('image.tif', 'w', **meta) as dst:
    dst.write(array)

This raises the following error for an array with a shape of (256, 256, 3).

ValueError: Source shape (256, 256, 3) is inconsistent with given indexes 3


texts: Comparable quantification across different maps

texts: I would like the color ramp / bubble size to be comparable over different maps; right now it auto selects the max value as the max in the data set. How can I set the max and min values so that I can compare across maps?


texts: What is the ideal equal-area projection to map Germany, Switzerland and Austria (DACH)?

texts:
I am wondering if there is a standard projection if you want to map the neighboring countries of Germany, Switzerland and Austria (known as "DACH" in equal-area?

For Europe there is EPSG:3035 which seems centered on Germany. Would that be the best choice?


texts: Automating attribute completion of new polygons

texts:
I am currently outlining property boundaries and then copying and pasting in the necessary attributes on a daily basis. I would like to speed up this process and remove the chance of error by introducing some automation.

From my database I am able to export (to an Excel file), all of the records for the coming day. I would like some kind of drop down box or search function which allows me to select the correct record from the Excel file for the property boundary which I have just created. I would like to achieve this using a Python script.

However, I have searched online for related problems or tutorials and can't find any solutions or similar problems. I would like some pointer about where to start with this and how to approach it. Are there any specific ArcGIS or Python tools/modules which can help me achieve this?


texts:
I try use arcpy.JoinField_management(inFeatures, inField, joinTable, joinField, fieldList) to join SQL server database tables to shapefiles directly, using an ArcSDE connection. The database connection from ArcMap is successful.

Following the tip from
How to run “sdetable register” command in SQL Server? at
How to run "sdetable register" command in SQL Server?

I ran the command:

sdetable -o register -t {table} -c {column} -C {database}

on the desktop that ArcGIS 10.3.1 is installed.

However, it gives error message:

  'sdetable' is not recognized as an internal or external command..

I searched "sdetable" in the desktop and didn't find some file. There are a lot of sde*** in the ArcGIS folder, but no "sdetable".

If "sdetable" command is not available for 10.3.1, is there other alternative option to join the tables? I don't prefer to export database tables to csv from the SQL Server, then import each one into ArcGIS to make it have ObjectID.


texts: Merging partially overlapping polygons into a table in PostGIS

texts:
I have a couple of tables (three) with polygons (landuse coverage in this case). Within each table, there is no overlap, but between tables there might be a lot of overlap.

I would like to combine the polygons from the tables into a single table with polygons with no overlap. To resolve situations where the source tables overlap, I would like to put a priority on each layer, so that polygons from the layer with highest priority is used for the overlap.

One way to visualize what I want to do is drawing the tables on top of each other, discarding already drawn polygons where there is overlap. You could also say I want to merge the three tables into one. (Please inform me if there is some more appropriate terminology to describe what I want to do.)

A big issue in this case is performance: each of the tables contain millions of polygons.


texts: Creating histogram for continuous raster dataset in Create Graphs wizard of ArcGIS Desktop?

texts:
I would like to create some histograms of continuous raster data (in TIFF format) to display in a report.

When I enter the Create Graphs wizard, I can only display histograms for feature classes, not rasters. I understand that I can also use 'Create Histogram' function but this does not give me all of the options for visual edits such as colour.

Is there a way I can create and edit a histogram for this raster data? It is a massive amount of data to export and I have had no success trying to make the histograms in Excel or SPSS (they just crash!).