texts: Set layer rendering 'refresh layer at interval (seconds)' using PyQGIS from python console?

texts:
I am trying to use the python console to achieve the following manual procedure of setting each layer to refresh automatically at a specified time interval.

The manual process:
Layer Properties > Rendering > [tick box] Refresh layer at interval (seconds) > [enter seconds].

How do I do this from the python console?

This specific rendering information is not stored in a style, and manually setting it for each new layer is time consuming. I often reload all of my layers into a new QGIS project, so being able to set the refresh interval rate with python will save a lot of time.


texts: Navigator for ArcGIS: search by coordinates

texts: I just installed Navigator for ArcGIS and start testing it. I tried to search using coordinates and it was not possible. It seems to search only addresses. Search by coordinates is not an included functionality?


texts: How to write GML and GeoJSON features while retaining attributes in .NET?

texts:
I am developing a .NET project [C#] that needs to export/write features as GML and GeoJSON, including, of course, the attributes; data could queried by a database or manually entered.

I checked NetTopologySuite and GeoApi, I have also tried DotSpatial and SharpMap, but none seems to be able to include attribute fields and values in the output.

Am I wrong? Is there anybody who can help me or suggest the right place to start from?


texts:
I'm trying to convert a map I found here into lat/lon compatible with another data set I created from Google Maps lat/lons.

From the directory with that .shp file, I run:

library(rgdal)
SG = readOGR('../data', 'MP14_PLNG_AREA_NO_SEA_PL')

A subset of the data I'm trying to co-plot is:

dim = c(5L, 2L)
# polygon coordinate matrices
poly_mats = list(structure(c(103.843, 103.843, 103.854, 103.854, 103.843,
                             1.318, 1.324, 1.324, 1.318, 1.318), .Dim = dim),
                 structure(c(103.843, 103.843, 103.854, 103.854, 103.843,
                             1.324, 1.329, 1.329, 1.324, 1.324), .Dim = dim),
                 structure(c(103.854, 103.854, 103.865, 103.865, 103.854,
                             1.318, 1.324, 1.324, 1.318, 1.318), .Dim = dim),
                 structure(c(103.854, 103.854, 103.865, 103.865, 103.854,
                             1.324, 1.329, 1.329, 1.324, 1.324), .Dim = dim))
# nest properly to create SP object
polySP = SpatialPolygons(lapply(seq_along(poly_mats), function(ii) {
  Polygons(list(Polygon(poly_mats[[ii]])), ID = ii)
  # assign lat/lon CRS
}), proj4string = CRS('+init=epsg:3857'))

Both of these objects look fine when plotted separately:

png('separate.png')
par(mfrow = c(1, 2))
plot(SG)
plot(polySP)
dev.off()

Not that it should matter, but for certain, the objects are currently in different CRS:

identical(proj4string(SG), proj4string(polySP))
# [1] FALSE

So we should be able to transform SG like so:

SG = spTransform(SG, proj4string(polySP))

But this fails:

png('together_no.png')
plot(SG)
plot(polySP, col = 'red', add = TRUE)
dev.off()

It seems spTransform has failed, as the coordinates in the output are non-sensical as lat/lons:

coordinates(SG)[1:5, ]
#       [,1]     [,2]
# 0 11559649 153646.0
# 1 11569258 147405.4
# 2 11559462 150848.3
# 3 11544079 146374.9
# 4 11549925 150925.5

What's going wrong/how can this be rectified?


texts: .getDownloadURL() does not respect image clipping (Earth Engine)

texts:
I want to download a clipped version of image with getDownloadURL(), but the result only maintains the bounds of clip polygon, not the shape itself.

This is an example script:

var no2_coll = ee.ImageCollection("COPERNICUS/S5P/OFFL/L3_NO2");

// irregularly shaped polygon
var aoi = ee.Geometry.Polygon([[[36.19717602142744, 36.235397497040246], [35.03262524017744, 34.38876424418324],
[39.17447094330244, 34.37062991132883], [39.41617016205244, 36.36820725445314], [38.17471508392744, 35.3890351383359]]]);

var image = no2_coll.select("tropospheric_NO2_column_number_density")
  .filterBounds(aoi).first();

// this clips an image, but has no effect on result of .getDownloadURL()
var img_clip = image.clip(aoi);

// this has effect on result of .getDownloadURL()
var img_clip_double_mask = img_clip.mask(img_clip.mask());

var downloadArgs = {
  name: 'expected_clipped_image',
  crs: image.projection(),
  scale: 7000,
  region: aoi
};

var url = img_clip_double_mask.getDownloadURL(downloadArgs);

print(url);

So while both img_clip and img_clip_double_mask are rendered as a clipped image on Map, it is only when you change  img_clip to img_clip_double_mask (essentially self masking - but seems redundant) when you get clipped image from .getDownloadURL().

Why is that?


texts: ArcPy Runtime error: Invalid SQL statement used

texts:
I am trying to use a search cursor to read fields from a table, assign their values to variables, find the difference and print it. So far I keep getting a runtime error. I'll list my code and the error below. I'm using a virtual desktop which may make the path look odd.

import arcpy

fc = r'\\directory\userprofiles$\me\Documents\ArcGIS\Default.gdb\usPOPS'
arcpy.env.overwriteOutput = True

with arcpy.da.SearchCursor(fc, ['Resident population 1990', 'Resident population 2010']) as cursor:
    for row in cursor:
        pop1990 = row[0]
        pop2010 = row[1]

del cursor

popDiff = (pop1990 - pop2010)

print (popDiff)

Runtime error
Traceback (most recent call last):
  File "<string>", line 12, in <module>
RuntimeError: An invalid SQL statement was used. [SELECT Resident population 1990, Resident population 2010, OBJECTID FROM usPOPS]

Any ideas? I've just been hitting dead ends.


texts: Custom Map Projection from a set of coordinates (no transformation/parameters available)

texts:
I want to create a Custom Map Projection in QGIS without knowing any parameters for the transformation.
My basis is a georeferenced raster image that contains a grid with coordinates from a local system that has been used in some area in the 1960s. However, there is no known transformation for this local coordinate system. I could generate an image/dxf or even a list with more coordinates in both a known projection and the desired local system from AutoCAD, if it helps. The system only covers an area of about 3 by 3 kilometers and I only need any kind of accuracy inside this small area, because the system has never been used nor was it meant to be used outside it.
I cannot find much on the internet about this topic. I am neither a geography nor a GIS specialist, and so I have no idea how to approach this, though I figure there must be a way of calculating the transformation. Is there any way to achieve this?


texts: QGIS v.transform output extent not lining up with input

texts:
I am using the updated Affine tool (v.transform) to rescale my vector layer (containing multiple polygons). I want to shrink each polygon to half of its size. So I enter in 0.5 for both the x and y scale. I leave everything else as default and then hit run. The output is in another extent and is located somewhere weird. I tried setting the extent argument to my original vector, but I still get the same output problem. I've also tried SAGA's transform tool. Same problem.
v.trasnform tool inputs:

original vector files extent:

transformed vector files extent:

It's also weird that the feature count increased. I went through the file and noticed that overlapping sections were being split into their own polygon. Not sure why that is happening.


texts: Spatial join features from different geodatabases

texts:
I have two geodatabases with several feature classes in both of them. I want to spatial join the feature classes in gdb1 to those in gdb2. I tried to use "arcpy.ListFeatureClasses()" function, but it needs to set the workspace to one gdb to get its list. Here is what I wrote:
import os
import pandas as pd
import arcpy

#set workspace as gdb1 and get its list
arcpy.env.workspace = r'P:\Jung\Brownfield\=2020\Parcels_poly_scag.gdb'
FC_list = arcpy.ListFeatureClasses()

#set workspace as gdb2 and get its list
arcpy.env.workspace = r'P:\Jung\Brownfield\=2020\Superfund_Brownfield_Database_2020.gdb'
ds_list = arcpy.ListFeatureClasses()

for i in ds_list:
    target_ft = i
    for j in FC_list:
        join_ft = j
        arcpy.SpatialJoin_analysis(target_ft, join_ft, target_ft + '_' + join_ft, 'JOIN_ONE_TO_ONE')

And I got the error -

ERROR 000732: Target Features: Dataset SF_Active_Site_Inventory does not exist or is not supported.

It seems when I changed the workspace, the feature list of the former workspace expired.
How can I solve this problem?


texts: Geonode importlayers 504 Gateway Time-out

texts:
I'm trying to upload layers to Geonode using importlayers. Some of the shapefiles were successfully uploaded. However, we've encountered some errors in other shapefiles during uploading. Below is the error message.

================
Traceback (most recent call last):
  File "/home/geonode/geonode/geonode/layers/utils.py", line 588, in upload
    title=title
  File "/home/geonode/geonode/geonode/layers/utils.py", line 470, in file_upload
    defaults=defaults
  File "/home/geonode/.virtualenvs/geonode/lib/python2.7/site-packages/django/db/models/manager.py", line 154, in get_or_create
    return self.get_queryset().get_or_create(**kwargs)
  File "/home/geonode/.virtualenvs/geonode/lib/python2.7/site-packages/modeltranslation/manager.py", line 341, in get_or_create
    return super(MultilingualQuerySet, self).get_or_create(**kwargs)
  File "/home/geonode/.virtualenvs/geonode/lib/python2.7/site-packages/django/db/models/query.py", line 383, in get_or_create
    obj.save(force_insert=True, using=self.db)
  File "/home/geonode/.virtualenvs/geonode/lib/python2.7/site-packages/polymorphic/polymorphic_model.py", line 90, in save
    return super(PolymorphicModel, self).save(*args, **kwargs)
  File "/home/geonode/.virtualenvs/geonode/lib/python2.7/site-packages/django/db/models/base.py", line 545, in save
    force_update=force_update, update_fields=update_fields)
  File "/home/geonode/.virtualenvs/geonode/lib/python2.7/site-packages/django/db/models/base.py", line 569, in save_base
    update_fields=update_fields)
  File "/home/geonode/.virtualenvs/geonode/lib/python2.7/site-packages/django/dispatch/dispatcher.py", line 185, in send
    response = receiver(signal=self, sender=sender, **named)
  File "/home/geonode/geonode/geonode/geoserver/signals.py", line 74, in geoserver_pre_save
    charset=instance.charset)
  File "/home/geonode/geonode/geonode/geoserver/helpers.py", line 1125, in geoserver_upload
    workspace=workspace)
  File "/home/geonode/geonode/geonode/geoserver/helpers.py", line 998, in _create_db_featurestore
    raise GeoNodeException(msg)
GeoNodeException: An exception occurred loading data to PostGIS-
504 Gateway Time-out

The server didn't respond in time.


texts:
I can make a vector tiles request with valid response against Geoserver as WMS but in OL you seem to have to use a ol.source.VectorTile source which is some kind of xyz-url.
Perhaps WMS is possible as a tileloadfunction?


texts: Is there a way to execute SQL queries on file geodatabase from Python? I'm looking for the same functionality provided by arcpy.ArcSDESQLExecute for SDE databases.


texts: Use lat / long to create polyline that can be passed into DrawPolyline

texts:
public void DrawPolyline(ESRI.ArcGIS.Carto.IActiveView activeView)
{
  if(activeView == null) return;
  ESRI.ArcGIS.Display.IScreenDisplay screenDisplay = activeView.ScreenDisplay;

  screenDisplay.StartDrawing(screenDisplay.hDC, (System.Int16)ESRI.ArcGIS.Display.esriScreenCache.esriNoScreenCache); // Explicit Cast
  ESRI.ArcGIS.Display.IRgbColor rgbColor = new ESRI.ArcGIS.Display.RgbColorClass();
  rgbColor.Red = 255;

  ESRI.ArcGIS.Display.IColor color = rgbColor; // Implicit Cast
  ESRI.ArcGIS.Display.ISimpleLineSymbol simpleLineSymbol = new ESRI.ArcGIS.Display.SimpleLineSymbolClass();
  simpleLineSymbol.Color = color;

  ESRI.ArcGIS.Display.ISymbol symbol = (ESRI.ArcGIS.Display.ISymbol)simpleLineSymbol; // Explicit Cast
  ESRI.ArcGIS.Display.IRubberBand rubberBand = new ESRI.ArcGIS.Display.RubberLineClass();
  ESRI.ArcGIS.Geometry.IGeometry geometry = rubberBand.TrackNew(screenDisplay, symbol);
  screenDisplay.SetSymbol(symbol);
  screenDisplay.DrawPolyline(geometry); //I WANT TO CREATE THIS FROM LAT / LONGs!
  screenDisplay.FinishDrawing();
}

The above code snippet works but the problem is I want to create the IGeometry from lat \ longs not from clicking the screen.  How can I create an polyline IGeometry object for the drawPolyline function using standard lat \ long coordinates?


texts: How to remove outlying islands from a DEM raster?

texts:
I have a raster image that looks like this (see below)

How do I remove the outlying islands and keep only the central main island in the image. Is it a geoprocessing clip function?
I am using QGIS 2.8.1, I added a .tif file from a downloaded DEM and then used the raster>analysis>DEM(terrain models)


texts: How to construct a gs:PointStacker parameters using openlayers against geoserver

texts:
I want use openlayers 2.13 dev and geoserver to achive the cluster of points.
I user openlayers 2.13 dev 's "OpenLayers.WPSClient".
here is my code:

<!DOCTYPE html>

<meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
<meta name="viewport" content="width=device-width, initial-scale=1.0, maximum-scale=1.0, user-scalable=0">
<meta name="apple-mobile-web-app-capable" content="yes">
<title>wps</title>
<link rel="stylesheet" type="text/css" href="./theme/default/style.css">
<link rel="stylesheet" type="text/css" href="style.css" >

    <button onclick="cluster();">cluster</button>

but I don't know how to construct the pointStacker.execute especially outputBBOX, chrome gives me the error Uncaught TypeError: Cannot read property 'left' of undefined, meaning the source code's options.bounds.left
could somebody give some help!


texts: How can i automatically display the output of my toolbox?

texts:
i have made a python toolbox (.pyt file) and run it in Arcgis 10.1. The Output is a raster and I want it to automatically show up in the Arcmap viewer. The option for displaying the output of geoprocessing under Geoprocessing Options is turned on, but the raster layer still does not show up. Of course I can add it the usual way but normally when running a toolbox, the results show up without doing anything. I also tried this: arcpy.env.addOutputsToMap = True but there´s no effect out of it.

Can somebody help me?


texts: I'm trying to save an Image in QGIS Composer by polygon borders as seen in the right figure. But when I try to save it, I get a white rectangle box around my figure as seen in the left figure. Is there a way to save the figure like the right one?


texts: Find all lines first point in radius 500 meter

texts:
I need to find all lines startPoint and lines attribut in a given radius of (0.5 km) from a point(X,Y) Then i need to find the disnace from the point(X,Y) to all of the startpoint in the given radius (0.5 km). Can someone give me a hint how to do this.

(You have to forgive me for my bad english)


texts:
I have built my tile server on AWS EC2 instance and want to use it as a base map in QGIS.
When I try to add it as a WMS layer on QGIS, I get an error saying The server you are trying to connect to does not seem to be a WMS server. Please check the URL.

Do I need to convert my AWS url to a WMS URL?
Are there any tutorial on how to do this?


texts: What is the difference between a multipoint shapefile and a point shapefile. in the creation of a shapefile In ArcGIS it has the choice of point shapefile and multipoint shapefile. But I do not know what their difference.


texts: I have set my time interval using CuRL tool for automatic Geowebcache.. How to set bounding box in curl command?


texts: Anyone used Informatica Powercenter and FME together?

texts:
FME apparently integrates with Informatica Powercenter. I use Informatica quite a bit to move data around, push and pull to ArcSDE, and to slurp mass quantities of XML into SQL Server. I'm particularly interested in extracting DWG attributes and features and pushing them to featureclasses in SDE and into tabular SQL Server databases as well. Does anyone have any experience with anything like this or just Informatica and FME together, in general?

I'm currently extracting data from DWGs through ArcObjects, and it can be painful. So, I guess my specific question to anyone with experience would be just how well does FME integrate with Informatica?


texts: Model converted to Script - No editing mode?

texts:
I have a model that an ESRI associate sent me as an example. When I run the model everything is great. I click on the model and because the seleting features is a parameter, it automatically starts an edit session for me to digitize the selecting features. However when I convert this model to a script and try to run as a script, when I run it, the selecting feature parameter is not in edit mode, a polygon does not pop up in the table of contents as a feature class to be digitized like it was with the model. Is there a way to change the model so that the selecting features can be digitized ? Thanks

The first picture shows the model, the second shows the model when I run it, and the third shows the scrip when I run it.

# -*- coding: utf-8 -*-
# ---------------------------------------------------------------------------
# toolbox_scrip3t.py
# Created on: 2013-09-30 13:13:52.00000
#   (generated by ArcGIS/ModelBuilder)
# Usage: toolbox_scrip3t <Selecting_Features> <value>
# Description:
# test
# ---------------------------------------------------------------------------

# Import arcpy module
import arcpy

# Script arguments
Selecting_Features = arcpy.GetParameterAsText(0)
if Selecting_Features == '#' or not Selecting_Features:
    Selecting_Features = "in_memory\\{44DC94EF-CABC-4233-A4AF-CEFD34618B2E}" # provide a default value if unspecified

value = arcpy.GetParameterAsText(1)

# Local variables:
Input_Points = "FirePerimeter"
Final_Output = value
Selected_Features = Selecting_Features
Input_Points_Layer = "SamplePoints_Layer"

# Process: Make Feature Layer
arcpy.MakeFeatureLayer_management(Input_Points, Input_Points_Layer, "", "", "OBJECTID OBJECTID VISIBLE NONE;SymbologyField SymbologyField VISIBLE NONE;DomainField DomainField VISIBLE NONE;Shape Shape VISIBLE NONE;Subtype Subtype VISIBLE NONE;EditField EditField VISIBLE NONE")

# Process: Select Layer By Location
arcpy.SelectLayerByLocation_management(Input_Points_Layer, "INTERSECT", Selecting_Features, "", "NEW_SELECTION")

# Process: Calculate Field
arcpy.CalculateField_management(Selected_Features, "EditField", "'%value%'", "PYTHON_9.3", "")


texts: Refering in one layer to another layer in ModelBuilder or ArcPy?

texts:
I have a problem in using Python in ModelBuilder.

My model is going to calculate pace of lateral erosion.

I subtracted two polygon layers by Symmetrical Difference and add field "tempo" to new layer with expression in "Wyrażenie na tempo" "Shape Area"/55.5 (55.5 is length of my segment). There is no problem when pace is positive because Symmetrical Difference's output is an absolute value. The problem is when pace should be below zero.

I should use some kind of if statement but I don't know how.

My statement should be like this:

If Shape Area of Dane 2 > Shape Area of Dane1: tempo * (-1)
else: tempo *1

I'm not sure where in model should I use it (maybe add next Calculate Field Tool or write a script) and how to refer in one layer to other layers.


texts: Defining circular neighborhood for custom raster calculator function?

texts:
I have a 1 meter resolution raster containing values from 0-3 which represent increasing levels of habitat suitability.  I need to identify circular regions that have contiguous pixel values.  For example, I define a 50 meter radius, and I am looking for regions with values of 3.  I need to know if every cell within a 50 meter radius of each target cell has a value of three.  If true, then I need to return a true value (1) to the output target cell, if false, then return a false value (0).

This is similar to a focal statistics function using a circular neighborhood.

Region group will identify regions of contiguous cells, but I need the regions of contigous cells to be circular and of a minimum diameter.

I'm using ArcMap 10.2.


texts: Viewing Google Earth or Google Maps satellite imagery on GPS handhelds?

texts:
I have managed to upload BirdsEye imagery to a Garmin GPS64s handheld using BaseCamp ( v 4.2.5 ) by following the instructions, however the imagery provided is useless to me.

I have determined Google Earth or even Google Maps to be a better alternative ( for my purpouses ) and I would like to upload selected segments of Google Earth imagery to my handheld unit.

I am technically literate ( IT background ), however I am new to the technology of uploading satellite imagery to handheld GPS units ( the terminology, file formats and technology are confusing me .... I can't see the forest for the trees ).

Where can I inform myself to gain an overview of what is involved in uploading imagery from Google Earth or Google Maps to Garmin handheld ( ... or any handheld unit for that matter )?

I am looking for resources and starting points ( public domain / on-line material ) which provide enough overview to understand what to ask and where to ask it ( if necessary ).


texts: Computing the elevation angle of a GPS satellite with respect to the given receiver location

texts:
can anyone source an explicit formula for calculating the elevation angle for a GPS satellite as seen from a receiver. I can't seem to find one. Even a text reference would be great. I couldn't quite find anything posted already that does the same thing.

I'm also looking for something similar for calculating elevation but for Glonass satellites as well.


texts: Difference between nodes and waypoints OpenStreetMap

texts: What is the main difference between a node and a waypoint in OpenStreetMap? I know both represent a point on the map. Also what is the significance of node 0?


texts: Map production checklist for organisation

texts:
I would like to create a map production workflow list for my GIS colleagues at a local public authority. When someone orders a map, I would like to have a kind of check list where we can just run through the standard questions and write them into an excel table. This would help us keep a record of what maps we have acually made, which datasets were used, who does the mapping, and a production history which is good for time management purposes.
I was thinking around the lines of creating a mask in an excel table including "Purpose","target Audience", "Output medium", "size", "Scale/resolution requirements", "Title", "north arrow"..etc...

Has anyone created something similar? Does this seem useful? Does anyone have similar documents which I could look at and compare?


texts: How to fix an error while registering PostgreSQL database to ArcGIS Enterprise Server

texts:
Software in use is:
ArcMap 10.8.1
ArcGIS Server (Enterprise) 10.8
PostgreSQL 13 + PostGIS 3.1
Windows 10 64bit
I followed a guide to register a PostgreSQL database (PostGIS) to ArcGIS Server. When I try to register the database, it gives an error as in the picture I provided. Connecting to the PostgreSQL database through ArcMap client works flawlessly (I can query, make edits etc.). What can be the problem? I googled for almost an hour without finding anything similar to my error.


texts: Scale and aggregation function in Google Earth Engine export

texts:
Given a 30 m resolution Image in Google Earth Engine, I would like to export it to Google Drive at a 1000 m resolution such that within each 1000 m pixel I obtain the sum of the underlying 30 m pixels.
I know the reduceResolution function, which however is suitable only for small areas (up to 65536 pixels) and does not allow to specify the output resolution.
In other words, I would like to achieve what in R can be done with the function aggregate of the raster package.
Take the following example:
var image = image.clip(geometry)
Map.addLayer(image);
Export.image.toDrive({
  image: image,
  description: 'image',
  scale: 30,
  maxPixels: 10000000000000,
  fileFormat: 'GeoTIFF',
  crs : 'EPSG:4326'
});

If I change the scale parameter in the export function to 1000, which is my objective, GEE will interpolate pixels to decrease the resolution, but not sum their value, which is my objective.


texts: Submit empty string via POST wfs in GeoExt.form.FormPanel

texts:
I have a two textfields in FormPanel. When one of them is empty, after pressing the search button the form sends the empty string in WFS query. How can I disable the empty field which sends the WFS request after I pressed the search button?

    formPanel = new GeoExt.form.FormPanel({
    ref: "formPanel",
    width: 300,
    height: 200,
    region: "west",
    protocol: protocol,
    buttons:[{
     text: "search",
     handler: function() {
         formPanel.getForm().search();
     },
     //scope: formPanel
 }],
    items: [{
        xtype: "textfield",
        name: "trans_id__eq",
        value: "",
        fieldLabel: "ID",
        disabled: true,
        ref: 'myTextField'
    },

    {
        xtype: "textfield",
        name: "nabywca__eq",
        value: "",
        fieldLabel: "nabywca",
        submitEmptyText: false
    },
    ],
    listeners: {
        actioncomplete: function(form, action) {
        // this listener triggers when the search request
        // is complete, the OpenLayers.Protocol.Response
        // resulting from the request is available
        // through "action.response"
            features = action.response.features;
            app.featureGrid.store.loadData(features);   //załadowanie zwróconego wyniku wyszukiwania do "app.featureGrid.store"
            vm=app.mapPanel.map.getLayersByName("Results");
            if(vm.length==0){
                vecLayer = new OpenLayers.Layer.Vector("Results");
               app.mapPanel.map.addLayer(vecLayer);
                app.featureGrid.store.bind(vecLayer);
                app.featureGrid.getSelectionModel().bind(vecLayer);
            //app.featureGrid.getSelectionModel().bind(vectorLayer);
        }
    }
    }

})

    );

When I leave blank first textfield and put the value to the second my wfs request send empty string and I get the error.

My protocol:

var protocol= OpenLayers.Protocol.WFS({

            version: "1.1.0",
            url: "http://localhost:8080/geoserver/wfs",
            featureType: "transakcje_punkty",
            featureNS: "http://geoserver.org/transakcje",
             srsName: "EPSG:900913"

        });

Wfs is working properly when I have one field, but I would like to add 7 "textfields".

The user fills the fields with the attributes he want to find and  rest leaves empty.

But empty field send some empty value to WFS and I get an error... instead results.


texts: (Arcpy) read values from a txt file and convert the data into a shapefile

texts:
I´m a beginner at python programming and wanted to start with a litte and hopefully simple project. I tried to read float values out of a textfile and write them into a shapefile.

The program is working fine until the writing in the shapefile. Following error occurs:

TypeError: value #0 - unsupported type: list
TypeError: sequence size must match size of the row

I tried different styles (last lines in the code) but i could´t work it out.
So my question is: What is my error and what could i do to solve this Problem?

Example of the Textfile ("test.txt"):

2500000.00 5600000.00 10
2500001.00 5600010.00 20

Script

import fileinput
import string
import os

# Declaration
env.overwriteOutput=True
env.workspace= outpath

x = []
y = []
z = []

newfc="newfc.shp"
outpath="c:/tmp/"

# Create new Shapefile and add FIELDS
arcpy.CreateFeatureclass_management(outpath, newfc, "Point")
arcpy.AddField_management(newfc, "X", "FLOAT", field_length = 50)
arcpy.AddField_management(newfc, "Y", "FLOAT", field_length = 50)
arcpy.AddField_management(newfc, "Z", "FLOAT", field_length = 50)

# Reference Cursors
cursorX=arcpy.da.InsertCursor(newfc,"X")
cursorY=arcpy.da.InsertCursor(newfc,"Y")
cursorZ=arcpy.da.InsertCursor(newfc,"Z")

# Read File
a = open("c:/tmp/test.txt","r")
inputF = a.readlines()

for line in inputF:
        xCoordinate, yCoordinate, zValue = line.split()
        x.append(float(xCoordinate))
        y.append(float(yCoodrinate))
        z.append(float(zValue))

print([x]) # [[2500000.0, 2500001.0]]
print([y]) # [[5600000.0, 5600010.0]]
print([z]) # [[10.0, 20.0]]

# Copy Arrays into the Field
cursorX.insertRow(x) # TypeError: value #0 - unsupported type: list
cursorY.insertRow((y)) # TypeError: sequence size must match size of the row
cursorZ.insertRow(z)


texts: How to remove layers by name - openlayers 3

texts:
I was add vector layer with wkt source to map with below code:

var SelectVector = null;

for (var i = 0; i < wktarray.length; i++) {
   var wkt = wktarray[i];
   var format = new ol.format.WKT();
   var feature = format.readFeature(wkt, {
       dataProjection: 'EPSG:4326',
        featureProjection: 'EPSG:4326'
    });

SelectVector = new ol.layer.Vector({
        source: new ol.source.Vector({
            features: [feature]
        }),style: new ol.style.Style({
            fill: new ol.style.Fill({
                color: 'rgba(255, 255, 255, 0.2)'
            }),
 stroke: new ol.style.Stroke({
                color: '#ffcc33',
                width: 2
            }),
image: new ol.style.Circle({
                radius: 7,
                fill: new ol.style.Fill({
                    color: '#ffcc33'
                })
            })
        })
    });

map.addLayer(SelectVector);
SelectVector.set('name', 'selectvector');
}

Now, I want remove this vector layer from map. wrote below code but, it does not remove all layer with selectvector name.

map.getLayers().forEach(function (layer) {
    if (layer.get('name') != undefined & layer.get('name') === 'selectvector') {
        map.removeLayer(layer);
    }
});

what is wrong?


texts: Use ModelBuilder to copy/paste attributes from Table to Feature?

texts:
I am using ArcMap 10.3. This question involves the use of tables and personal geodatabase features.

I am looking for an efficient way to populate Water Meter fields using data from a Table.

The table contains all of the billing information for the water meters.
The water meters layer only contains the physical location and blank fields awaiting the billing information.

I am seeking a way to highlight a single water meter feature from a layer (output/destination), highlight a single row from the table (input), run the Model and have the data from the highlighted row in the table be appended to the water meter feature.

The data I am looking to copy is three fields, all in String/Text format.

(Yes, the point is to do this one at a time with hand-selected data)
(Yes, there are better ways to accomplish my end-goal, but this is the method I need right now)


texts: How to set the transparency for multiple layers or add a global transparancy preference?

texts:
How can I set the transparency to 27 % for all 245 layers in QGIS? The only way I know of is by right clicking on the layer name, selecting Properties, then Transparency and then moving the transparency slider left or right.

This is simple enough. But this is only good for up to 10 layers maybe. What if you have 245 layers like I do? Do you just keep on repeating the process? Now surely, there there must be a way to apply this to all 245 layers at once!?

Alternatively, is there a global transparency preference setting I can add in so that when I add new layers they automatically get 27 % transparency?

The QGIS online documentation mentions something about exporting your transparency setting to a file for a later use.

  As you can see this is quite easy to set custom transparency, but it
  can be quite a lot of work. Therefore you can use the button Export to
  file to save your transparency list to a file. The button Import from
  file loads your transparency settings and applies them to the current
  raster layer.

This seems like a useful feature. But I don't think this is what I'm looking for.

I tried selecting multiple layers in the table of contents and then right click and select Properties and set the transparency level, apply changes and click OK. It applied the changes, but only to the last layer in the selection, the one I right clicked on. None of the other layers in the selection were affected. (This could be a bug actually.)


texts: pgRouting TSP with eventual order constraints?

texts:
Is-it possible to run Travelling Sales Person (TSP) in pgRouting (I've got 2.3... installed) with constraints?

A general TSP algorythm or those I've identified in pgRounting solve the TSP
problem visiting each nodes (cities) to create the shortest cost path at the end.

I need a more sophisticated solution that enables adding constraints
before the calculated related to node (city) order.

Eg.: Say we have 5 points to reach like A B C D and E.

In a regular TSP when starting from A it would answer A -> D -> E -> C -> B [-> A].

Let's assume that we have a package to carry meanwhile from C to E, so C must precede E on that day of traveling.

So, a constraint comes up like: C -> E

The full definition of the call would be:

Do TSP with

A B C D E

Starting from A (and also finishing here)
CONSTRAINT: visit C before E (and maybe more...)

Is there a solution, algorithm, program like this?


texts:
In my program:

I am selecting only schools of the specified type based on the "FACILITY" field and creating new feature class of the result.
Deleting records from the resulting feature class from "Select" where the "Verified" field does not equal "Y" or "y" (I need to use: Update cursor and deleteRow method)
For the records that do not equal "y" print the name of each deleted school

My code is not deleting the features that are not equal to "Y". And I am not sure how to print out the values that do not equal "Y".

Current error: Traceback (most recent call last):
  File "C:\Program Files (x86)\Wing IDE 101 4.1\src\debug\tserver_sandbox.py", line 26, in
  File "C:\Program Files (x86)\Wing IDE 101 4.1\src\debug\tserver_sandbox.py", line 19, in selectSchools
TypeError: 'Row' object does not support indexing

import arcpy

work = arcpy.env.workspace = "c:\Scripts\Lab 6 Data"
sType = "JUNIOR HIGH"
outFile = "output.shp"
arcpy.env.overwriteOutput = True

def selectSchools(work,sType, outFile):
    fc = "Schools.shp"
    arcpy.MakeFeatureLayer_management(fc, "schoolslyr")
    query = "\"FACILITY\" = '{0}'".format(sType)
    arcpy.Select_analysis(fc, outFile, query)
    arcpy.AddMessage("Select_Analysis layer tool completed!")

    where_c = "\"VERIFIED\" = 'y%'"
    update = arcpy.UpdateCursor(outFile)
    row = update.next()
    for row in update:
        row.setValue(where_c)
        updates.updateRow(row)
    if (row[0].upper() != 'y'):
        row.deleteRow(row)
        print ("NAME")
        row = update.next()

selectSchools(work,sType,outFile)


texts: Fulfilling removed cloud area with free cloud data in Google Earth Engine

texts:
I have completed doing cloud removal with Landsat image (in this case is Lansat 8) in the Google Earth Engine. Of course the result of this process is the cloud pixels become no data, consequently the area which was as the cloud is perforated.
I want to process this image to the nest step for classification and image transformation. So I want to fill the holey areas with no cloud data from other time period so that my study area has complete data within.

How can I solve this problem as to fill the holey image?
I use this cloud masking scirpt:
var roi = ee.Geometry.Point([97.90305, 3.98383]);

function maskL8sr(image) {
  // Bits 3 and 5 are cloud shadow and cloud, respectively.
  var cloudShadowBitMask = 1 << 3;
  var cloudsBitMask = 1 << 5;

  // Get the pixel QA band.
  var qa = image.select('pixel_qa');

  // Both flags should be set to zero, indicating clear conditions.
  var mask1 = qa.bitwiseAnd(cloudShadowBitMask).eq(0)
      .and(qa.bitwiseAnd(cloudsBitMask).eq(0));

  var mask2 = image.select('B.*').gt(0).reduce('min');

  // Return the masked image, scaled to TOA reflectance, without the QA bands.
  return image.updateMask(mask1.and(mask2)).divide(10000)
      .select("B[0-9]*")
      .copyProperties(image, ["system:time_start"]);
}

//Clip
var clipToCol = function(image){
  return image.clip(roi);
};

// Map the function over one year of data.
var collection = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR')
    .filterDate('2019-01-01', '2019-12-31')
    .map(maskL8sr)
    .map(clipToCol);

var composite = collection.reduce(ee.Reducer.percentile([25]));

// Display the results.
Map.addLayer(composite, {bands: ['B7_p25', 'B6_p25', 'B4_p25'], min: 0, max: 0.2});

How can I overcome this issue?


texts: Navigating to Lat/Long in QGIS with a relative coordinate system

texts:
It seems really fundamental, but how do I pan the map center to a decimal lat/long in QGIS? I would type in the lat/long at the bottom, but I use WGS 84/Pseudo Mercator to use the streetmap from openlayers which shows the x/y distance instead of absolute coordinates. If I change the projection, the streetmap gets unaligned even with OTF CRS transformation enabled.

I do not want to search for the coordinate with my mouse (i.e. using coordinate capture), I just need to pan the center of the map.
It's easy to do in google earth (just search for the lat/long), ArcMap (Go to XY in the toolbar) and MapInfo (Map->Change View). I can even geocode addresses (geosearch plugin) in QGIS easier.


texts:
I had plotted a Javascript gMap.Using the marker,the values of latitude and longitude are displayed.I need to get that values and save it to database.What is the best way to do this and ,How to accomplish this.

In model,i created both field as CharField.I want to know how Python get the values from HTML and save.


texts:
I wrote a toy plugin for QGIS-3 (2.99 Ubuntu 16.04), which worked for quite a while. But after a recent upgrade, the plugin disappeared from the QGIS toolbar while the rest of QGIS seem to be working fine.

I have checked that the ~/.qgis3/python/plugins/ folder is still present. I tried to use pb_tool deploy to re-install the plugin. But the plugin still wouldn't show in the toolbar.

Does someone know what could have gone wrong?


texts: Using SelectByAttribute_management after making XY Event layer in ArcPy

texts:
I imported a CSV file using the MakeXYEventLayer_management and trying to use the SelectlayerByAttribute_management function, but I am getting error message saying on "Object ID(OID)" is present..
Is there a way to work around this?

Below is my code so far:

import arcpy
arcpy.env.overwriteOutput = True

xFieldName = 'LON'
yFieldName = 'LAT'
outFolder = "C:\\GEOG485\\FinalProject\\Output"
asos = "ASOS Stations"
spatialRef = arcpy.SpatialReference(4326)
csvFilePath = "C:\\GEOG485\\FinalProject\\isd-history.csv"
targetState = 'NJ'

arcpy.MakeXYEventLayer_management(csvFilePath, xFieldName, yFieldName, asos, spatialRef)
arcpy.FeatureClassToShapefile_conversion(asos, outFolder)

arcpy.SelectLayerByAttribute_management(asos, "NEW_SELECTION", '"STATE_NAME" = targetState')


texts:
I am in the process of writing my two first QGIS plugins. (https://github.com/sickel/qgisSpectre and https://github.com/sickel/qgisstripchart) both of them contains two drop down boxes, one to select a layer, and then the other is populated with the fields from the layer when a layer is selected.

For this, I want to just list the vector layers in the first dropdown. After some fooling around, I ended up doing it like this, which feels like going around in circles:

layers = QgsProject.instance().layerTreeRoot().children()
self.dlg.cbLayer.clear()
for layeritem in layers:
    layername=layeritem.name()
    layers = QgsProject.instance().mapLayersByName(layername)
    if len(layers)==0: # Invalid layer name
        return
    layer = layers[0] # first layer
    if "fields" in dir(layer):
        self.dlg.cbLayer.addItems([layer.name()])

So I am getting the layer names from layerTreeRoot().children(), then I am using that name to pick up the layer from mapLayersByName.

Is there a way to either Get the layer directly from layerTreeRoot().children()
or Iterate over the map layers in another way, e.g. through QgsProject.instance().mapLayers()and get each layer's type and name?

(Another problem my current code has, is that it will not list layers inside a group, I think this may work fine if I can get what I need from QgsProject.instance().mapLayers())


texts: non alignment of base layer and other vector layers

texts:
i'm sure you've encountered this problem repeatedly. i've tried all different methods.

Enable on the fly crs transformation.
go to settings and options to use crs data.
go to properties and change the coordinate system to the same one.

I'm using google street as my base map.

CRS of vector data WGS84/48N

metadata is -

                                            General
Storage type of this layer
ESRI Shapefile
Description of this provider
OGR data provider (compiled against GDAL/OGR library version 1.11.2, running against GDAL/OGR library version 1.11.2)
Source for this layer
C:\Users\hp 1.000\Desktop\GIS\weatherdata.shp
Geometry type of the features in this layer
Point
The number of features in this layer
1253
Editing capabilities of this layer
Add Features, Delete Features, Change Attribute Values, Add Attributes, Delete Attributes, Create Spatial Index, Fast Access to Features at ID, Change Geometries, Simplify Geometries, Simplify Geometries with topological validation
Extents
In layer spatial reference system units
xMin,yMin 11303.3,25898.2 : xMax,yMax 44627.4,46712.2
Layer Spatial Reference System
+proj=utm +zone=48 +datum=WGS84 +units=m +no_defs


texts: (Float Type - Modis Data) Tile error: Expected a homogeneous image collection, but an image with an incompatible band was encountered

texts:
I'm stuck trying to combine two images due to the following error. The type for the bands are both floats, but as the error shows, are not quite the same:
Scaled LST: Tile error: Expected a homogeneous image collection, but an image with an incompatible band was encountered. Mismatched type for band 'LST_Day_1km':
Expected type: Float<0.0, 2184.5>.
  Actual type: Float<-3274.25, 2.5>.
     Image ID: 1
This band might require an explicit cast.

I've attached my script below, I've seen other threads seemingly addressing similar issues. I'm very new to Java.
//Loading in the Modis EVI data
var evi = ee.ImageCollection('MODIS/MOD09GA_006_EVI')
                  .filter(ee.Filter.date('2020-08-01', '2020-08-30'));
var colorized = evi.select('EVI');
var colorizedVis = {
  min: 0.0,
  max: 1.0,
  palette: [
    'FFFFFF', 'CE7E45', 'DF923D', 'F1B555', 'FCD163', '99B718', '74A901',
    '66A000', '529400', '3E8601', '207401', '056201', '004C00', '023B01',
    '012E01', '011D01', '011301'
  ],
};
Map.setCenter(6.746, 46.529, 2);
Map.addLayer(colorized, colorizedVis, 'EVI');

// Use the expression function to generate scaled EVI.

var scaledevi = evi.map(function (evi) {
  return evi.expression(
    '(evi - 0.1)',{
      'evi': evi.select('EVI')
    });
});

Map.addLayer(scaledevi, {bands: 'EVI', max: 1 , min: 0}, 'scaledevi',false);

//Loading in the Modis LST data
var lst = ee.ImageCollection('MODIS/006/MOD11A1')
                  .filter(ee.Filter.date('2020-08-01', '2020-08-30'));
var landSurfaceTemperature = lst.select('LST_Day_1km');
var landSurfaceTemperatureVis = {
  min: 13000.0,
  max: 16500.0,
  palette: [
    '040274', '040281', '0502a3', '0502b8', '0502ce', '0502e6',
    '0602ff', '235cb1', '307ef3', '269db1', '30c8e2', '32d3ef',
    '3be285', '3ff38f', '86e26f', '3ae237', 'b5e22e', 'd6e21f',
    'fff705', 'ffd611', 'ffb613', 'ff8b13', 'ff6e08', 'ff500d',
    'ff0000', 'de0101', 'c21301', 'a71001', '911003'
  ],
};
Map.setCenter(6.746, 46.529, 2);
Map.addLayer(landSurfaceTemperature, landSurfaceTemperatureVis, 'Land Surface Temperature');

// Use the expression function to generate scaled LST.

var scaledlst1 = lst.map(function (lst) {
  return lst.expression(
    '(lst/30)',{
      'lst': lst.select('LST_Day_1km')
    });

});

Map.addLayer(scaledlst1, {band:'LST_Day_1km' ,max: 1 , min: 0}, 'scaledlst1');

//

var scaledlst2 = lst.map(function (lst) {
  return lst.expression(
    '(lst * (-0.05)) + 2.5',{
      'lst': lst.select('LST_Day_1km')
    });

});

Map.addLayer(scaledlst2, {band:'LST_Day_1km', max: 1 , min: 0}, 'scaledlst2');

//select images

var ListOfImages1 = scaledlst1.toList(scaledlst1.size())
var ListOfImages2 = scaledlst2.toList(scaledlst2.size())

var LST = ee.Image(ListOfImages1.get(0))

var LST2 = ee.Image(ListOfImages2.get(0))

//Take min value
var scaledLST = ee.ImageCollection.fromImages([LST,LST2]).min()

Map.addLayer(scaledLST, {max: 1 , min: 0}, 'Scaled LST');


texts: Opening ArcGIS layer package in open source software?

texts:
It seems that ESRI introduced the layer package in ArcGIS 9.3.1.

I have been provided with layer packages containing essential data, but do not have access to this software.

Ideally, I'd like to be able to open them using some open source tool.  I mainly use R, but Quantum GIS or GRASS would be okay at a stretch.

Suggestions about the best way to access this data are appreciated.


texts: Selecting specific range pixel values using Google Earth Engine

texts:
I have calculated NDWI on my study area. After this, I want to select a specific range of values that I will use to make a mask. My specific values range  is between 0.2 and -0.1.
This is my code:
var sentinel = collezione.filterBounds(area)
                         .filterDate('2020-10-10','2020-10-28')
                         .select('B4','B3','B2','B8','TCI_B','TCI_R','TCI_G')
                         .sort('CLOUDY_PIXEL_PERCENTAGE')
                         .first()
                         .clip(area);

var ndwi = sentinel.normalizedDifference(['B3','B8']);
var palette_NDWI = ['de0000','ff4000','ff7a20','ffc522','fffd30'];
var maschera = ndwi.lt(0.2).and(ndwi.gt(-0.1))
var filtro = sentinel.updateMask(maschera2)

Map.addLayer(filtro)

I don't have my expected mask. What should I do? Is there some operator that can extract my range?


texts: white page border when converting geoPDF to geoTiff with gdal_translate

texts:
I am using gdal_translate to convert some USFS geoPDF maps to geotiff images.
I actually set more layers off than is listed in the code snippet below but I've pulled them for clarity.
The resulting geotiff contains the desired PDF content but also has a white page border. Is it possible to for gdal_translate to do the rendering without the white border?
I've seen a site with a reference to using esri ArcGIS to remove the border after gdal_translate processing but that is very cost prohibitive.
If it is not possible to do with any of the gdal utilities, what is the best way to remove the border? Crop with photoshop and add the geospatial references back?

gdal_translate -of GTiff -co "TILED=YES" -co "TFW=YES" D:\USFSTopo.pdf D:\USFSTopo.tif --config GDAL_PDF_LAYERS_OFF Other_4,Adjacent_Quadrangles,Quad_Location_Diagram,Map.Measured_Grid --config GDAL_PDF_DPI 300


texts: Populate SQL Expression in Toolbox based on hard coded parameter

texts:
I have an ArcGIS toolbox script that applies a user-defined definition query to multiple layers within the current MXD.

The only user variable currently required is the definition query itself, which I have listed as an "SQL expression". The script works fine, but it would be nice if the SQL validation box on the toolbox GUI was populated with the details of the first layer in the MXD, allowing the user to construct and validate a definition query against the first layer.

How can I achieve this using a hard-coded variable pointing to the first layer in the TOC?

It's fairly straightforward in a one-to-one scenario, by asking the user to input the feature layer as a parameter and then making the SQL expression parameter "Obtained from", but I can't figure out how to do it by defining it myself.


texts: Extracting or masking of very large raster to reduce file size on disk?

texts:
I have a very large 8-band raster (339806 x 61704) that I cannot load to memory to manipulate (R, raster:canProcessInMemory indicates I need 4999.16 GB of ram to load this file). However, I only need about 10% of this total area and I have a shapefile of the areas I need to keep.

What is the most efficient way to load my areas of interest to memory?

Masking of the original raster with my shapefile, or extracting my areas of interest and writing them to another raster (I need to retain the spatial reference of all of the sites)? Or some other process? I am new to this and prefer to use R if possible.


texts: Changing color of pin icon marker on CartoDB?

texts:
I am interested to use the pin icon (see attached image) in my cartodb map (cluster type).

Can someone tell me how do I change the color of the pin to reflect categorical data? I tried to click img fill but it the button is not active when I select the pin icon as my marker type.


texts: How to detect geometric network issues (exhaustive check within extent), and repair them using ArcObjects

texts:
I am trying to use this tool which can be found in the ArcToolBox under data management, geometric network

I need to access this tool using ArcObjects because I have loaded features into their respective feature classes, and I would like to run an exauhstive check on their geometric network status using the extent.

I came across IGeometricNetworkConnectivity2, and I can use it on each feature class separately, but it does not allow me to do an exhaustive search or use an extent for that matter, it will run on the entire .
Here a sample of where I am at now

var ftrClass =LoadFeatureClass(shapeFileObject)
var netWorkGeo = ftrClass as INetworkClass;
var nWCheck=(IGeometricNetworkConnectivity2)netWorkGeo.GeometricNetwork
string logger= string.Empty;
ISet set = new Set();
bool consistincy;
nWCheck.CheckAndRepairConnectivity(true, logger, out set,out consistincy,repairConnectivityProg);

I tried using the IGeoProcessor Interface to run the tool, but I am getting an exception eventhough I am passing correct parameters according to the ArcToolBox Help provided.
The exception is raised on Execute.

   try {
            var array = new VarArrayClass();
            var gp = new GeoProcessorClass();
            gp.AddToolbox(@"c:\program files (x86)\arcgis\desktop10.2\ArcToolbox\Toolboxes\Data Management Tools.tbx");
            array.Add("ElectCoilNet");
            array.Add(@"c:\OutLog.txt");
            array.Add("VERIFY_ONLY");

            var result = gp.Execute("VerifyAndRepairGeometricNetworkConnectivity_management", array, null);

        } catch (Exception ex) { Logger.LogException(ex);}


texts:
I am trying to extract zonal statistics from raster imagery using polygons in QGIS 3.2.1, but I only want to use values from those pixels that are completely included in the polygons, ¿how to do this?

Answer for the comment:
I'm agree but I have a problem when I polygonize my raster. This is because similar values form bigger polygons and then I have problems to select by location because the polygon is more big than a pixel –


texts: Mapinfo wor file not working when file paths are updated using python

texts:
We did a data move and so I scripted the updating of the wor, so the old file path is changed to the new file path. It looks like it works perfectly but when I open the wor file it doesn't. When I manually go through and re-assign the file paths and then look at the wor, it is exactly the same as the auto-updated one.

Why is this?

See the wor files in https://drive.google.com/folderview?id=0B1XGDTNNfbMJSjRkcTBrVzdkbEE&usp=sharing

_original : before the script below was run
_script_updated : updated by the script
_fixed : updated manually by assigning the new file name when you open the wor in mapinfo

Script follows

    #!/usr/bin/python

import glob, os, csv

def inplace_change(filename, old_string, new_string):
        s=open(filename).read()
        if old_string in s:
                print 'Changing "{old_string}" to "{new_string}"'.format(**locals())
                s=s.replace(old_string, new_string)
                f=open(filename, 'w')
                f.write(s)
                f.flush()
                f.close()
                text_file = open("Process_worked.txt", "a")
                text_file.write("%s;%s;%s\n" %(worfile,old_string,new_string))
                text_file.close()
        else:
            print 'No occurances of "{old_string}" found.'.format(**locals())
            text_file = open("Process_notfound.txt", "a")
            text_file.write("%s;%s;%s\n" %(worfile,old_string,new_string))
            text_file.close()

#worfile='Bruce Workspace_13-06-2008.wor'
csvfile='File_Locations.csv'
special='Open Table \"'
wor_list=glob.glob('*.wor')

for worfile in wor_list:
    print worfile
    r = csv.DictReader(open(csvfile, "rt"), dialect="excel")
    for row in r:
        old_string=row['Existing File Path']+row['File Name']
        new_string=row['New File Path']+'\\'+row['New Name']
        print old_string, new_string
        inplace_change(worfile,old_string,new_string)
        inplace_change(worfile,row['File Name'],row['New Name'])


texts: Where can I get road data that includes max speed?

texts:
Does anyone know how to detect the road speed

I. E.  If one of our vehicles was driving on a certain road above the speed limit,  we would want to know

But no idea where to get this data from
Any help appreciated
Thanks


texts: Importing a shapefile with no specifed projection to Postgres?

texts:
I am working with a shapefile that does not come with a .prj file - only .shp, .shx and .dbf files. The file's accompanying README says:

ALL DATA IS IN A GEOGRAPHIC COORDINATE SYSTEM (NOT PROJECTED)
USING DECIMAL DEGREES

I want to use shp2pgsql to import the data into Postgres, but I know I should use -s to specify the source SRID.

But what SRID should I be using in this case?

I tried importing it into Postgres without a source SRID, but looking at it in QGIS it seemed to have imported badly.

UPDATE: I tried shp2pgsql -s 4326:5070 and that seems to work OK. I think "geographic coordinate system" might actually mean EPSG:4326 - is that correct?


texts:
I have a polygonal layer, with adjacent and non adjacent territories. Each territorie has a population field. I would like to create a new polygonal layer based on location and attribute : adjacent territories with a population sum under eg 600.000 inhabitants.

I can first select objects with population over eg 20.000 inhabitants. But then?

I tried the "Select by attribute sum" tool, but it didn't work (failure message).

(The colours depend on the number of inhabitants.)


texts: QGIS spatial join - join all features, not just one

texts:
I'm trying to do a spatial join of a feature on a polygon layer (p1) that intersects with two features of another polygon layer (p2)) via QGIS.

What I want: get the attributes of both features of p2 (name as string) into the attribute table of p1, not just the first one. Sum feature just works with int/double and no with strings.

Any suggestions?


texts: Solar potential or irradiation GIS data for India

texts:
I am looking for solar irradiation data for India in a gis format. Like the ones on the link below:
http://www.arcgis.com/apps/OnePane/basicviewer/index.html?appid=bf886ca471a54e2c9127a54f60d4ca27


texts: How do I choose a specific position for ArcMap labels

texts: ArcMap 10.1: Labels on tiles I have set up for data driven pages randomly locate around the page instead of in the centre or in the corner where I want them to be despite allowing overlapping labels and setting labels weight to high. How do I choose a position for these labels pls?


texts: Equalize number of columns and rows number with different raster data which was taken from different sources

texts: I am interested in creating a flood map using machine learning methods. I need the same number of columns and rows for each raster data or layer. However, although I used the same boundaries for cutting raster datasets, I have a different number of columns and rows. By the way, cell sizes are the same for all data.


texts:
I'm using osm2pgsql to import an extract of OpenStreetMap's data so I can geocode an address to a lat/lng. But I'm new to OSM and mapping so I'm having trouble figuring out what data is relevant and what is unnecessary.

What is the bare minimum data I would need to import from an OSM data file in order to geocode an address?

Thanks!


texts: What projection appeared in the NYTimes this morning?

texts:
Perhaps not as GIS related as some, but definitely about maps.
Can anyone identify the map projection that appeared in the New York Times this morning? It's not one I've seen before, and a quick scan through Map Projections: A Working Manual did not yield the result I was looking for.


texts:
In my Openlayers app, since the base map projection is not epsg900913 or 4326, I need to include proj4js in order to use geolocation feature. but I got such an error notice:

 GET http://svn.osgeo.org/metacrs/proj4js/trunk/lib/projCode/null.js 404 (Not Found)

I am wondering how come this null.js is needed here. the funny thing is I can't find this file from proj4js lib files..

anybody can help me out?

Code:

<!DOCTYPE html>

<link rel="stylesheet" href="theme/default/style.css" type="text/css">
<link rel="stylesheet" href="style.css" type="text/css">


texts: Updating a Tool Validator Class from a CSV File

texts:
I have a CSV file which contains a list of clients and different project numbers and project names.  This list is generated periodically.  The list can contain multiple values for client name but unique values in the project number/name fields.  For example, the CSV might look like this:

COMPANY NAME, PROJECT NUMBER, PROJECT NAME

ABC Corp, P123, Mapping Project 1
ABC Corp, P124, Mapping Project 2
Acme, P38, Some Kind of Project 1
BBB Corp, P188, Another Kind of Project
BBB Corp, P189, Another Kind of Project 2
BBB Corp, P190, Mapping Project 1

What I would like to do is create a python tool in ArcGIS 10.2.1 which has two parameters:

Company Name
Project

The Project parameter would be a concatenation of PROJECT NUMBER, PROJECT NAME in the format: PROJECT NUMBER - PROJECT NAME.

I'm trying to program a validation class such that the user is presented with a drop down list of UNIQUE values of the client name.  Then, the second parameter would present only the project numbers associated with that client.  These values would then populate values in a map layout.

So, for example, if the user selects ABC Corp from the first drop down, the second drop down would only have two choices:

P123 - Mapping Project 1
P124 - Mapping Project 2

Similarly, if the user selects BBB Corp from the first drop down, the second drop down would have three choices:

P188 - Another Kind of Project
P189 - Another Kind of Project 2
P190 - Mapping Project 1

I've successfully got the first drop down to create a unique list of companies.  However, I'm having trouble getting the validator class to display properly.  I can generate a list in IDLE correctly, I'm just not sure how I can get it to display in ArcGIS.

My validator code looks like this:

import arcpy, csv
class ToolValidator(object):
  """Class for validating a tool's parameter values and controlling
  the behavior of the tool's dialog."""

  def __init__(self):
    """Setup arcpy and the list of tool parameters."""
    self.params = arcpy.GetParameterInfo()

    csvpath = "S:\projects\P000000 Pre-Sales and Internal GIS Activities\ArcGIS Toolbox Development\Get XIC Information\CSV\ProjectData.csv"
    xicData = []
    with open(csvpath,'rb') as csvfile:
        reader = csv.reader(csvfile, delimiter = ',', quotechar='|')
        for row in reader:
            xicData.append(', '.join(row))

    xicData.pop(0)
    self.xicData = xicData
    pass

  def initializeParameters(self):
    """Refine the properties of a tool's parameters.  This method is
    called when the tool is opened."""

    companyset=set()
    for i in self.xicData:
        companyset.add(i.split(', ')[0])
    companylist=list(companyset)
    companylist.sort()

    self.params[0].filter.list = companylist
    return

  def updateParameters(self):
    """Modify the values and properties of parameters before internal
    validation is performed.  This method is called whenever a parameter
    has been changed."""

    if self.params[0].value:
        projectList = self.params[1].filter
        companyName = self.params[0]
        projectList = []
        xicSplit = []
        for i in self.xicData:
            xicSplit.append(i.split(', '))
        for i in xicSplit:
            if i[0] == companyName:
                projectList.append(i[1] + ' - ' + i[2])

        projectList.sort()

        self.params[1].filter.list = projectList
    return

  def updateMessages(self):
    """Modify the messages created by internal validation for each tool
    parameter.  This method is called after internal validation."""
    return

The problem seems to be in the updateParameters function.  I think my logic is sound, but I don't know where I'm going wrong.


texts:
I have a very simple script which takes in table of X,Y --> Convert to point --> Convert to Shape
centerpoint = r'C:/Users/YYY/pointfile.xlsx|layername=Sheet1'
pointshape = r'C:/Users/YYY/pointshape.shp'

processing.run("native:createpointslayerfromtable", {'INPUT':centerpoint,
                                                    'XFIELD':'POINT_X','YFIELD':'POINT_Y','ZFIELD':'','MFIELD':'',
                                                    'TARGET_CRS': QgsCoordinateReferenceSystem('ESRI:102080'),
                                                    'OUTPUT':pointshape})

processing.run("native:rectanglesovalsdiamonds", {'INPUT':pointshape,
                                                  'SHAPE':0,'WIDTH':Width_of_Bld,'HEIGHT':Height_of_Bld,'ROTATION':Rotation_of_Bld,'SEGMENTS':2,
                                                  'OUTPUT': bldshape})

This script is in a simple for loop over a list. It works fine when it's processing the first item in the list, but when it reaches the second item, I get error:
Application path not initialized
ERROR 1: C:/Users/YYY/pointshape.shp is not a directory.

Which is weird because it's obvious that is name of file to be generated and not a directory. I don't need pointshape so I was hoping it would just simply overwrite it.


texts: Schema problem during tool creation using Append

texts:
I have a model that creates a census tract segment layer by intersecting with a .5 mile buffer around every park in the city. It iterates through each park buffer individually on the GISPROPNUM field, and then appends the census tract segments onto a master shapefile. This model works great but I would like to turn it into a tool that can iterate through any layer and intersect it with any other selected layer in any mxd it opens in.

My problem is that after the first iteration, it stalls at the Append step because the schema does not match the target shapefile (because there is no target shapefile yet). Ideally, I would like the schema of the output to match the source (in this case, the park layer).

SO in a nutshell, at the append step it needs a target shapefile whose schema matches the source (the park buffer layer), but how can there be a target dataset if the tool has never been run within an mxd? How do I create a shapefile after the first iteration that can be used as a target dataset?


texts: How do I update the postgresql database after I've already imported all the data?

texts:
I ran osm2pgsql and it worked great but then I spent a few days adding a whole other PBF only to find that I had forgotten the append flag and I can't render any of the new data.

Does anyone know what command I have to run to fix this without wasting all the CPU cycles and disk space of reimporting?


texts: Can BFAST algorithm be applied to to a whole study area/ ROI?

texts:
I am working on identification of breakpoints in Landsat based NDVI time series (1990-2018). I am using The Breaks for Additive Seasonal and Trend (BFAST) algorithm which is a pixel based time series decomposition method to detect and characterise abrupt changes. I wonder would it be possible to apply it for a region of interest or study area as a whole instead of applying it on selected pixels?

ndvi_stack <- stack(ndvi_list)
ndvi_ROI_1 <- crop(ndvi_stack, extent(ROI_1))
ndvi_ROI_1

So, ndvi_ROI_1 is a study area or region of interest. Now I select the pixel:

selected_pixel <- 100

I created an irregular time series for the selected pixel (100 in this case) using bfastts() from BFAST package.

(s <- bfastts(as.vector(ndvi_ROI_1[selected_pixel]), dates, type = c("irregular")))

The time series is converted to a regular time series using advice from https://philippgaertner.github.io/2018/04/bfast-preparation/.

I also tried not to select to a specific pixel and run bfastts() on whole study region and it seems to work (although I am not sure about the validity of results). I used the following code:

(s <- bfastts(as.vector(ndvi_ROI_1), dates, type = c("irregular")))

I am still not sure if this is technically right to apply bfast method to a whole study region instead of pixel by pixel basis? And if it is, is the coding correct or I have to introduce any other steps?

Looking forward to the kind advice.


texts:
In GeoNetwork 3.4, the default format for Metadata information "Contact" is the American Format i.e. street, nummer, city, Postal Code, Country

example: Rue Vautier 29, Bruxelles, 1000, Belgium

However in Europe (in Belgium at least) it should be:

Rue Vautier 29, 1000 Bruxelles, Belgium

Otherwise the address is not valid. Do you know how to modify that?

The metadata files are harvested on a regular basis from a WAF so I can't (and don't want to) hard-code it.

Here is an example of a MD file

http://geonetwork.bmdc.be/geonetwork/srv/eng/catalog.search#/metadata/bmdc.be:dataset:2023


texts: Trying to publish a feature layer with a very basic symbology; no fill, just a basic hatch to match the outline color. Adding the layer to a new web map in AGOL keeps the outline color but changes the hatch color to black. I have other polygons that used the same symbology conventions and their hatch colors were all changed to black, too. Is there a way around this, so that the hatch color and outline color match, like they do in ArcMap before publishing? I've tried publishing to Portal and out of Pro with the same results.]2


texts: gdal2tiles tile origin differences using Windows/OSGeo4W and Linux/GDAL

texts:
I am getting different X and Y values for tiles from gdal2tiles.py when in a PC
(OSGeo4W - GDAL 1.9.2) environment and Linux (GDAL
1.10.1).  It seems to me that the newer 1.10.2 x's and y's are wrong (not respecting the specified origin).

Is this a bug, which result is correct, and how do I get them to
behave the same?

===

Details:
I've created a sample tileset on a PC in windows using the
tools from GDAL 1.9.2 in OSGeo4W for a map in projection
NZTM2000 (EPSG:2193).  The resulting tiles are
generated as-per my expectations, and the origin from which
the tiles are numbered is the origin of the projection
(+lat_0=0 +lon_0=173).

Having perfected this with a small map area, I switched to a
cloud server running linux and gdal1.10.1.  The process
appears to run successfully - tiles are generated, BUT the
tile x and y numbers are different.

After much stuffing around, I've realised that the x and y
values used match those used for the same geographic area
when using the default EPSG:4326.

During troubleshooting I have run the gdal2tiles on the same
single map sheet with the same command line options in both
environments, manually specifying the SRS in case the
definitions differ, and got two different results.
I've also tested it using the built-in projection -s
EPSG:2193.

python gdal2tiles.py test.tif -s '+proj=tmerc +lat_0=0
+lon_0=173 +k=0.9996 +x_0=1600000 +y_0=10000000 +ellps=GRS80
+towgs84=0,0,0,0,0,0,0 +units=m +no_defs' -z 5-10

PC / Windows / Cygwin:
mbriggs@forestero
/d/OSGeo4W/maps_past $ ls test/11/
1116/  1117/  1118/
mbriggs@forestero
/d/OSGeo4W/maps_past $ ls test/11/1116
1303.png*  1304.png*
mbriggs@forestero
/d/OSGeo4W/maps_past $

Linux:
mbriggs@dev:~$ ls
test/11
2022  2023  2024  2025
mbriggs@dev:~$ ls
test/11/2022
767.png  768.png  769.png
mbriggs@dev:~$

mbriggs@dev:~$ gdalinfo
test.tif
Driver: GTiff/GeoTIFF
Files: test.tif
Size is 4193, 2861
Coordinate System is:
PROJCS["unnamed",
    GEOGCS["GRS 1980(IUGG, 1980)",
        DATUM["unknown",

SPHEROID["GRS80",6378137,298.257222101]],
        PRIMEM["Greenwich",0],

UNIT["degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",173],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",1600000],
    PARAMETER["false_northing",10000000],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]]]
Origin = (1810357.092181787360460,5488467.918922254815698)
Pixel Size = (10.000000000000000,-10.000000000000000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  ( 1810357.092, 5488467.919)
(175d29'27.21"E, 40d43'40.29"S)
Lower Left  ( 1810357.092, 5459857.919) (175d30'2.02"E, 40d59' 7.25"S)
Upper Right ( 1852287.092, 5488467.919) (175d59'12.51"E, 40d42'57.90"S)
Lower Right ( 1852287.092, 5459857.919) (175d59'54.22"E, 40d58'24.47"S)
Center      ( 1831322.092, 5474162.919) (175d44'39.06"E, 40d51' 3.45"S)
Band 1 Block=4193x1 Type=Byte, ColorInterp=Red
  Mask Flags: PER_DATASET ALPHA
Band 2 Block=4193x1 Type=Byte, ColorInterp=Green
  Mask Flags: PER_DATASET ALPHA
Band 3 Block=4193x1 Type=Byte, ColorInterp=Blue
  Mask Flags: PER_DATASET ALPHA
Band 4 Block=4193x1 Type=Byte, ColorInterp=Alpha
mbriggs@dev:~$


texts: Is it possible to to have higher variances for predictions in Universal Kriging?

texts:
I am using gstat to do Universal kriging.

I am finding that the variances for the  predictions of a hold out sample (i.e., new spatial locations for which I want to determine the response variable) are extremely large.

Infact, they are larger than my actual predictions. How is this possible?
Unless I am missing something , aren't the variances a representation of the confidence levels of the predictions ?


texts: I'm using the v.to.lines algorithm (GRASS algorithm) to convert a layer from polygons to lines in QGIS. The algorithm works without problems but the output layer only has an attribute called "cat", I would like to know if it is possible to inherit the attributes of the input layer.


texts: Deauthorizing ArcGIS for Desktop 60 Day Trial

texts:
I installed ArcGIS For Desktop Advanced 60 Day Trial. I'm not sure what went wrong, but I didn't get a clean install. I uninstalled it.

I later installed ArcGIS For Desktop, and went through the authorization wizard to configure it as a Single User license under my company's organization.

Everything looked like it worked fine, except that when I try to run ArcMap, it loads the splash screen, sits for 30 seconds, then disappears.

When I run ArcGIS Administrator, I see Desktop. I can run through the Authorization Wizard, and it tells me everything is authorized. But when I look at the availability list, everything is set to expire in June, 60 days after I installed the free trial, not in November when the company's license expires.

If I go in to Deauthorize, I see two items:

ArcGIS for Desktop Advanced 60 Day Trial - expires 27-Jun-2016
EDN ArcGIS for Desktop Basic - expires 01-Dec-2016

If I try to deauthorize the free trial, I get an error message: "7466: The return of the fulfillment is denied by the return policy because max return exceeded [Incident# 8414-03497135]"

I've uninstalled all the ArcGIS products, walked through and removed all the ArcGIS registry settings, reinstalled the EDN Desktop Basic, and I still see the free trial showing up in the ArcGIS Adminstrator lists.

Am I right in supposing that it's the presence of this aborted free trial install that's keeping my fresh install from working?

And how can I get rid of it?


texts: arcpy.Exists not working within standalone script

texts:
I have the following code. I want to see if a layer name exists within a map document, not if the shapefile exists in the directory. If I copy paste into the python window within ArcMap works fine. IE. Prints: True

If I run the program as a stand-alone script it does not work?  Prints False when should print True
I have code the sets up the mxd and dataframe.

if arcpy.Exists("project_aoa") == True:
    print("TRUE")
elif arcpy.Exists("project_aoa") == False:
    print("FALSE")

These layers are withing a group layer should this matter?

project_aoa is a shapefile.

EDITED:

I ended up doing this:

for item in mxds:
    mxd = arcpy.mapping.MapDocument(r"{0}".format(item))
    df=arcpy.mapping.ListDataFrames(mxd,"Project Area")[0]

    for lyr in arcpy.mapping.ListLayers(mxd, "project_aoa", df):
        if lyr.name == "project_aoa":
            print("True")
        else:
            print ("False")


texts:
We're trying to enable enterprise geodatabase functionality in an existing PostgreSQL database. We've already created the sde user (set as superuser) & sde schema.

We get the following error:

  "Could not create geodatabase tables and stored procedures. Failed to
  execute (EnableEnterpriseGeodatabase).

the database is: "PostgreSQL 9.2.9 on x86_64-unknown-linux-gnu"

client: ArcGIS 10.2


texts:
I'm quite new to geotools and I'm trying to get to grips with it.
Trying to load some shapefiles (see attachment) into a MapContent which in turn
is associated with a SwtMapPane I always get an exception that suggests that
there's something wrong with the current CRS of the map pane:

java.lang.IllegalArgumentException: Argument "sourceCRS" should not be null.
        at
org.geotools.referencing.operation.AbstractCoordinateOperationFactory.ensureNonNull(AbstractCoordinateOperationFactory.java:734)
        at
org.geotools.referencing.operation.BufferedCoordinateOperationFactory.createOperation(BufferedCoordinateOperationFactory.java:249)
        at org.geotools.referencing.CRS.findMathTransform(CRS.java:1205)
        at org.geotools.referencing.CRS.findMathTransform(CRS.java:1173)
        at org.geotools.swt.SwtMapPane.setCrs(SwtMapPane.java:436)
        at org.geotools.swt.SwtMapPane.layerAdded(SwtMapPane.java:812)
        at org.geotools.map.MapContent.fireLayerAdded(MapContent.java:468)
        at
org.geotools.map.MapContent$LayerList.addIfAbsent(MapContent.java:1063)
        at org.geotools.map.MapContent.addLayer(MapContent.java:369)
        at org.geotools.JSurveyMvn.App.fileOpen(App.java:87)
        at org.geotools.JSurveyMvn.App$1.widgetSelected(App.java:113)
        at
org.eclipse.swt.widgets.TypedListener.handleEvent(TypedListener.java:234)
        at org.eclipse.swt.widgets.EventTable.sendEvent(EventTable.java:84)
        at org.eclipse.swt.widgets.Widget.sendEvent(Widget.java:1258)
        at org.eclipse.swt.widgets.Display.runDeferredEvents(Display.java:3540)
        at org.eclipse.swt.widgets.Display.readAndDispatch(Display.java:3161)
        at org.eclipse.jface.window.Window.runEventLoop(Window.java:825)
        at org.eclipse.jface.window.Window.open(Window.java:801)
        at org.geotools.JSurveyMvn.App.<init>(App.java:47)
        at org.geotools.JSurveyMvn.App.main(App.java:58)

When trying to explicitly set a CRS directly after associating the (empty)
MapContent to the map pane (see line 149), I'll get a NullPointerException like that:

java.lang.NullPointerException
        at org.geotools.swt.SwtMapPane.setCrs(SwtMapPane.java:433)
        at org.geotools.JSurveyMvn.App.createContents(App.java:152)
        at org.eclipse.jface.window.Window.create(Window.java:431)
        at org.eclipse.jface.window.Window.open(Window.java:790)
        at org.geotools.JSurveyMvn.App.<init>(App.java:47)
        at org.geotools.JSurveyMvn.App.main(App.java:58)

I'm quite desperate of being unable to solve that problem. Can you please give  me a hand to make it work?

Obviously it's not possible to attach files here. So here's my original post to the Geotools-gt2-users list which seems to be inactive but allows adding code and map samples: http://sourceforge.net/p/geotools/mailman/message/34635944/ Attached is the source code of an example application that demonstrates the  problem (run under linux (amd64) using getools v. 13.3 or 14.0) as well as the shapefile.


texts: Display image from attribute table in infotemplate

texts:
I have added image in attribute table. It looks like this:

Now I want to display this image into InfoTemplate.
so far I am doing like this:
var content = "Name: ${name}" +
                  "
Picture: ${Pic}" +
                  "
Last Name: ${lname}" ;
    var infoTemplate = new InfoTemplate("${FIELD_NAME}", content);

But it does not display image. It is empty. Image is showing in identity window in ArcGIS, but not in the InfoTemplate.


texts: How to down load modis burned area extraction product . I found a link show the product but not found where can download the product geotiff image http://modis.gsfc.nasa.gov//gallery/individual.php?db_date=2015-07-30


texts: Continuously changing values in the same band of a multiband raster

texts:
I am having issues with continuously changing the values of pixels within a single band of a multiband raster using python 2.7. Essentially, what I am trying to do is continuously change (increment) the pixel values of a single band based on variables drawn from 2 lists that are used in the nested loop. I've included some code to help understand what I'm doing.

import arcpy
import os
from arcpy.sa import *
import math

arcpy.CheckOutExtension('spatial')

ConstantList = [2.1, 1.6]
TimeStep = ["4h", "8h", "12h", "16h", "20h", "24h"]

multiband = "Path:\\to\\my\\multibandraster.tif\\"

arcpy.env.workspace = multiband

bands = [Raster(os.path.join(multiband, b)) for b in arcpy.ListRasters()]
count = 1
for constant in ConstantList:
    for time in TimeStep:
        temperature = Raster("Path:\\to\\temperatureRasters\\" + time + "_1_1971.tif"

        bands[1] = Con(((temperature > 2.5) & (temperature < 35)), bands[1] + some equation, bands[1])

    bands[1].save("Path:\\to\\output" + str(count) + ".tif")
    count = count + 1

arcpy.CheckInExtension('spatial')

Using this code, I get a raster with all pixels = NoData. If I remove the 'time' for loop (for time in TimeStep:....) I can get the correct values, but only for one timestep. Why would including the 'time' for loop give me a raster with NoData, and how can I fix it?

Also, I should mention that when I remove the 'time' for loop, the values of bands[1] during the execution of the second constant in the constant list, are wrongly added to bands[1] values calculated during the execution of the first constant. At the start of the second constant, the values of bands[1] should be reset to the original values of the multiband raster.


texts: Creating new PostGIS table with records + geometries from another table?

texts:
I have a bunch of PostGIS tables (one per county) containing polygons for all of the property parcels in that county. Each layer also includes several attributes like street address, property value, and acreage that I'm interested in keeping. Each layer also has several dozen columns that I don't need and want to discard. I want to merge all of the data I'm keeping from these tables into a single PostGIS table parcels with an additional county column, that looks like this:

I understand how to use INSERT/SELECT statements to take values from one table and insert into another, but I am confused about exactly what I need to do to copy the geometries from each county's layers along with the attribute data.
Am I supposed to just copy the geom column along with the columns I'm interested in keeping? In addition to the geom field, several of the county layers have columns called x_coord, y_coord, and others have columns with names like shape_area, shape_len, shape_star, etc ... I'm wondering if I need to copy any of these over along with the geom field? Are any of these fields related to how PostGIS stores geometry data, or are these all just useless remnants from the ArcGIS tools that the counties used to generate these tables?
What exactly is the process for copying attributes + geometry data from one PostGIS table to another? I'm using QGIS and would like to know how to do it with DB manager, if possible, but I'd also like to know how to do if from the psql command line if anyone knows how.
Just to have a concrete example, let's say I have a table called jefferson_parcels that contains 40+ attribute fields, but I only want to copy the fields situs_addr, situs_city, cfmv_total, ttl_acres, land_use from jefferson_parcels and INSERT them into the parcels table in the fields parcel_addr, parcel_city, total_value, total_acres, land_use respectively, and set the county field equal to jefferson. Along with all of this attribute data, I want to copy over the geometry data for these parcels. How would I do this?


texts: How to remove a parcel from a larger area?

texts: New QGIS user here.  I'm trying to clip a large (2000-ac) property boundary with a smaller (500+/-) acre property.  Assuming the larger boundary is the input layer and smaller is the clip layer, my output layer then shows up the exact same 500 acres boundary I'm trying to clip with - which is clearly wrong (the larger tract is not clipped).  Seeing the error, I reversed the input and clip files - switched the 500 acre tract to input and the 2000 acre tract as the clip - and the clip file then shows up as the exact same 500 acre boundary as the first failed clip, again wrong.  WTH is wrong here?  How does the clip boundary not cut the input boundary?


texts:
I have a shp with size 673'476 KB, with shx size 192'422 KB (as shown in the picture below). This file is with points.
Is there any way that can I reduce the shx file size using QGIS?


texts: Determining if point on boundary using shapely?

texts:
I'm new to Python, and I want to test if a certain point lies ON the (boundary) of a polygon - not inside, not outside, just if it's on the boundary.
The points and polygon already defined

polygon = [(3, 2), (5, 1), (7, 2), (8, 6), (9, 7), (11, 6), (10, 4), (12, 2), (15, 2), (16, 5), (19, 7), (18, 11), (14, 12), (11, 9), (5, 9), (2, 6)]
Point_X = 14
Point_Y = 12

this point should be a boundary.
How can I do that?


texts: Failed to update point geometry through WFS-T request when sending it with python

texts:
I am trying to update some data on postgis through a WFS request to geoserver.
More specifically I want to update the geometry of a point.
My XML looks like this and works correctly when I execute it through geoserver Demos section:

        <wfs:Transaction service="WFS" version="1.1.0"
          xmlns:ogc="http://www.opengis.net/ogc"
          xmlns:wfs="http://www.opengis.net/wfs"
          xmlns:gml="http://www.opengis.net/gml">
          <wfs:Update typeName="geonode:idn_poi_tradecenter_sumba_rbi25k_big_a">
            <wfs:Property>
              <wfs:Name>the_geom</wfs:Name>
              <wfs:Value>
                  <gml:Point>
                    <gml:pos>-9.37861241467984 120.9814453125</gml:pos>
                 </gml:Point>
                </wfs:Value>
            </wfs:Property>
            <ogc:Filter>
              <ogc:FeatureId fid='idn_poi_tradecenter_sumba_rbi25k_big_a.11'/>
            </ogc:Filter>
          </wfs:Update>
        </wfs:Transaction>

When I execute exactly the same xml using python requests, the point is not updated. Although I do get a status_code = 200 back when I print the response of the request I get this:

    <?xml version="1.0" encoding="UTF-8"?><wfs:TransactionResponse xmlns:xs="http://www.w3.org/2001/XMLSchema" xmlns:ogc="http://www.opengis.net/ogc" xmlns:wfs="http://www.opengis.net/wfs" xmlns:geonode="http://www.geonode.org/" xmlns:ows="http://www.opengis.net/ows" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:gml="http://www.opengis.net/gml" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" version="1.1.0" xsi:schemaLocation="http://www.opengis.net/wfs http://localhost:8080/geoserver/schemas/wfs/1.1.0/wfs.xsd">
<wfs:TransactionSummary><wfs:totalInserted>0</wfs:totalInserted>
<wfs:totalUpdated>1</wfs:totalUpdated>
<wfs:totalDeleted>0</wfs:totalDeleted></wfs:TransactionSummary>
<wfs:TransactionResults/>
<wfs:InsertResults>
<wfs:Feature>
<ogc:FeatureId fid="none"/>
</wfs:Feature>
</wfs:InsertResults>
</wfs:TransactionResponse>

If you notice in the response its mentioned:

<ogc:FeatureId fid="none"/>

Does anyone has an idea what's going on?

EDIT

I found out what was the problem. It was related with encoding issues and hidden whitespaces or something like this.
The issue (and major mistake) was that I have copied my XML structure to my email and then I took it from there. Sometimes this cause trouble such as adding hidden spaces etc.

I had to retype the XML string and then it worked.
What made me suspect this issue was this line in geoserver log:

Content is not allowed in trailing section


texts:
Could anyone suggest a script/tool/software, that would:

recognize and
quantify the similarity/dissimilarity between the four attached images?

Meaning, it would show that there is a greater similarity between the top two images (images B1 and B2) and between the bottom two images (images A1 and A2) compared to the similarity between either of the top two images and either of the bottom two images. Meaning, for example, that B1 is more similar to B2 than it is to either A1 or A2. The similarity that is referred to here is the overall visual similarity between the images. When I look at those images, it is obvious to me that the top two are more similar between themselves  than they are to the bottom two images and I would assume that the tool that would recognize similarities would deal with the spectral, texture, pattern, etc., characteristics of the images. Presumably, it wouldn't work with entire images at once but would allow for a window of a certain size in pixels to be passed over the images. For example, depending on the pixel size, the window in these examples (these are trees, black spruce in the top two images and mostly poplar in the bottom two images), would be large enough to capture the 'minimal unit' of difference, i.e., several trees together.  Is there an approach/tool that can quantitatively express such a similarity?

I can use Python and I have tried SSIM (Structural Similarity Index Measurement) but that didn't produce satisfactory results.


texts: How to add tiger schema to the database search path?

texts:
The tiger schema must be added to the database search path for the functions to work properly.

I just started recently to use tiger dataset and postgres database and I have basic understanding of both. I would like to fix the following:

loader_generate_census_script: command not found


texts: What do cartocss attachments along filters mean?

texts:
I'm exploring openstreetmap-carto style and I wonder what does such syntax mean:

  [natural = 'water']::natural,
  [landuse = 'reservoir']::landuse,
  [waterway = 'riverbank']::waterway

Is this declaring multiple attachments at once, or is it redirecting filters somehow?
I'm clueless.


texts: Getting print statements in Python script to print through a batch file

texts:
I'm starting to learn Python and batch files. I have a Python script which selects data by location and exports the selection to a new feature class. The script runs just fine - included in this script I have the following lines:

matchCount = int(arcpy.GetCount_management(featureLayer).getOutput(0))

print str(matchCount) + " rows exported"

When I run this script in PythonWin the print statement above is produced in the Interactive Window. When I run this script through a batch file the print statement isn't produced. At the end of my batch file I have written PAUSE so the screen doesn't disappear after running my script. Anyone have ideas on how to produce my print statement in command line?


texts: Python script to update attributes crashing during loop

texts:
I have a script to update a field in one feature class from a selected record in another feature class. This is done through looping through records in one feature class selecting appropriate the record in the other then reading the ID and updating the first feature class. I ran some tests on small subsets of the data and it worked fine. When i run across the full data set the MXD crashes (there are 30,000+/- records in the first feature class and 80,000 in the second). I am assuming this has something to do with memory, but I am not sure. Any help in solving this issue or improving the efficiency of the script (included below) is much appreciated.
thanks.

import arcpy
cur = arcpy.UpdateCursor("layer_1")

for i in cur:
    st=i.getValue("Street")
    if st[-4:] == 'PKWY' or st[-4:] == 'EXPY':
        pass
    else:
        qry = "Street = '"+st+"'"
        arcpy.SelectLayerByAttribute_management("TAMED_BUF","NEW_SELECTION",qry)
        arcpy.SelectLayerByLocation_management("TAMED_BUF","CONTAINS","layer_1",1,"SUBSET_SELECTION")
        cur2=arcpy.SearchCursor("TAMED_BUF")
        for row in cur2:
            seg=row.getValue("SegmentID")
            i.T_SegmentID = seg
            cur.updateRow(i)


texts: Programmatically export OpenLayers map with maximum zoom/highest resolution

texts:
I have a problem to export a OpenLayers map to a png or jpeg file.

What I want is the full map and with maximum resolution, what I get is a part of the map with minimum resolution.

The map is hosted on a server on the Internet.
Here is the link to call the page:
http://bvbm1.bib-bvb.de/webclient/DeliveryManager?custom_att_2=simple_viewer&pid=233679

Here is The code (I found it in another thread, but I don't know how to zoom or how to set the size of the map.

// modify these parameters accordingly!
var url = 'http://digital.bib-bvb.de/webclient/DeliveryManager?application=&owner=&custom_att_1=&custom_att_3=&user_mode=&pds_handle=GUEST&search_terms=&adjacency=&frameId=1&user=&pid=233679&locale=de&zoom=100'
var mapID           = 'map';
var outputFilename  = 'hurrz.png';

var resourceWait  = 300,
    maxRenderWait = 10000;

var page          = require('webpage').create(),
    count         = 0,
    forcedRenderTimeout,
    renderTimeout;

// function to execute map clipping and saving when page has fully loaded
function doRender() {

    console.log("Page has fully loaded");
    console.log('clipping map...')

    var clipRect = page.evaluate(function(mapID){
        console.log('Map Id');
        console.log(mapID)
        return document.getElementById(mapID).getBoundingClientRect();
        }, mapID);

    page.clipRect = {
        top:    clipRect.top,
        left:   clipRect.left,
        width:  clipRect.width,
        height: clipRect.height
        };

    console.log('rendering image...')
    console.log('top, left, width, height:')
    console.log(clipRect.top)
    console.log(clipRect.left)
    console.log(clipRect.width)
    console.log(clipRect.height)

    page.render(outputFilename);
    phantom.exit();

}

page.onResourceRequested = function (req) {
    count += 1;
    //console.log('REQ > ' + req.id + ' - ' + req.url);
    //req.url = req.url.replace("zoom=0&", "zoom=1&");
    //req.url = req.url.replace("zoom=1&", "zoom=100&");
    console.log('REQ > ' + req.id + ' - ' + req.url);
    console.log(' ');
    clearTimeout(renderTimeout);
};

page.onResourceReceived = function (res) {

    //res.url = res.url.replace("zoom=0&", "zoom=1&");
    //res.url = res.url.replace("zoom=1&", "zoom=2&");

    if (!res.stage || res.stage === 'end') {
        count -= 1;
        console.log('REC  ' + res.id + ' ' + res.status + ' - ' + res.url);
        if (count === 0) {
            renderTimeout = setTimeout(doRender, resourceWait);
        }
    }
};

page.open(url, function (status) {

    if (status !== "success") {
        console.log('Unable to load url');
        phantom.exit();
    } else {
        forcedRenderTimeout = setTimeout(function () {
            console.log('Count');
            console.log(count);
            doRender();
        }, maxRenderWait);
    }
});


texts:
If my computer has only MS Visual Studio installed but not ArcGIS, can I test and debug the addin?

I know MS Visual Studio 2010 Pro can be used to create ArcGIS addin but it is not free. Can I use the free version of MS Visual Studio to build addin? If yes, may you suggest one?


texts: How to make rules so two QGIS contour maps don't cross

texts:
I'm creating multiple contour maps stacked on one another over one geographical area based on the depths of various attributes. In some cases, given the mass of data and occasional inaccurate data, the contour maps will 'cross' each other. For example, contour A should always have lower depth than contour B, but because of inaccurate data point A1, contour A has a greater depth than B at point A1.

It becomes too much work to just exclude certain data points by hand, as there are so many. Any ideas on how to set a rule where if contour A crosses contour B, the depth for A is automatically set to a certain distance above B? Or perhaps a way in QGIS/SQL to systematically exclude data? Note that only one depth reading accompanies a data point, which accompany a predetermined categorization for contour A or contour B.


texts: Is it possible to try Esri Developer Network (EDN) without purchasing license?

texts:
Albeit using QGIS for most of my GIS tasks, I still use ArcMap for final map production. I am curious and would like to test out EDN. According to Esri the annual subscription is $1,500. I'm not sure I'm ready to invest that much money into it as I'm not sure I'll even use it.

Do I have any options to develop on the ArcGIS platform without subscribing to EDN?


texts: Configuring WMS from a MSSQL data source with MapGuide Open Source

texts:
I've having trouble configuring WMS out of MapGuide Open Source. Here's what I've done so far:

Installed MapGuide
Downloaded MapGuide Maestro
Created a "feature" in Maestro that connects to MSSQL2008
Added the "feature" to a new layer
Added the new layer to a new map
Combined those three things into a "package", imported that into Site Administrator.

Thus, MapGuide has access to the data, and WMS in enabled. Yet, doing a GetCapabilities request returns a 404 out of my browser and a result but no information about the data from the MapGuide UI.

How do I "point" the WMS service to my MSSQL data?


texts: Extracting raster pixel values by polyline

texts: I am using ArcMap 10.3.1 to try and extract the values of raster pixels that individual polylines intersect. I have tried using extract by mask but the wide mask created for each polyline intersects a greater number of raster pixels than just the thinner polyline does, which is unfortunately not very useful for me. Are there any ways to achieve the desired result?


texts: QGIS Centre Pivot Wheel Vector Height Intersection With Raster

texts: I am trying to visually display a centre pivot irrigator against a DTM Raster file to check terrain clearances in QGIS.  I have a 1m DTM Raster and polygon/and line vector files of the pivot wheels.  Each pivot wheel runs approximately 45m apart there is a 2.5m tower at each wheel point with a straight arm connecting the two.  So at each wheel point as the pivot moves around the wheels will follow the contour/DTM but I want to check that between the two (straight line) there is no intersection with the ground/DTM.  I am struggling to figure out what to convert my vector data to and how to process but I feel there must be a way.


texts: How to display QGIS plugin panel "OSM place search"?

texts: It is possible to hide the panel for the experimental QGIS plug-in 'OSM place search' - but the panel is not listed (to re-display) under the menu item View | Panels. How can it be re-displayed?


texts:
My QGis Difference tool doesn't work really well, it gives me lots of little squigly lines. i read in the link below there is such a thing as a vector diffence tool from SAGA, but i can't find it.

Whats wrong with the difference tool QGIS 2.14.1?

Does anyone know where i can find the SAGA vector difference tool, or other ways to calculate the difference between 2 polygons?

i'm using QGis 2.6.0


texts:
My GeoServer for geOrchestra SDI has been recently upgraded from v2.12 to 2.17.4-georchestra by the instance admin.
I had a weird bug now: any text search ends with a false statement in the search bar and either nothing is filtered (workspaces) or nothing remains (stores, layers, and styles).
Did someone already face this?


texts: Geocode using layer of latitude and layer of longitude using R

texts:
I am working in R.
I need to geocode a raster image having a layer of latitude and a layer of longitude (that i’m importing with “raster” package).
To be more precise, I have the first layer with an image (for example 2000x2000) that is not geocoded.
I also have two more layers (raster layers) with the same dimensions (2000x2000) of the the first layer (the image) but one of these just cointains latitude and the other one just contains longitude.
I would like to use latitude and longitude layers to georeference the image layer.
In few words I need to do in R what in ENVI is called GLT (Geographic Lookup Table).
Can please someone help me?