texts: node.js GeoLib fails on polygon that straddles the anti-meridian

texts:
geolib.isPointInside does not appear to work correctly when the polygon longitudes straddle the anti-meridian.
I maintain a system with lots of polygons that do straddle the antimeridian (I'm in Australia :-) so it's important to me.

see the test code below that illustrates the problem when it's run in node:

var geolib = require('geolib');

// -------------
// Test 1 This succeeds
// -------------
var polygon = [
    {latitude: -46, longitude: 1},
    {latitude: -46, longitude: 2},
    {latitude: -47, longitude: 2},
    {latitude: -47, longitude: 1}
];
console.log('Test 1: Polygon:',polygon);
console.log('should be true', geolib.isPointInside({ latitude:-46.5, longitude:1.5}, polygon))
console.log('should be false', geolib.isPointInside({ latitude:-47.9, longitude:1.2}, polygon))
console.log('should be true', geolib.isPointInside({ latitude:-46.8, longitude:1.1}, polygon))
console.log('------------');

// -------------
// Test 2 This succeeds
// -------------
var polygon = [
    {latitude: 10, longitude: -10},
    {latitude: 10, longitude: 10},
    {latitude: 20, longitude: 10},
    {latitude: 20, longitude: -10}
];
console.log('Test 2: Polygon:',polygon);
console.log('should be true', geolib.isPointInside({ latitude:15, longitude:-5}, polygon))
console.log('should be false', geolib.isPointInside({ latitude:15, longitude:-15}, polygon))
console.log('should be true', geolib.isPointInside({ latitude:15, longitude:-5}, polygon))

// ----------
// Test 3 : AntiMeridian 1: This fails
// ----------
var polygon = [
    {latitude: 10, longitude: 175},
    {latitude: 10, longitude: -175},
    {latitude: 20, longitude: -175},
    {latitude: 20, longitude: 175}
];
console.log('Test 3: Antimeridian 1: Polygon:',polygon);
console.log('should be true', geolib.isPointInside({ latitude:15, longitude:177}, polygon))
console.log('should be false', geolib.isPointInside({ latitude:15, longitude:-170}, polygon))
console.log('should be true', geolib.isPointInside({ latitude:15, longitude:-176}, polygon))
console.log('------------');

// ----------
// Test 4: Antimeridian 2: This fails
// ----------
var polygon = [
    {latitude: 10, longitude: -175},
    {latitude: 10, longitude: 175},
    {latitude: 20, longitude: 175},
    {latitude: 20, longitude: -175}
];
console.log('Test 4: Antimeridian 2: Polygon:',polygon);
console.log('should be true', geolib.isPointInside({ latitude:15, longitude:177}, polygon))
console.log('should be false', geolib.isPointInside({ latitude:15, longitude:-170}, polygon))
console.log('should be true', geolib.isPointInside({ latitude:15, longitude:-176}, polygon))
console.log('------------');


texts: Startup/Shutdown script for ArcGIS Enterprise?

texts:
I am after a startup/shutdown script for the newly named ArcGIS Enterprise at 10.5.

This is my home machine, where I have installed ArcGIS Server, ArcGIS Data Store, Portal for ArcGIS and also GeoEvent Server.  Therefore I just want to startup the Web GIS when I want to do some development.  I know a Virtual Machine might be a good option here, but can we rule that out as an answer.

I remember back at 9.3.1, I used to make use of some shell scripts that I would have as shortcuts on my desktop to Start SQL Server, SDE, ArcGIS Server, but have lost them.  I remember they had some logic to wait for other service to finish starting up before the next one would startup.

Will accept an answer that includes everything but GeoEvent.


texts: Point within a Polygon within another Polygon

texts:
I have a table with towns (polygons) and another with company assets (points). I am trying to get the town these assets are located in spatially.

I have some point geometry that is contained within two polygons but I want to return only the inner polygon.  I am using PostgreSQL and all of the functions (i.e. St_Intersects, St_contains) are giving me the same duplicate records where the point is showing up in both Smithville and Smithville - Village.  How do I get only the inner polygon?  .

 select distinct c.id,t.division,t.city,t.state from compequip c, a.geom g,
---geometry table a.geom tg, ----geometry table a.towns t, a.station s
 where c.id = g.id and tg.id = t.id and c.stationid = s.stationid
 and st_intersects(st_makevalid(st_curvetoline(g.geom)), st_curvetoline(tg.geom))


texts: OpenLayers map shows partial map on page load, unless dev-tools are opened

texts:
I have a OpenLayers map with a few KML layers being added dynamically on top of it. On page load I see only a part of the map in top left corner of designated html div. I can't figure out the source of this issue. Neither of projection, zoom level seems or initial position is affecting the load.
It is strange that when I open the chrome dev-tools (or Inspector window for Firefox), map behaves exactly as expected and complete map is displayed.

Note: above issue is only found with base layer from Google such as "Google Streets" while with base layer such as OSM it works fine.

Browsers tested on:

Chrome Version 30.0.1599.101 m
Firefox 23.0.1
Firefox 24.0


texts:
I need to do some raster/vector processing on some non-overlapping tiles.  My python script needs to find the correct tiles and then extract the raster values and count them up.  My script needs to be flexible enough such that sometimes there are no tiles to process while othertimes there could be up to four (since they are tiles, and my vector is so small, it is impossible for there to be more than four)

There is no switch/case in python so I am struggling at an efficient way of doing this.  So basically I need to construct the variable inRaster out of 1,2,3,or 4 different rasters.

I can do it when there is only one but the problem is that with all these search cursors and loops, (a) my script is really slowing down and (b) it is not clear to me how to add 2, 3, or 4 rasters.  I'd like to do something like inRaster=Plus(tile1,tile2,tile3) or inRaster=Plus(tile1,tile2,tile3,tile4).

I figured I would see if the intersect between my polygon and tiles has any rows in it (it should have 0,1,2,3,or 4 rows) and if there is some data, then construct inRaster

While I can do it for one row, I have trouble understanding how to loop thru all the rows.

        arcpy.MakeQueryTable_management(tempIntersect, tempTable,"USE_KEY_FIELDS")

        # Print the total rows
        numberOfRows=int(arcpy.GetCount_management(tempTable).getOutput(0))
        arcpy.AddMessage("       numberOfRows for distance = " + str(dist) +"m is=" + str(numberOfRows) + " ...")
        if not (numberOfRows is None):
            if numberOfRows   == 1:
                for myTileName in arcpy.da.SearchCursor(tempIntersect, "TILENAME"):
                    myTileFullName=str(myTileName[0])+".tif"
                    inRaster= os.path.join(TSpace,myTileFullName)
                    arcpy.AddMessage("       numberOfRows " + str(inRaster) + " ...")

                    outExtCircle=ExtractByCircle(inRaster, arcpy.Point(x,y), dist, "INSIDE")
                    rtemp          = os.path.join(WSpace,"myRastTemp")
                    if arcpy.Exists(rtemp):
                        #arcpy.AddMessage("       Deleting " + rtemp + " ...")
                        arcpy.Delete_management(rtemp)
                    outExtCircle.save(rtemp)
            else:
                #  Here I am trying to loop through and get all the tiles (and failing)
                myTestString=""
                for myTileName in arcpy.da.SearchCursor(tempIntersect, "TILENAME"):
                    myTestString+=myTileName[0]+".tif"
                    arcpy.AddMessage("       tilename is " + str(myTestString) + " ...")


texts: Is it possible to generate a vector polygon from a polygon on a raster image?

texts:
I have a raster map which I have extracted from a government report providing flood extent return frequencies. I want to 'snatch' these flood extent polygons (i.e. convert to vector) and I was wondering if there is a simple way to do this in QGIS. Simply digitising these layers by tracing new polygons would be too time consuming and potentially result in human error.

Is there a tool or methodology someone could recommend to me which would recognise these polygons in raster format (there are clearly outlined and coloured distinctly differently to the background) so I can convert them to vector?

Any help with this is very much appreciated.

Camilla


texts: How to store binary data referenced by 3D coordinates (xyz)?

texts:
I would like to store huge amounts of voxels (3D cubes with color and some other values) in a gis database.

Does someone can recommend a database (indexing) system, where I can store binary data (a sparse voxel octree) and find it by a 3D coordinate query?


texts: Adding a column to polygon table with correct corresponding address

texts:
I have 2 shapefiles:

one with polygons (buildings) and their properties
other one is a point address shapefile - it contains points that are inside polygons from the other SHP, but it's still another layer.

I need to add 3 columns to Polygons layer from Address layer. How to do it so the records will match correctly?

So far I got to know I can use ST_Intersect in SQL, but I'm not sure it's the easiest way to do it.


texts: Aligning pixels for using GSW as mask layer with Landsat 8 and Global Surface Water (JRC) Layer

texts:
I want to isolate some surface water areas on Landsat 8 images by using the Global Surface Water dataset in GEE as mask. After mapping both in the Editor I realized the positions of pixels of GSW and Landsat8 do not align, instead there is an offset in both dimensions (see image).

I did simply load and display the data. Didn't reproject any, as GEE states to harmonize Layers automatically. Looking in the properties of the images there are different projections assigned to them (LS: EPSG:32624 and GSW: EPSG:4326).

Is it necessary to reproject one of the layers?

I did try the function .register(), which should find displacements and correct them. However the output layer didn't have any changes in pixel position.

var LC08_CE1_2019= ee.Image('LANDSAT/LC08/C01/T1_SR/LC08_217063_20190819');
var gsw = ee.Image('JRC/GSW1_0/GlobalSurfaceWater');
var extentGSW = gsw.select('max_extent')
  .clip(Ceara)
  .unmask(0);
var regGSW= extentGSW.register({referenceImage: LC08_CE1_2019.select('B1'), maxOffset:15.0, patchWidth:60.0, stiffness:10});


texts:
I'm having some difficulties using the hyperlink tool (ControlsMapHyperlinkTool) in ArcEngine. Does anyone have any experience with it?

It looks like it behaves correctly except that every time a user clicks a hyperlinked feature it shows a message saying "Unable to launch hyperlink". I'm using a "script" hyperlink and it performs fine in ArcMap. It doesn't seem to be trying to execute the script on click. I checked the IHotlinkExpressionProperties object and the Expression property is set correctly. Any thoughts? Thanks!


texts: PostGis - st_union polygons, keep sum of attribute

texts:
My current query is:
    SELECT
        (ST_Dump( ST_Union( way ) )).geom AS way,
        10000 AS way_pixels
      FROM (
        SELECT
          way,
          way_pixels,
        FROM polygons
      ) AS test

it works fine to obtain and dissolve and split into indivisual polygons. But is there a reasonable way to compute sum of way_pixels attribute over merged polygons?
For example polygons A, B, C were merged. A has way_pixels value of 20, B of 120 and C of 19. So I want merged polygon to have way_pixels set to 159, not 10000 as it is done currently.
Computing sum has a known problem. With overlapping rather than joining polygons it is going to produce misleading values. But it is preferable to assuming that size is always 10000.
ST_Area is not solution as way_pixels values are not real area but projected area (in Web Mercator).
Merge two polygons, attributes are added - Tool? and Sum of attributes in Spatial join when polygons cross boundaries in ArcMap? Sum attribute in overlapping polygons are similar but about ArcGis


texts: Fill color and outline color with the same value

texts:
I made a categorized symbology on a shapefile, but the color value of outline (stroke color) is different than the color of fill color option. How can I ensure that both colors are the same for fill color and outline color using categorized symbology?
See an example below.


texts: Creating 3D polyline in ArcGIS for Desktop?

texts:
I have a csv file containing a breakline points' coordinates and their elevation.

How can I create a 3d polyline from these points?

The CSV file looks like this:

  X              Y         z
384282.6    244807.4    22.1807
384376.4    244791.2    12.6306
384467.5    244773.5    12.634
384535.7    244745      12.5274
384593.4    244701.9    12.6257


texts: Converting from Washington State Plane Coordinates to Web Mercator?

texts:
I received coordinate data in Washington State Plane: EPSG:14602. I need to convert this to EPSG: 3857, aka Web Mercator.

I'm using Javascript. Here's what I've tried so far:

Using OpenLayers: ol.proj.transform. This doesn't work since it does not contain the definition for EPSG: 14602.

To add definition, I used ol.proj4js. Issue with that is that it also needs a definition to begin with.

Looking for that definition, I found this website: https://epsg.io/14602 If you try https://epsg.io/14602.js, you get an empty definition. But if you try https://epsg.io/3857.js, you get a full proper definition.

How/where do I get the definition as formatted in https://epsg.io/3857.js?

Do you know of any tools that convert or might be close?

Perhaps, guide me towards the math and I can do it manually?


texts: How to get 'IGraphicsContainer' from the map?

texts:
I tried to get IGraphicsContainer from Data View (active map) using FocusMap, ActiveView

IGraphicsContainer graphicsContainer = (IGraphicsContainer)neApp.Document.MxDocument.FocusMap;

but the variable always null, is there any different way to do this?

I already get from page layout view in this question right here:
Get IGraphicsContainer from page layout
but i need it from data view (map).


texts: Point/Area location of data in reference to geotransform of raster

texts:
I'm trying to understand how to interpret the position of pixel data within a GDAL raster, specifically as it applies to files containing height data.  I'm create Triangular Regular Network (TRN) for terrain from various GDAL sources and need to generate the terrain as precisely as possible.

I've read the follow sources several times, but I have failed to understand whether they answer my questions:

  Is there a standard for the coordinates of pixels in georeferenced rasters ?
  Seeking tool to generate Mesh from DTM?
  Does a value in srtm dem file refer to the center or corner of a cell?

The GDAL documentation has this to say about the GeoTransform:

    Xgeo = GT(0) + Xpixel*GT(1) + Yline*GT(2)
    Ygeo = GT(3) + Xpixel*GT(4) + Yline*GT(5)

  In case of north up images, the GT(2) and GT(4) coefficients are zero,
  and the GT(1) is pixel width, and GT(5) is pixel height. The
  (GT(0),GT(3)) position is the top left corner of the top left pixel of
  the raster.

  Note that the pixel/line coordinates in the above are from (0.0,0.0)
  at the top left corner of the top left pixel to
  (width_in_pixels,height_in_pixels) at the bottom right corner of the
  bottom right pixel. The pixel/line location of the center of the top
  left pixel would therefore be (0.5,0.5).

I understand the last paragraph in the quote above. However, my confusion begins when trying to understand the geo-referenced location of a pixels value within the pixel itself (height in my case) and the interaction of PixelIsPoint vs PixelIsArea.

So, by way example, if we consider a very simple raster:

float RasterData[4][4] = {
    { 1.0f, 2.0f, 1.0f, 3.0f },
    { 1.5f, 1.5f, 2.0f, 2.5f },
    { 2.0f, 2.0f, 2.5f, 2.0f },
    { 1.5f, 2.0f, 2.0f, 1.5f }
};

float GeoTransform[6] = {
    0.0f, // X origin
    1.0f, // Xaxis.x of transform
    0.0f, // Xaxis.y of transform

    0.0f, // Y origin
    0.0f, // Yaxis.x of transform
    1.0f  // Yaxis.y of transform
};

By the GDAL documentation,

the top-left corner of RasterData[0][0] has a geo-referenced position of (0.0, 0,0)
the bottom-right corner of RasterData[3][3] has a geo-referenced position of (4.0, 4.0)
the center of RasterData[0][0] has a geo-referenced position of (0.5, 0.5)

This part I believe I understand well enough.

However, I do not understand how to determine the (geo-referenced) location of a pixels data within the pixel itself, nor do I understand whether PixelIsPoint / PixelIsArea affects the location of the pixel corners, or the location of a pixels data within the pixel itself.

For example, I'm not entirely certain whether this is the correct way to create a vertex from this dataset:

// loop pixel_y 0 -> height (integer)
// loop pixel_x 0 -> width  (integer)
vertex.x = GeoTransform[0] + GeoTransform[1] * pixel_x + GeoTransform[2] * pixel_y;
vertex.y = GeoTransform[3] + GeoTransform[4] * pixel_x + GeoTransform[5] * pixel_y;
vertex.z = RasterData[y][x];

The method above assumes the value of RasterData[y][x] is geo-referenced at the top-left corner point of the pixel -- this is where my doubt begins.

So my question is:
How do I determine the correct geo-referenced position of a given pixel value RasterData[y][x] when interpreting the data under PixelIsPoint versus PixelIsArea?

IE: I want to take RasterData[y][x] and GeoTransform and generate (X,Y,Z) vertex positions where:

The X coordinate is the appropriate geo-referenced position for RasterData[y][x]
The Y coordinate is the appropriate geo-referenced position for RasterData[y][x]
The Z coordinate is the height value contained in RasterData[y][x]


texts: Deleting points outside polygon with PostGIS

texts: How can I delete all points (rows) in one table that fall outside a polygon stored in a different table? Assume the tables follow the common convention of storing the geometry in a column named geom.


texts: Merging edited data back into ArcSDE master geodatabase?

texts:
I need some suggestions/inspiration for a problem we have.

The problem:

We have an ArcSDE master GIS database. We will edit this data almost exclusively outside of ArcGIS and outside of ArcSDE. In other words, we will make an export of some small area and edit this. Then we will push this data back into an other SDE table (diff database).
We will also maintain an extra attribute that contains: add, delete, edit.

What we want to do now is write this diff table to the master. My question is what is the best tool for doing this?
Using Arcpy with cursor, using FME, or are there other Arcgis tools to help?
I was thinking using ArcPy and doing something like this:

simply add all objects with status add.
search for every id with status delete and delete it
search for all ids with status edit and copy the edited ones

What do you think would be the best tool to do this?

Arcpy?
FME
direct SQL (we are using MS SQL)
others...

The master database will be in the range of about 500.000 objects, so performance is maybe an issue. At the moment we are thinking about a sync operation where an operator starts the sync and can wait for it to finish so he can do some checks after it is completed.

Bonus questions:
We will need to sync dimensionings and labels too. The problem is that, in the diff database, we probably cannot store them as arcgis dimensioning feature classes but they will be line and point symbols with the necessary data attachted to it. Can we use the same technique for transforming and migrating this data so it fits in the dimensioning classes?


texts:
I am having problems with ArcGIS Pro and ArcGIS Online.

I am working in an already existing database and I am creating polygons to symbolise areas in which types of animals have been seen.

I am missing a value (a fish sort) but I cannot figure out how to add a value. When I do it in ArcGIS Pro it will not register in ArcGIS Online.

Should it be done via ArcGIS Online or ArcGIS Pro?

I need to add Arctic Char as a type But I am unsure of whether that should be done via ArcGis Online or ArcGIS Pro. My logic says that it should be doable through ArcGIS Pro. I have tried to add value in attribute table, without any luck, I have also tried out designs, without any luck.


texts: Set ZIndex for Google Fusion Tables in Google Maps API

texts:
Is there any way to set the ZIndex of a Google Fusion Table layer in the Google Maps API?  I have several Fusion layers added as overlays to a map with the ability to turn them on and off.  I want certain layers to always be ontop if they are on, however Fusion Tables/Google Maps just renders them in the order that they are drawn.

For example:

Layer 1 starts On.
Layer 2 starts Off. Layer 2 should Always be on top if it is being shown.

Click Layer 2 On, and it is on top.
Now, click Layer 1 Off and the On again, and now it is on top as it is the most recent layer drawn.

I'd ideally like to set the ZIndex of Layer 2 to higher than Layer 1 so that it is always on top, but there doesn't seem to be a way to do this.  Any ideas?


texts:
I'm a basic user of ArcMap and currently trying to figure out a solution for an issue i'm having with my data driven pages. The data frame layer that I'm using for the data driven page has a field with simple numbering 1,2 and 3 which correlate to a specific type of soil that it represents. I'm using this numbering as the page name field and using page definition queries on my other layers to only show attributes that correlate to that specific soil type. The issue is, a couple of other layers which I also need to change with each driven page do not contain that specific numbering as it is not related to soil types but instead it is the location ID, which is basically a redline boundary which would ideally change as driven page changes.

Is there a way I can create a custom page definition query which uses another attribute within the driven page attribute table and not the page name field?

At the moment my only solution to this is to simply use a simple definition query and change the location ID with each map but ideally i'd like this to be a little more automated.


texts: Does ArcGIS for Server support others using own credentials (Single Sign On) for Federated Services?

texts:
I have a number of secured ArcGIS 10.2 for Server (AGS) services which a number of other companies use as well. At the moment we have created domain accounts on our domain which they need to access the services.

However, does ArcGIS 10.2 for Server allow me to setup my site as a federated resource partner so that the other companies can use their own credentials instead (Single Sign On).

I know that arcgis.com and the Portal will allow you to do it but haven't found anywhere about doing it on your own AGS installation?


texts: What is the going rate for GIS freelancers?

texts: There appears to be very little information available about competitive hourly rates for GIS freelancers.  Additionally, GIS freelancers never post rates on their websites.  Upwork has some information on what people are asking, although I believe the rates are well below  the industry standard due to the structure of their bidding system (example).  The US Bureau of Labor Statistics states the median pay for Cartographers and Photogrammetrists is $27.62 USD, although this is presumably low because it does not taking into account fringe benefits. What is the going rate for GIS freelancers in the USA?  I am interested in hearing from both freelancers and businesses who contract out work.


texts: Constructing a weighted graph a.k.a. a road network with distances for a specific area

texts:
Pardon my ignorance in GIS/Routing tools as I'm new to it:

I have a map in OSM format and wanted to build a routable graph out of it ( to experiment a searching algorithm I have in my PhD thesis).

I went to www.openstreetmap.org and selected an area and exported it. I got an map.os.xml file.
I downloaded osm2po as I read somewhere here that it gives you a routing graph.
Running java -jar osm2po-core-5.1.0-signed.jar prefix="osm2poConversionResults" map.osm.xml in the command line.

osm2po is poorly documented after this step; It just give me a nice map on http://localhost:8888/Osm2poService . I need a graph where every edge is annotated with its distance i.e. a weighted graph. Looking at the directory osm2poConversionResults I found file .gph which seems to be the one I am looking for, but opening it through TextEdit in Mac OS, get unreadable symbols.

So my question is: How to construct a weighted graph out of .osm.xml file in a very basic way with minimal reliance on complex software? Or at least a graph (I can look for the distances later). Another small question: Does the way of converting .osm to a graph differ substantially from software to software or is there a standard for the conversion?


texts: QGIS 2.0 Automatic generation of a vector layer out of a point layer

texts:
I am an environmental engineer and the aim of my QGIS 2.0 project is to visualize the results of a sewage water system inspection. In a first step I included the hoppers from a csv file and added a google streets layer.
Now I have to draw the pipes between the hoppers.
One way would be to create a vctor layer and draw one pipe after the other.
However, I have 522 pipes to draw. So I wonder whether QGIS is able to do it automatically (draw the line of shortest distance between two hoppers?!)

I tried the "nearest neighbour analysis" After this tutorial (http://www.qgistutorials.com/en/docs/nearest_neighbor_analysis.html) it seems to work for different layers. I would like to find the nearest neighbour points inside of my hoppers layer.
I also installed the mmqgis plugin in order to make a hub diagram. I loaded a table with start and end hopper number of the different pipes into qgis and tried to link it with the hopper layer. I got this error message:

  Traceback (most recent call last): File
  "C:\Users\stutte/.qgis2/python/plugins\mmqgis\mmqgis_dialogs.py", line
  1427, in run message = mmqgis_hub_lines(self.iface, hubname, hubattr,
  spokename, spokeattr, savename, 1) File
  "C:\Users\stutte/.qgis2/python/plugins\mmqgis\mmqgis_library.py", line
  3377, in mmqgis_hub_lines spokex =
  spokepoint.geometry().boundingBox().center().x() AttributeError:
  'NoneType' object has no attribute 'boundingBox'

Do you have any idea?


texts: Experiencing transformation anomalies when performing a Save As to convert coordinates between NAD83-Z14 and NAD27-Z13

texts:
I have recently discovered QGIS and have been hoping to find a solution for converting bacthes of xy coordinates between two projection systems for our mapping system.  I have known good coordinates for both datums that I am trying to convert between, but even though I believe I'm doing the right steps in QGIS, my coordinates are still out by a small margin.
Prior to finding QGIS, for single coordinate conversion I have been using the free Java-based Coordinate and Data Transformations tool that is available on the internet from Charles L. Taylor. I've attached an image of the two known coordinates from this tool.
I've provided the text data of my coordinates.
Original coords in NAD83/WGS84 Zone 14N (EPSG:4326)
Lat/Long: 49.362977,-101-742026
UTMs: x 300914.226, y 5471423.864

Converted coords in NAD27 Zone 13N (EPSG:26713)
Lat/Long: 49.362924,-101.741634
UTMs: 736577.110,5472689.933

I had imported my list of coordinates (original datum is NAD83/WGS84) into QGIS to get my list of points.  After adding the corresponding UTM coordinates using the Export/Add Geometry function, I can see that the UTMs in the attribute table are the same as my original .csv file.
After conversion to what I thought was NAD27-Z13 using QGIS, my first coordinate comes back with the following coordinates:
Lat/Long: 49.36297733,-101.74205733
UTMs: 736548.431,5472694.612

When I plot the above coordinate on the map I have with the known good points already on the map, obviously they don't line up.  What I'm wondering is if I'm missing a step or if there are options that must be set during conversion.  From what I can tell of QGIS, if I can find a solution to this, it will make batch conversion between datums much quicker.  Many of the coordinate files that I receive are in .csv format and QGIS seems like the perfect application for what I need, however, I'm hoping that someone that has more experience can point me in the right direction.


texts: Installing Cartaro using Acquia Dev Desktop and link it with GeoServer and PostGIS

texts: How do I configure GeoServer and PostGIS so that I don't get the following error when installing Cartaro using the Acquia Dev Desktop. I am just installing Cartaro locally for testing, or else how can OpenGeoSuite be installed with Cartaro?


texts:
I show you my code, but basically, I have been following the steps in the official bibliography about Map Composition The code works, but I don't get the items I have been adding to the map (scalebar, labels, and legend). Sorry, some words are in Spanish

QGIS 2.14 (Essen)

I'm using Text Editor of QGIS

from qgis.core import *
from qgis.utils import iface
from PyQt4.QtGui import *
from PyQt4.QtCore import *
from qgis.core import *

lyr = QgsVectorLayer("\\ROUTE\\MUNICIPIOS.shp", "Municipios", "ogr")
reg = QgsMapLayerRegistry.instance()
reg.addMapLayer(lyr)
##etiquetamos
#etiquetas = QgsPalLayerSettings()
#etiquetas.readFromLayer(municipios)
#etiquetas.enabled = True
#etiquetas.fieldName = 'Municipios'
#etiquetas.placement= QgsPalLayerSettings.Horizontal
#etiquetas.setDataDefinedProperty(QgsPalLayerSettings.Size,True,True,'8','')
#etiquetas.writeToLayer(municipios)
##fin de etiquetas

#EMPEZAMOS CON LA IMAGEN

i = QImage(QSize(800,800), QImage.Format_ARGB32_Premultiplied)
c = QColor('#cda56d')
i.fill(c.rgb())
p = QPainter()
p.begin(i)

lyrs = reg.mapLayers().keys()
#AQUI EL MAP RENDERER
el_renderer = iface.mapCanvas().mapRenderer()
#el_renderer = QgsMapRenderer()
el_renderer.setLayerSet(lyrs)
rect = QgsRectangle(el_renderer.fullExtent())
rect.scale(1.2)
el_renderer.setExtent(rect)
el_renderer.setOutputSize(i.size(), i.logicalDpiX())
el_renderer.render(p)
p.end()

# Aqui tengo mi COMPOSICION DE MAPA ----> . QgsComposition
composition = QgsComposition(el_renderer)

composition.setPlotStyle(QgsComposition.Print)
#composition.setPaperSize(215.9, 279.4)
# Mi QgsComposerMap
w, h = composition.paperWidth() * .50, composition.paperHeight() * .50
x = (composition.paperWidth() - w) / 2
y = ((composition.paperHeight() - h)) / 2
componer_mapa = QgsComposerMap(composition,x,y,w,h)
composition.addItem(componer_mapa)

#Aqui termina el ComposerMap

#Aqui esta mi etiquetado. QgsComposerLabel
etiqueta = QgsComposerLabel(composition)
#QgsComposerLabel::adjustSizeToText
etiqueta.adjustSizeToText()
#QgsComposerLabel::setText
etiqueta.setText('Municipios')
#QgsComposerLabel::marginX()
etiqueta.marginX()
#QgsComposerLabel::marginY()
etiqueta.marginY()
#QgsComposerLabel::displayText()
#etiqueta.displayText()
#QgsComposerLabel::displayName()
#etiqueta.displayName()
composition.addItem(etiqueta) #DOESN'T WORK

legend = QgsComposerLegend(composition)
legend.model().setLayerSet(el_renderer.layerSet())
composition.addItem(legend)

scBar = QgsComposerScaleBar(composition)
scBar.setStyle('Numeric')
scBar.setComposerMap(componer_mapa)
scBar.applyDefaultSize()
composition.addItem(scBar) #DOESN'T WORK

i.save("RUTA\\mapa_con_todo.jpg","jpg")


texts:
Some months ago I was asking for ideas in order to do a displacement model of a lost person with GIS. This model would represent the displacement of a person by a terrain with some costs paramethers like slope, vegetation… So, I already have my model almost done. I have managed to do a final Raster where each pixel have a time value that represents the time in seconds that is needed to cross it.

Do you know if there is any tool in ArcGIS that accumulates these pixels from a starting point until a máximum of accumulated time? Let me explain:

In the model there would be a starting point, for example in .shp format. From this starting point, the tool could be able to accumulate the value of each pixel until it reaches 3600 seconds in order to simúlate the displacement of one person in one hour. The result should be a “circle” or isochrone, that represents the area until the accumulation of 3600 seconds.


texts: Copernicus Sentinel - National Mirrors and Collaborative Ground Segments

texts:
Is there an official list of all national mirrors of the collaborative ground segment members for data from the Copernicus Sentinel satellites?

So far I found the following:

Copernicus Sentinels Scientific Datahub (SciHub)
Australia National Mirror
Austria National Mirror
Finnish National Mirror
French National Mirror (PEPS)
German National Mirror (CODE-DE)
Greek National Mirror
Italian National Mirror

I am trying to piece together a list of all national or partial mirrors but could not find any official list or document about the available national mirrors. Since each mirror requires a dedicated access to the central archive I was hoping that there is some kind of documentation about the existing or planned national mirrors.


texts:
I want a cheaper solution than ESRI ArcGIS server for my online web GIS applications.
I want to:

display multiple SHP files and an aerial image,
retrieve the SHP attributes,
all in a web page inside my PHP applications.

I tried Geoserver, and got stuck in tiling the aerial images, plus the performance is not that good.
When I used ImagePyramid for example, I had the universal problem

  Could not list layers for this store, an error occurred retrieving them:
   This plugin accepts only File, URL and String pointing to a file

I tried to google it, found a lot of other users faced the same problem, no straight online tutorial describes the required steps or how to troubleshoot the mistakes.

Any help or suggestion for alternatives?

Thank


texts: Geometry Generator - default value not saved on insertion

texts:
I'm using QGIS 3.12.1 on Win10 x64.
Overview
I have a point layer representing tree positions. I have a calculated field ("rpaR") based on other attributes. I am using the geometry generator in styling to create a buffer around the point.
When I create a point, the geometry generator styling is not displayed. If I then edit that point, changing any attribute, it works. I thought this might be because the calculated field "rpaR" on which the buffer is based was not calculated before saving the new feature, but: (1) I can see the "rpaR" value is there before I save. If I then edit and change an underlying feature that changes the value of the "rpaR", the amended geometry generator styling is implemented immediately.
DETAIL
First, the drag and drop designer, showing the settings for the calculated field rpaR. The expression is as follows:
( if("dbh1" is null,0,("dbh1" / 1000)^2 * "stems1") + if("dbh2" is null,0,("dbh2" / 1000)^2 * "stems2") + if("dbh3" is null,0,("dbh3" / 1000)^2 * "stems3") + if("dbh4" is null,0,("dbh4" / 1000)^2 * "stems4") + if("dbh5" is null,0,("dbh5" / 1000)^2 * "stems5") )^0.5

I note the text tip below the default value field reads "column 'dbh1' not found". I think this is sympomatic of the same issue. After saving the feature and amending any field, 'dbh1' is found and this tip reports the correct value.

Second, the geometry generator field.

Third, plotting the feature (ref arrow) before saving (pressing OK). rpaR has been calculated correctly (blue arrow).

Four. Feature has been plotted - red arrow. But neither the simple marker, nor the geometry generator are displayed.

Five. I edit the feature and change any field. Press OK. Now the geometry generator has created the dashed blue line (red arrow). Also, at the feasture centre, the simple marker is also displayed.

Is this a bug, or some unseen error by me.
Any help appreciated.


texts: Inserting new data in PostGIS are not displayed automatically in Geoserver layer

texts:
I have build an application using PostGIS, geoserver and openlayers2.12 and I am facing the following problem:
When I add some data in my Postgis tables, geoserver doesn't get refreshed properly and the new data (lines and points) are not shown on the map.
If I zoom in the map I get some parts of my new data (parts of lines and all the points). But when I zoom out I can not see anything.
Also when I do "Layer Preview" I don't see the new geometries.

I tried to delete the cache of the layer and my browser..but no luck.
What solves the problem is if I delete and then recreate the layer.
Then everything its shown. But doing this I loose all my cache and the application becomes slow.

Any idea what It can cause this behavior?


texts: function ST_Line_Interpolate_Point with multilinestring

texts: Is there a function like ST_Line_Interpolate_Point but using a multilinestring?


texts: Showing density of overlapping features via transparency in ArcMap?

texts:
Below is a famous map produced by the Facebook team, which shows the connections between their users. Where there are multiple overlapping features, the lines show more brightly, giving a clear indication of feature density.

How can this effect be emulated in ArcMap? For example, is it possible to show individual features with a set transparency, and stack these features with the transparency values accumulating to show where multiple features coincide?

The standard > Display > Transparency option applies to the whole layer equally.

The Facebook example uses light features on a dark background, but I'm also interested in dark features on a light background.


texts: Mosaicking raster tiles with different color palettes


texts: QGIS 3.14.15: is there a build in option for storing additonal quality data from trackpoints captured with DGPS antenna

texts:
I'd like to use QGIS for collecting datapoints in field. Due to the 'GPS information module' it is no big deal to connect the R2-Antenna via Blutooth-link and read the NMEA stream.
As far as I figured out till now, QGIS- stores the XY(Z) data according to layer properties where I store my points.
I'd like to store information about acquisition quality e.g. fix-time, accuracy, satelites, etc. as attributes.
I could find some stuff in qgis.core.QgsGps...but obviously need to be pushed in the right direction.
Even the storing of the entire $G[PNA]GGA sentence would be a solution.
Can someone give me some hints or point to the right plugin to look for or to use the QGIS-internal Python classes?
Using Python with pymmea2 allows reading all that data on plain Python level but it's not my aim writing my own GUI for that purpose.


texts: What is the best method to snap features to road network?

texts: I have a road network with some of my tract centroids and facilities that fall off of it, however I would like to execute my analyses so that all features are part of the road network?  My facilities are hospitals and incidents are census tract centroids.  I am executing a closest facility analysis for drive time and distance.  I would like my output to adjust these features to the network with an output in an attribute table which shows the distance the feature was adjusted to the network.  Much like the near analysis, but actually adjusting features so that they fall within the network.


texts: Environment variables for QGIS 3 plugin development on Windows

texts:
I want to start learning QGIS plugin development, and I'm trying to get set up for Python 3, PyQt5 and QGIS 3.

I've installed QGIS 2.99, qt5-devel and setuptools via OSGeo4W. Now I get stuck on the environment config, following and modifying this quick guide to PyQGIS on Windows. I've modified it to:

@echo off
SET OSGEO4W_ROOT=C:\OSGeo4W64
call "%OSGEO4W_ROOT%"\bin\o4w_env.bat
call "%OSGEO4W_ROOT%"\apps\grass\grass-7.2.2\etc\env.bat
@echo off
path %PATH%;%OSGEO4W_ROOT%\apps\qgis-dev\bin
path %PATH%;%OSGEO4W_ROOT%\apps\grass\grass-7.2.2\lib
path %OSGEO4W_ROOT%\apps\Python36;%OSGEO4W_ROOT%\apps\Python36\Scripts;%PATH%
set PYTHONHOME=%OSGEO4W_ROOT%\apps\Python36
set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\qgis-dev\python
set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\Python36\Lib\site-packages
set QGIS_PREFIX_PATH=%OSGEO4W_ROOT%\apps\qgis-dev
cmd.exe

I've changed references to %OSGEO4W_ROOT%\apps\qgis to %OSGEO4W_ROOT%\apps\qgis-dev, and %OSGEO4W_ROOT%\apps\Python27 to %OSGEO4W_ROOT%\apps\Python36.

Since "%OSGEO4W_ROOT%"\bin\o4w_env.bat includes etc\ini\python-core.bat, which does SET PYTHONHOME=%OSGEO4W_ROOT%\apps\Python27 and PATH %OSGEO4W_ROOT%\apps\Python27\Scripts;%PATH%, I've just added Python 3.6 explicitly to the path and set PYTHONHOME to 3.6 as well, but that's ugly. Is there a better way?

It results in these env variables:

echo %PATH%
C:\OSGEO4~1\apps\Python36;C:\OSGEO4~1\apps\Python36\Scripts;C:\OSGEO4~1\apps\Python27\Scripts;C:\OSGEO4~1\bin;C:\Windows\system32;C:\Windows;C:\Windows\system32\WBem;C:\OSGEO4~1\apps\qgis-dev\bin;C:\OSGEO4~1\apps\grass\grass-7.2.2\lib
echo %PYTHONHOME%
C:\OSGEO4~1\apps\Python36
echo %PYTHONPATH%
;C:\OSGEO4~1\apps\qgis-dev\python;C:\OSGEO4~1\apps\Python36\Lib\site-packages

With the above I can run the correct python version, but I run into this:

E:\>python
Python 3.6.0 (v3.6.0:41df79263a11, Dec 23 2016, 08:06:12) [MSC v.1900 64 bit (AMD64)] on win32
Type "help", "copyright", "credits" or "license" for more information.
>>> import qgis.core
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\OSGEO4~1\apps\qgis-dev\python\qgis\__init__.py", line 71, in <module>
    from qgis.PyQt import QtCore
  File "C:\OSGEO4~1\apps\qgis-dev\python\qgis\PyQt\QtCore.py", line 1, in <module>
    from PyQt5.QtCore import *
ImportError: DLL load failed: The specified module could not be found.


texts: function, return, 'instancemethod' object has no attribute '__getitem__'

texts:
I'm still trying to write a python plugin to create an input file for an urban climate model out of three shape files. Now I stumbeled into another problem, that I think is an easy one, but I just can't see whats wrong.
I wrote a little function to get me the path and file name of the current project file (as I want the result file of my calculations to be called like the project file by default).

def projectfilepathname():
    # function to extract the file path and name of the current project file
    qgis.mainWindow().statusBar().showMessage("Getting project name and path")
    prjfi = QFileInfo(QgsProject.instance().fileName())
    path_name = []
    if prjfi:
        path_name.append(str(prjfi.absolutePath()))
            filename = str(prjfi.fileName())
            if filename.count(".") > 0:
              path_name.append(filename.split(".")[0])
            else:
              path_name.append(filename)
            path_name.append(prjfi.fileName())
    else:
            path_name.append(r'C:\Users')
            path_name.append(r'test')
    return path_name

I call the function by

currentprofile = self.projectfilepathname

but that doesn't seem to work. At least I get an error for the next line

currentprofile_dir = currentprofile[1]

telling me 'instancemethod' object has no attribute 'getitem'.
Maybe I'm wrong, but I have the impression that the function is not called, but its name is saved to currentprofile. Whats wrong?


texts:
I have a Geotiff map and I want to overlay a Shapefile onto it. Ideally I want to specify a 'nice' style for the polygons e.g. have them opaque (for now all polygons should be the same).

I want to print it so I don't want any resampling etc. going on.

How can I accomplish this?

edit: To add a little background: I generated Maps from Openstreetmap using Maperitive. I want to overlay Polygons from Shapefiles to visualize buildings in contrast to the osm buildings. I do not want to resize or resample the maps themselves as the workflow is adjusted exactly for these image sizes. I already converted coordinate systems of raster- and vector files so both are of the same type.


texts:
I've got a map composer in QGIS with a polygon layer (attributes: code, area).

I want to add a legend with the codes and the total area for each code. I used the following expression (wrong, as it happens) in the Expression string builder:
@symbol_label ||  ' ' || '(' ||  "area"  || ' ha' || ')'

And in the legend, this is what I get

How can I get what I want?


texts: Deleting all Features of a Vector Layer in pyQGIS

texts: I am developing a python plugin for QGIS. In QGIS map window, certain features of a vector layer are selected and those features are highlighted. Now I need to delete all the existing features from another vector layer without disturbing current selection in the map window. Is it possible to delete all the features of a vector layer without selecting them?


texts:
I must have a Linux library out of wack since the following image has "smears" on the colors.  Is there a better term?  On my Windows 7 RStudio it works fine.  this image is from a Ubuntu 16.04 LTS and the following are the libraries:

library(ggplot2)
library(ggmap)
library(maps)
library(mapdata)

the code that generates this image is:

ggplot(data = world2) +
  geom_polygon(aes(x = long, y = lat, fill = region, group = group), color
    = "white") +
  coord_fixed(1.3) +
  guides(fill=FALSE)


texts:
I am trying to identify lakes and other water bodies in a DEM file. The source file is located here. I am using gdaldem aspect to identify contiguous pixels with zero slope change.

gdaldem aspect -co compress=lzw infile_dem.tif outfile_aspect.tif

When I plot the output file, it looks like below (please ignore the white border):

Obviously, I am doing the exact opposite of what I need to be doing. Is there any way of inverting the operation of gdaldem aspect? The output of gdalinfo is give below:

Size is 2928, 2285
Coordinate System is:
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"]]
Origin = (81.381111111111110,22.564444444444447)
Pixel Size = (0.002083333333333,-0.002083333333333)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  81.3811111,  22.5644444) ( 81d22'52.00"E, 22d33'52.00"N)
Lower Left  (  81.3811111,  17.8040278) ( 81d22'52.00"E, 17d48'14.50"N)
Upper Right (  87.4811111,  22.5644444) ( 87d28'52.00"E, 22d33'52.00"N)
Lower Right (  87.4811111,  17.8040278) ( 87d28'52.00"E, 17d48'14.50"N)
Center      (  84.4311111,  20.1842361) ( 84d25'52.00"E, 20d11' 3.25"N)
Band 1 Block=2928x1 Type=Int16, ColorInterp=Gray
  NoData Value=-32768


texts: How to keep the DN value range (data type) after projecting MODIS using gdalwarp with python

texts:
I downloaded the MODIS images with HDF format from NASA EOSDIS Reverb website. And I am trying to project MOD09GQ from Sinusoidal to WGS84 using gdalwarp with Python. However, the DN value range (data type) in the output Geotiff after projecting MODIS using gdalwarp with python has been changed.

gdalwarp -of GTIFF -r near -t_srs EPSG:4326 -ot Float32 HDF4_EOS:EOS_GRID:"MOD09GQ.A2014279.h11v02.006.2015291113528.hdf":MODIS_Grid_2D:sur_refl_b01_1 20142792.tif
0...10...20...30...40...50...60...70...80...90...100 - done.

I used the above code. And I also projected the image using the MCTK of ENVI. The DN value range of the image which has processed by the MCTK is same with the original image. The DN value of the image which was processed by GDALwarp has been changed from [0,1.15] to [112,11422]. But, I set up the data type of the output Geotiff to float 32-ot Float32. Moreover, I have tried to set up -ot Int16 (is same with original data type of Modis). Its DN value range still was changed.

How do I keep the original DN value range (data type) after projecting MODIS using gadlwarp?


texts: Since VBA support is no longer supported in ArcGIS 10, what is being used out there for custom data entry forms.  I'd like to develop something that is very intuitive (for non-GIS users).


texts:
I have lat and long coordinates of Premier League football stadiums, which I want to map onto a NUTS-2 shapefile.

They are saved in a CSV document, which when I go to upload, not one of my coordinates are read in successfully. Instead the following error occurs:

  The following lines were not loaded into QGIS due to errors

Perhaps someone may know what my issue is?

I have attached a picture for reference.


texts:
I have a fairly extensive workspace which among other things has an output of shapefiles in several needed projections.
I have those going to the respective folders for use in different business needs.

There seem to be several problems with one of my processes (I am using the clipper to split datasets by county lines.)

First the configuration of my clipper is with the counties (polygon)[esri geodatabase sqlserver] as the clipper and the dataset (polygon)[esri goedatabase sqlserver] as the clippee.  The clippee has several polygons (20) that overlay several counties.
Several of those would "cut" a clippee polygon.

I have tried several scenarios within the parameters dialog (including single clipper, multiple clippers, and clippers first); but first don't see a change in the results and, second don't see what I want.

The results are that

some of the data is passed through the "inside" port of the
transformer, but is not clipped by the boundaries in the clipper dataset.
the attributes I want to use to fanout my dataset on are not in the data.
All of the data is passed through the "outside" port but is again not clipped, and has no attributes from the clipper dataset.

What I ultimately want is a separate dataset for each county with each respective clippee element having the attributes of the county that clipped them, and a few of their own attributes from the original data.
Edit:
I am confused.
I created a new workspace with only my clipper and one clippee.
I do get results in my inspector view but all values show clipped - no.
I have done as you said the clipper is on top and the multiple clippers is chosen.
I still can't make the 2 attributes from the clipper turn "green" but they are in the output (blank values).
HMM! Not sure what I am doing wrong at this point.

I see these values in the log.  REJECTING BELOW FEATURE:
Error - Expected an aggregate or area geometry.
REJECTING BELOW FEATURE:
Polygon feature must have at least 4 coordinates...rejecting
REJECTING BELOW FEATURE:
Error - Expected an aggregate or area geometry.
REJECTING BELOW FEATURE:
Polygon feature must have at least 4 coordinates...rejecting
Apperantly clipper doesn't take line features as clippee.  I added a bufferer (because of the confusing message) and not It is creating 3300 fan out datasets.
Now I have to figure out how to make it just create a dataset for clippers that have clippees in them.


texts: How to convert an 2D array image to table - Google Earth Engine

texts:
I'm trying to get a table from a 2d-array image. Basically, the Image (called segmentationInfo) looks like this:
The 2-D ‘LandTrendr’ annual segmentation array looks like this:

[
  [1985, 1986, 1987, 1988, 1989, 1990, 1991, 1992, 1993, 1994, 1995, 1996, 1997, ...] // Year list
  [ 811,  809,  821,  813,  836,  834,  833,  818,  826,  820,  765,  827,  775, ...] // Source list
  [ 827,  825,  823,  821,  819,  817,  814,  812,  810,  808,  806,  804,  802, ...] // Fitted list
  [   1,    0,    0,    0,    0,    0,    0,    0,    0,    0,    0,    0,    1, ...] // Is Vertex list
]

The question is: How to convert this to exportable table?
If someone can give me a hint, it will be much appreciated. Please find the commented code below.
// START INPUTS

var area = ee.Geometry.Polygon(
        [[[-64.77293019260956, -2.4770762036671767],
          [-64.77293019260956, -2.5250956006529575],
          [-64.71902852024628, -2.5250956006529575],
          [-64.71902852024628, -2.4770762036671767]]], null, false);

// define collection parameters
var startYear = 1986;
var endYear = 2020;
var startDay = '06-01';
var endDay = '12-31';
var index = 'B5';
var maskThese = ['cloud', 'shadow', 'water'];

// define landtrendr parameters
var runParams = {
  maxSegments:            6,
  spikeThreshold:         0.6,
  vertexCountOvershoot:   3,
  preventOneYearRecovery: true,
  recoveryThreshold:      0.25,
  pvalThreshold:          0.05,
  bestModelProportion:    0.75,
  minObservationsNeeded:  6
};

// define change parameters
var changeParams = {
  delta:  'loss',
  sort:   'greatest', //option: newtest, least, oldest, fastest, slowest
  year:   {checked:true, start:1986, end:2020},
  mag:    {checked:true, value:200,  operator:'>'}, //default 200
  dur:    {checked:true, value:4,    operator:'<'}, // default 4
  preval: {checked:true, value:300,  operator:'>'},
  mmu:    {checked:true, value:11},

};

// END INPUTS

// load the LandTrendr.js module
var ltgee = require('users/emaprlab/public:Modules/LandTrendr.js');

// add index to changeParams object
changeParams.index = index;

// run landtrendr
var lt = ltgee.runLT(startYear, endYear, startDay, endDay, area, index, [], runParams, maskThese);

print(lt, 'LandTrendr');

// selecting the LandTrendr band

var segmentationInfo = lt.select(['LandTrendr']); // subset the LandTrendr segmentation info
print(segmentationInfo, 'segmentInfo'); //It contains 4 rows and as many columns as there are annual observations for a given pixel through the time series.


texts: How can I sample points from images mosaic in Google Earth Engine

texts:
Rodrigo E. Principe did a great job for Mosaicking a Image Collection by Date (day) in Google Earth Engine (https://code.earthengine.google.com/20ad3c83a17ca27b28640fb922819208)
based the code, I mosaic a whole state temporal images, then I have tried to sample some points on a stacked images. While  an error  said: "The default WGS84 projection is invalid for aggregations. Specify a scale or crs & crs_transform." The tutorial says "An image which is a composite or mosaic of input images with different projections will have the default projection, which is WGS84 with 1-degree scale, generally, an aggregation at 1-degree scale is not desired or intended, so Earth Engine gives this friendly reminder to provide a complete specification for the output." Is there a better  method  except for reproject() which is not advocated in tutorial.


texts: Muhammad al-Idrisi lake Victoria map 12th Century

texts:
Wikipedia and many other sources state along the lines of, "An excellent map, known as the Muhammad al-Idrisi map from the calligrapher (sic) who developed it and dated from the 1160s, clearly depicts an accurate representation of Lake Victoria, and attributes it as the source of the Nile."

But when I look at the maps online I cannot depict Lake Victoria on his maps. Can anyone show me where Muhammad al-Idrisi depicts Lake Victoria on his 12th Century Maps?

I got the images below from http://www.alrahalah.com/2015/01/al-idrisi-12th-century-ce-5th-century-ah-geographer-and-cartographer/

On the top map I cannot depict sub-Saharan Africa as this map appears to be mostly North Africa, Asia, and Europe and just misses Lake Victoria.

On the lower map (it says north is down) I suspect the yellow circle is Lake Victoria but this is just a guess.


texts: Calculating land use area in specific area using ArcGIS Desktop?

texts:
I'm digitizing land use and I want to be able to calculate the land use percentage of area in three different areas (circles). I.e. I want to know the land use percentage in each of the three different circles.

How could I do that effectively?

One problem might be that the circles overlap, but that is intentional so that's part of the question. The image shows the three different circles that are my areas of interest and I will start to digitize every land use category in them.

I will then divide the area of each circle with the land use categories area, but is there some way of delimitate it by each circle?


texts: Create right to left legend in GeoServer

texts: Can I create right to left legend like this?


texts:
I'd like to change the symbology of my vector layer and can't figure it out how to read the documentation. I would like to have an polygone with an orange outline, a dashed fill style and an orange stroke color. However that would be nice to have an answer that explains how to access most of the symbology change! (if possible ;-) )

Here is what I managed so far:

import os
import processing
from qgis.core import *
from qgis.gui import *
import qgis.utils

# initialize project
project = QgsProject.instance()
project.setCrs(QgsCoordinateReferenceSystem(21781))
project.setFileName("C:\\test_symbology.qgs")

my_shp = "Y:\\shp_collection\\A300.shp"
my_shp_name = "A300"

# load polygon
shplayer = iface.addVectorLayer(my_shp , my_shp_name , "ogr")
if not shplayer.isValid():
    print(".shp failed to load!")
shplayer.setCrs( QgsCoordinateReferenceSystem(21781,QgsCoordinateReferenceSystem.EpsgCrsId))

layer = iface.activeLayer()

# make my polygone orange
single_symbol_renderer = layer.renderer()
symbol = single_symbol_renderer.symbol()
symbol.setColor(QColor.fromRgb(255,128,0))

How can I set the other symbology properties?


texts: QGIS transparency layer export to pdf not working

texts: I have in a drawing polygon layer with a degree of transparency setting. However, after exporting to PDF a layer is not transparent anymore. Where is the problem? I use QGIS 3.16.


texts: Displaying standalone attribute table as a feature layer

texts:
Using ArcGIS Pro, I have a feature class, which I converted the selected attributes into a standalone attribute table.
Having cleared my selection, I want to then know how to convert that table into a feature layer?


texts: Using QGIS with National Hydrography Dataset (NHD)?

texts:
I'm a trying to learn QGIS. I am trying to download and open with QGIS 2.0.1 state-specific data from the USGS NHD Viewer, and i think i have at least managed to do the downloading correctly. For example, in efforts to gather data for the state of Florida, i ended up finding, downloading & extracting a compressed folder called NHDH_FL_931v210. extracted from it were 244 items of various file-types (ATX, GDBTABLE, GDBINDEXES, FREELIST, SPX, GDBTABLX) which altogether i assume is an ESRI.gdb and would open in ARC like i am used to. I have QGIS 2.0.1 now and am as green as it gets... i have been trying for a few days to open this data but am not sure really where to begin. When i try to open the files, nothing shows up in the QGIS opener because, i imagine, the files are not recognized by the program.

Do i have the right file downloaded from NHD to display USGS rivers and lakes with QGIS?

Do these files need to be converted? if so, how?

Do i have to connect to the source folder within QGIS like i did in ARC? if so, how?

Do i need to create a GDB in QGIS or is this thing ive downloaded gonna work?

In my research so far i have seen posts from a few years ago (about Wroclaw, an older version of QGIS) that discussed certain steps to take during the installation process in order for ESRI .gdb files to be readable, and if this needs to be done could someone explain this process for 2.0.1 in steps?

Is the new version coming out going to be able to handle this problem i have easier or better?


texts: QGIS moving the layer/all features using Vector Affine Transformation

texts:
I'm using QGIS to analyse some flight tracks which were produced using a different software. Now when they are on QGIS the reference point of everything is on 0,0 which is the aerodrome reference point.

But in real world this aerodrome reference point is a value in Easting and Northing. i.e. Easting=286499.025 Northing=6230965.204, so the everything else (runway ends, tracks points should be exactly as they are in the real world)

I did some research and found the Vector Affine Transformation plugin. Now using this how can I move everything to the reference point of Easting=286499.025 Northing=6230965.204
Easting = X and Northing = Y.

I'm having trouble understanding the plugin's matrix.


texts: Select buildings in a html map and get osm

texts:
I have my own database of buildings with the architects and some metadata. In the past I tied them manually to the ID of openstreetmap.

But now I need a map where I can click on a shape of building or area and getting the id into Javascript for processing more. (its functionality is near to the edit mode of openstreetmap). I already read about leaflet and openlayers but as I recognized there is always the query of lat/lon and not getting the real id from the object underneath the cursors. That is not so nice because I have to very often select a specific shape or id.

Is there anything which might help me or already exists to build something like that? I can not imagine that it might not possible!


texts: Predict user's future location - location awareness mobile computing

texts:
This question is a part of an online course (Mobility Data Management) that I am currently auditing. A part of the project is to implement a system that can predict user's future location. This kind of application has recently got some interests from database and mobile application developers.

This is the problem statement:
Assume that we are given a road network G (V,E), the recent motion history (trajectory) T, and current instant velocity v of a user, a temporal duration dt, and (optional) a historical trajectory database D of other users’ movements. Our goal is to estimate user’s location after dt.

I found that the following paper has a good description about what a road network and the trajectory data is: http://www.sis.pitt.edu/~xliu/papers/ceus.pdf

Does anyone has an idea on how to formulate a simple solution of this mobility awareness problem (predicting user's future location), and how to start the implementation (how you formulate the problem in an OOP way, such as: what classes are needed to represent the problem correctly).


texts: Locating unknown destination with Trimble?

texts:
I would like to know if it is possible to input coordinates in Trimble Geo7x and be able to locate an unknown destination.

If so how do I input the coordinates and also to navigate to the input coordinate points destination?


texts: Copy all attributes of a polyline in a layer to a different polyline in the same layer

texts:
How do I copy all the attributes of a polyline in the attribute table represented by a row to a different polyline in a different row in the attribute table?

In the above figure, all the attributes of ID 1301 has to be copied and pasted to ID 1388.


texts: QGIS Cloud - partial display of polygons only

texts: I'm attempting to upload a layer of polygons (land zoning). The output in QGIS cloud has significant geographical gaps in some areas. The source data is not sound however. I'm using a free account but are well under the MB limit. Any thoughts or suggestions?


texts: How to perform unsupervised classification of LANDSAT time series using Grass?

texts:
I'm using GRASS to perform an unsupervised classification on a time series of LANDSAT images (monthly temporal res, over 30 years), in order to compare land cover change in a particular area.

I'm concerned that if I run i.cluster for each map in the time series the signatures used for the classes may be slightly different for each map, affecting the validity of the comparison. Would this be the case? If so, is there a way to ensure that the class signatures used on each map in the time series are identical throughout?

My initial idea was to run i.cluster on a particular map where all classes are well represented, then use the sigfile generated from this as the input sigfile to i.maxlik for all maps.


texts:
I need save selected attributes of layer as CSV with separator. This is my code:
def initAlgorithm(self, config=None):
        self.addParameter(QgsProcessingParameterVectorLayer('Input', 'Input', types=[QgsProcessing.TypeFile], defaultValue=None))
        self.addParameter(QgsProcessingParameterFeatureSink('Output', 'Output', type=QgsProcessing.TypeFile, createByDefault=True, defaultValue= ''))

def processAlgorithm(self, parameters, context, model_feedback):
        feedback = QgsProcessingMultiStepFeedback(1, model_feedback)
        results = {}
        outputs = {}

        alg_params = {
            'FIELDS_MAPPING': [{'expression': 'FID',
                                'length': 10,
                                'name': 'FID',
                                'precision': 0,
                                'type': 10},
                              ],

            'INPUT': parameters['Input'],
            'OUTPUT':  parameters['Output']

           }

        results['Result'] = processing.run('qgis:refactorfields', alg_params, context=context, feedback=feedback, is_child_algorithm=True)
        outputs['Outputs'] = results['Result']['OUTPUT']
        crs = QgsCoordinateReferenceSystem("EPSG:5514")
        uri = outputs['Outputs']
        layer = QgsVectorLayer(uri, 'name', 'ogr' )

        QgsVectorFileWriter.writeAsVectorFormat(layer, 'C:/Users/xxxx/OneDrive/xxxx', "utf-8", crs, "CSV", layerOptions=["GEOMETRY=AS_YX", "SEPARATOR=SEMICOLON"])

        return outputs

Problem is that the CSV is exported, but with comma not separator.
In layerOptions I have I want separator as delimiter, but it isn´t working. I can´t find mistake, maybe in Output?


texts: How to compile pgRouting 1.05 binaries for Windows?

texts: Anybody with a compiled 1.05 binary for (Windows / PostgreSQL 8.4 / PostGIS 1.5.3) or a detailed description on how to compile?


texts: I've been trying to get 2 raster layers to match up, so the cells are in the same locations for both layers (they have the same resolution already). But I can't seem to figure it out, already tried the SAGA and GRASS resampling tools with no result.


texts: Changing ArcGIS Viewer for Flex page title?

texts:
When working with the ArcGIS Viewer for Flex application is it possible to change to title page tab to something else e.g. sometitlepage.com or will the Internet police come for me if I try to alter this setting?

was going to post an image but cannot because I'm 'new'.

link here - http://i1100.photobucket.com/albums/g417/jdsegars/junkphotos/pagetitletoremore.png?t=1311342928


texts: Applying async/await feature to ArcGIS API for JavaScript code?

texts:
How do I apply the async/await feature to this code?
I need to wait for the reltionshipQuery and forEach to finish before I return contentSpecies.

        contentSpecies = function(value, key, data) {
          var contentSpecies = [];

            objectID = data.OBJECTID;

            var queryTask = new QueryTask({
                url: "https://webmaps.com/MapServer/1"
            });

            var relationQuery = new RelationshipQuery({
                objectIds: [objectID],
                outFields: ["Species", "Notes", "Status"],
                returnGeometry: true,
                relationshipId: 0
            });

            queryTask.executeRelationshipQuery(relationQuery)
            .then(function(rslts) {
                var features = rslts[objectID].features;
                features.forEach(function(ftr) {
                    var t = ftr.attributes;
                    var species = t.Species;
                    contentSpecies += "Species: " + species + "
";
                    var notes = t.Notes;
                    contentSpecies += "Notes: " + notes + "
";
                    var status = t.Status;
                    contentSpecies += "Status: " + status + "
";

                });

                //return contentSpecies;
            })
            return contentSpecies
        }


texts: Get attributes with Latitude and Longitude in shapefile as a CSV file

texts: I have some shapefiles of geodata. I need to add that into a database and access attributes with Latitude and Longitude in shapefile. So how can I solve using any type of programming language and database?


texts: Creating one-sided buffers or parallel lines in PostGIS?

texts:
I'm looking for a way to create so-called one-sided buffers or parallel lines in PostGIS. Example: 1. buffer, 2. one-side, 3. both-sides

I've found some discussion on the mailing list from 2009 and information that it has been implemented in GEOS, but nothing about the current PostGIS status.

In case the function is not implemented yet, do you know any workarounds? Is it possible to cut off one side of a normal buffer?


texts: Performance of Picture Marker Vs Character Marker (font)?

texts:
We have a point data type to symbolise.  The feature classes can often contain more than 6000 points up to 20000.  There are about 30 different symbol catergories to these points.

My question is what's more performant - creating EMF pictures and symbolising using picture markers or create a font and use character marker symbolisation?

We would be using this dataset in both ArcGIS Desktop and ArcGIS Server (silverlight API).


texts: What does either of the appreviations b and u refer to

texts:
In
https://wiki.openstreetmap.org/wiki/Overpass_API/Overpass_API_by_Example
there is an example as shown below in the code section. I would like to know what is
b
u
and b.u?

Code:
[timeout:600];
rel[admin_level=2][wikidata]->.a;
foreach .a -> .b (
  rel.a(if:t[wikidata] == b.u(t[wikidata]))->.d;
  rel.d(if: id() == d.min(id()) && d.count(relations) > 1 );
  convert rel ::id = _.u(id()),
             wikidata = _.u(t[wikidata]),
             all_relation_ids = d.set( id() );
  (._ ; .result; ) -> .result;
);

.result out;


texts: ogr2ogr problem converting sosi fil to kml (or anything else)

texts:
I am trying to convert a SOSI file to KML, but I am getting the following error.

Error message

ERROR 4: Unrecognized geometry of type 5.
ERROR 4: Unrecognized geometry of type 20.
ERROR 1: Feature 168 referenced by 1474, but it was not initialized.
Aborted

I have tried different output formats but that do not seam to effect the error.

Are there any settings I can pass ogr2ogr to ignore or get around the errors?


texts: Intersection of all layers of different feature classes using ArcMap

texts:
I have 20 different layers where a few of the layers are lines ( rivers) some are points and a few are polygons. I need to find all the single polygons where all 20 layers intersect as describe:

within distance 10 meters of layer #1
within distance 20 meters of layer #2.
... ...
within distance 100 meters of layer #20.

I have ArcMap 10.7 with a license that allows doing multiple-layer intersections. I use two methods:

Buffer (with dissolved) each layer by the distance and feed all in a intersect in my model.
I buffer the first one and I use "select by distance" within the distance specified by the next layer and keep doing this to layer 20!

But none of the methods returns a correct result.
Is this because the layers are from different feature classes?


texts: GEE Landsat SR year composite // mask cloud / shadow out w/ other method than quality pixel

texts:
Has anyone found a good solution to mask out cloud/shadow cover not using quality assessesment pixels or combining it with something else?

I am trying to generate cloud/shadow free composite for a large area [code with 1 point to make it easier] using Landsat 8 SR. I wrote this code based on qa_pixels. It works, but there is still cloud pixels which are described as clear land in the quality assessment [easier to see them when plotting an image w/ the max NDVI value per pixel].
Here is the code: https://code.earthengine.google.com/8b42d2a7c3c31ee3606a3994373ea59a

var UXRB = /* color: #0b4a8b */ee.Geometry.Point([-52.207, -14.349])
var colL8_SR = ee.ImageCollection("LANDSAT/LC08/C01/T1_SR")
    .filterBounds(UXRB)
    .filterDate('2015-01-01', '2015-12-31')
    .mosaic();
print("colL8_SR", colL8_SR);
Map.addLayer(colL8_SR, {bands: ['B7', 'B5', 'B4'], min:0, max:10000}, 'colL8_SR');

var colL8_SRnc = ee.ImageCollection("LANDSAT/LC08/C01/T1_SR")
    .filterBounds(UXRB)
    .filterDate('2015-01-01', '2015-12-31')
    .filterMetadata('CLOUD_COVER', 'equals', 0)
    .mosaic();
print("colL8_SRnc", colL8_SRnc);
Map.addLayer(colL8_SRnc, {bands: ['B7', 'B5', 'B4'], min:0, max:10000}, 'colL8_SRnc');

////////////////////////////////////////////////
// I can get cloud-shadow free images w/previous section - static scene;
// However, instead of static scene, I need a composed of the maximum number
// of cloud-shadow free images in a year as possible to calculate temporal indexes
////////////////////////////////////////////////

// Add additional bands to Landsat 8 collection
print('Creates add NDVI function');
var addBANDS = function(image) {
    var ndvi = image.normalizedDifference(['B5', 'B4']).rename('NDVI');
    var evi = image.expression('2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))', {
        'NIR': image.select('B5'),
        'RED': image.select('B4'),
        'BLUE': image.select('B2')}).rename('EVI');
    var time = image.metadata ('system:time_start').rename('Time');
    return image.addBands(ndvi).addBands(evi).addBands(time);
};

// Call L8 REFLEC collection T1 & filter by the bounds (UXRB), date (2015), pixel & add NDVI & Mosaic
print('L8 SR, Spatial/Temporal filter, greenest');
var colL8_SR15 = ee.ImageCollection("LANDSAT/LC08/C01/T1_SR")
    .filterBounds(UXRB)
    .filterDate('2015-01-01', '2015-12-31')
    .map(addBANDS);
print("colL8_SR15", colL8_SR15);

//A function to mask out cloudy pixels based on L8 SR qa_pixel
var cloud_shadows1 = function(image) {
    // Select the QA band.
    var QA = image.select(['pixel_qa']);
    // Get the values for clear & water & fill
    var mask = QA.eq(322) //Clear land: 322
              .or(QA.eq(324)) //Clear water: 324
              .or(QA.eq(336)) //Snow-Ice: 336
              .or(QA.eq(1)) //Fill: 1
    // Get updated image
    return image.updateMask(mask);
};

var colL8_SR15ncs1 = colL8_SR15.map(cloud_shadows1);
Map.addLayer(colL8_SR15ncs1, {bands: ['B7', 'B5', 'B4'], min:0, max:10000}, 'colL8_SR15ncs1');
print("colL8_SR15ncs1", colL8_SR15ncs1);

var L8_SR15ncs1m = colL8_SR15ncs1.mean();
Map.addLayer(L8_SR15ncs1m, {bands: ['B7', 'B5', 'B4'], min:0, max:10000}, 'L8_SR15ncs1mean');
print("L8_SR15ncs1m", L8_SR15ncs1m);

var L8_SR15gr11 = colL8_SR15ncs1.qualityMosaic('EVI');
Map.addLayer(L8_SR15gr11, {bands: ['B7', 'B5', 'B4'], min:0, max:10000}, 'L8_SR15gEVI1');
print("L8_SR15gr11", L8_SR15gr11);

var L8_SR15gr21 = colL8_SR15ncs1.qualityMosaic('NDVI');
Map.addLayer(L8_SR15gr21, {bands: ['B7', 'B5', 'B4'], min:0, max:10000}, 'L8_SR15gNDVI1');
print("L8_SR15gr21", L8_SR15gr21);
// A lot of bad pixels 'cause of cloud/shadow which are described in qa as clear pixels

I also run the code described here [Cloud mask for Landsat8 on Google Earth Engine ], but retrieves the same result - some bad pixels due to cloud-shadow which are identified as clear land in qa.

Has anyone found a better solution for this issue? Or provide any advise on the way to go?


texts: CartoDB: Can I use some "if..then" structure or "loop" structure in popup custom HTML?

texts:
In CartoDB (https://carto.com/)'s popup CustomHTML, can I use some logic like "If..then", "loop" or so?

I want to use it for usecase something like below:

Avoiding representation if attributes are null.
Showing all images if image url list is given.

Would you give me the information for such function in CustomHTML?


texts: How to convert Lat/Long to UTM using Proj(4j) similar to JScience UTM.latLongToUTM?

texts:
I'd like to use Proj(4j) to convert between UTM and Lat/Long coordinates.

I used to use JScience for UTM <-> LatLong conversions like this:

UTM.latLongToUtm(latlong, ReferenceEllipsoid.WGS84);
UTM.utmToLatLong(coordinate, ReferenceEllipsoid.WGS84);

But I would like to accomplish this using Proj4j. So I'm currently having set up two CoordinateReferenceSystems:

CRS_FACTORY.createFromParameters("EPSG:4326", "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs");
CRS_FACTORY.createFromParameters("?", "+proj=utm +? +zone=?);

Where I don't know how to fill up the ? gaps.

Why do I have to specify a zone-parameter? Is there a way to accomplish the same above in JScience using proj4j?


texts: I am trying to do principle component analysis, my file size is 1,5GB. It has been not responding past 30 mins,  but I decided not to close the program, is it running or I am just waiting for nothing?


texts:
I often come across code examples when searching for something like
GeoTiff size - what to do?.

Sometimes, the presented code examples seem to do exactly what I was searching for. There is just one big drawback which is that I just don't know which type of console is to be used.

In the example, it says:

  You may try out something like the following, using the GDAL
  utilities:

gdal_translate -co "TILED=YES" -co "COMPRESS=JPEG" -co "PHOTOMETRIC=YCBCR" input.tiff output.tiff

gdaladdo -r average --config COMPRESS_OVERVIEW JPEG --config PHOTOMETRIC_OVERVIEW YCBCR output

So I googled GDAL utilities, but that didn't help. I am used to R, but apart from that, I am quite unfamiliar with coding.

Which program do I use with such code examples?


texts: How to set order of operations for Field Calculation in Network Analyst Evaluator Dialogue?

texts:
I desire to have a ratio calculated each time I create a new route.

To get this ratio, I created and populated a field that multiplied the length of each segment by a designated number from 1-5, where 1 was most desirable, we'll call this [FieldValue]. Theoretically, when divided by the length, this would produce a value that, the closer to 1 would indicate a more desirable route.

In the field evaluator value dialogue box, I was able to divide my calculated value by the length of the segment ([FieldValue]/[Length]). I was hoping that they would sum each field first and then divide the two values at the end.

However, when I ran the route analysis, the values were dividing each row, and summing all of the individually divided values at the end.

I would like each field, [FieldValue] and [Length] to sum first, then divide the two summed fields to produce the ratio. It would appear that I can do this with a VB script or Prelogic within the network analyst evaluator value field, but I lack the programming skills to do so.

What would be the best route to take in order to have them calculate in my desired order?


texts: PyQGIS - How to import WKT string from a CSV to a polygon shapefile?

texts:
In QGIS, I can add WKT data stored in a CSV by clicking the 'Add Delimited Text Layer' tool and setting it to 'WKT'. This works fine:

But how can I do it as a python script?

Ideally I'd like a script which-
1. Imports the CSV into QGIS as a polygon (drawn from the WKT field)
2. Exports the new (temporary?) polygon layer to a shapefile on my disk

I have minimal experience writing scripts in QGIS. Using help guides, I've written a script to import a CSV to a point layer, with X and Y fields defined. (using this link)

But obviously a WKT is different as you only need to define the WKT field.

Here's what my source CSV looks like:


texts:
I have a folder on my local drive with 1,500+ JSON files (lines) that I need converted to shapefiles. (JSON to feature tool in ArcGIS works great on a single JSON file)

I've tried getting a batch model in Python to work, but this code seems to give me an ERROR 000206 after 31 JSON files...

Any ideas?

(I'm using ArcGIS 10.2 and PyScripter)

import arcpy
import os

dirOfJSON = 'c:/temp/jsons'
outputDIrectory = 'c:/temp/jsonSHPs'

arcpy.env.workspace = dirOfJSON

for f in arcpy.ListFiles('*.json'):
    j = os.path.join(dirOfJSON, f)
    outputShp = os.path.join(outputDIrectory, os.path.splitext(f)[0] + '.shp')
    print("Processing : {}".format(outputShp))
    arcpy.JSONToFeatures_conversion(j, outputShp)


texts: what is the difference between the Geographic coordinate, transformed coordinate and screen coordinate

texts:
What is the difference between the Geographic coordinate, transformed coordinate(projected) ,cartesian coordinate, and screen coordinate? Is there a live example?

And what is the association between them? For example, on geographic coordinate will have one or more transformed coordinate?

The following is a .prj file from this post:

PROJCS["Lambert Azimuthal Equal-Area",
       GEOGCS["GCS_WGS_1984",
              DATUM["D_WGS_1984", SPHEROID["WGS_1984", 6378137, 298.257223563]],
              PRIMEM["Greenwich", 0],
              UNIT["Degree", 0.017453292519943295]],
       PROJECTION["Lambert_Azimuthal_Equal_Area"],
       PARAMETER["false_easting", 0],
       PARAMETER["false_northing", 0],
       PARAMETER["latitude_of_origin", 0],
       PARAMETER["central_meridian", 20],
       PARAMETER["xy_plane_rotation", 0],
       UNIT["Meter", 1]]

Then what confused me is that Lambert Azimuthal Equal-Area is the name of the projection coordinate system, and it specified the geographic coordinate system here(GCS_WGS_1984),why? Does it mean that this projection is only fit for the GCS_WGS_1984 geographic coordinate? As far as I understand, the projection have nothing to do with the geographics coordinate.

For example,given a geographics coordinate system GCS_WGS_1984.

Now I can use different type of projection(like UTM, and Lambert Azimuthal Equal-Area in the post) to transfrom the geographics coordinate to projected coordinate. So the projection is not bounded with the geographics coordinate system. Isn't it?


texts: Skipping 'file already exists' error in ArcPy?

texts:
Doing an intersection and some files already exist.Is it possible to tackle this error specifically and not with try except in general?

try:
    for i in cur:
        Intersect_analysis (in_features, out_feature_class)
except arcpy.ExecuteError as e:
    print(e)
    pass

  file already exists.

This may pass any execute error instead of just this with a file being already there.

OR

Since the above error stopped the whole program and now has to restart,
another plan would be to search in the gdb which feature layers aren't there and proceed to do the intersect only to the new ones that haven't been processed due to the program stop.


texts: Batch field-calculator with a decimal result

texts:
Is it possible to carry out a batch field calculator process in QGIS across multiple layers resulting in a decimal number? I have manage to do it and get an integer but not a decimal. There are only options for an integer, float, date and string in the field type drop down.

I am trying to calculate the area of roughly 3900 layers so don't want to have to calculate it individually.


texts: How to reproject all featureclasses in a database to another database?

texts:
I am very new to arcpy and have ArcGIS 10.0

I would like to make a script where I would only have to change two variables in a script (database_input and database_output) so that I can reproject all vector featureclasses within the database and save them into another database(which exists already).

Would anyone be able to put me on the right track?

My code so far. It is basically an exported modelbuilder model using the batch reproject tool.

# Import arcpy module
import arcpy

# Local variables:
shp1 = "Path_to_Geodatabase_old\\old.gdb\\shp1"
shp2 = "Path_to_Geodatabase_new\\new.gdb\\shp2"
test_etrs_gdb__2_ = "G:\\Tifosy_Serverinstall\\test_etrs.gdb"

# Process: Batch-Projizierung
arcpy.BatchProject_management("Path_to_Geodatabase_old\\old.gdb\\shp1;Path_to_Geodatabase_new\\new.gdb\\shp2", test_etrs_gdb__2_, "PROJCS['ETRS_1989_UTM_Zone_N32',GEOGCS['GCS_ETRS_1989',DATUM['D_ETRS_1989',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Transverse_Mercator'],PARAMETER['False_Easting',32500000.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',9.0],PARAMETER['Scale_Factor',0.9996],PARAMETER['Latitude_Of_Origin',0.0],UNIT['Meter',1.0]]", "", "DHDN_ETRS_1989_8_NTv2")


texts: How to use python to set a variable to a value from a vector attribute table for GRASS

texts:
I'm adding to a python script for GRASS v.transects.py. I want to set a variable to equal a value from a field in the vector's attribute table using a for loop to cycle through each row. This is one of my first GRASS scripts, I'm much more comfortable with arc. For arcpy I would use a search cursor and getValue:

vector = "river.shp"
lines = arcpy.SearchCursor(vector)
for line in lines:
     attributeValue = str(line.getValue("FIELD"))

The existing script has the right for loop already (for line in vector), its just a matter of figuring out the last line, how to have it read the field in the attribute table of the line segment it is currently cycling through.


texts: QGIS - Can't see rows data but I can see the total number of rows in the attribute table

texts:
I tried to download some shapefiles data from here and the layer is shown correctly, however when I open the attribute table I can see just the number of rows but not the rows themselves and when I click organize columns no fields are detected, any tip on how to see the data in the attribute table?
The version of QGIS I am using is from the arch Linux repositories and my OS is Manjaro Linux if this can help.


texts: Open Source Alternatives to ERDAS IMAGINE ATCOR?

texts:
I have been using ERDAS IMAGINE's ATCOR extension for ATmospheric CORrection of ASTER imagery.

What opensource options exist for atmospheric correction of ASTER Level 1B products?

I'm primarily interested in converting registered radiance at the sensor to true reflectance values.


texts: Implementation of multiple-point statistics in GIS

texts: Is there any implementation of multiple-point statistical algorithms in GIS? at least for spatial interpolation. code for instance.


texts:
How do I detect a click on a layer and find out information about the layer that was clicked?

Currently I can do the usual binding a popup and adding a click handler on the layer:

L.geoJson(data.streets, {
    onEachFeature: function(feature, featureLayer) {

        featureLayer.bindPopup(feature.properties.name);

        featureLayer.on('click', function(e) {

            console.log('Layer clicked!', e);
        });
    ...

However the click handler on the layer doesn't give me any information about the layer, just the following:

containerPoint: o.Point
latlng: o.LatLng {
    lat: 50.804443085898185,
    lng: -1.089920997619629
}
__proto__: Object
layerPoint: o.Point
originalEvent: MouseEvent
target: e
type: "click"
__proto__: Object

How can I get a reference to the layer it's clicking?


texts: Converting vector layer to raster using QGIS?

texts:
In qgis in the raster menu there is an item "convert vector to raster".

This behaves strangely - I do not understand how it is supposed to work.

I can select a vector file, an attribute and an output raster. But, for some reason the output file is only accepted if it exists and even if I select an existing raster nothing seems to happen.

Also I do not understand how the program knows what kind of raster to create, what cell-size, projection etc.

How can I convert vector to raster?


texts: Do I need to define extent when I define a new projection?

texts:
GGRS87 (EPSG:2100) is the projection system for Greece. But there is also another projection system which caters for a complex of islands (Kastellorizo, etc) which are located east form Greece (near Turkey). I have created this new projection with the following parameters
(+proj=tmerc +lat_0=0 +lon_0=27 +k=0.9996 +x_0=500000 +y_0=0 +ellps=GRS80 +towgs84=-199.87,74.79,246.62,0,0,0,0 +units=m +no_defs).

I have digitized some polygons on these islands with the help of some orthophotos. However, these orthophotos were not on the correct location. They are inside the extent of the area for GGRS87 and not near Turkey. I saved this shapefile to the new projection. The map of the extent of the CRS for the shapefile shows that the area of the projection system is located on the west side of Africa.
 On the contary, the map of the extent of the CRS for the whole project shows that the area of the projection system is near Turkey (which is the correct location).

How is this possible? Has the fact that I did not determine the extent when I created the new projection played a role in it?
In addition you can download the polygon shapefile here.


texts: How to use r.patch with two maps having different colortables?

texts: I am using r.patch to combine two raster maps in order to fill in the null values of the other. Is there a way to use r.patch but still maintaining the colortables of the separate raster maps? I do not want to change the colortable of input1 raster map A and input2 raster map B.


texts: Ordering data through SQL Select query in MapInfo

texts:
I have a All India Village tab file(6Lac records) and an excel which has information on most of those villages(4Lac records).

I want to select 100 villages from a particular state whose area is the highest.

I used the below shown query:

But in this query the Area is sorted after selecting first 100 records. I want to select 100 records which have the highest area.

Is there a way the villages are ordered first as per the area and then the selection is carried out?


texts:
I have PRISM data which is in the following format.

Are there any ways (software and methods) to extract the DEM or elevation information from these data files?


texts: Create FileGDB using Geoprocessor on ArcGIS 10.3 ArcObject

texts:
Exact same code working fine with 10.0. But, 10.3 usual Esri "Error HRESULT E_FAIL has been returned from a call to a COM component".

                Geoprocessor geoprocessor = new Geoprocessor();

                //create fGDB
                CreateFileGDB createFileGDB = new CreateFileGDB();
                createFileGDB.out_folder_path = path;
                createFileGDB.out_name = name;
                createFileGDB.out_version = "10.3";


texts: AJAX WMS request from Linux to Windows ArcGIS Server

texts:
I am trying to figure out a CORS header or cross domain issue but to no avail.  We have an OpenLayers application hosted on a Linux server.
However, it calls a WMS service that is hosted on a Windows Server created through ArcGIS Server.  The request goes out via an AJAX call but never comes back.

The outgoing request looks like this:

http://.../WMSServer?TRANSPARENT=true&LAYERS=6&QUERY_LAYERS=6&STYLES=,&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetFeatureInfo&BBOX=-84.82432751264871,34.890657004827005,-79.56736950410722,36.72589197674798&FEATURE_COUNT=100&HEIGHT=824&WIDTH=1914&FORMAT=image/png&INFO_FORMAT=application/vnd.ogc.gml&SRS=EPSG:4326&X=1088&Y=271

I have verified this is something to do with a CORS header or cross domain I think because I took the same OpenLayers site and put in on the Windows server and it works.

I have tried a couple of thing like putting a crossdomain.xml in a couple of places on the Windows server but had no luck.  Any help would be greatly appreciated.


texts: Using nested If statement in Python parser of ArcGIS Field Calculator?

texts:
I'm trying to use field calculator to assign values to a new field I created. I wan't to assign the value from the first field unless that field is empty, if it is empty use field 2, if that is empty use field 3. Here is what I have thus far. I feel it is way off but not sure where to go from here.

New_Size: My new field

wmDeviceEq: First field to check for new field name

wmDevice_4: Second ""

wmDevice_5: Third ""


texts:
I am using QGis 2.18 on Windows 10. I want to execute an external executable file that is writen in python3 and show the output in the dialog provided by QgsRunProcess. The code I am using is:

from qgis.core import QgsRunProcess
my_cmdLine = r'"C:\Python36\Scripts\myExecutable.exe" -h'
my_process = QgsRunProcess.create(my_cmdLine, True)

I am not able to see any output, nor any error message. The executable works if executed on a cmd console, showing the help text of the executable in this case with the -h option.

How can I see the output that the program using QgsRunProcess or if it is the case see any possible error.

Moreover, once the output will be vissible, can I get it as is being generated by the executable?


texts: What does "CLEAR" mean on a land survey map?

texts:
I paid for a land survey of my house lot, and there was one term on the survey map that I wasn't 100% sure about: "CLEAR".  It's used on the map to refer to a short distance (e.g. 0.1' CLEAR or 0.2' CLEAR).

Is CLEAR short for "clearance"?  If yes, does that mean that my house is 0.1 feet and 0.2 feet, respectively, away from the property line?  Or does it more likely mean that my house is on the property line and those small distances are measurement errors?

The house is in California if this term is locale-specific.


texts: PyQGIS: How to send signal to call unsetMapTool()

texts:
I've created a map tool to add some functionality to adding point features in QGIS 2.0. I'd like to unset this map tool if editing has manually been terminated by user. I've attempted as shown below, but my map tool remains active even after editing is terminated.

if not self.myLayer.isEditable():
    self.canvas.unsetMapTool(self.clickTool)


texts: OpenLayers: How to add permalink long/lat anchor to links outside of map on page

texts:
I would like to extend the funtionality of OpenLayers that appends, on-the-fly as the user zooms and pans the map, the longitude/latitude parameters as a anchor to the URL in the browser's address bar, e.g.

http://mywebsite.com/#zoom=12&lat=50.90306&lon=-1.38329&layers=BT

I would like certain other links on the same page, outside of the map, to have the the anchor #zoom=12&lat=50.90306&lon=-1.38329&layers=BT appended to them.

The reason is so that when these links are clicked on, the map will be loaded at the same pan and zoom position. These links contain other 'arguments' that determine what data appears in the overlay. The site runs on the Drupal 7 content management system, with the OpenLayers module and these other arguments control a Drupal View to control what is displayed in the overlay.

So my approach is to implement some self-contained javascript/jquery that doesn't change openlayers javascript (to avoid overheads of burdening the authors with testing, commits to the source etc) but instead is loaded after openlayers javascript, this self-contained code would at run time 'register'/'bind' as necessary to be able to listen to pan and zoom events and obtain the longitude and latitude anchor data and add these to the certain other links. I would identify those certain other links with a class attribute in the html <a /> element and use the JQuery DOM functions to search for the class and prepend the anchor.

What I need is to know how the OpenLayers code is updating the location field/address bar with the # anchor in order to figure out how to write some separate code to listen for the pan zoom events.

Please share your thoughts. Thanks.


texts: Creating Thiessen Polygons (Analysis) by defining external feature boundary?

texts:
Although this issue has been addressed in this post some months ago: Thiessen Polygon delineation within a Feature, the appropriate answer has not been given. So, new ideas are welcome.

How can we use the Create Thiessen Polygons (Analysis) tool available with ArcInfo license by defining an external feature boundary (e.g. study area boundaries defined by a polygon) for the created thiessen polygons?

As I have seen this is not possible via the Environments...options in that tool.

Is there any available code in VBA for doing this? or How can I run this tool within a VBA code?

The comments indicate that ArcGIS 9.2 and 9.3 is being used, the tag indicates version 10, the first post is for that version not previous ones.


texts: Calculating nearest distance between two Linestrings using R

texts:
I have two simulated linestrings, one that represents a satellite ground-track and the other representing a linear surface feature (e.g., a road or river). I am trying to calculate the minimum distance the two features.
From searching stack exchange, I have only found minimum distance methods for polygons to points, points to points, lines to points, etc. but cannot seem to find line to line examples.
My current code approach is as follows:

library(sf)

# Create 'coordinate'  dataframes -----------------------------------------------

#Simulate satellite transect coordinates
x0 <- runif(1, 1,10)
y0 <- runif(1, 1,10)
x1 <- runif(1, 1,10)
y1 <- runif(1, 1,10)

#Simulate road coordinates
l_x0 <- runif(1, 1,10)
l_y0 <- runif(1, 1,10)
l_x1 <- runif(1, 1,10)
l_y1 <- runif(1, 1,10)

# Set up a dataframe for the satellite transect
transect.df <- data.frame(
  'x' = c(x0, x1),
  'y' = c(y0, y1)
)

# Set up a dataframe for the simulated road info
road.df <- data.frame(
  'x' = c(l_x0, l_x1),
  'y' = c(l_y0, l_y1)
)

# Create 'sf'  objects -----------------------------------------------

road.sf <- st_sf(id = 'L01', st_sfc(st_linestring(as.matrix(road.df)))) #create a formal sf feature
st_crs(lead.sf) <- st_crs(4326) # assign crs
road.sf <- st_transform(road.sf, crs = 4326) # transform

transect.sf <- st_sf(id = 't01', st_sfc(st_linestring(as.matrix(transect.df)))) #create a formal sf feature
st_crs(transect.sf) <- st_crs(4326) # assign crs
transect.sf <- st_transform(transect.sf, crs = 4326) # transform


texts: How do I mosaic several WorldView-2 images into single image in ArcGIS Pro?

texts: I have downloaded multiple WorldView-2 images that were outlined by an AOI, but it is unclear to me if it is possible to combine two overlapping images into a single image for downstream data analysis.  I have been unable to find a workflow for doing this. I would preferably like to perform this task using ArcGIS Pro.


texts: Why are split-lines saved in the attribute table? (QGIS 1.7.4)

texts:
A few months ago I splitted some lines in a shapefile with official roads because they are not accessible anymore. Now I loaded that file again into my project and the splitlines appears. Because I didn't made a legend yet, they show up. (see the picture below)

edited: I first selected the layer and turned on the toggle editing. Then I activated split features . I drew a line across the point where the road had to be split. After doing all the splitting, I turned off the toggle editing mode.

Is it normal that splitlines are saved in the attributetable? Can this be avoided? Because it is messing up my data.


texts: Can QGIS be coaxed into exporting text as text in an SVG?

texts: Having resolved my squished projection I'm now facing the fact that QGIS exports all my labels as vectors, which means I can't take the map over to Inkscape and fine tune the fonts without retyping them all. Is there any way around this?