texts:
After switching from 2.14.1 to 2.16.1 adding layers from a certain WFS results in an error, example is with version 1.0.0 forced:

  Layer is not valid: The layer srsname='EPSG:28992'
  typename='app:Bouwvlak'
  url='http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?'
  version='1.0.0' table="" sql= is not a valid layer and can not be
  added to the map

QGIS WFS log:

  Analysis of DescribeFeatureType response failed for url
  srsname='EPSG:28992' typename='app:Bouwvlak'
  url='http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?'
  version='1.0.0' table="" sql=: Cannot find element 'Bouwvlak'

I guess the problem is in parsing the DescribeFeatureType response, but I do find the element 'Bouwvlak':

...
<element name="Bouwvlak" substitutionGroup="gml:_Feature">
    <complexType>
        <complexContent>
            <extension base="gml:AbstractFeatureType">
                <sequence>
                    <element name="fid" type="integer"/>
                    <element name="bestemmingsvlak" minOccurs="0" type="string"/>
                    <element name="historisch" minOccurs="0" type="boolean"/>
                    <element name="identificatie" minOccurs="0" type="string"/>
                    <element name="naam" minOccurs="0" type="string"/>
                    <element name="plangebied" minOccurs="0" type="string"/>
                    <element name="planstatus" minOccurs="0" type="string"/>
                    <element name="typePlan" minOccurs="0" type="string"/>
                    <element name="versieIMRO" minOccurs="0" type="string"/>
                    <element name="verwijzingNaarObjectgerichteTekst" minOccurs="0" type="string"/>
                    <element name="geometrie" minOccurs="0" type="gml:GeometryPropertyType"/>
                    <element name="dossierID" minOccurs="0" type="string"/>
                    <element name="dossierStatus" minOccurs="0" type="string"/>
                    <element name="datum" minOccurs="0" type="string"/>
                    <element name="verwijzingNaarExternPlan" minOccurs="0" type="string"/>
                </sequence>
            </extension>
        </complexContent>
    </complexType>
</element>
...

What I tried:

loading the features with FME - works
manual request in browser - works
force version to 1.0 / 1.1 / 2.0 in QGIS - no luck
ignore axis orientation in QGIS - no luck
Using Wireshark I checked the differences in requests from 2.14.1 and 2.16.1:

2.14.1 requests:

GetCapabilities

http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?SERVICE=WFS&REQUEST=GetCapabilities&VERSION=1.0.0

DescribeFeatureType

http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?SERVICE=WFS&VERSION=1.0.0&TYPENAME=app:Bouwvlak&&REQUEST=DescribeFeatureType

GetFeature (1)

http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?SERVICE=WFS&VERSION=1.0.0&REQUEST=GetFeature&TYPENAME=app:Bouwvlak&SRSNAME=urn:ogc:def:crs:EPSG::28992&MAXFEATURES=1

GetFeature (all)

http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?SERVICE=WFS&VERSION=1.0.0&REQUEST=GetFeature&TYPENAME=app:Bouwvlak&SRSNAME=urn:ogc:def:crs:EPSG::28992

GetCapabilities

http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?SERVICE=WFS&VERSION=1.0.0&REQUEST=GetCapabilities&TYPENAME=app:Bouwvlak&SRSNAME=urn:ogc:def:crs:EPSG::28992

2.16.1 request

GetCapabilities

http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?SERVICE=WFS&REQUEST=GetCapabilities&ACCEPTVERSIONS=2.0.0,1.1.0,1.0.0

DescribeFeatureType

http://afnemers.ruimtelijkeplannen.nl/afnemers2012/services?SERVICE=WFS&REQUEST=DescribeFeatureType&VERSION=2.0.0&TYPENAME=app:Bouwvlak

No differences up to this step, so now I am stuck. Switching back to 2.14 works but I prefer using 2.16.

Can someone tell me this is a bug or a problem in the service?


texts: Qgis2web (Leaflet) - Search Function: Zooming in on search results

texts:
So I just made a nice Leaflet webmap in Qgis2web for work. The map shows the locations of all the projects we have worked on in the last few years as points, including information about which team worked on it, the clients and so on.
I included the search function where you can search different project-numbers, it works fine but I want the search function to zoom in on the search result, so its clear which point on the map is the result of your search (Sometimes two or more projects are in really close proximity). Does anyone have an idea on how to accomplish that?


texts: Detect the point that lies inside line containing area

texts:
I have Two layers  one is point another is line as indicated in picture.

I want to know they way that i can detect that point is contains/intersect/near to the line so that I can add it's attribute data into line attribute data. As you can see that neither line contains the point nor line intersect the point. So how to get that point? Like point under the bounding box of line etc


texts: Merging a .csv file with a .shp file in R

texts:
I am trying to merge a .csv file with a .shp file. The shapefile has more information that I need so I also want to limit it by only adding information for the counties in my csv file. I tried my code based on a response to a similar question posted here but received an error. Once I merge them, I would like to make it a shapefile. Here is my code:

library(raster)

county_info <- read.csv("/path/county_info.csv")
county_shape <- readShapePoly("/path/county_shape.shp")
data <- merge(county_info, county_shape, by="GEOID", all.x=TRUE)
data <- shapefile(data, "/path")

The error I receive is:

Error in (function (classes, fdef, mtable) : unable to find an inherited method for function ‘shapefile’ for signature ‘"data.frame"’


texts: How to convert a double field to a string?

texts: I would like to convert a column in a spreadsheet from a double to a string so I can join the data with a shp file. What is the easiest way to do it?


texts: How to insert images from geodatabase to layout in ArcMap?

texts:
I have a geodatabase downloaded from Survey123 for ArcGIS with 1.783 points with attachments (images) and I want to generate reports from the attribute table of the feature class.

I've designed a report in the layout tab using the Data Driven Pages dynamic text feature, but I have to insert the images in the report now. The images are stored in the geodatabase in a blob field format, and I don't know how to insert the graphical reference of the blob field in the layout.


texts: Expression in field calculator to pick out multiples of x

texts:
I have an attribute table of contours at intervals of 10m but I wish to label the layer at intervals of 50m. My method was to create a new field with an expression which would copy the entries in the parent field (PROP_VALUE) whose value would be a multiple of 50. i.e 50,100,150,200 etc leaving the entry empty (NULL) if not a multiple of 50.

I got somewhere with the CASE expression but QGIS wouldn't let me repeat the expression for 100m, 150, and 200m. I suppose I could update the existing field but there must be a simple/faster way.


texts: Can't create/save changes in shapefiles on QGIS

texts:
I can't create or edit shapes in QGIS, the files are created on my system but cannot "fill" it with the features. I'm using QGIS 2.14.2 64bit Windows 7

This is the error I get when trying to save a new feature to an existing shapefile:

  alriscin_ETRS89: Error de OGR creating the object -2: Error in
  psSHP->sHooks.FSeek() while writing object to .shp file.


texts: PostGIS table shows twice in QGIS Browser panel

texts:
I've observed for a while that for many (but not all) of my tables, QGIS 3 (2.99) shows the same table name twice in the Browser Panel on the left, when I expand my PostGIS data source.

The first entry in the Browser view can be loaded into QGIS normally, but the second instance gives a Coordinate Reference System Selector dialogbox and once loaded shows no data at all.

Does someone know why this happens and is it possible to not show the 2nd entry?

Related question: Duplicate entries in table list when adding PostGIS layer in QGIS said that PostGIS tables are duplicated in Table list, but

  everything's fine in browser panel


texts:
I have an input POINT-layer with lots of features in it. What I am trying to do is to calculate the azimuth for two features. And the way I am doing this is the following:

geomFrom = (...).geometry()     # this variable holds the geometry of the first point (the starting point of the imaginary line which connects both points.)
geomTo = (...).geometry()       # this variable holds the geometry of the second point (the end point of the imaginary line which connects both points.)

# Now I am trying to get the Point out of the geometry for calculating the azimuth.
# the method asPoint() causes the Minidump in C!
# Although I store it in a separate variable!!
fromPoint = geomFrom.asPoint()
toPoint = geomTo.asPoint()

# If the asPoint() method wouldn't cause the minidump crash, I would now calculate the azimuth like this:
azi = fromPoint.azimuth(toPoint)

In order to solve the problem, one could probably say: Upgrade to QGIS 3!

But I have to solve this in QGIS 2.18. Otherwise, I will have to rewrite the whole algorithm in QGIS 3 again as the python scripting has changed fundamentally.

Therefore, I am looking for a workaround. How can I calculate the azimuth without calling the asPoint() method?
What came in my mind: Is there a way to get the x and y-coordinates of a geometry variable? Maybe I can then create a new point variable like this:

from_X = geomFrom.getXCoord()
from_Y = geomFrom.getYCoord()
fromPoint = Point(from_X, from_Y)

to_X = geomTo.getXCoord()
to_Y = geomTo.getYCoord()
toPoint = Point(to_X, to_Y)

azi = fromPoint.azimuth(toPoint)

Any suggestions?


texts: Adding Vertical Reference in Python to ArcGIS feature dataset

texts:
I am creating a python script to generate a terrain. In the process I am creating two feature datasets. However, in order to successfully create this FD, i need both horizontal and vertical datums. I am able to code the horizontal datum successfully, but unable to input the vertical datum. Although the code does not error out, it does not create the requested feature dataset. I used the ESRI WKID's to reference horizontal and vertical datum as shown below.

Here is the code line for spatial reference that I used:

sr = arcpy.SpatialReference(horizontal_wkid)

I tried doing this sr = arcpy.SpatialReference(horizontal_wkid,vertical_wkid) but to no avail.  I am using ESRI ArcGIS 10.2 and coding in Python 2.7


texts: Generating polygons with shared boundaries using QGIS?

texts:
Essentially the problem is that if you digitise each polygon separately the boundaries will not match exactly. My instinct is to have a table of line segments and generate the polygons from that using postgis but then all the line segment end points have to join exactly which is painful as there are a lot of line segments.

I have thought about writing a script that uses heuristics to match up lines and then add extra segments to bridge the gaps as suggested in a previous question I asked: How to automatically convert line features into polygons using QGIS/PostGIS?

I feel this is such a common task that there must be some "best practice". If I were starting from scratch then how would I do it?


texts: MetadataImporter_conversion lags on first iteration of a loop

texts:
I have a Python script that is looping through a list of feature classes stored in a config file.  It is exporting the metadata for each feature class using the methodology in the answer posted to Creating table containing all filenames (and possibly metadata) in File Geodatabase?

The script runs as expected, however, when performance testing I'm finding that the first time the MetadataImporter_conversion command is executed, it is taking approximately 60 seconds to complete.  On subsequent iterations the same command takes less than a second to complete.

I have tried re-ordering the feature classes within the config file to test if the issue is being caused by a particular feature class, but it had no impact.  The first feature class in the list takes 60 times longer to process than the ones that follow.

Ultimately this script will be incorporated into a Geocortex workflow, so I'm trying to make it as efficient as possible so that the user doesn't have to hang around waiting for the results to be returned.

This is a snippit of the code I'm using:

# ---------------------------------------------------------------------------
def CreateDummyXMLFile():
    filepath = r"e:\temp\tmpXML.xml" # This path is for testing purposes only
    with open(filepath, "w") as f:
        f.write("<metadata />")

    return filepath
# ---------------------------------------------------------------------------
# ---------------------------------------------------------------------------
    def GetElementText(tree, elementPath):
    """Returns the specified element's text if it exists or an empty
    string if not."""
        element = tree.find(elementPath)
        return element.text if element != None else ""
# ---------------------------------------------------------------------------
# ---------------------------------------------------------------------------
def GetMetadataElementTree(dataset):
    """Creates and returns an ElementTree object from the specified
    dataset's metadata"""
    xmlfile = CreateDummyXMLFile()
    arcpy.MetadataImporter_conversion(dataset, xmlfile)
    tree = ElementTree()
    tree.parse(xmlfile)
    os.remove(xmlfile)
    return tree
# ---------------------------------------------------------------------------

for aIntersect in ListOfFeatureClasses:
    xmlfile = CreateDummyXMLFile()
    tree = GetMetadataElementTree(aIntersect)
    aTitle = GetElementText(tree, "dataIdInfo/idCitation/resTitle") # Title
    aPubDate = GetElementText(tree, "dataIdInfo/idCitation/date/pubDate") # Publication Date
    print "      Title: {}".format(aTitle)
    print "      Date: {}".format(aPubDate)
    del tree

Processing times with baam feature class as the first item in the list:

Processing times with esa_c feature class as the first item in the list:

Can someone explain this behaviour or replicate it?

Testing on my personal PC running arcgis-10.6 has not been able to replicate the issue.  Would be interested to know if this issue is related to arcgis-10.3 or perhaps just isolated to my work PC.


texts:
I'm quite new to OpenLayers and much of the documentation is just strings of jargon, so please forgive me if this is a noobish question.

I've intialised an OpenLayers map with Google physical and satellite layers plus a vector layer for users to add polygons. The init code is as follows:

//define extent
var miny=-2000000; var minx=6000000; var maxy=600000; var maxx=10003750;

        // allow testing of specific renderers via "?renderer=Canvas", etc
var renderer = OpenLayers.Util.getParameters(window.location.href).renderer;
renderer = (renderer) ? [renderer] : OpenLayers.Layer.Vector.prototype.renderers;

//define the layers needed
    g_phys  =   new OpenLayers.Layer.Google("Google Physical",{type: google.maps.MapTypeId.TERRAIN, numZoomLevels:20}, {isBaseLayer:true});
    g_sat   =   new OpenLayers.Layer.Google("Google Satellite",{type: google.maps.MapTypeId.SATELLITE, numZoomLevels:20},{isBaseLayer:true});
    boxes  =    new OpenLayers.Layer.Vector("Boxes", {rederers:renderer});

//create new set of option for the map
var options = {
            projection: new OpenLayers.Projection("EPSG:3857"),
            units: "m",
            maxResolution: "auto",
            maxExtent:  new OpenLayers.Bounds(miny, minx, maxy, maxx),
            theme: null,
            layers: [
                //define the base layers for the map
                g_phys,
                g_sat,
                boxes

                    ]

        };

//create new map object based on "map" div
map=new OpenLayers.Map('map',options);
map.zoomToMaxExtent();

However, the map allows a maximum zoom such that the window shows about the length of Canada and one of my clients would like to easily add polygons that span larger than this distance, so needs to zoom out further.

I've tried changing the maxextent values and this doesn't seem to help. Could anyone help explain to someone not overly familiar with Geography how to allow the user to zoom out further - ideally to view from the Arctic to Antarctica at maximum zoom.


texts: Display two different types of attributes, in same dataset, as pie charts in different colours?

texts:
Is there a python script or something (maybe there is a simpler way) to display data of the same shapefile as a pie chart but as different coloured pie charts?

For example I have road links and intersections in the dataset, with the same 7 different attributes columns which I use to make the pie charts, but instead of creating two separate shapefiles I want them to display the data in the same way but in different colours. So blue colour palette pie chart for intersections and yellow colour palette pie chart for links and no fiddling around with "why is a 5000 count pie chart smaller in links but the 2000 count pie chart in intersection is bigger" - would like the consistency.


texts:
So I have a shapefile of a vector for Philadelphia's city lines and a raster image of Philadelphia county that I got from the lovely people from USDA:NRCS. I converted the raster from a sid file to TIFF (using gdal_translate) and then reprojected the newly converted raster TIFF to 3857 (using gdalwarp).

This new raster now has the same native projection as the vector for the city limits: EPSG:3857. So I load both the shapefile for the city limits and this new raster tif that I prepped into QGIS and I fire up that useful clipper tool (Raster->Extraction->Clipper) and I go about my business of setting the raster as the source and selecting the vector as my mask layer and press "OK". #thatwassimple. The computer starts chugging along doing its thing and I'm getting pumped about my new raster coming in.

Ding! computer finished... only South Philly and a tiny part of Central Philly are extracted out--West, North, and East Philly seem to have ducked out of the party and it's all black.

Hm... that's weird, maybe I can use GDAL command line to force this, because you know, if it's done on the command line, I am a real hacker and it's not like QGIS is already doing this in the background, but I digress.

So I cmd those mofos (using gdalwarp again!) 0...10....100! DING! Alright, alright, alright, as I McConaughey myself over to the directory that I outputted the clipped TIFF. Open that sucker up and I'm with the exact same output of Southside and Central only left to represent.

I've attached two images, the first shows a rose colored polygon that I want to cut out of the raster, and the second shows the output both times.

Photo 1:

Photo 2:

SOURCE FILES

Shapefile used (in zip)

TIFF used (in zip)

UPDATE:
I noticed that the clipped areas are actually in fact regions outside of the city and just cuts north of the city. (Note south philly no longer has Schuylkill River going through it). I tried making a polygon of this so I could upload an image for clarity, but I couldn't figure it out. Please just replicate the raster clipping with the files provided above and you too can see what I'm talking about.

This leads me to believe that there must be a projection problem, however, after I gdalinfo-ed the raster and consulting the .prj of the city_limits file, I've confirmed that they have the same aforementioned projection. #sayhwatt


texts: Displaying OpenStreetMap in locally hosted web page?

texts: I'd like to set up my own tile server but I'm very new to it all and do not know what to do next. I followed the guide Switch2OSM on how to install OSM locally. Surprisingly I successfully got through it after a few wrenches thrown in. Now I have no idea how to display my map in a locally hosted web page. I saw the tutorial on leaflet through Switch2OSM but its not for a local server and is very uninformative to a beginner. Any help on where to go next would be great.


texts:
How can I have two sublayers on at once in carto.js? I'm very much struggling on how to accomplish this. My code is here, I have gotten so far as to hide and display each sublayer with a checkbox, but cannot get them to appear together.

Relevant code:

cartodb.createLayer(map, vizjson)
                .addTo(map)
                .on('done', function(layer) {

                    var subLayerOptions = {
                        sql: "SELECT * FROM final_1",
                        cartocss: "#final_1{polygon-fill: #F84F40; line-color: #000;}",
                        //interactivity: "name"
                        }

                    var sublayer = layer.getSubLayer(0);

                    sublayer.set(subLayerOptions);

                    sublayers.push(sublayer);

                    //alert("this is a map");

                });

            $('#bhp').on('change', function() {
                if($(this).is(":not(:checked)")) {
                    $('p1').html("something") } else {
                    $('p1').html("else")};
                    });

            var LayerActions = {
                bhp: function(){
                    sublayers[0].setSQL("SELECT * FROM final_1 WHERE mms_num= 2277 ");
                    sublayers[0].setCartoCSS('#final_1{polygon-fill: #e65400; line-color: white; }');
                    return true;
                },
                noLayer: function() {
                    sublayers[0].setSQL("SELECT * FROM final_1 WHERE 1=0");
                    sublayers[0].setCartoCSS('#final_1{polygon-fill: #e65400; line-color: white; }');
                    return true;
                },

                bp: function(){
                    sublayers[0].setSQL("SELECT * FROM final_1 WHERE mms_num= 2481");
                    sublayers[0].setCartoCSS('#final_1{polygon-fill: #009400; line-color: white; }');
                    return true;
                }
            };

            $('#bhp').on('change', function() {
                if($(this).is(":not(:checked)")) {
                    LayerActions['noLayer']()} else {
                    LayerActions[$(this).attr('id')]()}
                });

            $('#bp').on('change', function() {
                if($(this).is(":not(:checked)")) {
                    LayerActions['noLayer']()} else {
                    LayerActions[$(this).attr('id')]()}
                });


texts:
With a lot of help from gis.stackexchange.com I've created leader lines for point symbols. The lines are based on a geometry generator symbology:

Shortening line geometry using expression dialog in QGIS?

Now, I want to implement leader lines for moved labels. It's very importent to me, that the lines start at the center of the labels and end at the point coordinates (or vice versa).

-- Geometry generator:
make_line(
    make_point(label_x, label_y),
    $geometry)

Once again, the lines are covered by the labels:

Is there a way to address the label boxes and shorten the leader lines? I suppose this is quite a tricky question :/


texts: How to install Fiona to read Shapefile attributes with OSGeo4W?

texts:
I´m developing a plugin in Python and I need to read the shapefile attributes, so I found the Fiona API and dbfpy also but I don't understand how to install these packages because I install the QGIS through osgeo4W.
How do I do that?


texts: How is IValidate.GetInvalidRules supposed to work?

texts:
we are currently trying to programmatically  validate features against geodatabase-defined rules. This includes attribute rules as well as connectivity and relationship rules.

We have consulted the documentation for IValidate, IValidation and the Validating features article on EDN to learn how it's done.

After setting up the appropriate connectivity rules, invocation of the Validate Features in ArcMap correctly yields the number of invalid features in a selection. If I select a single feature, invocation of this command correctly displays what is wrong about that particular feature in a message box.

Now, calling IValidate.Validate on that feature works as well. Its return value indicates that the feature is indeed invalid, and the error message is identical to that displayed by the Validate Features command in ArcMap. However, calling this method has two disadvantages which are problematic for our purposes:

It performs the validation in a short-circuit manner. As soon as one of the assigned rules is invalid, it stops evaluating the remaining rules. We need to evaluate all the rules so that we can display the violations at once.
It provides the error information as a string. For our purposes, we would ideally like to examine the information in a more structured way, which means we need to access the IRules that turned out to be invalid.

So, this is where IValidate.GetInvalidRules comes into play. The documentation for this method states:

  The GetInvalidRules method returns an enumerator of rules that are currently violated by the row.

which is exactly what we are after. However, the enumeration returned by this method does not seem to ever contain anything. What is particularly weird is that:

Calling IValidate.Validate on a feature yields the appropriate validation error message as expected.
Calling IValidate.GetInvalidRules on the very same feature always returns an empty enumeration as if no connectivity or relationship rules were invalid. (It does not matter whether the call is made within an edit session or not.)

How can I validate a row against all the assigned geodatabase rules at once? Has anyone used IValidate.GetInvalidRules successfully?

UPDATE: IValidate.GetInvalidRules seems to return invalid rules which pertain to attributes. However, invalid network connectivity rules are still never returned from this method.


texts: ArcGIS Enterprise weblayer sharing automation with WMS/WFS capabilities enabled using Python

texts:
I have to create web services on ArcGIS Enterprise (via ArcGIS Pro) that offer WMS/WFS capabilities from several hundred datasets coming from Enterprise geodatabases. While this is not a problem when being done manually (I simply check the respective check boxes when using "Share as a Web Layer" - see image below), I am not able to find a way how to do that with Python.

In Python I have been using getWebLayerSharingDraft, exportToSDDraft, arcpy.StageService_server and arcpy.UploadServiceDefinition_server. This works pretty well, however, I only end up with a Map Image Service and Feature Service and have no idea how to get WMS and WFS. There seems to be no option to also create those services.


texts: ArcGIS Network Analyst Service Area Default breaks in Meters

texts:
I'm trying to use a meter value for the default break. The default value in length is Kilometers, and I can't figure out how to change that.

Additionally, the field doesn't allow a dot or a comma.

I need to create an catchment area for 600m, 1200m, 1800m.

How would I do this?


texts: PostGIS Geometric V Geography measurements when using ST_Length

texts:
I was measuring geometry using some PostGIS functions and noticed that I was getting different results from the two methods I used.

sum(st_length(st_transform(geom, 3857))) # Gave 3674

sum(st_length(geom::geography)) # Gave 1671

Why am I seeing different results?  I was expecting a slight deviation but this seems to be quite significant.

Is there something wrong with how am I using these queries?


texts: I have created webapplication in arcgisonline where I built a form fill widget. My question is how I can add this widget or activate this widget on popup via link i.e if user click on feature then popup appears and user can able to fill the form in popup only. No need to go form fill widget on the map.


texts: Legend Not Refreshing When Renderer Changes

texts:
I have the following code to grab a value from a drop down menu and use it to change the renderer of the featureLayer. However, it does not update the legend. The legend just stays on the initial, default renderer at page load.

I have tried adding a callback function after the switch that would legend.refresh(). That did indeed refresh the legend every time I changed the renderer, but oddly it was "one behind," always changing the legend to reflect the previous renderer instead of the current.

getRenderer = function(rendererFromSelect) {
 renderer = rendererFromSelect.value;
 switch(renderer) {
  case "approval_rate":
  map.removeLayer(featureLayer);
  var generateRenderer = new GenerateRendererTask(featureLayerURL);
  var params = new GenerateRendererParameters();
  var classDef = new ClassBreaksDefinition();
  classDef.breakCount = 5;
  classDef.classificationField = "approval_rate";
  classDef.classificationMethod = "quantile";
  params.classificationDefinition = classDef;
  generateRenderer.execute(params, applyRenderer);
  function applyRenderer(renderer) {
   featureLayer.setRenderer(renderer);
  };
...other cases...
};

I have also tried this:

(1) Create a function to build the legend
(2) Call this function when the body loads
(3) Every time a new renderer is selected via drop down, legend.destroy()
(4) Then, after the renderer ids changed and the feature layer loads, create legend again.

So here is my function to create the legend. (I declare it as a global variable and define it inside the scope of the AMD so I call access it from the HTML.) This works great.

createLegend = function(titleStr){
  legend = new Legend({
    map: map,
    layerInfos : [ {
      layer : featureLayer,
      title : titleStr
    } ]
  }, "legend");
 legend.startup();
};

And here is where I change the renderer, load feature layer again, destroy legend and recreate after feature layer loads.

getRenderer = function(rendererFromSelect) {
 renderer = rendererFromSelect.value;
 switch(renderer) {
  case "approval_rate":
   map.removeLayer(featureLayer);
   var generateRenderer = new GenerateRendererTask(featureLayerURL);
   var params = new GenerateRendererParameters();
   var classDef = new ClassBreaksDefinition();
   classDef.breakCount = 5;
   classDef.classificationField = "approval_rate";
   classDef.classificationMethod = "quantile";
   classDef.colorRamp = multiRamp;
   params.classificationDefinition = classDef;
   generateRenderer.execute(params, applyRenderer);
   function applyRenderer(renderer) {
     featureLayer.setRenderer(renderer);
   };
   legend.destroy()
   map.addLayer(featureLayer);
   map.removeLayer(highlightLayer);
   map.addLayer(highlightLayer);

   break;
   ....other cases...
 }
  featureLayer.on("update", createLegend(renderer));
}

When I try this approach, the legend is created successfully when the page loads initially. The legend.destroy() is also working properly. However, when createLegend is called again to rebuild the legend with the new renderer, I get init.js:148 Uncaught TypeError: Cannot read property 'firstChild' of null

So I'm wondering if legend.destroy() destroys something that I need to create the legend again later?

EDIT: Indeed, legend.destroy() was eliminating the destination div for my legend. So I now recreate that div before I build the legend. Now the legend refreshes almost correctly, but it is still refreshing based on the PREVIOUS renderer.

I tried calling createLegend on "update-end" instead of "update". That seemed to make no difference.


texts: How to create irregular shape grid in QGIS 3.0?

texts: As I mentioned in my question, I want to create a grid which I defined the row and column widths. However, difference between them are not equal. The shape example in my mind is:


texts: Loading WFS point layer with OpenLayers3?

texts:
I am trying to get a point layer from GeoServer to load in my web app as a WFS. I can get a sample set of polygon data to load, but my points will not. Here is the code I am using for this:

var geojsonFormat = new ol.format.GeoJSON();

var vectorSource = new ol.source.Vector({
  loader: function(extent, resolution, projection) {
    var url = 'http://........?service=WFS&' +
        'version=1.1.0&request=GetFeature&typename=Store:layer&' +
        'outputFormat=text/javascript&format_options=callback:loadFeatures' +
        '&srsname=EPSG:3857&bbox=' + extent.join(',') + ',EPSG:3857';
    $.ajax({url: url, dataType: 'jsonp', jsonp: true});
  },
  strategy: ol.loadingstrategy.tile(ol.tilegrid.createXYZ({
    maxZoom: 19
  }))
});
window.loadFeatures = function(response) {
  vectorSource.addFeatures(geojsonFormat.readFeatures(response));
};

var centerlocations = new ol.layer.Vector({
  source: vectorSource,
  style: new ol.style.Style({
    stroke: new ol.style.Stroke({
      color: 'rgba(0, 0, 255, 1.0)',
      width: 5
   })
  })
});

Is there something different that needs to be done for loading point data?


texts: Convert polygons into spatial points to use R "redlistr" package

texts:
I want to use the the makeEOO function from the redlistr package to calculate the Extent of Occurrence of 2 species ("binomial" A and B),  starting from a shapefile with the corresponding distribution maps in form of polygons. To use the function, I need to convert the polygons into points.
Here the starting shapefile.

I am using the sf package to import the shapefile

 library(sf)
 dummy <- st_read(dsn = "Calculate_EOO_dummy/dummy_maps_dissolved.shp")
 plot(dummy)

I converted the polygons into multipoint features, and the multipoints into points

multi.p <- st_cast(dummy.diss, "MULTIPOINT")
my.points <- multi.p %>%  group_by(binomial) %>% st_cast("POINT")
plot(my.points)
my.points

> my.points
Simple feature collection with 33 features and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: -0.9471488 ymin: -0.9276773 xmax: 1.155772 ymax: 0.8080668
epsg (SRID):    4326
proj4string:    +proj=longlat +datum=WGS84 +no_defs
# A tibble: 33 x 2
# Groups:   binomial [2]
   binomial                 geometry
   <fct>                 <POINT [°]>
 1 A          (-0.5855355 0.7913769)
 2 A          (-0.5771905 0.7774687)
 3 A          (-0.2267038 0.5771905)
 4 A         (-0.3936022 0.05702364)
 5 A         (-0.1098748 -0.1015299)
 6 A          (0.02642559 -0.340751)
 7 A        (-0.02920723 -0.4464534)
 8 A         (-0.3824757 -0.4297636)
 9 A         (-0.8553547 -0.2878999)
10 A         (-0.9471488 -0.1182197)
# ... with 23 more rows

However, if I try to call the makeEOO function with my.points, I get an error:

eoo <- my.points %>% group_by(binomial) %>% makeEOO()

Error in makeEOO(my.points) : trying to get slot "coords" from an object (class "sf") that is not an S4 object.

How could I convert such polygons (separately for each species A and B) into points to use makeEOO?


texts: OpenLayers re-projection of vector source GeoJSON

texts:
OpenLayers insists on drawing them over null island, what am I getting wrong?

            new ol.layer.Vector({
                title: 'vector',
                visible: true,
                source: new ol.source.Vector({
                    format: new ol.format.GeoJSON({
                        defaultDataProjection: 'EPSG:27700',
                        featureProjection: 'EPSG:3857'
                    }),
                    url: function(extent) {
                        return url;
                    },
                    strategy: ol.loadingstrategy.bbox,
                    crossOrigin: 'anonymous'
                }),
                style: new ol.style.Style({
                    stroke: new ol.style.Stroke({
                        color: '#ff9900',
                        width: 1
                    })
                })
            })


texts: Interpolating LiDAR intensity images in R?

texts:
I'd like to generate 0.5 m intensity images using lidR, but the resulting raster is full of pits. Here's a reproducible example:

library(lidR)
LASfile <- system.file("extdata", "MixedConifer.laz", package="lidR")
las <- readLAS(LASfile)
int <- grid_metrics(las, list(imean = mean(Intensity)), 0.5)
plot(int, col = height.colors(50))
summary(values(int))

Can the pit-free algorithm, or some other interpolation solution in lidR, be used to improve the result? It works great in the generation of CHMs with grid_canopy().


texts: Seeking GIS information on worldwide cases of wildfires

texts:
I am developing some case studies on wildfire cases around the globe, some examples include California (2007), Australia (2009), Moscow (2010), Valparaiso in Chile (2014), and so on.
Where could be a good start for gathering layers regarding the urban areas affected, infrastructure, nature features, and specifics of the fires such as perimeters?
Any tips on GIS libraries or sites from these cases?


texts:
I have two ortho images with having gap due to inproper cutting of images its overlay.

How to remove this gap in ArcMap 10.1?


texts: How do I embed the a Leaflet map into my web page, using QGIS2Leaf?

texts: I am generally new to web developement, and web gis for that matter, but I want to embed a leaflet map that I made in QGIS using the QGIS2leaf pluggin, into the html of my webpage, but am having troulbe getting it to work.


texts: Line and point feature data is not returned by getGetFeatureInfoUrl method after some time

texts:
I have web based GIS application with multiple layers.
Layers gets rendered as image from GeoServer.

Now upon clicking anywhere on the map data of overlapping layers (including line and points) shows up, However after 1-2 minutes upon clicking on map only data of overlapping Polygons is returned and not of lines and points.

If I refresh the web page then again data of lines and points show up on click.

System behaves irregular in returning data for point and lines ,I couldn't understand why ! Please suggest.

Below is my code to get the data on click :

 map.on('singleclick', function (evt) {

    var isAccordion = !!$("#popup-content").data("ui-accordion");
    if (isAccordion) {
        $("#popup-content").accordion("destroy");
    }
    var viewResolution = (view.getResolution());
    var url = multipleSrc.getGetFeatureInfoUrl(
      evt.coordinate, viewResolution, 'EPSG:3857',
      { 'INFO_FORMAT': 'application/json', 'FEATURE_COUNT': 50 }); //text/html

    if (url) {
        var parser = new ol.format.GeoJSON();
        var lookup = {};
        $.ajax({
            url: url,
            dataType: 'json',
            method: 'GET',
            success: function (response) {

                var result = parser.readFeatures(response);
                //if ((result.length > 0) && (wmsLayer.getVisible() == true)) {
                if (result.length > 0) {

                    var coord = evt.coordinate;
                    var content = "";
                    //content = "
";
                    for (var i = 0; i < result.length; i++) {
                        var Feature = result[i];
                        //content = content + "
";
                        // content = content + "
";
                        content = content + "
" + Feature.getId() + "
";
                        content = content + "
";
                        content = content + '
';
                        content = content + visitResultProperties(Feature.getProperties(), 1, "");

                        content = content + '
';
                        content = content + "
"; // According tab div
                        // content = content + "
";
                    }

                    content_element.innerHTML = content;
                    $("#popup-content").accordion({ collapsible: true });

                    overlay.setPosition(coord);

                }
            },
            error: function (jqXHR, exception) {
                if (jqXHR.status === 0) {
                    alert('Not connect.\n Verify Network.');
                } else if (jqXHR.status == 404) {
                    alert('Requested page not found. [404]');
                } else if (jqXHR.status == 500) {
                    alert('Internal Server Error [500].');
                } else if (exception === 'parsererror') {
                    alert('Requested JSON parse failed.');
                } else if (exception === 'timeout') {
                    alert('Time out error.');
                } else if (exception === 'abort') {
                    alert('Ajax request aborted.');
                } else {
                    alert('Uncaught Error.\n' + jqXHR.responseText);
                }
            }
        });
    }
});


texts: Raster area calculation using python script?

texts:
I am a beginner to python programming in GIS and am trying to calculate the area of a raster.

When I ran the program, I got error message below.

I am also providing my code.

import arcpy
from arcpy import env

#To overwrite output
arcpy.env.overwriteOutput = True

#Set environment settings
env.workspace = "C:/Subhasis/Test/raster-read"

#set local variables-STI extraction
inRaster ="sti"

# Create a search cursor for raster attribute
read = arcpy.SearchCursor("inRaster","","","Value;Count","")

#read the counts in raster
for row in read:
    count = row.getValue("Count")

    print count

Error Message:

"Traceback (most recent call last):
  File "C:/Subhasis/Test/raster-read/read-raster.py", line 14, in <module>
    read = arcpy.SearchCursor("inRaster","","","Value;Count","")
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\__init__.py", line 1167, in SearchCursor
    return gp.searchCursor(dataset, where_clause, spatial_reference, fields, sort_fields)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\geoprocessing\_base.py", line 359, in searchCursor
    self._gp.SearchCursor(*gp_fixargs(args, True)))
IOError: "inRaster" does not exist"


texts: Reference Error: selectZIP is not defined

texts:
In Firebug I keep getting this error when I try to click a button I made.

<button type="button" onclick="selectZIP()">Execute</button>

I declared a function called selectZIP in a .js file and put it in my HTML file. Yet am still having this error. SelectZIP is a geoprocessor which tries to...select ZIP codes.

function selectZIP(){
    gp = new esri.tasks.Geoprocessor("******/arcgis/rest/services/EDDM/selectZIP/GPServer/selectZIP");
    myMap.selectZIP.setOutputSpatialReference({
        "wkid": 102100
    });

In my .js file I have this, before my button statement

<script src="js/EDDM_API - Copy.js"></script>

I don't know where to go from here.


texts: Constructing graph from geometries and tieing points to graph after graph creation using PyQGIS network analysis

texts:
I am solving some Euclidean shortest path problems, thus I create my lines/arcs myself. Following code from PyQGIS cookbook I need to import a vector layer to start building a graph:

vl = qgis.utils.iface.mapCanvas().currentLayer()
director = QgsLineVectorLayerDirector(vl, -1, '', '', '', 3)
properter = QgsDistanceArcProperter()
director.addProperter(properter)
crs = qgis.utils.iface.mapCanvas().mapRenderer().destinationCrs()
builder = QgsGraphBuilder(crs)
pStart = QgsPoint(-0.743804, 0.22954)
tiedPoint = director.makeGraph(builder, [pStart])

As far as I understand the code the director is the one actually initializing the graph and according to the API documentation there is no other way to do so. I was wondering if this is correct or if I miss something. (It would be kind of cumbersome to first create a layer to then load it, if I already have all geometries.)

Also as far as I understand it I can only tie points to the graph while initializing it with the director. This would be pretty expensive, if the start/end points are variable depending on my Euclidean shortest path end point (I want to start a normal shortest path analysis from the end point of the Euclidean shortest path with an end point on a line).

Do I miss something?

I found QgsGraphBuilder().addEdge/Vertex and QgsGraph().addEdge/Vertex but am quite unsure what they exactly do (and if they are able to do what I want) and what input they expect.

I am using QGIS 2.8.1.


texts: Using domain coded values to define arcpy SQL query variable based on field.type

texts:
I am attempting to export selected features from many feature classes within a series of geodatabases into new ones - maintaining schema and all that. The trouble is that each geodatabase has different "ownership" field types. Some are coded value integers, while the rest are strings. Of those that are coded values, the values aren't even consistent!

I am trying to create a script that will reliably search for the coded value descriptions, then apply the corresponding value to the SQL query in the FeatureClasstoFeatureClass_conversion function. It works on some feature classes, but I am getting an invalid SQL query error on others.

env.workspace = item
for ds in arcpy.ListDatasets():
    listfc1 = arcpy.ListFeatureClasses(feature_dataset=ds)
    print listfc1
    for fc in listfc1:

        fcname =os.path.split(fc)[1]

        ownfieldlist=[]
        for field in arcpy.ListFields(fc,"*owne*"):
            if field.name != "DOWNELEV":
                ownfieldlist.append(field)
            else:
                pass
        if len(ownfieldlist) == 0:
            naughtylist.append(fc)
        try:
            ownfield = ownfieldlist[0]
            print ownfield.name + " {}".format(ownfield.type)
        except:
            print (fcname+" has no ownership field")
        if fc not in naughtylist:
            if ownfield.type == "String":

                arcpy.FeatureClassToFeatureClass_conversion(fc, os.path.join(output,itemname), fc,
                                                            """{} = 'USU Eastern'""".format(arcpy.AddFieldDelimiters(fc,ownfield.name)) )
            elif ownfield.type == "SmallInteger":
                number = getthenumber()
                arcpy.FeatureClassToFeatureClass_conversion(fc, os.path.join(output,itemname), fc,"""{} = {}""".format(arcpy.AddFieldDelimiters(fc,ownfield.name), number) )

        else:
            pass

The script above is just adapted to work within a single GDB, but normally, it refers to a function that creates a list of databases and iterates over those.

The getnumber call refers to another function, which passes the coded value which matches the ownership description I want. Like I said, it works for the first geodatabase, but not the second.

Here is the getnumber function:

def getthenumber(item):
    import arcpy
    from arcpy import env
    import os
    from arcpy import da
    item = r"D:\Eastern02\Electric.gdb"
    env.workspace = item
    for domain in arcpy.da.ListDomains(item):
        print domain
        if domain.name == "AssetOwner":
            coded_values = domain.codedValues
            answer = "USU Eastern"
            for val, desc in coded_values.items():
                if desc == answer:
                    thenumber = int(val)
                    print ("{} is represented by {}".format(desc,val))
                    return thenumber
        else:
            pass


texts:
I have a shapefile and I want to share it via email. I noticed that ArcMap saves this file as a bunch of little files with various types (.cpg, .dbf, .prj, .shp, etc.).

Is the best way to share a shapefile to zip a folder containing all these little files and send it that way?


texts:
I want to add a query layer to an arcmap document so that whenever i have  a db data update, it is reflected in my ma and I can update my service through python script. This is my script so far...
import arcpy
from arcpy import env
arcpy.env.overWriteOutput = 1

# define local variables
wrkspc = 'C:/ArcGIS_Workspace/'
mxd = arcpy.mapping.MapDocument(wrkspc + 'scripts/Trial.mxd')

input_db_name = "Database Connections/xxx.sde"
outLayer = "DynamicQueryLayer"
out_layer_name = wrkspc +"/shapefiles/DynamicQueryLayer.lyr"
query = "SELECT * FROM TEST_TABLE"

#unique key for the table
oid_field = "schema.TEST_TABLE.TEST_ID"

try:

    # Create a query layer
    arcpy.MakeQueryLayer_management (input_db_name, outLayer, query, oid_field)

    # save the created query layer
    arcpy.SaveToLayerFile_management(outLayer, out_layer_name, "ABSOLUTE")

    # Get the created layer
    lyr = arcpy.mapping.Layer(out_layer_name)

    # get the data frame of the given mxd
    data_frame = arcpy.mapping.ListDataFrames(mxd)[0]

    # Switch to data view
    mxd.activeView = data_frame.name

    arcpy.mapping.AddLayer(data_frame, lyr, 'TOP')
    print "Added layer to the MXD ..."

    mxd.save()
    print "MXD saved..."
except Exception, e:
    # If an error occurred, print line number and error message
    import traceback, sys
    tb = sys.exc_info()[2]
    print "Line %i" % tb.tb_lineno
    print e.message

I get this error

ERROR 000732: Input Features: Dataset DynamicQueryLayer does not exist
or is not supported.

I am sure there is a way to save the created query layer as a layer file. Any help is greatly appreciated.


texts: using changeset for larger region to update OSM data in PostGIS

texts:
I'm using the OSM import and update workflow described here to maintain OSM data in a PostGIS database.

Import (osm2pgsql):

osm2pgsql -c -d gps -U gps --cache 8000 --number-processes 4 --slim --flat-nodes europe_nodes.bin europe-latest.osm.pbf

Update (osmosis):

osmosis --read-replication-interval workingDirectory=/osmosisworkingdir/ --simplify-change --write-xml-change - | osm2pgsql --append -d gps -U gps --cache 1000 --number-processes 4 --slim --flat-nodes europe_nodes.bin -

I want to gradually import new countries into my database while keeping all the data up to date. To avoid the hastle of creating osmosis working directories for each new country and managing more and more separate updates, I would like to use changesets for the entire region.

Example: Germany and France are imported into the database using osm2pgsql. Changesets for all of Europe are then applied daily using osmosis.

Will this work properly, considering the changeset contains irrelevant changes to data I don't have? Will this approach have any undesired consequences?


texts: Clipping DEM / DTED data selecting by extent

texts:
similar to a few queries but really struggling to clip a DEM file from say USGS.

Going to Raster > Extraction > Clipper and then "select the extent by drag on canvas", does not work.

The process will continue for up to 1 hour then give error messages implying not enough storage?

I'm on QGIS 2.8 Tien, running a MacbookPro on El Capitan... didn't work prior to OS upgrade though from Yosemite.

Any ideas? Tried out thinking this and tried to clip after extracting contours, and then trying to clip by extent but this isn't possible at moment I think unless a tool is developed. My project is to Just make a file at 10km square of contours to input into Google Earth.


texts:
I want to print a map to scale on 11x17 paper. I understand you can change the scale bar to 'numeric' in print composer, and then set the scale under properties for the map.

What I want to know is what units those are, and how to change them. I want to print at 1:4,000 (one inch = 4,000 feet), but when I do that it is wildly off, showing something closer to 1:300(?). Everything is in EPSG 2229, which I understand to be in feet. My page properties have the paper size set for A3.

What am I missing?


texts:
When i am trying to display the Marker Cluster on Google Maps with the below code it is trhrowing some errors. Please find below code and errors.

<apex:page standardController="Account">
    <apex:includeScript value="{!URLFOR($Resource.MarkerCluster)}"/>

    <meta name="viewport" content="initial-scale=1.0, user-scalable=no"/>
    <meta charset="utf-8"/>
    <title>Marker Clustering</title>

 </apex:page>

Errors:


texts: Looping for importing and interpolating a lot of LIDAR files?

texts:
I'm quite new to Grass Gis and I'm trying to make DEM (digital elevation model) for modelling flooding. I have 170 *.xyz files, which need to be scanned for extents (r.in.xyz -s input=...\filename.xyz output=filename fs= ' '), then region needs to be set with previously recieved data (g.region n=... s=... e=... w=... res=2), importing lidar data (r.in.xyz input=...\filename.xyz output=filename fs= ' ') and at last interpolating that data (r.surf.idw input=filename@Siim output=filenameidw).

I did it successfully with one file, but it would be waste of time to do it manually for 170 files. Can anyone recommend me a loop, which I can use?


texts: ESRI JSAPI: How to use GMap3 as 1st layer when using gmaps-utility-gis.googlecode.com?

texts:
This question relates to ESRI JSAPI and Google Maps V3 integration...

I've put together this fiddle http://jsfiddle.net/br4pk/1/ which demonstrates a mismatch between the LOD of the map and gmaps.GoogleMapsLayer instance.  The map operates great until you introduce layers which specify an lod[0] with a level:0.  I believe this is because the gmaps layer lod[0] looks like this:

{
    "level": 1,
    "resolution": 78271.5169639999,
    "scale": 295828763.795777
}

and the map lod[0] like this:

{
    "level": 0,
    "resolution": 78271.5169639999,
    "scale": 295828763.795777
}

Should the map lod contain a level:0 when the basemap does not specify a level:0 lod?
What is the solution to getting this to work with other layers which do specify a level:0?

I did find some success in modifying the lods of the gmap by adding a level 0 detail, and perhaps that is the answer (to modify http://gmaps-utility-gis.googlecode.com/svn/tags/gmapslayer/1.0/src/gmapslayer.js):

    { "level": 0, "resolution": 156543.033928, "scale": 591657527.591555 },
    { "level": 1, "resolution": 78271.5169639999, "scale": 295828763.795777 },
    { "level": 2, "resolution": 39135.7584820001, "scale": 147914381.897889 },
    { "level": 3, "resolution": 19567.8792409999, "scale": 73957190.948944 },
    { "level": 4, "resolution": 9783.93962049996, "scale": 36978595.474472 },
    { "level": 5, "resolution": 4891.96981024998, "scale": 18489297.737236 },
    { "level": 6, "resolution": 2445.98490512499, "scale": 9244648.868618 },
    { "level": 7, "resolution": 1222.99245256249, "scale": 4622324.434309 },
    { "level": 8, "resolution": 611.49622628138, "scale": 2311162.217155 },
    { "level": 9, "resolution": 305.748113140558, "scale": 1155581.108577 },
    { "level": 10, "resolution": 152.874056570411, "scale": 577790.554289 },
    { "level": 11, "resolution": 76.4370282850732, "scale": 288895.277144 },
    { "level": 12, "resolution": 38.2185141425366, "scale": 144447.638572 },
    { "level": 13, "resolution": 19.1092570712683, "scale": 72223.819286 },
    { "level": 14, "resolution": 9.55462853563415, "scale": 36111.909643 },
    { "level": 15, "resolution": 4.77731426794937, "scale": 18055.954822 },
    { "level": 16, "resolution": 2.38865713397468, "scale": 9027.977411 },
    { "level": 17, "resolution": 1.19432856685505, "scale": 4513.988705 },
    { "level": 18, "resolution": 0.597164283559817, "scale": 2256.994353 },
    { "level": 19, "resolution": 0.298582141647617, "scale": 1128.497176 }


texts: How to show statistically significant incident locations

texts:
First let me preface with I tried to search for this solution first but I don't really know how to word it and could not find any results so excuse me if this is a duplicate. I would like to know how to prove or test to show if this problem is statistically significant or not.

I am using Arc 10.3.1 and have a feature class of incidents and a feature class of roads.  I queried out the major roads and buffered them to find out how many incidents fall within that road buffer.

I have approximately 6000 miles of roads. Of that 6000 miles, 900 miles are a major road. I also have 2010 incidents. Of those incidents 725 occurred on a major road. (ie 15% of the roads account for 35% of the problems) I would like to know of tutorial or guide that would show me how to prove or not prove if the incidents on a major road are significant or not.


texts: Quantitative evidence supporting Solid State Drive (SSD) adoption for ArcGIS Desktop workstations?

texts:
Has anyone seen quantitative evidence supporting the use of a SSD (solid state drive) in an ESRI ArcGIS Desktop workstation? Or any GIS workstation for that matter?

I keep reading that SSDs can be a huge performance boost but I have not seen any numbers to support the claim. I have a feeling that the productivity increase could be comparable to a switch from single to dual monitor setup.


texts: Obtaining geographic points with CSV file

texts:
I have a xlsx file with a column called CoordX and another with CoordY, I save the file in CSV format (semi-colon) and use the delimited text plugin and follow all the information in the text box.

But QGIS doesn't find the data, which are "recognize" as belonging to the first line.

What can I do?


texts: Will ActiveSync 4.1 work with ArcPad 10.0.2?

texts: Because Microsoft's ActiveSync 4.5 will not work on machines with Windows XP (SP3), from what I've read I need to install ActiveSync 4.1.  Has anyone tested this on ArcPad 10.0.2?


texts: Buffers between center points of irregular polygons

texts:
I'm having some problems, I hope some of you have some experience with buffers :)

I have the road network of a city (40.000 polygons approx), and I am trying to create a buffer between the central point of the line between each segment of the road that crosses (perpendicularly). My goal is to create a buffer that has a different area on each side of it (having the polygon inside the buffer). This is because the road network is not uniform, so I cannot create a buffer with a specific radius that means the same for each and every road segment.

Do you have any suggestions?


texts: Performing unsupervised k means classification in QGIS?

texts:
I'm running QGIS 2.16 on Windows 10.

I'm looking to perform an unsupervised K means classification on a geoTiff image I have from a multispectral camera surveying sugar cane fields.

Alas, the documentation for the algorithm in QGIS needs work done to make the code user friendly to use in the Python console.

Parameters
Input Image [raster]

Available RAM (Mb) [number] Default: 128

Validity Mask [raster] Optional.  ## I don't know what this means.

Training set size [number] Default: 100  ## Or this

Number of classes [number]

Default: 5

Maximum number of iterations [number]Default: 1000 ##Don't know either

Convergence threshold [number] Default: 0.0001

OUTPUTS:

Output Image [raster]   ## I imagine it's the layer with the classification done
Centroid filename [file]   (Not exactly sure, I know what centroids are but                     ##those work for polygons in vector files not raster files).

I'm thinking maybe it's because I don't have the OrfeoToolbox installed but I´m positive I do because I downloaded and installed QGIS 2.16 along with GRASS and the OsGeo4w packages.


texts:
I'm very new to open street map so apologies if this is a daft question!
I'm trying to find all fast food locations in an area. When I pull the data using osmdata in R there are a lot of rows with no information at all. How are these classed as fast food locations?
Here's a reproducible example:
library(osmdata)

t = opq('Adur') %>%
  add_osm_feature(key = 'amenity', value = 'fast_food') %>%
  osmdata_sf()
View(head(t$osm_points))

The locations with OSM IDs 2439442254, 2439828308, and 2439828379 all have no information. Why are they classed as fast food amenities?
Thanks!


texts: PostGIS geography conversion not working

texts:
I am new to PostGIS and QGIS. I am trying to convert geom to geog and then cast it back to geom in order to buffer my streets by 1,000 ft with this query:

SELECT a.gid, (st_buffer(a.geom::geography,304.8)::geometry) as geom, a.street_ext
FROM public.orangeburgroads a
JOIN public.orangeburg b ON st_within(st_centroid(a.geom), b.geom)
LIMIT 1000;

And here is what I am getting:

Does anyone have a clue why I am getting this result?

FYI orangeburg is a municipal boundary and orangeburgroads is a street center-line file.


texts: Gdal_grid gives blank .tiff file after conversion

texts:
After using gdal_grid I get blank .tiff file. I use the following .vrt file :

<OGRVRTDataSource>
    <OGRVRTLayer name="p">
      <SrcDataSource>p.csv</SrcDataSource>
      <GeometryType>wkbPoint</GeometryType>
      <GeometryField encoding="PointFromColumns" x="field_1" y="field_2" z="field_3"/>
    </OGRVRTLayer>
</OGRVRTDataSource>

and my data is in following format :

318600.00,323600.00,640.00
318575.00,323525.00,639.73
318550.00,323475.00,644.10
318625.00,323575.00,645.10

I use the following command parameters for gdal_grid :

$ gdal_grid -l p p.vrt p.tiff
Grid data type is "Float64"
Grid size = (256 256).
Corner coordinates = (313814.744746 323793.349824)-(319085.945254 318715.120176).
Grid cell size = (20.510508 19.759648).
Source point count = 18922.
Algorithm name: "invdist".
Options are "power=2.000000:smoothing=0.000000:radius1=0.000000:radius2=0.000000:angle=0.000000:max_points=0:min_points=0:nodata=0.000000"

0...10...20...30...40...50...60...70...80...90...100 - done.

This is the resulting .tiff file

EDIT

This is the info returned from gdalinfo

**Driver: GTiff/GeoTIFF
Files: p.tiff
Size is 256, 256
Coordinate System is `'
Origin = (313825.000000000000000,318725.000000000000000)
Pixel Size = (20.510507812500009,19.759648437499891)
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  313825.000,  318725.000)
Lower Left  (  313825.000,  323783.470)
Upper Right (  319075.690,  318725.000)
Lower Right (  319075.690,  323783.470)
Center      (  316450.345,  321254.235)
Band 1 Block=256x4 Type=Float64, ColorInterp=Gray**


texts: Add time restriction in ArcGIS Network Analyst

texts:
I'm trying to use train tracks as a network in ArcGIS Network Analyst.

We are trying to create a model for security agents on board our trains and we'd like to optimize their daily "route". For now, it's all made by hand.

My approach is to take our gtfs data and create a network using stop times to determine when the train segment is "traversable". When there's a train on the segment, it's open, if not, it's restricted.

The problem is I can't set a restriction based on time in Network Analyst (or wasn't able to). I would like to be able to say by a formula or script than you can take that segment only from 7:03 AM to 7:10 AM for example.

Is it possible?

Thanks!


texts: Adding elevation data to existing vector layer using QGIS?

texts:
I'm new to qgis and I'm having trouble adding elevation data to my existing vector file.

I currently have two files:

a road network with the all attribute information and
a vector file with point elevation data.

What I want to do is add an elevation column to the attributes table of my road network.  How do I go about doing this?  I've looked for an answer, but I'm having trouble finding one.

The relative height difference between the start and end point of the line would be the most useful, but the average height would also work.


texts: Diffusion/Kernel Interpolation with barriers?

texts:
I need to interpolate water quality data in a rather small area with peninsulas, islands and developed coastline that will act as barriers between my testing/sampling points.

I understand I have to use non euclidean method to do so, but does anyone know the difference between "Diffusion Interpolation with Barriers" and "Kernel interpolation with barriers"?


texts: Is there a way to Alter Elevation Values of Pixles of a DEM based on BM and Control point elevation from actual Field survey?

texts: We have got a quick bird DEM of our site. After getting actual on ground elevation values with a differential GPS device, we have found that there is a difference of average 12 feet between DEM and Actual ground Points. nw i want to alter elevations of that DEM based on actual values.


texts: How do I write my equation with conditions into the QGIS raster calculator?

texts: I have two raster layers - one with slope(degrees), called theta, and one with slope(%). I need to calculate the S-factor by using attatched formula. How do I write this in the raster calculator? Is it possible?


texts: Create a reference number (id) based on x coordinate in QGIS

texts:
I would like to add a reference number to a point table. I would like this value to increase as I move from west to east across my data. So the values read 1 to n from west to east. I have an attribute table on my composer which will contain this reference number, and the name of each point.

The reason for this is so that it is easier to find the points by the reference number on a map. I want to use this reference number as a label for each point. At the moment the reference numbers follow no particular sequence.

Any tips on how I can achieve this in QGIS?


texts: Removing ZM From Geometry Column Polygon Using ArcPy?

texts:
Using ArcGIS 10.2 for Desktop, I have some Shapefiles which I have no idea why their Geometry column (shape) contain ZM! now I need to remove all those ZM and get a clean polygon instead.

I thought I can use the FeatureClassToShapefile_conversion() but not sure how to disabled both Output has Z Values and Output has M Values through following ArcPy

# Import system modules
import arcpy
from arcpy import env

# Set environment settings
env.workspace = "C:/data"

# Set local variables
inFeature = ["climate.shp"]
outLocation = "C:/output"

arcpy.FeatureClassToShapefile_conversion(inFeatures, outLocation)


texts: Google Earth Engine fails to export variable

texts:
I am exporting zonal statistics for 0.01-degree tiles in various locations in Africa. When I try to find the tile's 2010 percent forest cover from the GFCC data, I use a reducer function:

var gfcc_2010 = ee.ImageCollection('NASA/MEASURES/GFCC/TC/v3')
      .filter(ee.Filter.date('2010-01-01', '2010-12-31'))
      .select('tree_canopy_cover').mosaic();

//Mean pct FC per tile 2010
var mean_2010 = tiles.map(function(feature) {
  return feature.set(gfcc_2010.reduceRegion({
    reducer: 'mean',
    geometry: feature.geometry(),
    scale: 30
  }));
});

When I print mean_2010 in GEE, I can see the variable produced, called tree_canopy_cover. But when I export this feature collection as a table:
Export.table.toDrive(mean_2010, "gfccMeanGroup1_2010", "deforTilesGEE");

The resulting table only has the tile ID and geocoordinates. Why isn't tree_canopy_cover exporting? Full code here. (I use 0.5-degree tiles in this demo code so that I can print in the terminal).


texts: Using File Geodatabase coded domain with ArcGIS Online?

texts:
I created a File GDB and set up a domain, and the domain was assigned to a field in a feature class. The domain behaves as expected on ArcGIS Desktop (dropdown list when editing the content of that field).

But if I zip up the FGDB and upload to ArcGIS Online to created a feature layer and edit the layer online, there is no dropdown when editing.

I don't see any additional step in this tutorial, other than a different mechanism of getting data to ArcGIS Online. I don't think I can do it the same way in ArcGIS 10.0:
http://doc.arcgis.com/en/collector/android/create-maps/prepare-data-desktop.htm

Anyone has any idea what's going on?


texts: Stop rounding values in Qgis attribute table

texts: I am calculating the areas of my polygons in Qgis using $area and everything works fine until I click on "save". All my values are rounded to an integer. I tired setting the output field type to decimal but still see the same problem? Any solution to that? I am using Qgis 1.8


texts: How can I draw custom shapes on Graphics layer?

texts:
Hello everyone I need help. Any clue, guidance or help that would be appreciated, please bear with me as I am still struggling with Flex 4.

What is the simplest way to draw custom shapes using < s : path> on graphic layer. I had data for the  already created.
How can I pass my own styles in SimpleMarkerSymbol(SimpleMarkerSymbol.STYLE_CIRCLE, 22, 0x009933) in place of STYLE_CIRCLE?
I want to create my own styles so that I can explicitly use them there like STYLE_RECT.
Can I pass my path data in any symbol like fillsymbol, linesymbol or markersymbol so that they draw the shape?
How can I extend these symbol classes so that I can use member functions and properties of these symbols + my MXML code?

I am using SFV 2.4.
Thanks


texts: How to capture feature selection event using FeatureSelectedEventHandler in MapXtreme?

texts:
I am trying to implement functionality, where I need to run a piece of code when selection of a feature(point feature selection with mouse click) happens.
I am using PointSelectionCommand (i.e., javascript clientcommand) and PointSelection(i.e., C# class command) provided in MapXtreme SDK. The code I tried is below:

public event MapInfo.Tools.FeatureSelectedEventHandler OnFseh;
private void selectedEvents(object sender, FeatureSelectedEventArgs e)
{

}

private void Page_Load(object sender, System.EventArgs e)
{
    this.OnFseh += new FeatureSelectedEventHandler(selectedEvents);
}

the problem is control is not going into the selectedEvents method when I select a point feature on map.


texts: Using 3d Polyline in QGIS to calculate Volume?

texts:
This question may be completely off the mark. I have a survey drawing in DWG format, I use AutoCAD 2d but the survey was created in LSS.  The Countour lines in the DXF have XYZ coordinates.

Firstly how do i use this layer (i can save it out seperately) in QGIS and how do I display it?

Secondly, is it possible to run a plugin or something (i dont know python) to calculate the volume of this area....

thirdly, i know that the maximum depath of the area is 7m below the highest point. How to create a new layer with the a z coordinate 7 m below maximum height...eg. max height is 110 therefor the z would be 103.

Lastly is any of this possible and am i completely wrong in assuming this can be done?

Any pointers or advice is appriciated.


texts: arcpy script tool: creating a point feature from user input XY coordinate

texts:
This part of my tool is giving me issues.
I want to create a new feature class based on the user's input workspace, then generate a new point feature within that feature class based on the user's input XY coordinates.

So far this is what I have.
My script tool calls for the input DEM (Parameter 0) (this is needed for a later process).
an XY coordinate (Parameter 1)
and an output workspace (Parameter2)

The input XY coordinates are split, space delimited.
From here the Xcoord and Ycoord can be read in as type floating point and used to generate the point feature.

Problem is: Everytime, no matter what coordinate inputs used, the point is placed at

-81.511127  37.576795 Decimal Degrees

Even if I enter 0,0 for coordinates...

Is there something I am missing?
How can I create a new point feature using the Xcoord and Ycoord input from the tool?

Here is my code:

import arcpy
from arcpy import env
from arcpy.sa import *
from arcpy.da import *

env.overwriteOutput = True

arcpy.CheckOutExtension("3D")
arcpy.CheckOutExtension("Spatial")

DEM = arcpy.GetParameterAsText(0)

CoordsIn = arcpy.GetParameterAsText(1)

Workspace = arcpy.GetParameterAsText(2)

CoordSplit = CoordsIn.split(" ");

Xcoord = (CoordSplit[0])
Ycoord = (CoordSplit[1])

arcpy.AddMessage(Xcoord)
arcpy.AddMessage(Ycoord)

PointName = "pourpoint.shp"
GeomType = "POINT"
Template = ""
Has_m = "DISABLED"
Has_z = "DISABLED"

SpRef = arcpy.Describe(DEM).spatialReference

arcpy.CreateFeatureclass_management(Workspace, PointName, GeomType, Template, Has_m, Has_z, SpRef)

fc = Workspace + "/" + PointName

cursor = arcpy.da.InsertCursor(fc, ["SHAPE@XY"])
xy = (float(Xcoord), float(Ycoord))
cursor.insertRow([xy])
del cursor

arcpy.AddMessage(xy)


texts: Erdas Imagine DN to TOA reflectance for SPOT images

texts: I would like to convet DN to TOA reflectance in Erdas Imagine 2013. I see there is a way of doing this for Landsat7 in the Radiometric tab. How can I do this for SPOT2,4,5 images?


texts: Where is the Spatial reference ID stored in SpatiaLite?

texts: I have loaded a shapefile into a SpatiaLite Database using the spatialite-gui. In which table is the Spatial reference ID or EPSG code (of the stored data)stored in SpatiaLite?


texts: Why vegetation index, NDVI or SAVI is not working in the desert?

texts:
Could anyone explain why the vegetation indices are not working in the desert region?

I used Landsat 8 image which has been already radiometric/atmospheric corrected.

I want to see the fractional vegetation cover in super arid environment in part of middle east country. Most of the area is covered by sand sheet, but there has a few trees and I wish I could extract them.

So far, NDVI, SAVI, and EVI seems not working well. NDVI and SAVI output are quite similar, and somehow both have high value in some place which do not have any vegetation at all (blue circles in the image below). In contrast, the area contain a few trees has lower value (red circle). Only the value seems matching with real land cover types is agricultural area which shows in red colour.

The SAVI image below is displayed as a rainbow scale, from low (-1): dark blue to high (+1): red.


texts: Problem writing list of features to text file

texts:
I am attempting to write / run a simple script that list all the points within a directory, and the subdirectories contained within.

Here is the code:

import arcpy, os

workspace = r'c:\Users\athom\Documents\ArcGIS'

txtfile = r'c:\users\athom\desktop\PointList.txt'

featureClass = []
for dirpath, dirnames, filesnames in arcpy.da.Walk (workspace, datatype = 'FeatureClass', type = 'Point'):
    for filename in filesnames:
        fc = featureClass.append (os.path.join(dirpath, filename))
        print featureClass
        openfile = open (txtfile,'a')
        openfile.write(str(fc) + '\n')

openfile.close()

The script appears to work fine, as it lists the points in the Python Interpreter at the bottom of the PyScripter.

The issue is with the output to the .txt file as it comes out as the features are written as this: ਍潎敮਍潎敮਍潎敮਍潎敮਍潎敮਍ over and over and over.

I copied the characters and pasted them into a word document, but they remain as strange characters.

If the script points to a directory with only a few point features, it still outputs the same characters.

Does anyone have any suggestions? I think it may be related to this question, but am not 100% certain.


texts: Arcpy ListLayers function returns only the first sub layer in a group layer (as opposed to the group layer AND all sub layers of that group)

texts:
I thought the arcpy.mapping.ListLayers function in ArcGIS 10.1 was supposed to return group layers AND all of their sub-layers, but my script is only outputting the first sub-layer in a group layer. Can anyone see why by looking at my code? Thanks! (p.s. the purpose of my code is to output a text file - code omitted - that lists all of the layers in a specified directory and the layers' properties).

In ArcGIS 10.1 Desktop Help it says:

  Group layers are treated just like layers. The index values are simply generated from top to bottom as they appear in the table of contents or the way they would appear in a layer file.

Code:

layerList = arcpy.ListFiles("*.lyr")
for layer in layerList:
    path = os.path.join(directory, layer)
    lyrFile = arcpy.mapping.Layer(path)

    for lyr in arcpy.mapping.ListLayers(lyrFile):

        try:
            if lyr.supports("DATASETNAME") and lyr.supports("TRANSPARENCY") and lyr.supports("DEFINITIONQUERY") and lyr.supports("SHOWLABELS"):
                doc.write("{0}, {1}, {2}, {3}, {4}, {5}, {6}, {7}"
                          .format(lyr.name, lyr.datasetName, str(lyr.transparency), str(lyr.visible), lyr.definitionQuery, str(lyr.showLabels), str(lyr.minScale), str(lyr.maxScale)))
                doc.write("\n")
                break

            if lyr.supports("DATASETNAME") and lyr.supports("TRANSPARENCY") and lyr.supports("DEFINITIONQUERY") and not lyr.supports("SHOWLABELS"):
                doc.write("{0}, {1}, {2}, {3}, {4}, {5}, {6}, {7}"
                          .format(lyr.name, lyr.datasetName, str(lyr.transparency), str(lyr.visible), lyr.definitionQuery, " ", str(lyr.minScale), str(lyr.maxScale)))
                doc.write("\n")
                break

            if lyr.supports("DATASETNAME") and lyr.supports("TRANSPARENCY") and not lyr.supports("DEFINITIONQUERY") and not lyr.supports("SHOWLABELS"):
                doc.write("{0}, {1}, {2}, {3}, {4}, {5}, {6}, {7}"
                          .format(lyr.name, lyr.datasetName, " ", str(lyr.visible), " ", " ", str(lyr.minScale), str(lyr.maxScale)))
                doc.write("\n")
                break

            if lyr.supports("DATASETNAME") and not lyr.supports("TRANSPARENCY") and not lyr.supports("DEFINITIONQUERY") and not lyr.supports("SHOWLABELS"):
                doc.write("{0}, {1}, {2}, {3}, {4}, {5}, {6}, {7}"
                          .format(lyr.name, lyr.datasetName, " ", " ", str(lyr.visible), " ", " ", str(lyr.minScale), str(lyr.maxScale)))
                doc.write("\n")
                break

            if not lyr.supports("DATASETNAME") and not lyr.supports("TRANSPARENCY") and not lyr.supports("DEFINITIONQUERY") and not lyr.supports("SHOWLABELS"):
                doc.write("{0}, {1}, {2}, {3}, {4}, {5}, {6}, {7}"
                          .format(lyr.name, " ", " ", str(lyr.visible), " ", " ", str(lyr.minScale), str(lyr.maxScale)))
                doc.write("\n")
                break

        except KeyError:
            pass


texts: A method to match narrow bands from in-situ spectrometer with Landsat-8 central band

texts:
Is there a way to resample Landsat-8 bands to match different bands from a spectrophotometer. I'm doing atmospheric correction, comparing in-situ reflectance measurements performed using 10-nm wide spectral bands with reflectance from satellite images with much broader bands.

In-situ visible bands are for example: 441.8nm, 491.6nm, 530.5nm and 555.4nm; whereas Landsat-8's corresponding bands are 443nm, 482nm, 561nm and 655nm respectively.

I know ENVI can be used but don't know how to go about it.


texts: Random selection of polygon but with certain value

texts:
i'm using QGIS, and i have a multi polygon which has more then 6000 features with population data in it, so i want to select only 400 features randomly but the selected feature need to have a certain population sum, lets say around 20.000.000 people.

I've been using random selection but it only select 400 features randomly.

is there any other way?


texts: Changing Size of Labels in Survey123 for ArcGIS?

texts: I am facing a problem when changing the size of labels. The html codes are appearing in my fields.


texts: Sequentially distilling several raster data sets into final map based on first time pixel’s value was greater than zero?

texts:
I have 120 daily data sets that contain a value only when a threshold temperature is exceeded. If the threshold has been reached that pixel was assigned a value indicating the day it occurred (the Julian day), all the rest of the cells for that day were assigned a value of zero. I’d like to combine these raster data sets to create a final map indicating the Julian day when the threshold was first reached or exceeded.

I could use the raster calculator (in QGIS it would look something like this):

(“julian1” = 1) * “julian1” + (“julian2” = 2) * “julian2” + (“julian3” = 3)* “julian3”

up to 120 but  I have 30 years of data parsed into 120 day blocks and it would take a really long time. I’ll bet there is a way better solution in either QGIS or ArcMap or a python script.


texts:
I have 65 monthly MOD11A2 (LST) images. The images came from the average of the original 8 day products. I have made them a time series stack using the timeStackMODIS function of the bfastSpatial package in RStudio.
When I am running the bfmSpatial function with the following parameteres:

bfastmonitor(data,start=c(2016, 1), formula = response ~ trend + harmon, order = 1,history = c("ROC"), type = "OLS-MOSUM", h = 0.25, end = 10, level = 0.05, verbose = TRUE, plot = TRUE)

the first 2 out of the 3 output rasters are empty, and the 3rd raster it looks the same as my input time stack layer (my study area). I have the same problem no matter what the parameters.
Below it's the entire code:

library(devtools)
install_github('loicdtx/bfastSpatial')
library(bfastSpatial)
library(raster)
library(rgdal)

setwd("my_dir")

#create a list of the images
list = list.files(path =my_dir", full.names = FALSE, pattern = "\\.tif$" )

#time series layer stack of the list
m = timeStackMODIS(list, pattern = glob2rx("\\.tif$"))

#check the names of the timestack
names(m)

#check the number of layers and their dates
getMODISinfo(list)

#The number of available observations in modis timestack
obs = countObs(m8, navalues = c(NA), sensor = NULL, as.perc = TRUE)
plot(obs)
summary(obs)

#apply bfm to a ts pixel
plot(m, 4)
bfm = bfmPixel(m, start = c(2010, 1), interactive = TRUE) #monitoring period = c(2016,1)
#stats of bfm for the specific pixel
bfm$bfm
#breakpoint diagram
plot(bfm$bfm)

#apply bfm to the whole tsstack
t2 <- system.time(bfm <- bfmSpatial(m8, start=c(2016, 1),
formula = response ~ trend + harmon, order=1, history = c("ROC"), type = "OLS-MOSUM", h = 0.25, level = 0.05, verbose = TRUE, plot = TRUE))
# plot the result
plot(bfm)

Below is the result after the executing the algorithm. Only the 3rd layer has an image (error flag is the name of the layer). When I run the bfmPixel the algorithm finds a break. Any thoughts of what might caused that issue?
output raster layers by bfmSpatial]1
If someone wants, I can give him/her a set of my data.


texts: Argument 1 has unexpected type 'set' Error using processing in QGIS?

texts:
I am currently working on a code that adds a field to a shapefile and inputs its filepath. I understand that the processing.run has changed from qgis 2 to qgis 3 but I am not sure why I am getting this error. Below is the code that I am working on.

enter code here
from qgis.core import *
import sys
import glob, os, shutil, time, string
import qgis.utils
import processing
from qgis.analysis import *

Location = r"\shapefile"

for root,d_names,f_names in os.walk(Location):
    for f in f_names:

        if f.endswith('.shp'):
            if ".kml" in f:
            continue
            src = os.path.join(root,f)
            layer = QgsVectorLayer(src, f, "ogr")
            param = {layer, 'Filepath' , 2 , 150 ,0 , True, root, src}

This is the error I get:

TypeError: QgsProcessingAlgorithm.parameterAsBool(): argument 1 has unexpected type 'set'


texts: Assigning unique Asset IDs using ArcGIS Desktop?

texts:
I have a layer with street lights on it and they all have unique asset ID numbers assigned to them already.

Our data collectors are re-picking up all these street lights again but this time it will have all attribute information we need.

Is it possible to transfer that unique ID from the old layer to the new layer?

There will be 2 point features that are pretty much the same x,y coordinates.

Could it be possible to spatially relate them or something?

I am not a GIS expert. I just do a little data management and cartography.

At this point all I can think of is manually copying and pasting the ID numbers into the new attribute table.


texts: Replace several coordinates by its mean value grouped by multiindex Pandas

texts:
I have a table with some info that are getting by a GPS, in my table I have this structure:

        Date    Time    Latitude    Longitude   Accuracy    Weather Visibility  ViewPoint   Species Bearing Distance    Count
29/12/2018  10:20:45    39.85528    -5.36466    1.500000001 Soleado Buena   Otero 4 Anser anser 210 600 70
29/12/2018  10:23:32    39.85528    -5.36466    1.500000001 Soleado Buena   Otero 4 Anser anser 260 900 20
29/12/2018  10:24:17    39.85528    -5.36466    1.500000001 Soleado Buena   Otero 4 Larus fuscus    230 500 2
29/12/2018  10:05:16    39.85528    -5.36466    1.299999952 Soleado Buena   Otero 4 Fulica atra 100 450 60
29/12/2018  10:09:47    39.85528    -5.36466    1.299999952 Soleado Buena   Otero 4 Anser anser 160 560 110
29/12/2018  10:12:41    39.85528    -5.36466    1.299999952 Soleado Buena   Otero 4 Anas penelope   160 560 200
29/12/2018  10:17:33    39.85528    -5.36466    1.299999952 Soleado Buena   Otero 4 Anser anser 194 560 45
29/12/2018  10:20:01    39.85528    -5.36466    1.299999952 Soleado Buena   Otero 4 Anser anser 200 560 60
29/12/2018  08:50:07    39.86689    -5.35504    2.000000001 Soleado Buena   Otero 5 Anser anser 0   0
29/12/2018  08:55:49    39.866885   -5.355045   1.500000001 Soleado Buena   Otero 5 Anas penelope   200 550 60

I need to get the mean coordinates (to improve the GPS error) by the ViewPoint and Date. So, for the Date 1, I am going to have a mean coordinate for the ViewPoint1 (Otero1), and other mean coordinates for every ViewPoint. Later on, in the next Date, I am going to have the same structure, but the ViewPoint can be slightly displaced from one date to another, that's why I need group also by Date and not on ly by the ViewPoint.

I can apply the groupby and I get almost what I need:

df.groupby(['Date', 'ViewPoint'], as_index=True)['Latitude', 'Longitude'].mean()

Which throws this output:

                       Latitude Longitude
Date        ViewPoint
2016-12-29  Otero 0 39.859567   -5.399055
            Otero 1 39.871807   -5.385762
            Otero 3 39.853757   -5.381488
            Otero 4 39.855280   -5.364660
            Otero 5 39.866887   -5.354937
            Otero 9 39.830826   -5.356414
            OteroX  39.854657   -5.391452
2017-01-10  Otero 0 39.858376   -5.395767
            Otero 1 39.871402   -5.385639
            Otero 2 39.874892   -5.379403
            Otero 4 39.855324   -5.364828
            Otero 5 39.866883   -5.355035
            Otero 9 39.830876   -5.356457

So, what I am missing is how to replace the mean values in every row of the original table, with its mean values in the same combination  "Date" and "ViewPoint"


texts: Exporting kriging prediction map in .asc format?

texts:
I would like to save the information generated from an ordinary kriging (processed with the function autoKrige from the automap R package), into a file extension readable in ArcGis for desktop 9.3.1 (more specifically, an .asc file - aka ESRI ASCII raster format).

I tried to search for functions inside the raster (raster) and rgdal (writeOGR) packages but without success.

require(gstat)
require(automap)
require(sp)

#Reading data to provide a reproducible example.
data(meuse)
coordinates(meuse) =~ x+y
data(meuse.grid)
gridded(meuse.grid) =~ x+y

#Perform ordinary kriging and store results on object of type "autoKrige" "list"
kriging_result = autoKrige(zinc~1, meuse, meuse.grid)

plot(kriging_result) #I clipped one panel from the entire output plot.

#Extract coordinates and variable of interest to object of type "data.frame"
kriging_variable = data.frame(kriging_result[[1]])

head(kriging_variable)

#Example of observations for the variable of interest "var1.pred"
x      y       var1.pred  var1.var   var1.stdev
181180 333740  777.5909   96809.44   311.1421
181140 333700  839.7388   78629.08   280.4088
181180 333700  773.9795   84108.76   290.0151
181220 333700  710.5693   90272.17   300.4533
181100 333660  911.2036   61045.30   247.0735
181140 333660  836.4262   65746.42   256.4107

Does anyone have a clue?


texts: Repaint seems not working when called multiple times from the same method

texts:
I have the following scenario:

displaySelectedFeatures(IDs);
then some code... and then:
displaySelectedFeatures(features);

Display selected features looks like this:

public void displaySelectedFeatures(Map<Color, Set<FeatureId>> features) {

    Style style;

    style = createSelectedStyle(features);

    Layer layer = mapFrame.getMapContent().layers().get(0);
    ((FeatureLayer) layer).setStyle(style);

    mapFrame.getMapPane().repaint();
}

What I need is the following. First IDs object should be drawn and then after some time e.g. 5 seconds, object - features should be repainted.
What happens it that only the IDs is drawn while features(second repaint) is omitted.
How to solve this?


texts: How to clip mesh files with a mask, and calculate the (3D) area of triangles inside?

texts:
I have a georeferenced mesh produced by photogrammetry software such as Agisoft Photoscan. Within this mesh, I have some areas of interest, which I can define with 2D vector masks (on the x and y plane) in the form of a shapefile or a SpatialPolygonsDataFrame (R).

Can anyone reccomend a method to clip mesh files (typically in .obj format) using a mask?

This is like an extract-by-mask workflow, but on a mesh as opposed to a typical raster.

A typical mesh in this case looks something like this:

Ideally I wish to use R or python.


texts: Clipping multiparted shapefiles to multiparted shapefiles?

texts:
Using ArcGIS Desktop 10.4. I am trying to clip a multiparted (polygon) layer to a multiparted polygon layer.  Specifically, I have a statewide dataset consisting of multiple irregular polygons confined to the entire state boundary;  I want to clip those features to the individual state counties, which I have in both a complete shapefile consisting of constituent county shapes as well as individual, separate county shapes.  I don't want to have to run a clip on each of our 100 counties, so I figured there'd be a way of using the complete statewide county layer.

Can this be done or am I doomed to batch clipping?


texts: Pro vs Con on using export to JPEG or JP2000 files for large imagery files

texts:
What is the preference among GIS Specialist and Remote Sensing Specialist using the JPEG vs JP2000.

Yes, I am aware of other compressed files such as ECW or MrSID...

I am more of concern of how to add them to ArcMap or ENVI or ERDAS because it can take some time for these software to load them.

Currently, I have ArcGIS 10.1, ENVI 5.0 (plan on to upgrade to 5.1), and ERDAS 2010.


texts: What is this line in the middle of nowhere?

texts:
Coordinates:https://www.google.com/maps/place/Ramsey+Island/@73.3566979,177.7371323,1017153m/data=!3m1!1e3!4m5!3m4!1s0x30fa7771d2910233:0xbf0ed91f32e1737c!8m2!3d12.4930556!4d97.8388889

What is this, is it a rendering failure or a something real?


texts: Rasterize shapefile to raster based on a raster

texts:
I'm wondering if it's possible to rasterize a shapefile based on raster cells as shown below.
I've tried several algorithms on qGIS, but they haven't the option to create a raster based on another raster.


texts: OSRM server for public transport not a trip on roads by cars

texts:
I implemented the OSRM server for routing but I want the path of a person using the public transport lines (Bus) ; instead of getting the road path.

On building the OSRM from the pbf file I used the command

osrm-extract *.pbf -p profiles/car.lua

May be because I am extracting car paths?


texts:
Is it possible to conditionally fill in county polygon of a state based on a date field?

Example:

ABCD County has a plan that was approved by FEMA on 10/01/2015 and expires on 09/30/2020.
For ABCD county polygon I want to show "Green" until the current date is 2 years out.

Once the plan is two years from expiring I want the county polygon to turn "Yellow" automatically.

Then once the plan expires, I want the county polygon to turn "Red" until the date field for the county polygon is updated with a five (5) year out or greater.

Here is a sample of the project.[


texts: Convert 800 .hdf4 files to Tiff using GRASS/QGIS/SAGA/GDAL/R

texts:
I've have MODIS ndvi data (MODIS_Grid_16 days NDVI) from the iberian peninsula, from February of 2000 to August of 2012 (really big data).

How do I batch convert about 800 .hdf MODIS files to .tif ?

I can use QGIS/GRASS/SAGA/GDAL/R but I have no experience with bash scripts.


texts: How do I fix this QGIS sextante error? NameError: name 'sextante' is not defined

texts:
A colleague of mine posted this the other day. We haven't got any nibbles yet and I wanted to elaborate on the issue.

We are having trouble running the r.mapcalculator function and receive the error below

I have gone into the file "HistoryDialog.py" and changed a reference to Sextante on line 101 to 'sextante' (original text shown below). I think I am clutching at straws here, because after this change I receive another error altogether.

A good googling returned nothing related to the 'sextante is not defined error'. Surely someone else has come across this before. Also for us it is happening on Ubuntu and Windows. We installed the latest dev version through Osgeo4w installer but no luck!!


texts: How do I remove hyphen from field column

texts: I want to use the Python field calculator to remove the hyphen from a field column. For example, I would like to parse out the hyphen from an entry like "123-45-6789" and read "123456789".


texts: How to convert ArcMap *.shp to MapInfo *.tab file?

texts:
Possible Duplicate:
  Tool for converting .shp file to .tab file?

I have files in ArcGIS format, I want to export some layers to Mapinfo files as tab file


texts: File GeoDatabase Table is too large to export?

texts: We need to query data from the national soils database. The data is in a file geodatabase (which is about 20 GB), but we need to query the data in a SQL Server environment.  In ArcCatalog, I am able to import the tables I need into SQL Server Enterprise 2014 (not SQLServerExpress which has a 10GB limit!) from the file geodatabaes, except for one of them.  This one table is about 12GB in size (and does not contain spatial info).  When I try to import the table, ArcGIS Desktop (Advanced, ArcInfo license) crashes/closes. When I open ArcCatalog again  I can see that the table/structure  is there but it is empty. The file is too large to export to a dbf. How can I get this table OUT of the file geodatabase?


texts: FME tester/testfilter equivalent for ModelBuilder?

texts:
I've made an FME workbench with a series of  test filters that I want to replicate in ModelBuilder for Arc 10.2.

Anyone know whether this can be done easily?


texts: Incomplete Delaunay triangulation unsing pyqgis3


texts:
it would be absolutely great if someone could help me with the following problem:

I have a simple Python-Plugin that starts another executable in a subprocess. From within this application (which can't be modified on a large scale) I want to start another python-script that shows some pois in QGIS. But the problem is, that this script doesn't seem to have access to the current QGIS instance (no qgis.utils.iface object).

Is there ay way to remember the QGis-instance object or access it by looking through the currently running processes?

Thanks a lot in advance!
Kian


texts: Calculating distance between polylines from same feature class using ArcGIS Desktop?

texts:
I have a feature class with multiple lines outlining the extent of a glacier through time. I would like to calculate the distance between consecutive years (i.e two polylines from the same feature class).

Is there a way I can do that?


texts: How to subtract Geometry Point Polygons in OpenLayers?

texts:
I have two json arrays, one for background polygon and one for other polygon.
And from them, I got two polygon layers.

backgroundLayer = new OpenLayers.Layer.Vector("Background Layer", {styleMap: backgroundStyleMap} );
....
totalAreaPolygon = new OpenLayers.Feature.Vector(new OpenLayers.Geometry.Polygon([linearRing]))
backgroundLayer.addFeatures([totalAreaPolygon])

And same for other polygon.

How can I subtract the other polygon from background polygon as shown in the picture (right polygon).

This is one of the examples from OpenLayers.
http://openlayers.org/dev/examples/donut.html


texts: MapInfo TAB - Integer fields and null values

texts:
Can MapInfo TAB integer fields hold null values?

Is there a spec document on MapInfo TAB attributes available for download?


texts:
I am struggling with georeferencing over 300 aerial photographs. All of those photographs are geotagged and generally are rather easy to georeference but I am wandering if that can be done automatically?

What I know:
 - xy position for every photo (one point)
 - all of the photos were taken in similar distance (from 88 - 92 meters interval)
 - all of them has the same size, format and pixel size

I was thinking about georeferencing one picture and then use the same link points for the next photo but the problem is I don't know how to make it to get the enough precision.

When I tried to use the same link points for the next picture it worked well, so the photo was well georeferenced but in the wrong place (good shape but the location of previous photo) that's why I thought about moving that link points, maybe depending on loacation of those yellow points representing picture location?

Is there any way to joixn all the photos into one raster or anything which will let me make the georeferencing faster or automatically?


texts: Pixelated low res OSM XYZ tiles on retina display

texts: After installing the latest Qgis3 for Osx on my Macbookpro w Retina display the Qgis canvas is still pixelated and appears at a (seemingly) low resolution (when loading OSM XYZ tiles).  I knew this was an issue on previous versions of Qgis but I thought this was resolved in the latest version.  Is this just a settings adjustment issue?


texts: Creating equidistant sampling points in polygon in QGIS 3.6.0?

texts:
I have figured out how to add random points inside polygon boundaries. Although this works well, nevertheless I need a few different things to happen for my application.

What I would like to do is get equidistant points (5-acre square) within multiple polygon boundaries. I can somewhat replicate this result with random point generation but I can't buffer the points from the edge of the polygon which would be nice.

Additionally, when I use Vector -> Research Tools -> Random points in layer bounds, the points are generated with random numbered attributes, meaning there is no relation to the polygon boundary. Because these are used for sampling I assume that will get a bit frustrating when one area is not in order ie. "Field 1" has samples 1, 5, 8 but "Field 2" has 2, 4, 20 etc.

I have tried the regular point generation but it does not clip to each boundary which in a large data-set will take forever to do it with each boundary.

If anyone has any idea how to get regular points within each polygon on a multipart polygon layer?


texts: Launch Atlas composer generation from command line

texts:
I am using QGIS 2.2 and have configured a composer with atlas generation. I would like to launch this atlas generation every day.
The idea is to write a batch and launch it with a cron.

Is it possible to launch atlas generation from command line ?


texts:
I am getting trouble with creating a basemap with my own data. I have huge data. So I don't want to use any other basemap platform like leaflet or Mapbox or Google.

I have tried with ArcGIS server and successfully complete the process. But Now my challenge is to create it with an open source platform.

By which platform I can create Basemap with my data?