texts: Batch conversion of polygonz into polygon in QGIS?

texts:
I have a bunch of kmz files which I need to merge into one so I can convert it into one polygon shapefile using QGIS.

The problem is that I have some kmz that are polygon, and some that are polygonz and QGIS wont allow the merging of different types of geometry.

So, how can I convert all those polygonz at once (batch) into polygon? (without the z attribute)?

(I know I can save one by one as shapefile than change the geometry type, but I wonder if there is a faster way since I have many layers)


texts:
I use OpenLayers 2.13.1 along with GeoServer and PostgreSQL

I have a code for a vector layer in my OpenLayers like so

//style for layer
var style_d = new OpenLayers.Style({
 'fillColor': '${p_color}',
 'fillOpacity': .8,
 'fontColor': 'black',
 'pointRadius': 8,
 'strokeColor': 'black',
 'strokeWidth': 1,
 'fontFamily': 'Arial',
 'fontSize': 10,
 'label': '${p_name}',
 'labelAlign': 'tc',
 'labelSelect': true,
 'labelXOffset': 4,
 'labelYOffset':4,
 },
 {
    context: {
        label:function(feature) {
            if(feature.attributes['p_name']) {
                return feature.attributes.p_name;
            } else {
                return '';
            }
        }
 },
});

//style for layer - when click a geometry
var style_select = new OpenLayers.Style({
 'fillColor': '${p_color}',
 'fillOpacity': .4,
 'fontColor': 'black',
 'pointRadius': 5,
 'strokeColor': 'black',
 'strokeWidth': .5,
 'fontFamily': 'Arial',
 'fontSize': 10,
 'label': '${p_name}',
 'labelAlign': 'tc',
 'labelSelect': true,
 'labelXOffset': 4,
 'labelYOffset':4,
 },
 {
    context: {
        label:function(feature) {
            if(feature.attributes['p_name']) {
                return feature.attributes.p_name;
            } else {
                return '';
            }
        }
 }
});

//Add a select feature control for pins
var select_feature_control = new OpenLayers.Control.SelectFeature(mylayer, {multiple: false,toggle: true});

map.addControl(select_feature_control);

//Activate the control
select_feature_control.activate();

//do stuff when clicked
function selected_feature(event){
 document.getElementById("wesmes").innerHTML="Loading";
 sf = event.feature.fid;
 var sfs = sf.split(".");
 var sff = sfs[1];
 //call a function
}

Problem :

If there is a polygon or a line on the map I can click right on them and everything works.

If there is a point, I have to click outside the point, near at it, so it will be clicked and activate the function selected_feature(event) etc etc. If I click on the point, nothing happens

How can I fix this? Why is this happening? I think it is a problem of style. If you want I can provide more code


texts: ST_GeomFromText fails on syntax error inside tuple of coordinates?

texts:
I am receiving the following error when trying to execute ST_GeomFromText(wkt)

  psycopg2.ProgrammingError: syntax error at or near "1187780.76" LINE
  1: ...0822.05,(ST_GeomFromTEXT(POLYGON ((7458952.505715 1187780.76...

The project interpreter displays an arrow pointing at the 1 in 1187780.76 to show where it runs into the "error".

Here is the code

def import_shape(db_conn, db_cur, file_input, SRID, output_name):
     connection = db_conn
     cursor = db_cur
     cursor.execute("""DROP TABLE IF EXISTS {output_name}""".format(output_name=output_name))
     srcFile = file_input
     shapefile = osgeo.ogr.Open(srcFile)
     layer = shapefile.GetLayer(0)
     for i in range(layer.GetFeatureCount()):
         feature = layer.GetFeature(i)
         name = feature.GetField("CTUID").decode("Latin-1")
         wkt = feature.GetGeometryRef().ExportToWkt()
         cursor.execute("""INSERT INTO {output_name} (name, geom) VALUES ({name},(ST_GeomFromTEXT({outline})))""".format(output_name=output_name,name=name,outline=wkt,srid=SRID))

     connection.commit()

I have also tried ST_GeomFromEWTK and it returns the same error.

I am using PostGIs 2.4.4 Postgres 9.6 Geos 3.6.2


texts: ST_Intersects returning incorrect result on routing problem?

texts:
I have a table that contains the shortest path between destinations for survey respondents on an osm2po created network.  When calculating the shortest path result, I store the geometry in a column.  Now, using the same osm2po network table and each shortest path geometry, I am trying to determine whether the shortest path crosses or intersects a highway.  I created another table that just contains the osm2po network with clazz 11-16 and run an st_intersects(shortestpathgeom,highwaygeom) and get all false results.  The shortest path geometry is multilinestring and the highway geometry is linestring, both are SRID 4326 and have been created from the same table.  Why am I getting no intersections?  Here is a picture of the shortest path in purple, highway in green, osm2po in red.

here is my code where po_routing.route is the shortest path, po_highway.the_geom is highway only network, and obstacle is the intersection test.

update analysis.po_routing
set obstacle = st_intersects(route,b.the_geom)
from analysis.po_highway as b;

I don't understand why this fundamental function is not working as expected.


texts: Using OpenStreetMaps (*.pbf) data with Mapnik?

texts:
I am using mapnik to directly render a map using xml sheet.

I downloaded world boundaries, created the mapnik stylesheet from wiki page, and ran the python code which created a world image :

Now I want to create one of my college.

Would .osm.pbf format work to include its data?


texts:
I have created a table, dc_crimes, that has eight columns. I would like to fill the geom column with point data created from the lat and long column. The syntax I used to create the table is:

CREATE TABLE tmorr.dc_crimes
(
gid serial NOT NULL,
shift character varying(80),
method character varying(80),
offense character varying(80),
ward character varying(80),
latitude numeric,
longitude numeric,
geom geometry(Point,102285)

I then inserted data from a previous table.
Syntax used:

INSERT INTO tmorr.dc_crimes (gid,shift,method,ward,latitude,longitude)
SELECT gid,shift,method,ward,latitude,longitude FROM tmorr.dc_crimes_v

This is where I am a bit lost. I tried to update the dc_crimes table and use the ST_MakePoint syntax to create points. Syntax used:

UPDATE tmorr.dc_crimes SET geom = ST_MakePoint(longitude, latitude);

This returned an error:Geometry SRID (0) does not match column SRID (102285)
SQL state: 22023

So I tried to set the geometry SRID using this syntax:

SELECT UpdateGeometrySRID('dc_crimes','geom', 102285);

This also gave me an error: ERROR: column not found in geometry_columns table
CONTEXT: PL/pgSQL function updategeometrysrid(character varying,character varying,character varying,character varying,integer) line 36 at RAISE
SQL statement "SELECT UpdateGeometrySRID('','',$1,$2,$3)"
PL/pgSQL function updategeometrysrid(character varying,character varying,integer) line 5 at SQL statement

********** Error **********

ERROR: column not found in geometry_columns table
SQL state: P0001
Context: PL/pgSQL function updategeometrysrid(character varying,character varying,character varying,character varying,integer) line 36 at RAISE
SQL statement "SELECT UpdateGeometrySRID('','',$1,$2,$3)"
PL/pgSQL function updategeometrysrid(character varying,character varying,integer) line 5 at SQL statement

I believe I need to put values in the ST_MakePoint() but I am not sure which points.


texts: Why I cannot use mean() after merging two imageColletions?

texts:
I have two ImagesCollections. I can do mean() reducer before merging them. But after merging them I got error when I was trying to do mean() to the merged imagecollection. See below the code.

//The image collection of Chlorophyll from MODIS Aqua
var chla_Modis = ee.ImageCollection("NASA/OCEANDATA/MODIS-Aqua/L3SMI")
.filterDate('2010-7-04', '2011-3-14')
.select('chlor_a');
print(chla_Modis);

var seawifsData = ee.ImageCollection("NASA/OCEANDATA/SeaWiFS/L3SMI");
var chla_SeaWiFS = seawifsData
.filterDate('2010-1-04', '2010-12-11')
.select('chlor_a');
print(chla_SeaWiFS);

//merge the two imageCollections
var mergedChla=chla_Modis.merge(chla_SeaWiFS);

print(mergedChla);

var meanChla=mergedChla.mean();
print(meanChla);


texts: What is a reasonable time to get a gps fix with 10 meters precision?

texts:
I hope I'm asking this at the right part of Stack Exchange. Please bear with me if I'm wrong.
I'm developing some gps based applications. The demands I have for precision are not very high, but I need to know the possible errors.
I have learned that the fourth decimal gives you ~10 meters precision. That should be enough for me.
The real question is how fast will I get that precision realiably in different environments (indoors, outdoors free sky, forest, cloudy, city etc).
The applications I'm developing is for handheld devices so I prefer to have the gps active in as short intervals as possible.
As I do now, the intervals I use the gps are more governed by battery life than precision. Now I'm trying to balance the two.


texts:
I have ENVI 5.0 , ERDAS 2013 SP2, and ArcGIS 10.1 SP1.  Every time I save an Imagine file from ENVI 5.0 and open in ERDAS 2013 and it seems like ERDAS does not like it and it would not open and would close the program.

Does anyone out here have the same experience ? or am I the odd man out here ?


texts: Connecting multiple polygons with lines using shortest distance as parameter

texts: Is there any plugin in QGIS (or code) that is used to connect multiple polygons with lines using distance as parameter. Let say that I put 10 meters, than all polygons within that distance will be connected with lines where is shortest distance from polygons. Something like on images below:


texts:
What actually does the createObject(Point.getClsid()) do?Does it create a new Point Object,and return the class id of the Point.Does it mean that each point created will have the same class id.

I went through ESRI documentation but all it said was Identify class ID for an object.


texts: How do I get QGIS to draw the rotated symbols on my map?

texts:
I'm working on QGIS 3.2. In my'normal' window I've managed to rotate symbols, so that three attributes now make up three parts of a pie chart. But in the Print Layout composer these rotated symbols are put on top of each other. This results in the top third of the pie being visible, but the other two barely.

How do I get the Print Layout to display the three symbols as a pie chart?


texts: No Spatial reference in output TIFF file

texts:
The spatial reference is missing when I try to save tif file (name_ext1.tif). I am using the following code:
name=name_ext1.tif
savepath=os.path.join(r'H:\Data\imagery\test',name)
kwargs=src.meta
kwargs.update(dtype=rasterio.float32,count=1)
with rasterio.open(savepath,'w',**kwargs) as dst:
    dst.write_band(1,file.astype(rasterio.float32))

src.meta
{'affine': Affine(30.0, 0.0, 607485.0,
   0.0, -30.0, 3782115.0),
'count': 1,
'crs': CRS({'init': 'epsg:32642'}),
'driver': 'GTiff',
'dtype': 'uint16',
'height': 7281,
'nodata': None,
'transform': (607485.0, 30.0, 0.0, 3782115.0, 0.0, -30.0),
 'width': 7481}**

What is the possible workaround?


texts:
I'm trying to vectorize a raster image and when I'm running polygonize within QGIS 2.6.1. I'm getting the following message:

  Traceback (most recent call last):  File
  "C:\PROGRA~1\QGISBR~1\bin\gdal_polygonize.py", line 36, in  import
  gdal, ogr, osr  ImportError: No module named gdal

I know that gdal is installed.
I don't know what to do! Need help!


texts: Cut linestring segment at centroid in PostGIS

texts:
In a table of linestring three segments are connected.

first--second--third

could Postgis generate the centroid of second segment and merge first and third at center?
Is it possible?

first--centroid--third


texts: Inserting title block when making print templates

texts:
I am creating a print template for some figures for a report.

These figures will need a title block with all my company's details. Typically when I use it in surfer it is just a bunch of lines surrounding the image.

How would I insert a title block into the print composer? Or will I just have to make one from scratch.


texts: Performing action upon refresh command in TOC using ArcObjects?

texts:
We are using custom layers in our plugin used in ArcScene/ArcMap.

For quite some time we were trying to figure out how to hook into the event occuring when the user selects "Refresh" in the right click menu item of a layer in the TOC.

Any advice how to do that?


texts: Normalizing negative pixel values in a GeoTIFF with GDAL in Python

texts:
I have a GeoTIFF pulled from Google Earth Engine. I'm using the following code to see the red values of each pixel:

from osgeo import gdal

gdal.UseExceptions()
ds = gdal.Open(myfile)

band = ds.GetRasterBand(4)
arr = band.ReadAsArray().astype(int)

However, the array contains numerous pixel values below zero. Is there a way to convert these values such that they cannot be less than 0? Or are these pixels just simply error pixels?

Gdalinfo dump:

Driver: GTiff/GeoTIFF
Files: LC08_014030_20170220.tif
Size is 749, 1570
Coordinate System is:
PROJCS["WGS 84 / UTM zone 18N",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",-75],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AUTHORITY["EPSG","32618"]]
Origin = (604200.000000000000000,4855080.000000000000000)
Pixel Size = (30.000000000000000,-30.000000000000000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  604200.000, 4855080.000) ( 73d42'13.61"W, 43d50'29.23"N)
Lower Left  (  604200.000, 4807980.000) ( 73d42'46.30"W, 43d25' 2.88"N)
Upper Right (  626670.000, 4855080.000) ( 73d25'27.68"W, 43d50'16.59"N)
Lower Right (  626670.000, 4807980.000) ( 73d26' 7.41"W, 43d24'50.42"N)
Center      (  615435.000, 4831530.000) ( 73d34' 8.82"W, 43d37'40.10"N)
Band 1 Block=256x256 Type=Float32, ColorInterp=Gray
  Description = B1
Band 2 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B2
Band 3 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B3
Band 4 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B4
Band 5 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B5
Band 6 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B6
Band 7 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B7
Band 8 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B8
Band 9 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B9
Band 10 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B10
Band 11 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = B11
Band 12 Block=256x256 Type=Float32, ColorInterp=Undefined
  Description = BQA


texts: WMS Request not served to OpenLayers 2 by GeoWebCache off GeoServer?

texts:
I'm developing a GIS app that makes use of GeoServer and GeoWebCache for serving requests, and OpenLayers for frontend. WMS requests to the GeoServer are served with no issues. But when I make a request through GeoWebCache, nothing is returned to the client. I've enabled Direct WMS Integration in GeoWebCache configuration. I'm using GeoWebCache integrated in GeoServer. My request is as shown below:

var wmslayer = new OpenLayers.Layer.WMS(
   "Kenya Administrative Boundaries",
   "http://localhost:8081/geoserver/gwc/service/wms?",
   {layers: 'kenya_admin', format:'image/png', transparent: true},
   {opacity: 1.0, isBaseLayer: false, visibility: true}
);
map.addLayer(wmsLayer);

But the layer is not loaded.

What am I doing wrong?


texts: GDAL 1.10 ogr2ogr OSM "Unable to open datasource"

texts:
I have a problem with ogr2ogr and osm file.

I install from source (on Ubuntu 12.04.2 LTS) newest GDAL libary and OGR still can't read osm files.

gdal-config --version
1.10.0
ogr2ogr --version
GDAL 1.10.0, released 2013/04/24
ogrinfo
Supported Formats:
  -> "ESRI Shapefile" (read/write)
  -> "MapInfo File" (read/write)
  -> "UK .NTF" (readonly)
  -> "SDTS" (readonly)
  -> "TIGER" (read/write)
  -> "S57" (read/write)
  -> "DGN" (read/write)
  -> "VRT" (readonly)
  -> "REC" (readonly)
  -> "Memory" (read/write)
  -> "BNA" (read/write)
  -> "CSV" (read/write)
  -> "GML" (read/write)
  -> "GPX" (read/write)
  -> "KML" (read/write)
  -> "GeoJSON" (read/write)
  -> "GMT" (read/write)
  -> "PostgreSQL" (read/write)
  -> "PCIDSK" (read/write)
  -> "XPlane" (readonly)
  -> "AVCBin" (readonly)
  -> "AVCE00" (readonly)
  -> "DXF" (read/write)
  -> "Geoconcept" (read/write)
  -> "GeoRSS" (read/write)
  -> "GPSTrackMaker" (read/write)
  -> "PGDump" (read/write)
  -> "GPSBabel" (read/write)
  -> "SUA" (readonly)
  -> "OpenAir" (readonly)
  -> "PDS" (readonly)
  -> "HTF" (readonly)
  -> "AeronavFAA" (readonly)
  -> "EDIGEO" (readonly)
  -> "SVG" (readonly)
  -> "Idrisi" (readonly)
  -> "ARCGEN" (readonly)
  -> "SEGUKOOA" (readonly)
  -> "SEGY" (readonly)
  -> "ODS" (read/write)
  -> "XLSX" (read/write)
  -> "PDF" (read/write)

When I try this:

ogrinfo UMP-PL-CCBySA_20130717.osm

and I get this:

ERROR 0:
FAILURE:
Unable to open datasource `UMP-PL-CCBySA_20130717.osm' with the following drivers.
Supported Formats:
  -> "ESRI Shapefile" (read/write)
...
...

Anyone know what I do wrong?


texts:
I have a thousand raster tiles in un-projected decimal degrees from which to build a seamless mosaic in a metres-based coordinate system. The final image will be about 70,000px by 50,000px. Which is the better approach, thinking in terms of processing (time) efficiency and quality of result:

mosaic into a single image and then transform to the final
projection?
or transform each tile and then mosaic?
(or don't mosaic at all and use a virtual raster table in the
final projection?)

I can use either or any of ArcGIS Desktop, GDAL Utilities, QGIS and a have 64bit OS with 12GB RAM.


texts: Adjusting the map ratio in such a way that Time Manager videos are consistent with YouTube resolution formats

texts:
I've recently made a time lapse video depicting migrants deaths in the Mediterranean Sea and surrounding areas in 2020. It was made in QGIS by means of the Time Manager plugin. I've uploaded it here on youtube.
As one can see, the points that represent the fatalities are not circles, but ellipses. This is because the original images that comprised the video had to be adjusted so as to fit the appropriate aspect ratio of YouTube videos, which in this case is 640 x 360 pixels (360p).
However, I'd like to prevent these elongated circles / ellipsoids from being shown in future videos. I want them to be circles.
How do I change the ratio of maps in QGIS in such a way that Time Manager videos are consistent with YouTube's image resolution formats?


texts: Call field from zip() into Field Calculator expression ArcGIS

texts:
I've been trying to find a way to call in field names from a zip(). Basically I have 12 fields called "Flw_days1", "Flw_days2" etc and corresponding "Depth" fields that need to be calculated using the inputs of the "Flw_days" fields. I.e. "Depth1" is calculated using "Flw_Days1" and so on.

I've used the zip() function to match these together and then tried to call these into a codeblock and expression for use in the field calculator.

See example below:

depth_fields = ["Depth1", "Depth2", "Depth3", "Depth4", "Depth5", "Depth6", "Depth7", "Depth8", "Depth9", "Depth10", "Depth11", "Depth12"]
flow_fields = ["!Flw_Day1!", "!Flw_Day2!", "!Flw_Day3!", "!Flw_Day4!", "!Flw_Day5!", "!Flw_Day6!", "!Flw_Day7!", "!Flw_Day8!", "!Flw_Day9!", "!Flw_Day10!", "!Flw_Day11!", "!Flw_Day12!"]

zipped = zip(sorted(depth_fields), sorted(flow_fields))

for input in zipped:
  codeblock = """
def depth(flow, kfc, rf):
  if kfc == 0 or rf == 0:
    return ""
  else:
    return (flow * (kfc / rf) + 2)
"""
  expression = "depth(input[0], !K_FC!, !Ret_Fact!)"
  arcpy.CalculateField_management(CMPP_fc, input[1], expression, "PYTHON_9.3", codeblock)

As you can see I have tried to call the "Flw_Days" field into the expression and the "Depth" field into the field calculator.

The error I get is:

Runtime error <class 'arcgisscripting.ExecuteError'>: ERROR 000539: Error running expression: depth(input[0], 0, 0) <type 'exceptions.TypeError'>: 'builtin_function_or_method' object is unsubscriptable Failed to execute (CalculateField).

So it seems to be calling the input[0] into the expression that is causing the issue.

I am using ArcGIS Desktop 10.


texts:
I have tried examples from the GIS support website and many more, but I cannot figure out how to get a codeblock for a field calculation to work within ArcPy.

Part of my code is underneath, it is a continuation of a table being added (multiple calculate fields). However the earlier calculations do not use a codeblock, the last one does. As a side note to the codeblock section, wing.IDE will not let me indent the code further than this. I hope that that is not where it goes wrong

Codeblock:

Codeblock_for_geometry_10 = """
    def val:
    if [AFSPRAAK] = 'Huur' then
    val = 'D102'
    elseif [AFPSRAAK] <> 'Huur' then
    val = 'D101'
    end if
    """

Field calculation (RVB_geometry_9 being the previous calculation on the table):

RVB_geometry_10 = arcpy.CalculateField_management(RVB_geometry_9, "NIVEAU_3", "val", "VB", Codeblock_for_geometry_10)

The error I get:
A field name was not found or there were unbalanced quotation marks.

As far as I can tell the quatation marks are correct and the fieldnames are aswel. I would like to learn what I am doing wrong here for future calculations.


texts: Extricating and sorting data from attribute table in QGIS

texts:
I download shapefile of Alberta from GEOFABRIK in order to work with the roads shapefile.
The problem was that I get an attribute table with 94,644 rows that was slowing down the computer, additionally to the fact that I get every small road and just the main one as I wanted.

I looked in the attribute table and found that the road type are classified in the fclass column. As described here I use the categorized function to classify the different road types, so I was able to display only the primary roads.

Also, it is possible to change the line style, by clicking the "Symbel" bar.

However, it was only a partial solution, I still have a lot of unnecessary data seem to slowing down me computer.
I was looking for a way to extract the data in need, and to work with this kind of big database.


texts: Generating series of map for each field of attribute table in QGIS?

texts:
Using QGIS 2.8, I have a shapefile with more than 30 columns (or fields) in the attribute table. I would like to automate the creation of one map for each of them.

I gathered that the QGIS Atlas generator can do that by "line" of the field, how to do that by column?

I found Generate series of maps with Python script, but I am on MacOS ...


texts:
I have created a few GeoExt webmapping projects using Ubuntu10.04 / Tomcat6 / Geoserver /Postgis but have no good solution regarding printing the maps. Using GeowebCache I am not able to export to pdf using the mapfish module in geoserver due to tiling errors. If I do a simple html print of the mappanel then the toolbar in the mappanel is also printed and distorted....or I risk being sued by google if the maps are printed.

I need to be able to serve the maps with geoserver - ie not as tms or filebased tiling methods so that the tiles can be loaded as an ogc service.

Can anyone help?

yours,

Rob


texts: How to get the intersecting polygon of two polygons? (using JSTS intersection)

texts:
I want to intersect two polygons within an OpenLayers layer. I found the JSTS library, which should be able to intersect two polygons, and return the new created polygon of the intersection.
The function to use: intersection
polygon1.intersection(polygon2)

Using the OpenLayers Polygons directly, I get errors, resulting of the missing functions, added by JSTS. And if I use JSTS Polygons, there seems still something missing.
this.getFactory().createGeometryCollection is not a function
How do I get an intersection of those two polygons?


texts: I am having an issue with the PLSS boundaries and AIANNH boundaries aligning with each other. Both datasets are plotting in the correct region, but the difference between the two datasets varies roughly between 600 and 1200 feet. See the attached screenshot. The red hatched layer is the AIANNH dataset. Is there a citeable explanation for this? Both datasets are in GCS_North_American_1983. I'm working with the most recent PLSS dataset for New Mexico (2016) along with the TIGER 2017 AIANNH dataset and am focused in the western/north-western region of New Mexico.


texts: Converting X and Y coordinates Data to Latitude Longitude

texts:
I have a data file that contains locations of different places of Wilmington, North Carolina. The data file contains the X Coordinates, Y Coordinates, and Latitude. I want to calculate the Longitude value of those places? How can I do this? I have used the arctan2(y,x) formula but it is not giving the correct location.
Sample Data:
X Coord: 2579880
Y Coord: 12401479
Latitude: 34.1213976435541
Longitude: ?


texts: How to label one point bigger than all others in a layer?

texts:
I have a point dataset in one layer. I labelled all points but I would like to label one point larger size than other point labels. How can I do it?

Merzifon and Bahçecik are in the same point layer. I want "Merzifon" text larger size than "Bahçecik".


texts: Problem with OverviewMap control projection in Openlayers 3

texts:
I am using the OverviewMap control on my ol3-based website. Now, since I am showing a map of Greenland, I use a UTM 24N projection (EPSG:32624). When I see the overview map, however, it looks like it is using Openlayer's default projection instead of UTM 24N. Here is an example:

I am not using any "special" code to generate the OverviewMap control:

map.addControl(new ol.control.OverviewMap({
                collapsed: true
            }));

Any idea what I am doing wrong? Or is this maybe a bug? Could it have something to do with Geoserver (where the tiles come from)?


texts:
Does anybody know some Exonyms dataset?

What I would like to find is a list of Europe cities/toponyms in different languages where all the cities refer to the same.

Example:
London, Londra, Londres, Londino --> All refer to London

I don't know if this exists.
The idea is to find some database that I can download, in order to make as authomatic as possible, and not copy-paste it.

I've searched in the net, but have not found anything.


texts: Edit polygon borders in QGIS preserving topological correctness (planar partition)

texts:
Given a layer of polygons that is topologically perfect, i.e., has no overlaps or gaps between the polygons, is it possible to manually edit the edges of the polygons in QGIS in a way that preserves the topological correctness?

As illustrated In the image below, suppose I want to change:

the "border" between polygons A and B to the red dotted line
the "border" between polygons A and C to the yellow dotted  line

This is the POSTGIS code to generate the polygons:

CREATE TABLE a_toy_polygons AS
SELECT      1 AS gid, 'A' as name, ST_MakePolygon(ST_GeomFromText('LINESTRING(0 0,1 0  ,1 1,0.5 1,0 1, 0 0)')) AS geom
UNION   SELECT  2 AS gid, 'B' as name, ST_MakePolygon(ST_GeomFromText('LINESTRING(0 1,0.5 1,1 1,1 2  ,0 2, 0 1)')) AS geom
UNION   SELECT  3 AS gid, 'C' as name, ST_MakePolygon(ST_GeomFromText('LINESTRING(1 0,2 0  ,2 1      ,1 1, 1 0)')) AS geom

Notice that in the common "border" between A and B there is already a vertex in the middle. There are actually two superimposed vertexes, one for each polygon. Is there a way to move them simultaneously? Or should I move one, then move the other close while having "snap" activated?

In the "border" of A and C, however, there are no such intermediate vertexes to be moved.

This is, of course, a toy example, my end use would be to fix manually some features of a Postgis topology. I am aware of the "Topology Checker" QGIS plug-in, but that is not capable of moving nodes yet. So I figured it could be possible to cast the topology as a geometry (of polygons) and edit in QGIS.


texts: Embedding 2500 images into feature class using ArcGIS Desktop?

texts:
How can I automate embedding images into a feature class attribute table?

I am working with a point feature class within a geodatabase that has a total of 2500 rows of data that need two image attachments each. Each row has a unique ID that corresponds with the image (Labeled by ID) that needs to be tied to it. I am working within ArcMap 10.4.


texts: How to apply a Majority filter in GRASS GIS?

texts: I was wondering if someone could help me write a majority filter for GRASS GIS. I have this map (see below) which I need to run a filter over to make the red areas a bit more dominant.  I know I can use r.mapcalc but I'm just not sure of the formula.


texts: I would like to use osm2po as a routing engine/library. My problem is that in geojson response, there is no osm_id value returned, only id field which, if I assume correctly is not a fixed value. I mean if I made new graph file from newer osm data, id of the same feature/edge could be different. Is there a way to get osm_id from service response?


texts:
I often need to apply functions on every raster cell. Unfortunately I don't know how to use a efficient method like calc, so I am using very slow for-loops.

For example:

logistic_function <- function(x, x_null, L = 1, k = 0.5) {
  y = 1 / (L + exp(1)^(-k*(x - x_null)))
  return(y)
}

#random input raster
r <- raster(ncol = 370, nrow = 202)
values(r) <- runif(ncell(r))

#output raster
r_log <- raster(ncol = 370, nrow = 202)

#inefficient for-loop
for (this_col in 1:ncol(r)) {
  for (this_row in 1:nrow(r)) {

      r_log[this_row, this_col] <- logistic_function(x=r[this_row, this_col],
                                                 x_null = cellStats(r, median),
                                                 k=0.5)
  }
}

I would like to do something like this:

r_log <- calc(r, cellStats(r, median), fun=logistic_function)


texts:
I need to do something apparently simple: adding a new attribute to an existing simple feature type. Something like:

SimpleFeatureType sft = DataUtilities.createType("", "the_geom:Point:srid=3035,ex:String,nb:int");

sft.addAttribute("newAttribute", String.class);

There is unfortunatelly no such addAttribute method and I did not find any direct way, or even a work around. The class DataUtilities allows a number of operations on feature types, but not that one. Maybe it is not possible. Does anybody has an idea how such operation could be performed ?


texts: Merging lines with same name and intersect each other?

texts:
I have shp file with 50k lines (roads).
Many roads have same name like 'high road' . some of them are the same road (they intersect each other) and some are located in different places of the city.

I want to merge all lines with same name and that intersect to single road. How can I do that? I have available SQL 2012 spatial with geometries and geographic data of roads and QGIS with shp file,

Someone told me in other post to use dissolve, but it looks like dissolve works on polygons only (i don't have selection for line vector layer)


texts: pgrouting librouting.dll library is supposedly missing when it is not

texts:
Yesterday, I successfully installed pgrouting from the most recent x64 2.0.0 binaries, at least this is what I thought. I was able to perform first query - select pgr_createTopology() from the tutorial i was following (http://anitagraser.com/2013/07/06/pgrouting-2-0-for-windows-quick-guide/) . On the next step (step 4 from tutorial), however, I received an error:

ERROR: cannot load library "C:/Program Files/PostgreSQL/9.4/lib/librouting.dll": The specified module could not be found.
SQL state: XX000

I copied downloaded files to PSQL folder again, but it did not help. I found out they have been "blocked" by windows security module, so spent some time learning a new trick on how to unblock them. Finally, after making sure all the files are at the correct places, I did a DROP on pgrouting extension hoping to reinstall it with with replaced files. However, this time even create EXTENSION pgrouting; fails, resulting with the same error as above.

I assume my steps were not necessairly correct since select pgr_createTopology() actually worked, so librouting.dll was probably fine. But nothing explains why subsequent SQL query (with pgr_dijkstra())  prompted an error to librouting.dll file?

Do you maybe have any thoughts on that?
I have installed PSQL, PostGIS, pgrouting with the latest files for x64 Win7 as discrepancies between versions were mentioned as a probable cause of such error in Installing PgRouting.


texts: Creating a roughness ASCII file from a raster dataset

texts:
I'm trying to create a roughness file (derived from topography) in the ASCII format for a European raster dataset. This dataset is from the EEA.
I've already made an ASCII grid file of the raster data as a whole (euroascii.asc) using the command:
gdal_translate -of AAIGrid ~/elevation3x3/elevation3x3.tif euroascii.asc
The roughness file has to be continuous (AKA 0.1, 0.2, 0.3, and so on), and has to account for every value that is already in the euroascii.asc file.
There is a command that is meant to do just this (I think), but I've tried it in the QGIS Python Console and it doesn't work.
This is what I've done in the QGIS Python Console, taken direct from the manual:
import processing # this gets no errors or anything out
And also:
from qgis import processing
for alg in QgsApplication.processingRegistry().algorithms():
        print(alg.id(), "->", alg.displayName()

To test whether this has any effect, I tried the code from https://docs.qgis.org/2.8/en/docs/user_manual/processing_algs/gdalogr/gdal_analysis/roughness.html
processing.runalg('gdalogr:roughness', EEAdata, 1, compute_edges=F, roughnessFile)
Which gets the error:
AttributeError: module 'qgis.processing' has no attribute 'runalg'
I'm only using QGIS Hannover, with no access to ArcGIS or any other commercial software. I'm also running this on Ubuntu 18.04.


texts: Field calculater in Python console from QGIS

texts:
I have some doubts in syntax of processing.runalg(qgis:fieldcalculator) in python console of QGIS. I want to create a new field in a shapefile and calculate this formula: [D_PATRICIA]*1000. (D_PATRICIA= is another field form the table of shapefile that contain values).

I try to do this and retrieved me a shapefile with the new field (p_mm) but without values:

processing.runalg('qgis:fieldcalculator', str4, "p_mm", 1, 6, 6, True, '[D_PATRICIA]*1000', str5)

What should be the problem? I think the problem is in the formula parameter.
Any idea to help me?

My code is:

import sys
import processing
import os
indir="D:/PATRICIA/ECMWF/2014"

for root, dirs, files in os.walk(indir):
    for file in files:
        if file.endswith('.nc'):
            fullname = os.path.join(root, file).replace('\\', '/')
            filename = os.path.splitext(os.path.basename(fullname))[0]

            str1="D:/PATRICIA/ECMWF/2014/"+filename+".tif"
            str2="D:/PATRICIA/IG/ESTACOES METEOROLOGICAS/ESTACOES_SPITFIRE74.shp"
            str3="D:/PATRICIA/ECMWF/2014/"+filename+".vrt"
            str4="D:/PATRICIA/ECMWF/2014/"+filename+"_VERAO.shp"
            str5="D:/PATRICIA/ECMWF/2014/"+filename+"_VERAO1.shp"

            processing.runalg('gdalogr:buildvirtualraster', fullname, 0, True, False, str1)
            processing.runalg('saga:addgridvaluestopoints', str2, str3, 2, str4)
            processing.runalg('qgis:fieldcalculator', str4, "p_mm", 1, 6, 6, True, "D_PATRICIA" * 1000, str5)


texts:
I am writing a Python script, which can literate through all main folder and Subfolder. My main folder names are 2001, 2002.....and under main folder differents subfolder is there, Subfolder name : JANUARY, FEBRAURY so on,
So my existing script I can iterate all subfolder under one main folder, but I can't iterate through all the main and subfolder at the same time,
Below is my script:

Folder and subfolder structure like this:

import arcpy, os, calendar
from arcpy import env
from arcpy.sa import *
arcpy.CheckOutExtension("Spatial")
topWorkspace = r'D:\Workspace\2006'
env.workspace = topWorkspace
months = {calendar.month_name[i].upper(): str(i).zfill(3) for i in range(1, 13)}

 # Step through list of all folders
for folderPath in arcpy.ListWorkspaces():
    baseName = os.path.basename(folderPath).upper()
    if baseName in months: # Test that subfolder is a month name
        monthNumber = months[baseName] # Get month-number for use in output filename
    arcpy.env.workspace = folderPath

    #Loop every month folder ta take shape file from each month
    shapeFile = (arcpy.ListFeatureClasses('rain_*.shp')[0])

I have updated the code using arcpy.da.walk but its producing error

    import arcpy, os, calendar
    from arcpy import env
    from arcpy.sa import *
    arcpy.CheckOutExtension("Spatial")
    topWorkspace = r'D:\Workspace\2006'
    env.workspace = topWorkspace
    months = {calendar.month_name[i].upper(): str(i).zfill(3) for i in range(1, 13)}
    # Get list of all folders
   walk = arcpy.da.Walk(topWorkspace, topdown=True,datatype="FeatureClass")
     # Step through list of all folders
    for folderPath, folderNames, fileNames in walk:
        baseName = os.path.basename(folderPath).upper()
        if baseName in months: # Test that subfolder is a month name
            monthNumber = months[baseName] # Get month-number for use in output filename
        arcpy.env.workspace = folderPath

        #Loop every month folder ta take shape file from each month
        shapeFile = (arcpy.ListFeatureClasses('rain_*.shp')[0])

Error message:

File "<module2>", line 27, in <module>
IndexError: list index out of range

How to address here two  loop condition for main and sub folder


texts: I downloaded a shapefile of USGS DOQQ boundaries a few years ago and cannot seem to find one now.  The attached image shows an example of the shapefile I'm looking for (green polygons showing the quarter quad boundaries).  Where can I download a shapefile covering the contiguous USA that shows USGS DOQQ boundaries with relevant attributes?


texts: D8 algorithm to get flow direction in ArcGIS Desktop

texts:
I am using ArcGIS 10 for the get the flow direction from DEM. I wonder about the algorithm they have used while calculating flow direction. It is mentioned there that it is using D8 model. D8 model should have the value 1,4,8,16,32,64 and 128. But I am getting the value from 1 to 255. Please correct me if my understanding is not correct.
I want the actual values of D8 model ( 1, 4,8,16,32,64 and 128) in flow direction rather than 1 to 255.
Is it possible to get that value?


texts: Iterate Graphical Model over a set of numbers in a non geometric layer

texts:
I have built a Graphical Model which, among other actions, extracts attributes of a vector layer by using the 'Extract by Attribute' algorithm. The extracted features are those whose Time in the Date&Time field lies within one specific hour of a day (say between 08:00 and 09:00) and then certain operations are performed on them. I was trying to iterate the Model to perform this over every hour of the day in a loop. So it would start from 00:00 to 00:59, then 01:00 to 01:59...and so on till 23:00 to 23:59.
I tried to use Vector Features as an input since this can be iterated. The vector feature referred to a CSV file which had only one field of 'Hours'. This field has 23 rows which were integer values (epoch values of hours).
I wanted to give each feature of this CSV file as an iterative input to the Extract as Attribute function and was to use it in the 'Formula' to incrementally increase the hours to be picked up.
However, the values are not being picked up. I tried to experiment with a simple model which used 'Field Calculator' to insert a new field with the current value of the Feature. however, here again, the new field does not contain any value. Just not able to get it to work.
Model -

If I keep the formula as Value, I get 'Null' in the New field, if I keep it as 'pre-calculated value', I get an error 'No root node! Parsing failed?'. In another case, I got an error'Cannot convert to feature' I think when I was using it directly in a formula to modify time.
Where am I going wrong?


texts: Turning code into a script and parameterizing?

texts:
I have python code that I want to turn into a script.

I want the end user to be able to navigate to the feature class the code should be run on.

The code I have currently speaks to a specific feature class and doesn't allow for picking what feature class the code is to be run on.

How do I change my code so I can have the parameters in the script allow the end user to pick what feature class the code should be run on?

This is the piece of code I need to change to allow the end user pick the feature class..

What do I need to change in my code to allow for the peramiters I set up in the script to replace this line of code?  I need fc to equal the input the end user selects not a specific feature class as the code currently stands


texts:
Possible Duplicate:
  How to fix QGIS error “Entry Point could not be located”?

Why can I not start Quantum GIS?  What do I do about this following error message?

"The procedure entry point TIFFMergeFieldInfo could not be located in the dynamic link librarylibtiff.dll".


texts: Apply a different color ramp to subset of categorized symbols in QGIS?

texts:
The goal is applying a color gradient to a subset of symbols, instead of applying it to all the existing symbols. I used to do this in ArcMap. In QGIS I have to do it manually, and it's very time consuming. Also, the results are not as beautiful, because I am not as good in choosing right colors.

This is the final result of doing individual ramps for each symbol group or subset:

And this is the styled map:


texts: QGIS Essen 2.14.0 OpenLayers plugins running Windows 7 64 bit

texts: I cannot make any OpenLayers plugin visible for the part of the world for which I plot vector points. I set CRS to WGS 84 / Pseudo Mercator. After adding a layer of points it should be trivial to click on one of the opened layers and have a map, but the map does not show.


texts: I want to a buffer around road (dwg file) in Arcgis. I am merging road lines but there are separate buffers for every line of road.


texts: How to know which Tiles build a rendered Image obtained by a Map Query?

texts:
as Test Purpose for my app, I want to know exactly, which Tiles were used to render the result image I get by query.

I.e. I want to see a bounding box of 500px*300px with zoom 13 of Location 48.42,12.93 and I query that in api.mapbox.com like this http://api.tiles.mapbox.com/v3/examples.map-zr0njcqy/12.93,48.42,13/500x300.png, then I would get the next image

Then I would like to know how I can know the specific tiles*, which were used to render this image 500x300.


texts: Possible bug in attribute ordering for QGIS 2.16.0-Nødebo?

texts:
I recently downloaded QGIS 2.16.0-Nødebo and noticed the ascending/descending ordering of the attribute table has changed (tested with integer, integer64, real and string type fields). Couldn't find a bug report for this so thought I would ask here if anyone has experienced this before I add a report?

Ascending order:

Descending order:


texts:
I am trying to use Leaflet.Geodesic.js geodesic to draw a geodesic line between two points, but I don't see any line at all.

Based on the Leaflet.Geodesic.js example script, this is what I have tried...

var Geodesic = L.geodesic([[52.5, 13.35, 33.82, -118.38]], {
opacity: 1,
color: 'red',
steps: 50
}).addTo(map);

Can anyone help me see what I'm doing wrong?


texts:
Hello I'm trying to create a label in QGIS that is formatted similar to the "table" styled yellow labels in the photo below. The attached image was from a map created in AutoCAD. I would like the labels linked to my attribute table which includes all that info.

I tried to use an expression to create a custom label but it is quite limiting when it comes to formatting. Does anyone have any ideas on how I could do this?

Name || ': ' || Depth_1 || '\n' || 'CL  ' || cl ||'(mg/L)' || '\n' || Name || ': ' || Depth_2 || '\n' || 'S04  ' || s04 ||'(mg/L)'


texts: Installation error - PostGIS on Ubuntu 18

texts:
I have installed PostgreSQL 11 on Ubuntu 18 and while creating PostGIS extension it gives this error.

  could not load library "/usr/lib/postgresql/11/lib/postgis-2.5.so":
  /usr/lib/postgresql/11/lib/postgis-2.5.so: undefined symbol:
  AllocSetContextCreate

The command I used to download PostGIS

apt install postgis postgresql-11-postgis-2.5


texts:
I'm not familiar with coding.

In QGIS, I would like field attribute to auto populate X and Y value and combine in one field just like highlighted below,

I know that if auto calculate value in field attribute by using $x_at(0) for start X coordinate and $y_at(0) for start Y coordinate,

How to calculate if one field read combination 'Y'N OR 'X'E just like highlighted below?

Has try out this expression $y_at(0) ||  $x_at(0) but turn like below, still not sure how to space between two value with N and E back at each value.


texts:
Currently, I am planning a research on the relation between mobile signal and poverty.

I found that telecommunications companies provide maps of their signal coverage like the links. http://www.sensorly.com/map/2G-3G/ID/Indonesia/Telkomsel/gsm_51010#|coverage

Is there any way that I can import the points in the map into qgis as vector or raster?


texts: GRASS7.8 setup to work externally in Python

texts:
I am following the instructions in https://grasswiki.osgeo.org/wiki/Working_with_GRASS_without_starting_it_explicitly#Python:_GRASS_GIS_7_with_existing_location to setup GRASS 7.8 to work with Python 3.6 for Windows. GRASS 7.8 was installed with a QGIS 3.10 installation.
The script throws an error on the last line of my code and I am not able to get any further - TypeError: a bytes-like object is required, not 'str' .
Here is my code:
import os
import sys
import subprocess

grass7bin_win = '\QGIS 3.10\bin\grass78.bat'
gisdb = os.path.join(os.path.expanduser("~"), "grassdata")
location = "/newLocation"
mapset = "mapset_first"

if sys.platform.startswith('win'):
    grass7bin = grass7bin_win
else:
    raise OSError('Platform not configured.')

startcmd = [grass7bin, '--config', 'path']

p = subprocess.Popen(startcmd, shell=False,
                     stdout=subprocess.PIPE, stderr=subprocess.PIPE)
out, err = p.communicate()
if p.returncode != 0:
    print(sys.stderr, "ERROR: Cannot find GRASS GIS 7 start script (%s)" % startcmd)
    sys.exit(-1)
gisbase = out.strip('\n\r')

I have tried using gisbase = out.strip(b'\n\r'), and this prevents the error from occurring on this line but it causes an error on the next line  - TypeError: str expected, not bytes - for os.environ['GISBASE'] = gisbase


texts: problem with Bounding Box GetFeature WFS

texts:
I have a WFS with two different feature types and I am trying to work with BBox.

This KVP is working: http://inspire.maaamet.ee/arcgis/rest/services/public/au/MapServer/exts/InspireFeatureDownload/service?SERVICE=WFS&VERSION=2.0.0&REQUEST=GetFeature&TYPENAMES=au:AdministrativeUnit&BBOX=59.664,26.3415,58.6268,26.3677,urn:ogc:def:crs:EPSG::4258&COUNT=10 for the feature type called au:AdministrativeUnit.

But this one does not work:http://inspire.maaamet.ee/arcgis/rest/services/public/au/MapServer/exts/InspireFeatureDownload/service?SERVICE=WFS&VERSION=2.0.0&REQUEST=GetFeature&TYPENAMES=mu:MaritimeBoundary&BBOX=59.664,26.3415,58.6268,26.3677,urn:ogc:def:crs:EPSG::4258 for the feature type called mu:MaritimeBoundary

Why this is happening when I change only the feature type name? Could it be a problem from the server itself?


texts: Clipping raster by mask layer with rule that it only clips raster if shapefile overlaps by 10% using QGIS?

texts:
I am using QGIS.

Is it possible to clip a raster when the shapefile overlaps with it by 10%?

Clip raster by mask layer tool is fine for Ontario (left) but makes me lose a lot of data from smaller locations (right). Need a rule to ensure consistency but don't know how to go about it.

Does anybody know of a function / have a code which could allow me to only cut the raster if 10% of the shapefile overlaps it?


texts: Using findAndReplaceWorkspacePaths in one directory (divide into sub folders) using Python?

texts:
I try to use 'findAndReplaceWorkspacePaths' function on MXD files using ArcView 10.3.1. The MXD files located in one directory (which is divided into folders and sub folders). Also in the sub folders there are MXD files.

This is the directory content location:

and this is the sub folders location (there are more MXD files in the other folders):

I write in python 2.7:

import arcpy,os,sys,fnmatch
from arcpy import env

rootPath = r"G:\desktop\Project"
for root, dirs, files in os.walk(rootPath):
    env.workspace = r"G:\desktop\Project"
    for mxdname in arcpy.ListFiles("*.mxd"):
        print mxdname
        mxd = arcpy.mapping.MapDocument(r"G:\desktop\Project\\" + mxdname)
         mxd.findAndReplaceWorkspacePaths(r"D:\PROJECTS\green_line\gis\layers", r"G:\PROJECTS\green_line\gis\layers")
        print 'findAndReplaceWorkspacePaths'
        mxd.save()
del mxd

The result is:

>>>
landuse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePaths
landuse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePaths
landuse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePaths
landuse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePaths
landuse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePathslanduse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePaths
landuse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePaths
landuse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePaths
landuse__3___20013.mxd
findAndReplaceWorkspacePaths
landuse__4___20013.mxd
findAndReplaceWorkspacePaths
>>>

Now, there are two problems:

Why the code work only with the MXD files in the directory and not in the sub folders?
Why the code printed so many lines?


texts: I want to make an uncoupled Kalman filter. The input data are GPS and the output of my strapdown Algorithm from the IMU. The Strapdown Algorithm is used to calculate the actual position and velocity from the given daten from an Inertial Navigation System (Inertial measurement unit [received data are from the accelerometer and gyroscope]). The GPS data are given as phi, lambda and h with a time. Do I have to transform them into x,y and z? My Strapdown output is also in phi, lambda and h. How do I know the dimensions of my matrizes? The most important thing to know is: My Strapdown output do not provide a time. Therefore when do I know if I have to use the Strapdown output? Only if there is no data provided by GPS?


texts:
I am trying to build the newest qgis project on windows 7 in my 32 bit machine.

I followed the official instructions and during the cmake configure I get the following error.

  CMake Error at C:/Program
  Files/CMake/share/cmake-3.0/Modules/FindPackageHandleStandardArgs.cmake:136
  (message):   Could NOT find PythonInterp: Found unsuitable version
  "1.4", but required

The executable python is in cygwin/bin path. the problem is that it finds a lower version of pythoninterp.

Does anyone have ideas?


texts: I'm new to QGIS, and I'm trying to plot settlements on a map with their elevation data (using Openstreetmap). Problem is, I find most of the villages/ towns have their 'ele' data missing in the attribute table. Is there a work around for this?


texts: Why are z-values replaced with zero after in_memory use in Python?

texts:
My script takes separate 3D feature classes, calculates a field, appends to a target feature class, and then deletes the originals. This is a huge data processing task of over 4700 mass point feature classes appending to a single one. I know it was a big mistake to put so many feature classes into SDE, but given they're in there now, I've got to deal with them there (yes, the sheer fc count is bringing SDE performance to its knees--hence my urgency).

Early on it worked great, in that the merged result showed valid z-values. I was able to build a terrain and it looked great. I'd been running this in the background for weeks, and have just noticed z-values of 0 in the output, where the originals were not 0! Oh no!!! I'm certainly glad I still have the file-based originals!

Pseudocode:

  An in_memory feature class is created from a z-aware template

  For each fc in the feature dataset:

    It is appended to the empty in_memory fc

    The TILE_ID field is calculated using a portion of the name of the source fc

    It is then appended to the target SDE fc

    The in_memory feature class is emptied

  Each processed fc is then deleted

When I initially posted, I was curious if it was my use of in_memory that caused this trouble. I had not initially been using it, and only switched to it because performance was far too slow otherwise, and WAY faster in_memory. To find where trouble started, I ran the Add XY Coordinates tool, which proved the trouble started with in_memory use. At that time I had also switched from CalculateField_management to da.UpdateCursor, so that's likely the bigger reason for the performance enhancement. I will now change my code to avoid in_memory.


texts: Converting latitude and longitude to point in geodjango?

texts:
I am doing the geodjango tutorial and I went on google and chose a random longitude and latitude in Kansas, I am using srid=4326 for my polygons

So I got this random lat and long on google maps in Kansas USA

lat, long = 38.670049, -99.565798

Now I am trying to see where in the world it is

>>> from django.contrib.gis.geos import Point
>>> pnt = Point(38.670049, -99.565798)
>>> WorldBorder.objects.get(mpoly__intersects=pnt)
...
...
world.models.DoesNotExist: WorldBorder matching query does not exist.

Using the tutorial point I get

>>> pnt = Point(-95.3385, 29.7245)
>>> WorldBorder.objects.get(mpoly__intersects=pnt)
<WorldBorder: United States>

So my database is imported and working.

So clearly Points are not the same as latitude and longitude, how do I convert between the two

It seems google maps has the lat/long the other way around, see coordinates in image below at the bottom from google maps:

# Try enter the coords in the same order as in image above
# Get wrong intersection
>>> from django.contrib.gis.geos import Point
>>> from world.models import WorldBorder
>>> pnt = Point(38.907636 , -77.036905)
>>> WorldBorder.objects.get(mpoly__intersects=pnt)
<WorldBorder: Antarctica>

# Swap them around, and it works
>>> pnt = Point(-77.036905, 38.907636)
>>> WorldBorder.objects.get(mpoly__intersects=pnt)
<WorldBorder: United States>

Can anyone explain what is going on exactly. For a Point(x,y)   which is the longitude and which is the latitude, is it

Point(x,y)

x=longitude

y=latitude

And in which order does Google maps present them in the picture above

I am not a GIS specialist or even novice.


texts: Is it possible to change default encoding for QGIS?

texts:
I often transfer layers from Windows machine to Linux and vice versa and it causes me problems (incorrectly displayed labels and same with text in attribute table) when I use drag and drop to put vector layers in QGIS.

On linux there is System encoding set to UTF-8 while on Windows it is Windows-1250. Is it possible to change the setting of the Windows version of QGIS to be also set to UTF-8 by default?

Tried to google a solution but no success.


texts: Automate creating Rectified files from georeferenced images

texts:
A python 2.6.5 / ArcGIS 10 solution is preferred.

UPDATE:
I found a tool to do this "Create Ortho Corrected Raster Dataset" but the python script for it gives errors when I try.

import arcpy
import glob
import os

path = os.getcwd()
filename = path+"\Stage2_19"
#for pathname, directories, filenames in os.walk(path):
#    for filename in filenames:
#        if filename.lower().endswith(".TIF"):
arcpy.GetMessages()
arcpy.CreateOrthoCorrectedRasterDataset_management(filename+".TIF",\
                                                   filename+".jpg",\
                                                   "CONSTANT_ELEVATION", "0")

arcpy.GetMessages()
print "Created: " + filename + ".jpg"

Gives the following errors.

        Traceback (most recent call last):   File
        "Q:\Software\GIS\ESRI\ArcGIS10\CustomScripts\temp\Georef-Rectified.py",
        line 16, in
            "CONSTANT_ELEVATION", "0")   File "C:\Program
        Files\ArcGIS\Desktop10.0\arcpy\arcpy\management.py",
        line 7815, in
        CreateOrthoCorrectedRasterDataset
            raise e ExecuteError: ERROR 000425: Failed to create
        ortho-corrected raster dataset Failed
        to execute
        (CreateOrthoCorrectedRasterDataset).

      >

ORIGINAL:

I have georeferenced a few hundred survey plans using the Georeferencing toolbar in ArcGIS but did not Rectify the images (i.e create a new image) as I didn't want to create new files and take up more space. Now when I want to open them up in ERDAS or MapInfo they don't accept the prj created in the georeferencing process and require each using Rectify-->Save As

Can I use python to automate the process - open the georeferenced Tiff files and save them as Rectified Tiff or JPG files with the same name.

I noticed that there is an option for using gdal_translate but this gives me a syntax error and http://www.gdal.org/gdal_translate.html doesn't help to much (even when trying the example gdal_translate -of GTiff -co "TILED=YES" utm.tif utm_tiled.tif). Also not sure how it would handle the fact that the original file is only georeferenced and not rectified.

Please help...thanks,


texts: I want to browse a city map offline, what are my options?

texts:
I found exports on http://download.bbbike.org/ and http://metro.teczno.com/

But it's not tiles, so it's not really usable unless I render it.

Are there software where I can browse openstreetmap.org tiles, but completely offline, for a small region ? Not editing, just browsing and maybe searching streets.

Any software which allows me to generate those tiles I want using those exports, and another software to browse tiles ?

I don't have a smartphone, I don't have a GPS, but I have a computer and an internet access, but sometimes I don't have a fast enough internet access, so I wonder why I can't find any way to host a limited set of tile map data...

I guess you could create a search index on smaller areas on a local computer if you give it enough time...

Does JOSM has a viewing mode ? garmin and qlandkarte seems to just be tools to transfer data to GPS devices.


texts:
I'm trying to have a zoom to extents of WFS layer. I'm using this code:
var map = new ol.Map({
  controls: ol.control.defaults().extend([
    new ol.control.FullScreen({className: 'ol-full-screen', tipLabel: 'Toggle full-screen'}),
    new ol.control.ScaleLine({className: 'ol-scale-line', target: document.getElementById('scale-line')}),
  ]),
  target: 'map',
});

var view = new ol.View({
  center: ol.proj.fromLonLat([9.25040, 45.54318]),
  zoom: 12,
});
map.setView(view);

var osm = new ol.layer.Tile({
  source: new ol.source.OSM()
});
map.addLayer(osm);

var treeSource = new ol.source.Vector();

var treeLayer = new ol.layer.Vector({
  source: treeSource,
  style: customStyle,
});
map.addLayer(treeLayer);

var featureRequest = new ol.format.WFS().writeGetFeature({
  srsName: 'EPSG:4326',
  featurePrefix: 'workspacename',
  featureTypes: ['trees'],
  outputFormat: 'application/json',
});

fetch('http://localhost:8080/geoserver/wfs', {
  method: 'POST',
  body: new XMLSerializer().serializeToString(featureRequest)
  }).then(function(response) {
    return response.json();
  }).then(function(json) {
    var features = new ol.format.GeoJSON().readFeatures(json);
    console.log(features);
    treeSource.addFeatures(features);
    map.getView().fit(treeSource.getExtent());
});

What I see is a zooming on the extent of the WFS but the features are out of scale because I see 50mm on scale bar. I exclude that the problem is on the features because I can use this features using WMS service withouth problems.
Another strange thing is that for a few second I can see the background map from OSM and after this I see only the WFS and the map is putted on 0,0.
I've used almost the same code of the official example.
If I use this code for the same layer, but using WMS service, I've no problems:
var wmsTreesSource = new ol.source.TileWMS({
  url: 'http://localhost:8080/geoserver/workspacename/wms',
  params: {
    'LAYERS': 'trees',
    'TILED': true,
  },
});
var wmsTrees = new ol.layer.Tile({
  title: 'Trees',
  source: wmsTreesSource,
  visible: true,
});
map.addLayer(wmsTrees);


texts: rescale MODIS NDVI value- fill value and background value after clip study area

texts:
I want to use MODIS NDVI time series to extract phenological parameters using TIMESAT software.
I have question about data prepration, Is it alwaysd necessery to rescale NDVI data?
valid range in MODIS NVDVI is -2000 to 10000.
I have cliped NDVI images with my study area map and now I have NDVI images with background value 0. I think I must set background value to -3000, because zero is a valid value in NDVI data, and the fill value / no data value in modid NDVI is -3000.


texts: How to upload NetCDF and HDF files to Geonode?

texts:
What is the best way to upload NetCDF and HDF files to Geonode?

Are those formats supported natively, do I have to access them programatically, or do I need to convert them to another format first?


texts: Retrieve value by geographic coordinate from GeoTIFF with GDAL/C++

texts: Given a GeoTIFF height map image, how can an altitude value be read at a given latitude/longitude using C++ and GDAL?


texts: Using distance and transform from geomfromtext from postgis in spatialite

texts:
I have the following query which I want to apply on my from postgis to spatialite imported database. The problem is that the result of the distance is always null, while the other selected coloumns are containing the expected values. So how can I use this query in spatial lite? It should calculate the distance from the given point to all in the databases stored points and sort them by the distance.

"select x1, y1, source, st_distance(st_transform(st_geomfromtext('POINT(' || x1 || ' ' || y1 || ')',4326),32633),st_transform(st_geomfromtext('POINT(15.442161 47.068170)',4326),32633)) from tablex order by st_distance asc limit 1;"

I tried it in spatialite with the following query:

"select x1, y1, source, distance(transform(geomfromtext('POINT(' || x1 || ' ' || y1 || ')',4326),32633),transform(geomfromtext('POINT(15.442161 47.068170)',4326),32633)) from tablex order by distance asc limit 1;"

If I do not use the transform function I get the right result, but the distance is in degrees and not in metres, as I would need it.


texts:
I have a polygon layer in a PostGIS database and connected it with GeoServer. When I publish the layer using a WFS service, the polygons are deformed, vertices are lost.
Where can I change GeoServer settings so that the layer is published complete?

In the image, the blue outline layer is the layer that is generated with the WFS service and the green layer is the original layer.
What GeoServer settings should I change for the WFS service layer to be published as the original layer?


texts: QGIS-Server Plugin Process per client (WMS)

texts:
I have a working infrastructure containing Apache2 with FCGI and QGIS-Server bundled in a Docker Container to deliver a WMS. My written Plugin (Python) for QGIS-Server is working like expected. It takes information from a WMS request, performs some calculations and updates my WMS Layer. It has one endpoint to setup a new project (creating a new folder with a .qgs file) and one endpoint to send a json-string to configure some settings.

I want now to go one step further. My plan is, that every connected client has their own specific configuration and data in the QGIS-Server Plugin. As far as I know it is not possible with the standard configuration of Apache2 with fcgi. Every new request will override my variables in the QGIS-Server Plugin (Is my assumption right?).

To be more clear:

I have one Docker Container to deliver a WMS
Inside of the Container is a Python Plugin for the QGIS-Server
The Plugin has project specific configurations and data
Normally the project is specific to a client
If a second client creates a new project to get a WMS which fits his need, he will override the settings of the first client (my assumption)
To solve this I need your help

Here are a few questions:

Is it possible to serve one process per connected client (so that one QGIS-Server Plugin exists for one connection)
If that is not possible, what are valid alternatives to solve the problem (My idea is to outsource client specific configurations and data in separate files on the filesystem which will be loaded into the QGIS-Server Plugin on every request)

Hints for the configuration of Apache2 and fcgi are welcome.


texts:
I noticed that with at least QGIS 3.14 (maybe earlier too) all objects are exported to shapefile as multi-geometry ie. points are exported as multipoints, lines as multiLines, etc, even if objects are internally stored as single-part geometry (after running multipart to single part). I simply noticed that when some of my GDAL-based scripts stopped working with new data processed with this version of QGIS.
I consider it a bug. But I'm afraid it is a feature with some narrow-targeted explanations.


texts:
I need to select a layer with all attributes while the selection comes from a table. With my code I get only empty feature classes and emtey Excel files. The tab_list Elements are the selection value. If Name2 = Name1 the selection should starting.

 field = ["SHAPE_Length", "SHAPE_Area"]

arcpy.MakeFeatureLayer_management(workspace + "\\" + gdb + "\\" + "Verknuepfung", "selectlayer")

bet_tab_v_list = [row[0] for row in arcpy.da.SearchCursor(convert_table, ["Name1"])]
s_cur = arcpy.da.SearchCursor("selectlayer", ["Name2"])

tab_list = set(bet_tab_v_list)

counter = 1
for elm in tab_list:
    arcpy.SelectLayerByAttribute_management("selectlayer", "NEW_SELECTION", "Name2 =' {}  '")
    arcpy.CopyFeatures_management("selectlayer", speicherpfad + "\\" + gdb + "\\" + text_pruefung(elm))
    arcpy.TableToTable_conversion(workspace + "\\" + gdb + "\\" + text_pruefung(elm), workspace + "\\" + gdb,"var")
    arcpy.DeleteField_management(workspace + "\\" + gdb + "\\" + "var", field)
    arcpy.TableToExcel_conversion(workspace + "\\" + gdb + "\\" + "var", workspace + "\\" + elm + ".xls")
    arcpy.Delete_management(workspace + "\\" + gdb + "\\" + "var")
    arcpy.AddMessage("Round" + str(counter))
    counter+=1


texts: Programmatically defining schema not working for nulls

texts:
I am attempting to create a table with appropriate schema. I have created a data dictionary with the appropriate values stored in a list. See below

schema_dict = {"borough" : ["borough","TEXT","","",2,"borough","NON-NULLABLE"],
"block" : ["block","LONG",10,"","","block","NON-NULLABLE"],
"lot" : ["lot","SHORT",5,"","","lot","NON-NULLABLE"],
"cd" : ["cd","SHORT",5,"","","cd","NULLABLE"]}

Then I use this data dictionary to call the appropriate parameter value when looping through a list of my fields as such:

for item in schema_dict:
    arcpy.AddField_management(table_name, schema_dict[item][0],
                                          schema_dict[item][1],
                                          schema_dict[item][2],
                                          schema_dict[item][3],
                                          schema_dict[item][4],
                                          schema_dict[item][5],
                                          schema_dict[item][6])
    print("Converting to appropriate schema for field " + str(schema_dict[item[0]]))

For whatever reason I am receiving the following errors:

  ExecuteError: ERROR 000622: Failed to execute (Add Field). Parameters
  are not valid.

  ERROR 000800: The value is not a member of NULLABLE | NON_NULLABLE.

According to the examples provided in the add field management doc, a string of "NULLABLE" or "NON-NULLABLE" should be sufficient, no? I notice the type is Boolean, should I be giving the parameter a value of True or False rather than the strings I supplied?


texts: How to receive input name by script in the toolbox on QGIS 3.14?

texts:
In QGIS 3.4.7 I have a script that I use in the toolbox and one of the steps is to store the name of the input tiff in a variable. I accomplish this perfectly with this part of the code:
from qgis.core import QgsProcessing
from qgis.core import QgsProcessingAlgorithm
from qgis.core import QgsProcessingMultiStepFeedback
from qgis.core import QgsProcessingParameterRasterLayer

class tiffname(QgsProcessingAlgorithm):

    def initAlgorithm(self, config=None):
        self.addParameter(QgsProcessingParameterRasterLayer('tiff', 'TIFF', defaultValue=None))

    def processAlgorithm(self, parameters, context, model_feedback):
        feedback = QgsProcessingMultiStepFeedback(25, model_feedback)

        NAMETIFF = parameters['tiff'].name()

        feedback.pushInfo(NAMETIFF)

        return {}

    def name(self):
        return 'tiffname'

    def displayName(self):
        return 'tiffname'

    def tr(self, string):
        return QCoreApplication.translate('Processing2', string)

    def createInstance(self):
        return tiffname()

However in QGIS 3.14 I get this error:

NAMETIFF = parameters['tiff'].name()
AttributeError: 'str' object has no attribute 'name'

How to call name() in this version?


texts: Losing large portions of elevation contours from QGIS to PDF

texts:
Please forgive me as I'm new to GIS and all this might be a silly question, the wrong place for this question, or a hack way of going about what I'm getting at. So any suggestions on better practice is greatly welcome. This might not be a GIS question entirely, but I've been trouble shooting this for a bit and thought the expertise of this form would be worth a shot.

I'm trying to make vectors from QGIS editable in something like Illustrator. So, I'm currently taking elevation contours I created from a .dem and saving them into a .pdf from the print window, as to make them editable in some design software. The .pdf appears as expected, but when I bring it into Illustrator from there I lose large swaths of some elevations.

I'm going to give ArcGIS a shot to see if its something crooked in the .pdf export, but I would rather stay with a QGIS or GDAL solution, as open is better.

Thats roughly it, so thanks in advance for any help.


texts: Shortest route batch calculation in QGIS 3.2.2 with v.net.path GRASS GUI

texts:
I have three shapefiles: road network, postcodes and destinations. I`d like to calculate the shortest route between the postcodes and destinations (multiple points in each layer)

Could anyone describe the batch process for this? I have tried the v.net.path but it does not allow batch process (There is no way I can import a CSV file or add a layer).

I am working in QGIS 3.2.2

I do not think that this is a duplicate of QGIS Shortest Route Network Analysis using multiple features because that question is pretty outdated. The tutorial is for QGIS 2.18, Im using QGIS 3.2.2. Also, Im asking about batch process, therefore it`s not a duplicate.


texts:
I want to use topojson on the command line to simplify some files. npm will install topojson without any errors:

npm install -g topojson

But then running:

topojson

produces -bash: topojson: command not found. This article suggests that:

  topojson/topojson is no longer the command-line tools for manipulating
  TopoJSON — that’s all been put in the topojson/topojson-client
  repository.

So I am trying:

npm install -g topojson-client

But then trying topojson-client produces command not found, and so does topojson and I can't see any documentation in either repo explaining what to do. Can anyone suggest any ideas?


texts: WFS with authentication to PostGis (or to shape)

texts:
Is there a way to import data from a WFS layer to a Postgis table (or to a shape)?
there are blogs explaining how to do it in 3 steps (http://anitagraser.com/2012/09/26/wfs-to-postgis-in-3-steps/) or just in 1 (http://anitagraser.com/2012/09/29/wfs-to-postgis-in-one-step/) but in none of them is explained how to pass the USERNAME and the PASSWORD required by the WFS service.

If i could export the WFS layer to the shape, then i could convert the shape easily in a postgis table.
But, again, the problem is in the HTTP authentication.

Any suggestion?

I think that the problem is that "GDAL/OGR must be built with Curl support in order to the WFS driver to be compiled" as is written here: http://www.gdal.org/drv_wfs.html.

And i suppose i didn't have GDAL built with Curl support then i receive the FAILURE message (Unable to open datasource).

Now the problem is:
Is this hypotesis true?
and..
how to build GDAL with Curl on CentOS?


texts: Get measure for selected vector objects in Openlayers 4

texts:
I need to get area and edge size values from selected object on the map. I worked on this code below. But I can not find the how i get these values.

function(polygon) {
    var area = ol.Sphere.getArea(polygon);
    var output;
    if (area > 10000) {
      output = (Math.round(area / 1000000 * 100) / 100) +
          ' ' + 'km2';
    } else {
      output = (Math.round(area * 100) / 100) +
          ' ' + 'm2';
    }
    return output;
  };


texts: Choose the value of the new feature from an existing feature in QGIS

texts:
Suppose that you are adding another feature and suppose that you want to automatically fill all the fields with some values.
These values should be copied from an existing feature on the map by selecting it on the map.

Example: you want to copy all the values from the red feature (labeled with road) to the new one (bottom right corner). The new feature table should have type = road, subtype = one-way, width = 5.

I cannot use the Default Value of the Attributes Form because there is not fixed rule.

The button that I'm looking for is already available in the Attribute Form when you change to Widget Type to Relation Reference, but you need to have a Relation set up with 2 layers (so not this case).


texts:
I have to convert ellipsoidal height to Mean sea level on 5 million points which only possible on GPU segments.
Below is my code
import cudf
import pandas as pd
import numpy as np
from laspy.file import File
import gdal
import pyproj
from numba import cuda
from scipy.interpolate import RectBivariateSpline as Spline

infile=File("pointcloud.las",mode='r')
lasfile_df = cudf.DataFrame()
lasfile_df['X']=infile.x
lasfile_df['Y']=infile.y
lasfile_df['Z']=infile.z
# http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm96/binary/binarygeoid.html
# Download WW15MGH.DAC
egm=np.fromfile("/content/drive/MyDrive/WW15MGH.DAC",">i2").reshape(721,1440)/100

@cuda.jit
def ellipsoid_to_msl(lon,lat,height,out):
  for i,(lon,lat,height) in enumerate(zip(lon,lat,height)):
    longs=np.arange(0,360,0.25)
    lats=np.arange(-90,90.1,0.25)
    interp=Spline(lats,longs,egm)
    lat*=-1
    if lon<0:
        lon+=360
    z=interp.ev(lat,lon)
    out[i]=height-z
lasfile_df=lasfile_df.apply_rows(ellipsoid_to_msl,
                                 incols={'X':'lon','Y':'lat','Z':'height'},
                                 outcols={'out': np.float64},
                                 kwargs={})

In the lasfile_df, the apply_rows is not working properly.
Can anyone explain how to do it using a cuda dataframe because in a pandas dataframe it takes a lot of time?


texts: Duplicating a polygon and offsetting it by distance and bearing from the original polygon

texts:
I am trying to locate a historical distribution/denning point for a species with the data "35km NW of x", with x being a small town. I have created a rough polygon around the current extent of the town using a satellite image. To define my search area I therefore have to duplicate the polygon and then off-set it from the original according to the bearing and distance.

How can I do this?


texts: Using MakeFeatureLayer to select by location in arcpy

texts:
I'm trying to assign a state plane zone based on if a point intersects it.  I make a feature layer from a feature class.  It gets through the loop once and calculates everything correctly, but when it loops again it says that the dataset already exists.  Can I put the MakeFeatureLayer command at the top and have it loop through underneath that?

This is the line that it gets stuck on after it calculates the first group of points:
arcpy.MakeFeatureLayer_management(fc,"fcTemp")

    for fc in SMON:
        print fc + " Started Processing"
        arcpy.AddField_management(fc, 'LATDMS','text')
        arcpy.AddField_management(fc, 'LONGDMS','text')
        arcpy.AddField_management(fc, 'DESC_','text')
        arcpy.AddField_management(fc, 'StatePlane','text')
        arcpy.AddField_management(fc, 'State','text')
        arcpy.AddField_management(fc, 'Unique_ID','text')
        arcpy.AddField_management(fc, 'DateImport','text')
        print "Fields added"
        arcpy.CalculateField_management(fc, "LATDMS", "decdeg2dms ( !GlobalLati! )", "PYTHON", "def decdeg2dms(dd):\\n    negative = dd < 0\\n    dd = abs(dd)\\n    minutes,seconds =  divmod(dd*3600,60)\\n    degrees,minutes =   divmod(minutes,60)\\n    if negative:\\n        if degrees > 0:\\n            degrees = -degrees\\n        elif minutes > 0:\\n            minutes = -minutes\\n        else:\\n            seconds = -seconds\\n            seconds = round(seconds, 5) #my attempt to define it\\n   \\n    dms = '{0:.0f} {1:.0f} {2:.5f}'.format(degrees, minutes, seconds)\\n    \\n    return (dms)")
        arcpy.CalculateField_management(fc, "LONGDMS", "decdeg2dms ( !GlobalLong! )", "PYTHON", "def decdeg2dms(dd):\\n    negative = dd < 0\\n    dd = abs(dd)\\n    minutes,seconds =  divmod(dd*3600,60)\\n    degrees,minutes =   divmod(minutes,60)\\n    if negative:\\n        if degrees > 0:\\n            degrees = -degrees\\n        elif minutes > 0:\\n            minutes = -minutes\\n        else:\\n            seconds = -seconds\\n            seconds = round(seconds, 5) #my attempt to define it\\n   \\n    dms = '{0:.0f} {1:.0f} {2:.5f}'.format(degrees, minutes, seconds)\\n    \\n    return (dms)")

        #Assigning State Plane and Unique ID for SMON
        #fcTemp = "fcTemp"
        arcpy.MakeFeatureLayer_management(fc,"fcTemp")
        arcpy.SelectLayerByLocation_management("fcTemp", "INTERSECT", CON, "", "NEW_SELECTION")
        urows = arcpy.UpdateCursor("fcTemp")
        for urow in urows:
            urow.StatePlane = 'Colorado NorthTEST'
            urow.State = 'Colorado'
            urow.Unique_ID = str(urow.PointID) + str(Con)
            urows.updateRow(urow)

        arcpy.SelectLayerByLocation_management("fcTemp", "INTERSECT", NDN, "", "NEW_SELECTION")
        urows = arcpy.UpdateCursor("fcTemp")
        for urow in urows:
            urow.StatePlane = 'North Dakota North'
            urow.State = 'North Dakota'
            urow.Unique_ID = str(urow.PointID) + str(Con)
            urows.updateRow(urow)


texts: VB.net search for address in Geodatabase in ArcObjects

texts:
I am try to search for address in Geodatabase in ArcObjects using VB.NET

I tray this sample but my database is large so when I use this sample in TextBox1_TextChanged it is slow !
How can I search and refresh the query list when I type the address fast
https://www.youtube.com/watch?v=i7FLASb-ByY


texts: Split feature by polygon and do multipart to singlepart in QGIS

texts:
I am searching for a tool or a combination of tools in QGis to get done some line-cuttings. They have to follow some rules:

Lines should be cut at the borders of a raster. For example, the purple line in picture 1 should be cut into 4 pieces as can be seen in picture 2 (rasters 5-8).

Lines that leave a raster and enter the same raster afterwards should be cut into 2 separate lines (e.g. the red line in raster 1).

Lines that leave a raster and enter the same raster afterwards AND cross themselves inside this raster should not be cut by the line itself (e.g. the yellow line in raster 3) but only by the raster-borders.

picture 1 - start

picture 2 - what I try to get

So here's what I got so far:
I managed to cut lines (picture 1) by raster-borders to get line-pieces for each raster. Therefor I used the intersection-tool. The result is shown in picture 3. Using then "multipart features to single parts" the red line in raster 1 is cut into 3 pieces and the yellow line in raster 3 is cut into 6 pieces as shown in picture 4.
If you have a look at picture 2 again, the red line in raster 1 should only be cut into 2 pieces and the yellow line in raster 3 should be cut into 3 pieces.
picture 3

picture 4 - what I got


texts: Can qgis2web display rule-based layer style?

texts:
I’ve just made a map that has two layers: one that contains the polygons and  other that contains points. For the last one, I used rule-based symbology where I define the size and colour of the points by the following rule structure:

If  “X”(string) = 'a'  and   "Y"(real) <= 0.04   and   "Y" > 0.03

The problem is that when I try to export it using Qgis2web the layer points disappear no matter if a use Openlayers or Leaflet.

I’ve changed the rule-based symbology and used a simple one. That just worked so I deduce that rule-based symbology and qgis2web don’t get well togheter.

I found here that the rule-based renderer is tied heavily to Qgis expressions, and only a small subset of Qgis expression functions are so far supported in Qgis2web.

I wonder if is there any possible way that I can use my first style layer and not the simple one to export it.

By that I mean, a different rule expression allowed in qgis2web that gives the same result layer,  a different way in qgis to recreate my different point size and colour layer that’s also complatible with the qgis2web pluggin, or even any other pluggin that I can use to export my rule-based symbology layer.


texts: Extract administrative boundaries with physical borders at the coast line from OSM data

texts:
I successfully downloaded the OSM data from https://www.geofabrik.de/ and converted the .osm.pbf file to a .gpkg using og2ogr. Now I want to extract administrative boundaries using DB-Manager and SQL:
SELECT *
FROM multipolygons
WHERE admin_level = 4

However, I need the physical borders, which do not fall in a line with the administrative borders in coastal areas (see picture). Having a look at OSM, for example via QGIS, the information is there, but I can not figure out which key to use and how to extract them.

What would be the correct way?


texts: Adding layer to MXD from different SDE slows down ArcPy script

texts:
My python script looks at layers in an MXD, then does a couple of clip and measure operations and generates a PDF.

It works accurately and quickly (0.5-3 seconds) when all the layers in the MXD are from SDE database #1.

It works accurately and slowly (20-50 seconds) when one of the layers in the MXD is from SDE database #2 -- even though:

the analysis process is only happening on layers from SDE #1 -- the data from #2 is just there to look pretty
all the data is in the same projection

The only time that the script might be dealing with that particular layer is when the MXD and layer variables are first all defined. (Each of these layers -- Annexations, Customers, Pole -- are in SDE #1.)

mxd = arcpy.mapping.MapDocument("CURRENT")
df = arcpy.mapping.ListDataFrames(mxd, '*')[0]
lyrAnnexations = arcpy.mapping.ListLayers(mxd, "Annexations", df)[0]
lyrCustomers = arcpy.mapping.ListLayers(mxd, "Customers", df)[0]
lyrPoles = arcpy.mapping.ListLayers(mxd, "Pole", df)[0]

I rebuilt the MXD from scratch with the same layers, and saw the exact same issue. Add layer from SDE #2, the script is much slower. Take layer from SDE #2 out of the map, the script goes back to its original speed.

I asked the network/database guy if he knew what was going on, and he'd never seen anything like it. How can a layer that isn't involved in an analysis have any impact on performance, let alone slow it down 10x?


texts: Inserting a ENVI dataset with API into QGIS

texts:
I got an API with a ENVI Data source and an Authorization and I would like to add the data to my project in QGIS
Is there a way to do that?


texts: ArcGIS Online Free Account Transferred to Paid Creator Account Will Not Upload Data

texts:
I created a web map using a free ArcGIS Online (AGOL) plan. After trying to upload a vector file I was prompted to upgrade the map to a paid subscription to accommodate a larger shapefile.  I upgraded the account and purchased a Creator account as prompted.
I can log into my new AGOL Creator account and see the original map I made while still on the free plan. I would like to proceed forward to develop this map taking advantage of the hosted storage of the Creator plan.  After opening the map, and selecting my file to upload (from my computer) I am again prompted to purchase a paid plan.  Apparently, the Creator account's increased storage privileges did not get transferred to this map.  How do transfer/alter this existing map to my Creator account?


texts: Multiply raster by Numpy derived float value

texts:
So I have a python script that uses numpy to calculate the mean of a field, then I try to use that mean value in a raster multiplication as follows:

na = arcpy.da.TableToNumPyArray("%SCRATCHWORKSPACE%\\rc_instruments_sr_noNULL2", "ratio")
    dMeanRatio = numpy.mean(na["ratio"])
    outCorrectedSR = Raster(outSimulated) * dMeanRatio

This gives an error of:

TypeError: expected a raster or layer name

This error is not referring to "Raster(outSimulated)" because I've replaced dMeanRatio with a hard-coded value and the script runs just fine. So my question is, how do I convert the dMeanRatio into a value that is acceptable to the raster multiplication?


texts: What is QGIS equivalent of ArcGIS GetCount for use in Python script?

texts:
I have a number of layers in QGIS that I would like to select attributes on, and get the feature count for, within a Python script.

In ArcGIS for Desktop I have been using SelectLayerByAttribute and GetCount.

What is the QGIS equivalent?


texts: GEE: Extracting country names from FAO dataset

texts:
The FAO GAUl dataset linked below has an 'ADM0_NAME' attribute, and I would like to see/ export these country names. However, I cannot print or export these names because I get an error.
var countries = ee.FeatureCollection("FAO/GAUL/2015/level0");
var names = countries.select('ADM0_NAME')
print(names)

I then get the following error:
FeatureCollection (Error)
User memory limit exceeded.

(https://developers.google.com/earth-engine/datasets/catalog/FAO_GAUL_2015_level0#table-schema)


texts:
We have a little side project running where we need to create maps of different countries containing postcodes and postcode area boundaries, which are put in a database. Those maps are supposed to be usable for several reports and visual display - through the database entries.
It appears it needs to be done from scratch (meaning we won’t pay for the data).

I’m completely new to QGIS and georeferencing and searched through several tutorials, but some clues are missing.

So far I have a blank shapefile of a country, which is already georeferenced.
Now what i need to do (from my understanding):

Get the post code data and import it (which is hard for countries like    Romania, since there seems to be nothing like tables available for free), OR put only the post code areas into the map manually.
-> But how can you do that somewhat accurately? Just using polygons?
Make it usable for a database, e.g. using SpatiaLite.

Are there any good resources about the steps of such a specific task?

edit: Sorry for causing some confusion. I adjusted the question a little.


texts: What do different realizations of the same datum have in common?

texts:
Firstly, am I correct in assuming that the process of determining the shape and the size of the spheroid is independent from the process of determining the origin and the orientation?

For example, for NAD27, it seems that the following equation -- called the Clarke 1866 ellipsoid (reference 1) -- was fitted to the surveyed points:

where the parameters to be determined were x_0, y_0, and z_0 instead of a general equation of an oblate spheroid:

Suppose we have the complete equation (with parameters already determined), and the origin and the orientation relative to some arbitrary standard coordinate system, of each different realization of all the datums. Further suppose we do not know the name of each. For example, NAD83(CORS94) may be called "1" and WGS84 (G730) perhaps "2". So they are all jumbled up. Is it possible to sort them into different groups so that each group consists of the different realizations of the same datum?

References:

ARSITECH "Constants for Reference Ellipsoids used for Datum Transformations"
http://www.arsitech.com/mapping/geodetic_datum/


texts: Updating repository XML after deploying QGIS plugin

texts:
When I deploy my plugin to the qgis repository on my server, I usually do:

pb_tool zip
scp plugin.zip user@ip:/location
ssh user@ip
nano repository.xml # and change stuff

Last step is what this question is about - update repository XML on the server, so it matches metadata.xml of the plugin.

How to automate this? Keeping in mind, that one repository can have multiple plugins?

I thought of having a bash/Python script that reads metadata.xml of given plugin and substitutes appropriate strings in repository XML. But that seems to be really ugly. How do you handle this?


texts:
I am using a PostGIS function (ST_StraightSkeleton to be exact) and I am getting an error with my data that reads:

ERROR:  Intersection does not support polygon with connected rings

I understand what rings are, and my data does indeed have holes, but I don't think this is the problem because this process uses the same function on a dataset with holes.

I don't understand what a "connected ring" is and Google is failing me. Can someone explain it? I need to find a way to remove them from my dataset apparently.


texts: What is the UOM of Microsoft Footprint 2017 height data

texts:
It says on the website that height's unit of measurement of MSFP 2017 is in meters. However, as I looked into some real building examples on Google Map street view, they height looks very inaccurate (seems to make more sense if it's in feet). Here's are some examples:

Coordinates: 36.717682, -119.78708. Height according to MSFP 2017: 10.11 meters (doesn't make sense, 10.11 feet is more reasonable?)

Coordinates: 36.719514, -119.78763. Height according to MSFP2017: 10.51 m.

If you know more accurate building height data sources please let me know.


texts: How do I determine the hardware requirements for importing an OSM extract?

texts:
I have a 950MB OSM extract I'm trying to import into Postgres using osm2pgsql. All attempts thus far have failed as I've been running out of hard drive space. I have no problem processing metro or individual state extracts. I've encountered nothing but problems trying to import the US Midwest region extract.

I'm likely to build a new VM to process this file so I'd like to get this right. How do I determine the hardware requirements for this process?


texts:
I have a las catalog with overlapping tiles. LidR warns me of this, and I'm wondering how concerned I should be. The tiles are in las version 1.2 with point record format 1. None of the points are marked as Withheld so there is no filtering to be done.  I did a test on two of the tiles and found that the area of the intersection is 2200 square meters, and each of the tiles are about 225 hectares. This seems like a relatively small portion of overlap. Do I need to modify the tiles? Or can I process them like normal?

Las Check:

> lascheck(ctg)

 Checking headers consistency
 - Checking file version consistency... ✓
 - Checking scale consistency... ✓
 - Checking offset consistency... ✓
 - Checking point type consistency... ✓
 - Checking VLR consistency... ✓
 - Checking CRS consistency... ✓
 Checking the headers
 - Checking scale factor validity... ✓
 - Checking Point Data Format ID validity... ✓
 Checking preprocessing already done
 - Checking negative outliers... ✓
 - Checking normalization... no
 Checking the geometry
 - Checking overlapping tiles...
    ⚠ Some tiles seem to overlap each other
 - Checking point indexation... no
 - List item

Catalog Info:

> ctg
class       : LAScatalog
extent      : 554999.8 , 562500.1 , 4555500 , 4561501 (xmin, xmax, ymin, ymax)
coord. ref. : NA
area        : 33.61 kunits²
points      : 370.01 million points
density     : 11 points/units²
num. files  : 16

Catalog Plot:


texts: Export some specific pages from data driven pages to pdf file by python

texts: I need to export pages 14, 13, 23 and 24 of 26 pages of data driven pages to pdf file by python and then add the pdf file to the final pdf file by arcpy (python), any help please?


texts: How to mosaic multiple MODIS images of Multiple months using MRT tool?

texts:
I have monthly MODIS images for 2001-2014.
For few images mosaicing with reprojection can be done using MRT tool. However, for large number of images it's a cumbersome process.

There are multiple folder in a pattern, where each image is saved in a particular folder (representing MONTH) and all the months of a particular year are stored in a specific folder named after corresponding YEAR. Below is the example for reference

MODIS\Y2000\M01
Image#1
Image#2
...

MODIS\Y2000\M02
Image#1
Image#2
...
...
...
...
...

MODIS\Y2014\M12
Image#1
Image#2
...

Is there any way out to automate this process using MRT tool?


texts:
I converted DTM Data from LAZ to ASCII.
Now I got some artifacts with an obviously wrong cell values.
Is there a simple way to edit just a single cell?
I tried looking for the cell value in Notepad++, but there are several cells all over the area that have a similar value, so I can't find my cell easily.
Any ideas?


texts: Editing with projected coordinate system defined for data frame

texts:
I am looking for some confirmation on my understanding of the information message Spatial Reference Does Not Match that I receive when starting an edit session.  I have been told by so many that you can ignore it in most cases; however, I have done research that has indicated that there is a danger when editing data that is projected, and has a different coordinate system than the data frame, so I want to make sure that I have not been incorrect all along.  The coordinate system of my data frame is almost always set to the regional State Plane projected coordinate system, as well as the new shapefiles I create.  I have always thought that as long as the layer you are editing has the same geographic coordinate system as the data frame, your edits will be accurate because even though State Plane is a projected coordinate system, it is based on the same datum.  In regards to ignoring the message, I believe it can be accurately ignored as long as the layer you are editing is not in the list of layers whose spatial reference do not match (i.e. some statewide or nationwide datasets that likely have a different coordiante system than my project-level data).
Can anyone confirm or deny my thoughts so I can have greater confidence?


texts: Slow counting of points in polygon in QGIS 3.4

texts:
I am working on a plugin to count the point in polygons by PyQGIS in QGIS 3.4, and I found that compared to the same function in QGIS 2.18, it is very slow. I do not know the reason. My code is listed as follows.

processing.run("qgis:countpointsinpolygon",{'POLYGONS':mesh_layer,'POINTS':select_crash_layer,'WEIGHT':None,'CLASSFIELD':None,'FIELD':'Count','OUTPUT':category_mesh_layer})

I used this function to calculate approximate 25000 points to 3761 polygons to find how many points in each polygon. I realized it in QGIS 2.18 within 10 s. However, in QGIS 3.4 it needs more than 5 minutes.


texts: Buffer a single line, don't dissolve overlapping areas of buffer

texts:
Is it possible in ArcGIS to buffer one line and have the buffer not dissolve together the overlapping areas. Other questions here addressed dealing with overlapping areas of multiple features within the same feature class. But I'm curious about the same feature.

A buffer like what is drawn is what I'm looking for.


texts:
I am using Shapely's wkt.dumps function to get WKT from geometry objects.

How do I limit the precision of the coordinates?


texts: Using timeslider in QGIS to show census data per year?

texts:
I have several shapefiles with the same geographic extend (in my case: Berlin and all its administrative subareas) and statistical (census) attribute data for every subarea displayed trough a choropleth map.

Now I want to show the changing of the figures throughout the given period of time (2007-2015) with QGIS.

Is it possible with several shp or do I have to save all the attribute data for every year into one single shapefile? But if so, how to apply with only one time column (year) but repeating geographical reference (code of subarea) every year then?