texts:
When I execute this code:
list1 = [2500, 4600, 6500, 7500, 10500, 12500, 14500, 16000, 17500, 19000, 20500, 22500, 24500, 26500, 29500, 31500, 33500, 35900, 37500, 39500, 41500, 43200, 45500, 47500, 49200, 51500, 53500, 55500, 58000]
for i in list1:
    for row in arcpy.da.SearchCursor('shapefile.shp', ["SHAPE@"]):
        uu = (i/60910.0)
        print uu
        point = row[0].positionAlongLine(uu, True)
        cursor = arcpy.da.InsertCursor('pointsshp.shp', ("SHAPE@"))
        cursor.insertRow([point])
        del cursor
Then the 1st point is 2800 meters instead of 2500 meters from the start of the line. All points are not precise - I can understand a few meters but not 300 meters! The line is 6910 meters long for sure. Both shapefiles are in GCS_WGS_1984 and the data frame is projected to WGS_1984_UTM_Zone_33N.


texts:
I am a newbie to QGIS. Currently trying to augment a dted raster value file with vector shapefile(polygon) containing height of buildings.

Have converted the shapefile into raster via rasterize function. WheneverI mouse over the converted shapefile(now raster), I am able to get the value of the polygonal shape (from the band value) corresponding to the old shapefile polygon elevation

For the dted, the mouseover seem to give me the quantized color depth(band value). How could we get the real height value? Or display them?

Lastly, how do we add the two values together and merge them into a single raster data set. I have already projected them to use the same reference coordinate system

Thanks


texts: Converting raster to vector from multiple raster values?

texts:
I have a raster GeoTIFF with values from 0 to 255 representing intensities. I want to Vectorize this information, but in ranges of 25.
As an example, I want to create a layer in which polygon islands representing values grouped  between 230-255; 25-50; etc.

In another thread someone recommended using the "Vectorising Grid Classes" but that only works for 1 value!


texts: What's the meaning of numbers before an EPSG code, e.g. the 6.9 in EPSG:6.9:4326?

texts:
A WFS server I was using offered boundaries data in the following EPSG codes:

EPSG:6.9:27700
EPSG:6.9:4326

...defined in the XML giveCapabilities spec like this:

<DefaultCRS>urn:ogc:def:crs:EPSG:6.9:27700</DefaultCRS>
<OtherSRS>urn:ogc:def:crs:EPSG:6.9:4326</OtherSRS>

As I understand it, EPSG 4326 is the "standard" latitude and longitude ellipsoid used by GPS systems, and EPSG:27700 is the older British system behind the Ordnance Survey National Grid.

But what's the 6.9: mean?

I tried inputting it into ogr2ogr and it just got confused:

  {"errors":["ERROR 6: EPSG PCS/GCS code 6 not found in EPSG support files.  Is this a valid","EPSG coordinate system?","Failed to process SRS definition: EPSG:6.9:4326",""]}

...and I can't find any reference to what the component parts of an EPSG code mean, only lists of 'standard' EPSG codes that don't have these initial numbers.


texts: House numbers as multiple stops route with pgRouting

texts:
How do I create a routing with multiple stops?
Currently I have uploaded on PostGIS the nodes and edges of the road graph but I want the stops to match the housenumbers. The ideal work would be to considerate the time windows of each stop and the capacity of the vehicle (once reached the max capacity, the vehicle must return to a point 0, resume the route and then complete it).
Which is the best algorithm to solve this problem?
Here is a small area of housenumbers and nodes (orange points) on QGIS


texts: Creating DEM from cross sections in ArcMap?

texts:
I have a DEM in arcmap and cross section of a river which flows on the Dem in Hecras.

I now the cross sections are more accurate. I want to make a DEM from some of my cross sections and stitch it on my original Dem to analyze 2d flow in the other hydraulic soft ware.

How can I make Dem from cross sections and How can I stitch  it?


texts: Is it possible to create Non Earth map in openlayers?

texts:
I'm working on a little hobby project for a game, it has some information about locations that are non earth coordinates with extents of 99999999999999999 to -99999999999999999 might even be larger.

I have a few questions regarding such large coordinates.

Can openlayers handle such large coordinates
Would it be best to scale it down to fit lat lon coords to mimic earth coords
what would the maximum decimals be for lat long accuracy in openlayers if I cannot use the full non earth coordinates?


texts: Performance tuning Python code to query and update data in PostGIS stored on local machine?

texts:
I developed this code to query and update data in a postgres/postgis database that is stored on a local machine. The table being analyzed has only around 150K records. The code takes 6-7 hours to run, and I am wondering why? The data does have a spatial index. Is there a way the code could be written more efficiently, or could it be the way the database is set up, or my machine?

The code took about 6-7 hours to run on a Windows 2007 64-bit machine. I also tried to run on 2016 Macbook Pro and had to kill it; I set it to run overnight and after 9 hours or so I killed it.

How can I optimize this?

I am running PostgreSQL 9.5.

import psycopg2

FROM_ST_DICT = {}
TO_ST_DICT = {}

CONN = psycopg2.connect("dbname='test' user='postgres' host='localhost' password='pg*admin'")
print 'successfully connected to the database'

cur = CONN.cursor()
cur.execute('SELECT * FROM sf_citylots_unique')
for row in cur:
    GEOM = row[4]
    FROM_ADDRESS = row[0] + ' ' + row[2] + ' ' + row[3]
    TO_ADDRESS = row[1] + ' ' + row[2] + ' ' + row[3]

    if GEOM in FROM_ST_DICT:
        if FROM_ADDRESS not in FROM_ST_DICT[GEOM]:
            # append the new number to the existing array at this slot
            FROM_ST_DICT[GEOM].append(FROM_ADDRESS)
    else:
        # create a new array in this slot
        FROM_ST_DICT[GEOM] = [FROM_ADDRESS]
    if GEOM in TO_ST_DICT:
        if TO_ADDRESS not in TO_ST_DICT[GEOM]:
            TO_ST_DICT[GEOM].append(TO_ADDRESS)
    else:
        TO_ST_DICT[GEOM] = [TO_ADDRESS]
print 'adding FROM address data...'
for f in FROM_ST_DICT:
    cur.execute('UPDATE sf_citylots_unique set from_address = '+"'"+
                ','.join(FROM_ST_DICT[f])+ "'" + ' WHERE geom ='+"'" +f+"'")
print 'adding TO address data...'
for t in TO_ST_DICT:
    cur.execute('UPDATE sf_citylots_unique set to_address = ' + "'" +
                ','.join(TO_ST_DICT[t]) + "'" + ' WHERE geom =' + "'" + t + "'")

cur.execute("""select distinct a.from_address, a.to_address, a.geom into sf_citylots_unique_geom_with_addresses
               from sf_citylots_unique a""")

CONN.commit()
CONN.close()

Some background, the data represents city loth with associated address ranges. There are duplicate geometries where multiple addresses/ranges exist. I want a final layer that contains only unique geometry, and to retain all addresses associated with each lot. I set up 'to_address' and 'from_address' fields that should contain all full street addresses associated with each lot, like 123 Main St, 456 Park St for the from_address column, and 125 Main St, 458 Park St for the to_address column. In some cases there could be a large number or to/from addresses, and in some cases the to and from values are the same.


texts: How to check in Python if a model has been successfully executed from a toolbox?

texts:
I have created a Toolbox (alias: mytool) and a model (name: mymodel; and it takes 2 parameters which are 'this' and 'that') with ModelBuilder. I need to execute this model from a Python script.

arcpy.ImportToolbox(r"c:\temp\MyTool.tbx")
arcpy.mymodel_mytool(r"D:\temp\this", r"D:\temp\that")

The script runs well but didn't print anything. My question is how to check if this model has been successfully executed by printing out messages from the command prompt.


texts: Basic MapServer query of PostgreSQL/PostGIS database

texts:
I've just installed MapServer on a linux server and I am trying to understand how to display a basic map using a PostGIS database. I have displayed a shapefile in the browser, so I know that the mapserver setup works. My mapfile is as follows:

MAP

        NAME "sample"
        EXTENT -180 -90 180 90 # Geographic
        SIZE 800 400

        LAYER
          CONNECTIONTYPE postgis
          NAME "roads"
          CONNECTION "user=username password=**** dbname=databasename host=hostdb.com port=5432"
          DATA "geom from public."MLSOA""
          STATUS ON
          TYPE LINE
          CLASS
                STYLE
                  COLOR 0 0 0
                END
          END
        END

END # end of map file

I then call this by using a url string:

http://mysite.com/mapserv.cgi?map=/home/user/webapps/htdocs/mapfile_pg.map&layer=roads&mode=map

However, I get the error:

loadLayer(): Unknown identifier. Parsing error near (MLSOA_London):(line 11)

I have also tried calling my database just MLSOA and 'MLSOA'. If I use the line
DATA "geom from MLSOA" I get the following error:

msDrawMap(): Image handling error. Failed to draw layer named 'roads'. msPostGISLayerWhichShapes(): Query error. Error (ERROR: relation "mlsoa" does not exist ) executing query: select encode(ST_AsBinary(ST_Force_2D("geom"),'NDR'),'hex') as geom,"gid" from MLSOA where geom && GeomFromText('POLYGON((-180.225563909774 -90,-180.225563909774 90,180.225563909774 90,180.225563909774 -90,-180.225563909774 -90))',find_srid('','MLSOA','geom'))

I think that the problem is with the table name, but I'm not sure how it should be named. Using pgAdminIII, if I want to run a query on this table I use the following syntax:

select getsrid(the_geom) from public."MLSOA" limit 10

What is the correct syntax for my mapfile?


texts: Expression to show specific attribute of feature within (aggregate?)

texts:
I simply want to get a specific attribute (f.e. as label or default value in form) of a feature within another feature. Aggregate could do the job but it don't.

Hard for me to write it down, so let me explain on the screenshot below: All polygons are in the same layer. The green polygon is labeled with the expression:

aggregate('SameLayer','concatenate', to_string("Number"), intersects($geometry, geometry(@parent)))

That works fine so it tells me, I'm on the right way! All intersection features (101, 102, 103) are written to the label (or later to the attribute form). Since I just want the number of the the polygon, where the green one complete fits within, I changed 'intersects' to 'within' (or 'contains') and expect as result: 102  - but it don't!

aggregate('SameLayer','concatenate', to_string("Number"), within($geometry, geometry(@parent)))

The reason is, the green polygon touches with 3 of its sides, the parent polygon. If I draw the polygon complete within the parent, the expression works.

I wonder what could be the right geometry function in my case, and if Aggregate is the right way to solve that trivial task?


texts: RuntimeError: cannot open, when using arcpy outside ArcGIS

texts:
When I ran the below code, I got a runtimeerror: can't open the shapefile. The file was not occupied by ArcGIS. I even restarted the computer, but it still gave me the error. However, I can run the code in the built-in python with ArcGIS. Why?

import timeit
start = timeit.default_timer()

import arcpy

fc = r'C:\Users\TL\Desktop\1016test\CursorDelEmpRd\TigerMerge.shp'

with arcpy.da.UpdateCursor(fc, "FULLNAME") as cursor:
    for row in cursor:
        if row[0] == '':
            cursor.deleteRow()

stop = timeit.default_timer()

RunningTime = stop - start
print "The running time of the whole process is: ", RunningTime


texts: Easy way to Create Vector Grid polygons?

texts:
We're currently using MapInfo & MapBasic to analyze daily tests results for planning purposes.For each test,exported result may contain empty values (say 999,we know Mapinfo has no NULL values),so we decided to maintain each result (by point values) by averaging and cleaning empty results for each parameters inside a square(piece of grid),that way results may not repeat at the same coordinates and separation of empty values.(Or you can say get results by grid where each grid contain how many points and joining results for each cell by averaging)

In this scenario,firstly we thought about creating these grids via gridmaker.mbx which is free tool for Mapinfo,but when the area is too large or the cell size is too small or even if grid maker has no chance to imply our algoritm that like in the following

We approach this problem by creating rectangles by 10 to 5000 meter values which are square cell sizes,and insert these values onto these tests via MBR or Convex polygon objects by checking out whether cell intersects these generated object from test.

The problem starts here even if we process these cells,generating results may last weeks,months even.So my question is likely in the title,so i ask a different and easy way to approach this problem.

Solutions,scenarios,algorithms you answer may not be solved only in mapinfo,but other tools as well.


texts: Third Party Vector Tile Source in Mapbox Studio

texts:
I am currently styling a map in Mapbox Studio (the web interface, not desktop application) and would like to add a layer that draws on vector tiles from a third party.

When I select "New Layer" and am presented with "No tileset, click to select" for the source data, and then select "New tileset" the only options are to upload a file or create from an existing dataset.

Is there a way to add third party vector data in Mapbox Studio?


texts: How to add labelling to a polygon layer?

texts:
I'm trying to add labelling to a polygons vector layer called "myVectorLayer" with the following code:

palyr = QgsPalLayerSettings()
palyr.readFromLayer(myVectorLayer)

palyr.enabled = True
palyr.fieldName = 'name'
palyr.placement= QgsPalLayerSettings.OverPoint
palyr.setDataDefinedProperty(QgsPalLayerSettings.Size,True,True,'8','')

palyr.writeToLayer(myVectorLayer)

It does not give me any error but the labelling does not appear. Any idea how to fix this?


texts:
Is there a way to compute the wide function in QGIS?. I want to obtain this function to compare the topology of two stream networks.

According to Bras (1990) the wide geometrical function N(x) measures the number of channels at a distance from the outlet of the basin. The distance x can be geometrical (i.e. real distance) or topological (measured in terms of the number channels from the outlet of the basin)

References.

Bras, R. L.  1990.  Hydrology. An Introduction to the
Hydrologyc  Science.  Massachusetts  Institute  of
Technology.


texts: How do I convert ESRI X,Y to Lat, Long for Google Maps?

texts:
I've downloaded the data from NYC Parking regulation. I want to plot that on Google Maps. I use http://www.shpescape.com/ft/ to convert to Google Fusion Maps but the data is in X,Y coordinates instead of Lat, Long. How do I convert X,Y coordinates to Lat, Long? I'm using OSX.

Here's the sample data.

1005740.867110
243957.356623


texts:
I want to publish my .shp file using QGIS Cloud plugin.

I signed in to qgiscloud.com (free). On create new database it will auto create dbname='fbeudp_mfaphk' user='fbeudp_mfaphk' password='68456287' on host='spacialdb.com' port='9999'.

I am unable to change the dbname, user or password (no option is given). On uploading this (.shp)file, an error is showing like this:

Windows firewall is off. I am using QGIS 1.8.0.

My question is: How can I remove this error, so that I can publish my file or map on QGIS cloud (website). Is it possible to change the dbname?  If yes, How?, Why my account not shown any created db i.e. dbname='fbeudp_mfaphk'.


texts: "Merge Selected Features" Persistently Grayed Out

texts:
I am perplexed by a new occurrence in QGIS,  I seem unable to merge selected features e.g. Polygon.   Toggling Editing on a polygon layer,  and selecting two polygon features I would like to merge into one,  the merge polygon feature is always greyed out using the digitizing tools toolbar.  Now I have tried the following:

Saving the polygon in a different CRS (Save as)
Creating a simple polygon feature to compare the results (same results cannot merge any polygon features as the function is greyed out)
Restarting in a new project.
Updating QGIS to the latest version 2.18.7
Testing it on multiple polygon datasets.
Checking geometry.

Still the merge selected features is always greyed out no matter what I do.

Am I doing something wrong?   I have used this feature in the past but now it is not working at all. I am not a GIS expert; just using it for creating simple maps.


texts: Performing viewshed analysis from each pixel (now points) in raster and getting JSON out of it using QGIS?

texts:
My question can be divided in 3 steps and I guess I can answer 1.5 of them already.

To perform the viewshed analysis from every single pixel in a raster I will need:

A viewshed plugin (which I am willing to use either the one from GRASS or QGIS plugin for visibility analysis

Now the 3 steps are:

STEP 1: Get the points layer from the raster

This is my raster:

What I've done was used Raster Pixels to Points tool:

and then I got this layer:

which I think it's was I need for step 2.

STEP 2: Perform Viewshed Analysis for each Point

As I said I would use one of those plugins/tools, but since the one use in GRASS requires to select one pixel at a time, I was trying to used the other one (QGIS plugin for visibility analysis).
So here is the configuration of the tool that I used:

--> From this point on I have problem achieving what I wanted <--

After I run it I have an error associated with the fact that I am not passing some fields ( which I think are also fundamental, observ_hgt and radius) and I don't know either where to put them...

I even tried to add the fields observ_hgt with the value 1.5 and radius 2000 to each pixel (by clicking Toggle Editing > Open Field Calculator) but afterwards I got another error (as can be seen below):

From this point on I don't know how to proceed and even Step 3 I don't know how would it be possible to achieve. But what I had in mind was to get an output like this...

{
  "70.24,214.22": {
    "visible_coordinates" : [ "70.23,214.24", "70.27,216.24",.... ]
  },
  "70.24,215.22": {
    "visible_coordinates" : [ "70.23,214.24", "70.27,216.24", .... ]
  },
  "70.24,216.22": {
    "visible_coordinates" : [ "70.23,214.24", "70.27,216.24", "70,4...." ]
  },
  .........
}

I said JSON but what I need is a way to get the information easily using python of knowing this: given a pixel coordinate, from which pixels is it visible? (something like reverse viewsheding, hope this isn't that dumb to say)

How do I solve Step 2?


texts:
I have a GPS track in a FGDB that has an ID field with A or B. I converted these points into a line using Points To Line, grouping it by the ID field.

I then want to split the line at the same points I used to create the line. When I run the Split Line at Point tool I get an empty feature class as the return. Outputs are placed into a FGDB as well.


texts: Importing shapefiles in folders into geodatabase?

texts:
Using ArcGIS Desktop 10.3, I need to import thousands of folders containing shapefiles into a geodatabase. How can I do this?

Here is a part of the files that need to be exported into a geodatabase:


texts:
Followings is my code. I use it in QGIS 3.8.2.
It shows "Expected an indented block" error.
Please help me to sort the error out.

from PyQt5.QtCore import QVariant
for layer in QgsMapLayerRegistry.instance().mapLayers().values():
    if layer.name() == "Tree":
with edit(layer):
        layer.dataProvider().addAttributes( [ QgsField("Name", QVariant.String, ",", 255) ] )
        layer.dataProvider().addAttributes( [ QgsField("Height", QVariant.String, ",", 255) ] )

else:
    pass


texts: Button to add a new line in a non-spatial table without opening the table itself using QGIS?

texts:
While using QGIS and having loaded a non spatial table (from PostGIS), is it possible to add a new row to its table attributes without having to open the attribute table first as in here Is it possible to add a new line in a non-spatial table using Qgis?

Same principle as add a point but with no geometry, just information in a table.

I was able to create a plugin that chooses and starts editing the table but can't get the form window to open so I can fill the fields.


texts:
I've a custom WPS function gs:GetStyleSize, that I used to call from SLD like the following:

<FeatureTypeStyle>
    <Rule>
        <PointSymbolizer>
            <Graphic>
                <Mark>
                    <WellKnownName>circle</WellKnownName>
                    <Fill>
                        <CssParameter name="fill">#FF0000</CssParameter>
                    </Fill>
                </Mark>
                <Size>
                   <ogc:Function name="gs:GetStyleSize">
                    <ogc:Function name="parameter">
                        <ogc:Literal>style</ogc:Literal>
                        <ogc:PropertyName>style</ogc:PropertyName>
                    </ogc:Function>
                  </ogc:Function>
                </Size>
            </Graphic>
        </PointSymbolizer>
    </Rule>
</FeatureTypeStyle>

I need to convert this SLD to CSS, I tried the following:

/* @title gold point */
*{
    mark: symbol(circle);
    mark-size: gs:GetStyleSize(style);
    :mark {
        fill: #ffcc00;
    }
}

But when I click validate in the Geoserver's New Style page, I gives me the following error:

Invalid input ')', expected NameCharacter, OptionalWhiteSpace or Colon (line 4, column 37)

What's the correct way to call that function, or any other WPS function?


texts: Buffering around point using WFS Filter?

texts:
I'm trying to make a buffer around a geometry (point/s,linestring,etc) using WFS queries without luck.
This is what I've been doing:

   <ogc:Filter>
     ...
     <ogc:Function name="buffer">
       <gml:Point>
         <gml:coordinates cs="," ts=""> x,y </gml:coordinates>
       </gml:Point>
       <ogc:Literal>300</ogc:Literal>
     </ogc:Function>
     ...           </ogc:Filter>

I read I could do it using Dwithin, but is it possible with the buffer?


texts: Legal Issues Concerning Intellectual Property

texts:
I'd like your opinion:

When the GIS Analyst leaves his job (not that I'm in this position, thank god), is he allowed legally to take with him odds and bits he collected or made during his work (small stuff like e.g. small VBA/Python scripts, examples of complex labels), or Like all his bigger Projects it is copyrighted and belongs to the hiring Company? I mean, it's clear that he can't take programs, layers and big stuff (or secret data for that matter), but what about the small non-specific odds and bits?


texts: Visualizing serialized GeoJSON data in HTML table

texts:
I am using Django and I have serialized data as GeoJSON format. I can visualize that with HTTP response but I want to put these data in table.
views.py
from django.shortcuts import render
from django.http.response import HttpResponse
from geojson import Feature, Point, FeatureCollection, GeometryCollection
from django.template import loader

from crt_ast.models import Asset
from django.core.serializers import serialize

def cesiumAsset(request):
    geojs = serialize('geojson', Asset.objects.all(),
              geometry_field='geom',
              fields=('name', 'markercolor', 'markersymbol', ))

    return HttpResponse(str(geojs))

I  added serialization module on my settings.py
SERIALIZATION_MODULES = {
    "geojson": "django.contrib.gis.serializers.geojson",
 }

This is the models.py created for database
class Asset(models.Model):
    type = models.ForeignKey(AssetType, on_delete=models.CASCADE)
    lc_phase = models.ForeignKey("Asset_Life_Cycle.LifeCyclePhase", on_delete=models.CASCADE)
    name = models.CharField(max_length=30)
    # BaseSpatialField.srid()
    geom = models.GeometryField()
    elevation = models.IntegerField(null=True, blank=True)
    photo = models.ImageField(blank=True, null=True, upload_to='asset_photo/')
    comissioning_date = models.DateField(null=True, blank=True)
    decommission_date = models.DateField(null=True, blank=True)

    active_choices = (('Yes', "Active"), ('No', "Inactive"),
                      ('Maintenance', "Maintenance"))
    active = models.CharField(max_length=11,
                              choices=active_choices, null=True, blank=True)

    description = models.CharField(max_length=256, null=True, blank=True)

    markersize= models.CharField(max_length=7, null=True, blank=True)
    markercolor = models.CharField(max_length=7, null=True, blank=True)
    markersymbol = models.CharField(max_length=15, null=True, blank=True)

    def __str__(self, *args, **kwargs):
        return "{}, Type: {}".format(self.name, self.type.name)

I can visualize with http response like that:
{"type": "FeatureCollection", "crs": {"type": "name", "properties": {"name": "EPSG:4326"}}, "features": [{"type": "Feature", "properties": {"name": "Keban Barajı 315. Nolu Boru", "markercolor": null, "markersymbol": null}, "geometry": {"type": "Point", "coordinates": [38.75894165877253, 38.88087765578684]}}, {"type": "Feature", "properties": {"name": "deneme1", "markercolor": "#FF0000", "markersymbol": "1"}, "geometry": {"type": "Point", "coordinates": [62.578125, 43.58039085560785]}}]}

I don't know how to put them inside the table.


texts: Openlayers vector feature coordinates are not transforming,

texts:
I have a OpenLayers vector with a single feature [0]. It's coordinates are in the map projection (EPSG:3413) and I want them in WGS84 (EPSG:4326) decimal degrees. I have the following code which executes with no errors.

    // Transform the polygon to WGS84 (EPSG:4326)
        polybound.features[0].geometry.transform(new OpenLayers.Projection("EPSG:3413"),new OpenLayers.Projection("EPSG:4326"));

        // Get the extent/vertices of the polygon drawn
        var polyBounds = polybound.features[0].geometry.getVertices(); alert(polyBounds);
        var polyExt = polybound.features[0].geometry.getBounds(); alert(polyExt);

However both the alerts show the untransformed coordinates.

Eg:

alert(polyBounds);

POINT(-2486558.28125 -1947299.71875),POINT(-4121511.40625 -1737690.34375),POINT(-3555566.09375 -2722854.40625)

alert(polyExt);

-4121511.40625,-2722854.40625,-2486558.28125,-1737690.34375

How do I get the vector feature coordinates to transform correctly?


texts: I installed  OCHA icons on windows 8,But once that I try to open the styles with arcmap (10.2.2) I have some funny icons, do you have some idea how to fix this bug...


texts: Using multiple fields for random points in polygons QGIS?

texts:
I'm quite new to QGIS and I have a question about the 'random points in polygons' functionality.

Is it possible to use more than one field from the attibute table when using 'random points in polygons' (in the expressions field)?

Currently I can select a field from the attribute table to generate the number of points I want. But now I would like to select two fields from the attribute table and have their sum as the number of points. And of course giving them both a different color.

Is this possible to do in the expression field or elsewhere?

Current situation: 10 points in polygon X from one field in one color

=> Polygon X,
Field A = 10,
(Field B = 0)

Desired situation: 10 points in polygon X from two fields in two colors

=> polygon X,
Field A = 5,
Field B = 5


texts: WMS GetFeatureInfo and Google maps api V3

texts: Is there a way to request GetFeatureInfo from a WMS service with the Google maps api v3?


texts: QGIS: put vector values into separate layers?

texts: QGIS 1.8.0 Mac: I've got a polygon shape file with values 1-5 in it i.e. contours. How do I get the polygons for value 5 into one layer, 4 into another... so I can assign each a colour and see them?


texts:
I want to plot my RGB raster on a basemap in R.

So far I have a RasterStack object and I have downloaded the basemap with ggmap. I can plot them individually but not together.

stack_rrs <- stack(red4, green3, blue2)
plotRGB(stack_rrs,r=1,g=2,b=3, stretch = "lin")
bm <- ggmap(get_map(location = "Location", maptype = "satellite"))


texts:
How do I use GDAL tools to rasterise a shapefile, while also reprojecting to a different target coordinate system?

The gdal_rasterize documentation says on-the-fly re-projection is supported in recent versions, but doesn't explain how to perform it. It seems not to document any projection related arguments except -a_srs (which it says is for overriding not reprojecting). Do I need to first create an empty destination file with the target attributes, and if so then what tool do I use to do so?


texts:
I have trying to find way to get google map display route polygon, so that i can filter data from my database(using polygon):

Google directionsService.route

There is some option I found but it is not helpful because it is related with google map : drawing tools polygon


texts: I would like to change all the values in a attribute field to the same text

texts: So I recieved a great answer to this question, How to change all values in a field in one operation?, although I need some further information. I would like to change the values in the entire field to the font name that I would like it to be. When I type in the font name, it claims that it is an invalid expression. Another example of what I would like to do is change all values of another column to half, related to the vertical alignment.


texts: Calculate midpoint from a series of latitude and longitude coordinates

texts:
I've got a series of longitude and latitude coordinates which represent a building outline

e.g.

-0.5485381346101759,53.2285150736142
-0.5482220594232723,53.22842450827133
-0.5482298619861881,53.22841205254449

...(intermediate points not listed)...

-0.5483123769301657,53.22882101914848

How can I work out the midpoint? I've found tutorials that show how to do it if you've got three coordinates (e.g. http://mathforum.org/library/drmath/view/68373.html), but in many cases I've got more than three.

Thank you


texts: Appending characters before the filename of all layers in a GeoPackage in QGIS?

texts: I have a geopackage which contains around 30 files. I want to put a "mmbl_" before the actual file name for all the layers in the said geopackage. Is there a an SQL expression to automate this renaming process in QGIS?


texts: Filter features outside the viewport in Mapbox

texts:
I have a huge GeoJSON file having around 100K features and I want to filter features outside the viewport. I tried multiple approach given on:

https://github.com/mapbox/mapbox-gl-js/issues/8350
and
How to load geojson features outside view with Mapbox GL?

But nothing seems to work. My Code is:

//Here i have already added empty geojson source
 this.instance.addSource('seats_geojson', {
            type: 'geojson',
            data: {
                "type": "FeatureCollection",
                "features": []
            }
        });

//Here i have axios call with await and i am getting all the feature's in console
//console.log(e.data.featuresCollection.features) <--- 100K features

this.instance.once('idle', () => {

          //Once i have the reponse ready I am setting source using setData
          this.instance.getSource('seats_geojson').setData(e.data.featuresCollection);

          //Here i am creating empty layer with above source
                         this.instance.addLayer({
                         id: "rowSelectionDup",
                         type: 'circle',
                            source: 'seats_geojson',
                            paint: {
                                'circle-color': "#4400d9"
                            },
                            filter: [
                                "in", "s", ""
                            ]
                        });
                    });

but on doing something like:

this.instance.querySourceFeatures('seats_geojson', {
  sourceLayer: 'rowSelectionDup',
  filter: ['in', 's', "1C"] //This section "s" 1C is in viewport and i am getting 207 feature's
});

this.instance.querySourceFeatures('seats_geojson', {
  sourceLayer: 'rowSelectionDup',
  filter: ['in', 's', "7C"]. //This section is outside viewport and result is []
});

this.instance is new Map(config...) in mapbox-gl

Am I doing something wrong here or is there any other approach to get feature's from GeoJSON?


texts: Exporting a CSV file from QGIS into several CSV files based on a field attributes

texts:
I have loaded a CSV table (no coordinates, just data) into QGIS. The table consists of data with a field containing timestamps.

Is there a way to export the same table to CSV and creating separate files containing the same timestamp? Each filename would be the unique timestamp.

I have tried the Export ---> Save feature as..., I can export the entire table back into a single file but I cannot separate into several files containing the same timestamp.

I can use FME to do that using the "fanout" technique, but I need to use QGIS if possible.


texts: Add polygons to existing shapefile in QGIS

texts:
I'm completely new to GIS and this is a very basic question about editing existing shapefiles.

I am working off an existing shapefile which was digitized (by a colleague) so that polygons represent building cover. The shapefile contains information about the land use type and area. I am using QGIS to edit the shapefile.

Now, I wish to add/remove/edit polygons to improve accuracy as the extent of building area has changed a lot recently.

However, when I add new polygons into the shapefile (by digitizing it using an updated basemap), they are not added into the same polygon class as the earlier shapefiles, even though I input the right LandUse type (i.e. Buildings) in the attribute form. Each new polygon seems to create a new LandUse category, even though they're all namd "Buildings".

How do I ensure that any new polygons I add are added into the same LandUse class (all previously created polygons are in this "Buildings" class)?

Here's the attribute table

And here's an image of the map. The light pink polygons are the new ones I've added. And the dark pink ones with the red outlines are from the old dataset.


texts: Generate polygon automatically based on locations

texts:
I'm trying to achieve some similar to Draw a polygon automatically through markers (points) in google maps but not exactly.

Given a set of locations I'd like to know if it is possible to generate the polygon container for these locations automatically. Something like what Google does with fitBounds method. But if possible, I'd like the bounds snapped to streets. For example:

Given the current set of locations:

lat: 41.87368, lon: -87.65606
lat: 41.87388, lon: -87.65452
lat: 41.87327, lon: -87.65525
lat: 41.87241, lon: -87.65636
lat: 41.87317, lon: -87.64001
lat: 41.87397, lon: -87.64186
lat: 41.87397, lon: -87.64027

my desired polygons would be:

The one in red and also acceptable would be a variant with the one in green.

Link to map with locations: http://bl.ocks.org/anonymous/raw/f0e736dadb07a0634b228e9ff55fcbf9/

Is this possible? I'm doing this on my backend, so my question narrows it down to available APIs that could do it. I'm currently working with OSRM, Google Maps (regular and static APIs) and OpenStreetMap.

To summarize, given a set of locations are there any APIs that would automatically generate the bounds polygon for me?


texts: How to uninstall standalone QGIS & reinstall in OSgeo4w on Win7?

texts:
I need to apply a patch to QGis but installed as standalone instead of OSGeo4w: foolish.

I'd like to uninstall then reinstall via OSGeo4w so as to have linux-like  updating support, but before I do so, can anyone advise me whether my maps & WMS info will be saved, or if there's anything one specifically needs to do to avoid extra headaches/work during this process?

Many thanks in advance.


texts: What is point density when using the random points processing algorithm?

texts:
I'm assuming the density of points means the number of points within an area. So if I use the settings shown in the image, the tool would make 1 point per 100 meters?

I tested this on a single polygon which is around 1,000 squared meters but the tool just hangs. After about 30 mins, I just quit it. Is my understanding correct? Should it take this long to generate a few points?

The project CRS is the same as the layer CRS (EPSG:27700).


texts: georectification via SLD is not appearing at the correct place

texts:
I am trying to georectify layers on the fly (feeded by Geoserver) in OpenLayers.

To do so, we have created a style (SLD) in Geoserver that does the georectification provided that we feed it a projection and coordinates via the GET request (env variable):

<?xml version="1.0" encoding="ISO-8859-1"?>
<StyledLayerDescriptor version="1.0.0" xsi:schemaLocation="http://www.opengis.net/sld StyledLayerDescriptor.xsd" xmlns="http://www.opengis.net/sld" xmlns:ogc="http://www.opengis.net/ogc" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance">
    <NamedLayer>
      <Name>GeorectifyCoverage</Name>
      <UserStyle>
        <Title>Georectify Coverage</Title>
        <Abstract>Warps image to specified GCP</Abstract>
        <FeatureTypeStyle>
          <Transformation>
            <ogc:Function name="gs:GeorectifyCoverage">
              <ogc:Function name="parameter">
                <ogc:Literal>data</ogc:Literal>
              </ogc:Function>
              <ogc:Function name="parameter">
                <ogc:Literal>gcp</ogc:Literal>
                <ogc:Function name="env"><ogc:Literal>wms_gcp</ogc:Literal></ogc:Function>
              </ogc:Function>
              <ogc:Function name="parameter">
                <ogc:Literal>targetCRS</ogc:Literal>
                <ogc:Function name="env"><ogc:Literal>wms_crs</ogc:Literal></ogc:Function>
              </ogc:Function>
              <ogc:Function name="parameter">
                <ogc:Literal>transparent</ogc:Literal>
              <ogc:Function name="env"><ogc:Literal>wms_transparent</ogc:Literal></ogc:Function>
              </ogc:Function>
            </ogc:Function>
          </Transformation>

          <Rule>
           <RasterSymbolizer>
             <Opacity>1.0</Opacity>
           </RasterSymbolizer>
         </Rule>

        </FeatureTypeStyle>
      </UserStyle>
    </NamedLayer>
</StyledLayerDescriptor>

I then display the layer over an OSM map with and pass the coordinates of the georectification through the env argument:

var wms5 = new OpenLayers.Layer.WMS("Geotiff5", "http://mydomain.com/foobar",
                    {
                        layers: data[0].imagename+":"+data[0].imagename,
                        styles: "georectify2",
                        env: "wms_transparent:false;wms_crs:EPSG:900913;wms_gcp:"+getGCP(pixel,coord3),
                        transparent: true
                    },
                    {
                        singleTile: true,
                        isBaseLayer: false,
                        displayOutsideMaxExtent: true
                    }
                );

where

getGCP(pixel,coord3)="[[[0, 0], [725801.41007875, 5861382.4105691]], [[0, 450], [725725.04490808, 5858019.2685367]], [[600, 0], [725801.41007875, 5861382.4105691]], [[600, 450], [733521.52808368, 5859280.356713]]]"

I get a result where the layer is rectified, but the it doesn't appear at the right place, and the image is not cropped. I did a screenshot of the result. The polyon represents where the image should be, and the picture is just a random image I used for testing purposes. One can see that the upper left and lower right corner are missing.

I don't know if using a SLD is the most elegant solution, but it is the only one I found to do it on the fly.

Any Ideas?


texts: How to create a clickable hyperlink in a QGIS form?

texts:
I have a hyperlink field in my table. I wish I could go to those addresses without using actions. I want to go directly from the form by clicking on them as happens in a hyperlink. It is possible?

I want to load web explorer doing double-click in url.


texts: How can QGIS 3.X create atlas by a group field?

texts: I have a line layer with 8 features and a field TYPE (5 features is type_A and 5 features is type_B). I want to create a atlas which output is two single image (one is grouped by typ_A with 5 features, other is grouped by type_B with 3 features) by using Rule-based Symbology to control the features showing or not. But it seems the atlas create 8 single images not what I want. So is there a way to achieve my goal?I have more than 100 features group by more than 10 different types.


texts: Task Scheduler doesn't like Floating License

texts: While trying to run an FME workspace weekly it has become clear that Task Scheduler won't open FME. I think this is due to missing the floating license credentials. Has anyone else run into this? If so, have you found a workaround?


texts: PyQGIS: How to find features  another feature

texts:
Using a spatial query in QGIS Desktop to select features of layer MTA_ACQ_TE_polygon that are within features of layer nts_snrc_50k works, but I am attempting to do this programatically in PyQGIS, with no luck.

And this is my attempt:

from qgis.core import *

app=QgsApplication([], False)
app.setPrefixPath("C:/OSGeo4W/apps/qgis", True)
app.initQgis()

mapsheet = QgsVectorLayer("E:/MTOdata/archives/0/MTA_ACQ_TE/nts_snrc_50k.shp", \
    "mapsheet", "ogr")
tenures = QgsVectorLayer("E:/MTOdata/archives/0/MTA_ACQ_TE/MTA_ACQ_TE_polygon.shp", \
    "tenures", "ogr")

if not mapsheet.isValid() or not tenures.isValid():
    print("failed to open shapefiles")
    exit(0)

print("-> loaded mapsheet (%s), tenures (%s)" % (mapsheet.crs().authid(), tenures.crs().authid()))
print("-> finding mapsheet")

for m in mapsheet.getFeatures():
    if m["NTS_SNRC"] == "094E06":
        print("\tfound it!")
        sheet = m.geometry()
        break

print("-> finding claims within")

for t in tenures.getFeatures():
    tenure = t.geometry()
    if tenure.within(sheet):
        print("\tgot one: %s" % t["OBJECTID"])

print("-> done")

Which provides the output:

-> loaded mapsheet (EPSG:4140), tenures (EPSG:3005)
-> finding mapsheet
        found it!
-> finding claims within
-> done

It fails to find any features of MTA_ACQ_TE_polygon within nts_snrc_50k.  Why?


texts: Finding start and end junctions of road from coordinate on road

texts:
I have a list of co-ordinates where traffic counts were carried out.

I wish to associate these traffic counts to the link they are selected with. For that link I want to know the location of its two opposing junctions. For example, on a road that goes north south - I want to know the coordinates of the south junction and the north junction.

If I have a shapefile of points and another shapefile of the road graph (the nodes and links) is there a GIS method for

finding the road link that coordinate is on and
the coordinates of the closest junctions on that link?

In this case I have a gml file containing the Integrated Transport Network from OS and the traffic count coordinates come from the Department for Transport.

For illustration. Ive used Google maps to plot the coordinates of the traffic count. This is shown as the red flag.

View post on imgur.com

I would then like to find the coordinates of A and B.


texts: Leaflet MapBox Issue of poly line which crosses the International dateline

texts:
I am using  mapbox leaflet js in my application "https://api.mapbox.com/mapbox.js/v3.0.1/mapbox.js". Previously i was using google map,
I am facing an issue when i draw a polyline in map which crosses the international datelines, the Line which is drawn on map using MapBOx is not wrapped in one copy of world as in Google map.

I had done lot of Reaserach and found that the coordinates which are crossing the dateline add/subtract 360 from the longitude, which will draw a contineous line on map.

But as my client is habitual with google map services so he needs the polyline to be displayed as in google map.

I had attached a Zip file which contains three Screenshots

1) Original_functionality_mapBox.png : This image shows when coordinates of polyline which crosses the international date line are not changed.

2) change_360_Degree.PNG : This image shows when coordinates of polyline which crosses the international date line are changed by subtracting 360 from longitude.

3) google.PNG: when the original coordinates are passed to google map, it automatically adjust the polyline which crosses the international date line.

Please help me urjently to do the same as in google map.

                     1(Original_functionality_mapBox.png)

                     2(change_360_Degree.PNG)

                     3(google.PNG)


texts: Convert Longitude/Latitude to census block

texts: I would like to convert Longitude/Latitude to census block for the nation.  Is there an easy way to do it? If there is public data available, do I have to do it for each state, or I could do for the whole nation in one batch?


texts: Is it possible to get a bounding box for the rectangular region of unmasked pixels in an Image in Google Earth Engine?

texts:
I'm working with the 2019 Hansen Global Forest Change Image:
var gfcImg = ee.Image("UMD/hansen/global_forest_change_2019_v1_7");

Even though every band of gfcImg is masked everywhere above about 80° latitude and below about -57° latitude, gfcImg.geometry() covers the entire planet. I need a bounding box in the form of a Geometry which encloses only the rectangular region between these latitudes in which gfcImg isn't masked. Is there a way to do this?


texts: will ogr get coordinates in geographic coordinate system or projected coordinate system when using getPoint(int,OGRPoint*)

texts:
Consider the code that I have written in the stackoverflow question How to extract vertexes of geometries in ESRI shapefiles using ogr library with c++, I have used the code to show shapefiles on an OpenGL window and my problem is it seems the layer is stretched comparing to what is shown when we render it in ArcGIS, QGIS and etc that I think it is maybe because of not using the correct projection.

for a layer, like this:

are the vertex coordinates stored in Geographic coordinate system or projected coordinate system when we use functions like:

OGRLineString::getPoint(int,OGRPoint*);
OGRPoint::getX();
OGRPoint::getY();

when  I debug my code and see the vectors stored:

I come to the point that coordinates are stored in the PCS not GCS, but why my layer is stretched.
In fact I'm not sure what's the problem to start to solve. So just tell me:

Are the coordinates fetched when we use the above functions in
geographic coordinate system or in projected coordinate system?

Edited section based on the Andre Joost's comment
I have downloaded Country Regions dataset here, you see its coordinate system and bounding box here:

and the data stored is:

that you see is in decimal degrees. Now in ArcCatalog I change the coordinate system:
As you see the bounding box is not projected to the new coordinate system and they are the same number just their unit is in meters. Which is obviously not true:

and the data stored with for the same point in my program is:

So I think the projection is not done automatically and I should do the calculations with PROJ.4 library!


texts: Add KML file with associated PNG images as OpenLayers layer

texts:
I have private sonar data in KMZ format, and can easily produce more. I understand that OpenLayers cannot use KMZ files, so my question is if it can handle the uncompressed KML file with all the PNG images. Since all geographic information is in the KML file it should work, right? I've tried with GML and vector layers but just haven't managed to make it show up yet...

My solution so far has been to use MapTiler, create tiles and then use a xyz layer. But I would rather use my KMZ if possible.


texts:
I would like to use the Advanced Digitizing tool and the Move tool to move entire polygons by specified distances and angles.  For some reason the "d" (distance), and "a" (angle) fields are grayed out.

Both my project and my layer are in a projected CRS (EPSG:29611).

I would like to do what he does in this video.

How do I get this tool working?


texts: Calculating nearest point from a road, its name, distance and what road it intersects

texts:
I have a road and I need to calculate the nearest restaurant from it, its distance, its name and which road it is on(the restaurant). I have done all the things above except for the last one. Does anyone have any ideas?

 SELECT a.name, b.name, ST_Distance(a.geom, b.geom) AS distance
    FROM roads a, restaurants b
    WHERE a.name='Roehampton'
    ORDER BY distance ASC
    LIMIT 1


texts: Why can't I modify a python script tool or delete the toolbox it's in?

texts:
So I started getting this weird error for what I thought was no apparent reason:

Executing: saveTile
Start Time: Tue Apr 28 20:19:17 2015
Running script saveTile...
Failed script saveTile...

Traceback (most recent call last):
  File "C:\Image_Processing\Python_Scripts_Automation\Save_Tile_Automation.py", line 236, in <module>
    arcpy.createDashboardReport_ipauto(sReports3, sDatetime)
  File "C:\Image_Processing\Image_Processing_Automation.tbx", line 23, in createDashboardReport
    AttributeError: Object: Tool or environment <createDashboardReport_ipauto> not found

Failed to execute (saveTile).
Failed at Tue Apr 28 20:20:29 2015 (Elapsed Time: 1 minutes 12 seconds)

The tool had been working and it all of a sudden stopped working. Then I got looking, and it seems I no longer have access to the tool. When I open its properties, I can't change the path of the script. I also can't modify the parameters. And the delete option is grayed out when I right click on the toolbox:

Is there a lock on it somewhere? Has anyone seen this before, or know how to fix it?


texts:
I'm using the GDAL Raster>Translate to create JPGs from Georeferenced TIFs. The translation is working, but the resulting JPGs are significantly darker than the TIFs they are being created from.

Here is a screen shot showing the difference in brightness:

Anyone know why this is happening and (more importantly) how I can fix it?

Thanks!


texts: MapInfo Universal Translator turning integer columns to decimal

texts:
I've got a TAB file that I am using FME Universal Translator in MapInfo 12.5 to turn into a SHP.

The input TAB file has a variety of columns, some are floats and some are integers.

The output shapefile has the same fields, but it turns the integer fields into decimal(11,0). The float fields remain floats.

Is there a way to get the output to have the integer fields just like the input?


texts: Troubleshooting and fixing TypeLoadExceptions in ArcMap AddIn when executing compiled c# script?

texts:
I am attempting to compile a c# script (class) on the fly and run it from an Arcmap AddIn (the script contains no ArcObjects). My unit test, in Visual Studio 2013, works just correctly and I am able to compile the c# code, invoke the constructor, and call my method with params and return the expected string value. When the same code is invoked from the ArcMap AddIn, I receive this exception from the method of the compiled class:

Could not load type 'Oldow.Runt.Core.Models.ExpressionObject' from assembly 'Oldow.Runt.Desktop.AddIn

So from the AddIn , the compilation is successful and so is the invoking the constructor. The exception is thrown from the compiled assembly's Run method.

This is my compilation code (no errors or warnings during compilation):

    Microsoft.CSharp.CSharpCodeProvider csProvider = new Microsoft.CSharp.CSharpCodeProvider();
    CompilerParameters options = new CompilerParameters();
    options.GenerateExecutable = false;
    options.GenerateInMemory = true;
    options.ReferencedAssemblies.Add(Assembly.GetExecutingAssembly().Location);
    CompilerResults result;
    result = csProvider.CompileAssemblyFromSource(options, code);

And here's an excerpt of the execution code for running the compiled script (TypeLoadException occurs in scriptObject.Run():

    var exportedTypes = compiledScript.GetExportedTypes();

    foreach (Type type in exportedTypes)
    {
        foreach (Type iface in type.GetInterfaces())
        {
            if (iface == typeof(ICustomExpressionEvaluator))
            {
                ConstructorInfo constructor = type.GetConstructor(System.Type.EmptyTypes);
                if (constructor != null && constructor.IsPublic)
                {   var scriptObject = constructor.Invoke(null) as ICustomExpressionEvaluator;
                    if (scriptObject != null)
                    {
                        returnText = scriptObject.Run(inputs);

Has anyone encountered this before in the AddIn framework and figured out how to resolve?

For reference, here's the test method that successfully runs.

[TestMethod]
public void CompileScriptTest()
{
    ExpressionObject p = new ExpressionObject()
    {
        SKTNumber = "asdasd3312"
    };

    string scriptText = @"namespace Oldow.Runt.Scripting {
using Oldow.Runt.Core.Models;
public class XScript1 : ICustomExpressionEvaluator
{
public string Run(object[] expressionInputs)
{
    var plan = expressionInputs[0] as ExpressionObject;

    return string.Format(""Test SKT Number: {0}"",plan.SKTNumber);
}
}}";

    object[] inputs = new object[] { p };

    var text = ScriptCompiler.EvaluateExpression(
        scriptText,
        inputs);

    Assert.AreEqual("Test SKT Number: asdasd3312", text);
}


texts: MS SQL table with geometry field to feature class

texts:
I have a table in a native SQL database (no geodatabase, dc_geobank_public.sde in the code) with cca. 300 000 records and 18 fields. Two of the field contains projected X and Y coordinates, and one field contains geometry with SQL geometry type. The spatial index is rebuilt regularly with default parameters in MS SQL Server. The projection code is EPSG:23700. Based on this data I would like to create and update a feature class in a geodatabase (dc_geobank_map.sde in the code) on a daily basis.

So far I have written the following script, and it does the job, but clicking on the table in ArcCatalog and using the "Create table... from XY data" runs 100 times faster. (The script runs for 6 hours, and the "Create table... from XY data" tool runs for 3 minutes). Here is the script:

#importing arcpy library
import arcpy

#define workspace, this is the native SQL database
arcpy.env.workspace = r"C:\scripts\dc_geobank.sde"

print "changed workspace!"

#we let the script overwrite existing feature class
arcpy.env.overwriteOutput = True

print "output overwrite set!"

#creating xy event layer from table
arcpy.MakeXYEventLayer_management("geobank.dbo.FRS_fullos","EOV_X","EOV_Y","furasxy","23700")

print "eventlayer ready!"

#change workspace, that's where the feature class should be created, and THIS is a geodatabase
arcpy.env.workspace = r"C:\scripts\dc_geobank_map.sde"

print "changed workspace!"

#convert event layer to feature class
arcpy.CopyFeatures_management("furasxy","furas")

#project the feature class to Web Mercator projection
arcpy.Project_management("geobank_map.DBO.furas"","furasprojected","3857")

Anyone knows what the "Create feature class from XY data" tool might use that's so much faster? So far I am not making use of the geometry field, but I guess it should be useful.

EDIT:
For better understanding of the issue: later on I would like to make a map service and a web mapping application based on this feature class. So my question:

How to make a feature class from the table faster?

I am on the 10.2.1 ArcGIS stack.


texts: Lengths of lines within individual but overlapping buffers?

texts:
I have a set of points for which I have created a 12 mile buffer.

Eventually, I would like to know how many kilometers of road exist within each buffer zone individually?

The problem I see going forward is that the buffers have overlap.

I need each buffer to be individually selectable so that I can apply my roads layer and obtain the desired information.


texts: i am looking for a simple app that will record gps coordinates against photos and videos

texts: is there a simple app that will allow my Android phone to record not just a location name against photos and videos, but also the gps coordinates?


texts:
I am using QGIS 3.4 and I work with a shapefile and its Attribute Table with a few rows where some fields are just blank (no zero value or NULL).

When I export the shapefile file to ArcGIS Pro all the blank fields become zero but they should be NULL. How can I handle this problem?


texts: st_difference (erase overlapped feature) line from polygon

texts:
I have two layers.
One is MultiLineString and the other is MultiPolygon.

I am trying to erase lines that are on the polygon feature
With QGIS difference tool, it worked just fine

This is my code and I do not understand what is wrong with my query.
I do not see anything wrong.

select a.objectid, st_difference(a.geom, b.geom) as shape
from line a, mesh b
limit 500

The results return line segments  that are also in the mesh.


texts:
I have made a basic unsupervised k-means clusterer function in Earth Engine that creates two classes, which is what I want. While the function does apply consistent labels when run multiple times the labels are not consistent between images, i.e. sometimes clusters based on lower values are labeled 1, and sometimes they are labeled 0 (can be seen in code below, where features are black in one image and white in another). I need to run the clusterer on each image in the collection.

Is there any way to force the naming structure for the result to be linked to the image value it is based on? I know that sounds a lot like supervised classification, but I want to be able to run this without input training data and come out with the same naming scheme in each image.

Link to sample code

Reproduced here:

var input = ee.ImageCollection('COPERNICUS/S1_GRD')
 .filterBounds(geometry)
 .filterDate('2015-10-01', '2016-05-10')
 .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'HV'))
 ;

var cluster= function(img){
  var training = img.select('HH').sample({region: trainer, scale:10, numPixels:1000});
  var clusterer = ee.Clusterer.wekaKMeans(2).train(training);
  var result = img.cluster(clusterer, 'kMean');
  return result;
};

var clustered= input.map(cluster);
var clust = clustered.toBands();
print(clust)

Map.addLayer(clust, {bands: ['S1A_EW_GRDM_1SDH_20151024T120105_20151024T120209_008294_00BB14_F72B_kMean'], min:0, max:1}, 'One scheme')
Map.addLayer(clust, {bands: ['S1A_EW_GRDM_1SDH_20151216T120911_20151216T121016_009067_00D054_ADF1_kMean'], min:0, max:1}, 'Other scheme')


texts:
What is the unit of argument "tolerance" used in mosaic {raster}? Is it pixel or raster unit?

x <- srtm_v4.1_90m
names(x)[1:2] <- c('x', 'y')
x$fun <- mean
x$na.rm <- TRUE
x$tolerance = 0.5

I am trying to mosaic srtm tiles using following script but I am not sure about unit of argument "tolerance" used in the script.


texts: How to get bounds coords by center LatLng point, width and height in pixels and zoom?

texts:
I have to calculate bounds coordinates to query mongodb box on server side to get data before google maps is loaded.

How can I get four LatLng corners by center LatLng point, window width and height in pixels and given zoom level?

EDIT:
This answer might be helpful but I don't know how to calculate for google maps resolution or scale denominator.


texts: Preserving measurement between QGIS and Google Earth

texts:
Has anyone found a way around this?

I am having issues when importing shp files into Google Earth. I need reasonably accurate measurement data and my measurements change significantly when imported from EPSG 3857 to Google earth (EPSG 4326). I comprehend the difference between the two but am having issues overcoming the problem since my layers have to be based upon locations as projected in the google openlayers plugin, which, to my knowledge, requires EPSG 3857 to be used. Converting the shp file layer projection to EPSG 4326 when saving as a kml has not worked. Neither has changing the layer or project projection in a new project without adding the google openlayers base map.

Alternatively, if there is a time efficient way to set kilometer points along a >500km line in Google Earth without using QGIS I have open ears.


texts: Mapping data from related table using ArcGIS Desktop?

texts:
I have a polygon shapefile that contains multiple lakes.  Some of these lakes have multiple residuals.  I have all of the residuals for the lakes in an Excel file, and I want to show residuals for the lakes as dots proportional to the size of the residual to see if the size of the residual is influenced by geographic location.

I did a relate, but how would I show all the residuals as symbols on top of the map of lakes?

This is an image of when there is only one residual per lake.  However, I had multiple samples from the lake.


texts:
I want to count points (1000 of students home) in district polygons with Qgis 2.16. I didn't find the tool in the analysis tools plugin.

Is there another plugin or method I could use to perform this task ?


texts: How to calculate zonal statistics from multi-band raster with no attribute table

texts:
I have a raster image (.tif) that is multi bands (35 band), floating point data type that does not have attribute table, and a grid shapefile of 107 polygon created by the fishnet function in ArcGIS. Now, I need to overlay them to get the mean value of each polygon of that grid from underlying raster image values.

How can I do that using ArcGIS or ENVI?

I tried to use zonal statistics (tabulate) but nothing would be computed as there is no attribute table for the image!


texts: How to get the value of longitude/latitude if I already know one of them for a road on a map of roads contained in a shapefile using java?

texts: I've shapefiles containing a map of roads. I need to get the value of latitude on a given road if I already have the longitude in java. Similarly, the value of longitude if I already have it's latitude. I'm using org.geotools for interpreting the shapefiles in java. Please help.


texts: How to write Python code with QList class

texts:
I need to write an interpolation script for QGIS. I want use TINInterpolation method from Interpolation plugin. I found class in API documentation but I don't know how to write QList object with attributes for interpolation. This is constructor from QGIS API.

QgsTINInterpolator (const QList< LayerData > &inputData, TIN_INTERPOLATION interpolation=Linear, bool showProgressDialog=false)

I cannot find any code snippet only this API documentation - https://qgis.org/api/classQgsTINInterpolator.html#a81c71e2be210a4b082ea807b0626a1bd

and QList documentation from Qt Api - http://doc.qt.io/qt-4.8/qlist.html#details

I tried to find any other code with QList class but I found only QgsVectorLayer and QgsRasterLayer everytime.

In the documentation is not available in Python bindings it means I can't use it with python or what does it mean? Got this from here: https://qgis.org/api/classQgsInterpolator.html

I'll be very glad for any help.

EDIT

In QgsInterpolator.LayerData I don't understand why interpolationAttribute is integer. It is number of column? In Interpolation qgis plugin it is a name of column. It is just next question for discussion. Still don't know how to crete LayerData structure.

I got this error:

Traceback (most recent call last):
  File "", line 1, in
TypeError: arguments did not match any overloaded call:
  QgsInterpolator.LayerData(): too many arguments
  QgsInterpolator.LayerData(QgsInterpolator.LayerData): argument 1 has unexpected type 'int'

when tried this line of code:

lay = QgsInterpolator.LayerData(8, "BREAK_LINES", layer2, True)

I want interpolate vector line data, then use Contour and select only one value from features, this all in script. I started write code to console, but this was first problem which I have. I never use python object programming just procedural. Basics of object programming should I have from Java but I don't understand that documentation. I don't know what means struct reference in python and how should I use it.

EDIT

layer = QgsVectorLayer(pathToFile, 'test', "ogr")
QgsMapLayerRegistry.instance().addMapLayer(layer)
ld1 = qgis.analysis.QgsInterpolator.LayerData()
ld1.vectorLayer=layer              #variable with vector layer
ld1.zCoordInterpolation=False
ld1.interpolationAttribute=7       #number of column with elevation value
ld1.mInputType=2                   #BREAK_LINES

...next two layers and LayerData object (ld2 and ld3)...

#TIN Interpolation
itp=qgis.analysis.QgsTINInterpolator([ld1, ld2, ld3])

rect = layer.extent()
res = 10
ncol = int( ( rect.xMaximum() - rect.xMinimum() ) / res )
test = qgis.analysis.QgsGridFileWriter(itp,pathToSave,rect,ncol,ncol,res,res)
test.writeFile(False)

Until this part I understand that, but when I start QgsTINInterpolator should take a long time for execution but it is too fast. It is weird...

Next code I don't understand, concretely what is zi variable? (line 693 on github)


texts:
I have a csv file with some info and coordinates and I want to create a table in PostGIS with that CSV file. This is what I have so far:

import pandas as pd
import geopandas
from datetime import datetime
#import psycopg2
from geoalchemy2 import Geometry, WKTElement
from sqlalchemy import *

ruta=('csv_file')
df= pd.read_csv(ruta,sep=';')
df['Tiempo'] = df['Fecha'].map(str)+ " "+ df['Hora']

now = datetime.utcnow()
tiempos=[]
registro=[]
datetime_str = df["Tiempo"]

for i in datetime_str:
    datetime_object = datetime.strptime(i,'%Y%m%d %H:%M:%S.%f')
    check = now - datetime_object
    check.total_seconds()
    if (check.total_seconds() > 300):
        registro.append(1)
    else:
        registro.append(0)
df['Ocurrencia']=registro

gdf=geopandas.GeoDataFrame(df, geometry=geopandas.points_from_xy(df.lng, df.lat))

this is the new geodataframe

#database connection
engine = create_engine('postgresql://postgres:postgres@localhost:5432/Rayos')
#conn = psycopg2.connect("dbname='Rayos' user='postgres' host='localhost'
#password='postgres'")

gdf.to_sql('table_name', engine, if_exists='append', index=False,
                     dtype={'geometry': Geometry('POINT', srid= 4326)})

this is the error:

#old error
ValueError: geometry (geometry(POINT,4326)) not a string
#new error
ProgrammingError: (psycopg2.ProgrammingError) can't adapt type 'Point'

I made some changes according to this post.What I'm doing wrong?


texts: Trouble using "Clip Raster by Mask Layer" algorithm QGIS

texts:
I'm having issues using the "Clip Raster by Mask Layer" algorithm to clip a raster by a shapefile. When I run the algorithm, I get an error saying that the mask layer could not be loaded, even though it is already present in the legend. I've tried this after freshly opening it as well.

Though I can't go into the details about the process, I can tell you that a shapefile with the boundary of a field is used to make a 20ft inner buffer. This is then cut away from the original boundary file to make a new boundary that is 20ft smaller on every side than the original. This is then used to cut a very large raster with the Clip raster by mask layer. So it goes:

Boundary (shapefile) --[SAGA:Shapes Buffer (Fixed Distance)]-> Inner Buffer
Boundary + Inner Buffer  --[ QGIS:Symmetrical Difference]-> SymDif (Shapefile)
Raster + SymDif --[GDAL:ClipRasterByMask]-> clippedRaster

Step 3 is where I get the errors. Sometimes it just doesn't write any file (probably going to make another post about that), but in this case, it can't read the output from step 2 (SymDif file).


texts: Show in QGIS the number of overlapping lines in the Label

texts:
Is there a way in QGIS to show in linestring's labels the number of overlapping linestrings (with exact same geometries)?

I've got cables that have the same geometry (are in same pipeline) and I need to know how many of them are there as I only see one linestring...


texts: R fractional area calculations using the exactextractr package

texts: I have a question related to this question. I have a gridded shapefile, which has 256 grid polygons. This shapefile represents an area of interest, and that area of interest also has a raster. I want to know fraction area of the each of the polygons covered by the raster pixels. I have two classes in the raster o and 1. So based on the number of cells that have a value one, the function will return a value telling me that, for example, 60% of grid polygon 1 are covered with raster pixels that have a value one and 40% of grid polygon 1 are covered by raster pixels that have value 0. For such a task, can I use the exactextractr package in R?


texts: Loading QgsComposition from template without throwing "QgsComposition constructor is deprecated"?

texts:
I Googled the deprecation message and the only result was a GIS.se question I asked that referenced the error. Despite this constructor having been deprecated in QGIS 2.4, there's a recent question still using the deprecated construction, and one of the highest voted answers on the topic relies on it as well:

renderer = QgsMapRenderer()
myComposition = QgsComposition(myMapRenderer)
with open(template_composer, 'r') as template:
    new_composer_content = template.read()
    new_composer_document = QDomDocument()
    new_composer_document.setContent(new_composer_content)
composition.loadFromTemplate(new_composer_document)

How should I update this code for PyQGIS versions greater than 2.4?


texts:
I have this table (all row aren't in order):
time, id, latitude, longitude

time is character varying and the format is: dd/mm/year hh:mm
id, latitude and longitude are numeric
I would create lines from latitude and longitude but from id and time. In each row it is possible found same id, but I want associate to the same id (example 20 same id is like 1 id) multiple time (each time have a latitude and longitude). So the line that will be create must respect a specific order. After I create the lines I will simply import the table in QGIS from postgres.

example

line1:
id 12345 ---> time 01/01/2001 00:00 ---> (lat, long)
id 12345 ---> time 01/01/2001 00:02 ---> (lat, long)
id 12345 ---> time 01/01/2001 00:04 ---> (lat, long)

line2:
id 54663 ---> time 01/01/2001 00:24 ---> (lat, long)
id 54663 ---> time 01/01/2001 00:54 ---> (lat, long)

line3
id 234252 ---> time 01/01/2001 00:54 ---> (lat, long)

....

lineN:
id 3223423 ---> ........
id 3223423 ---> ........
id 3223423 ---> ........


texts: Clip on ArcMap not clipping out the correct layer

texts:
In ArcMap 10.8 I have two shapefiles, the yellow area and the green square.
I am trying to create a shapefile that is the yellow area without the green square.

I want to clip the small square out of the large yellow area (so I get a shapefile with a yellow area with a missing square), however- no matter how I input the data onto the "Input" and "Clip" features, it only clips the yellow area out of the square, leaving me with the square which I am trying to get rid of.
Here is the image pre-clip, and the post-clip below and when I clip it with the green box being clipped out of the yellow box, it just gives me the green box again because it's in the yellow area. The yellow area does not get clipped at all.


texts: Checking if geopoint is in shapefile containing world cities

texts:
I have approximately 13 million (lat,long) geocoordinates. I also have a shapefile containing 3460 cities.
I need to find the City in which this geocoordinate can be found.
This is what I am currently doing:
from shapely.geometry import shape,mapping, Point, Polygon, MultiPolygon
import shapefile

with shapefile.Reader(r"shapefile-name") as sf:
        point = Point(x,y)
        for shapeRec in sf:
            polygonalRegion = shape(shapeRec.shape.__geo_interface__)
            if polygonalRegion.contains(point):
                countryCityRegion = [shapeRec.record[2],shapeRec.record[3],shapeRec.record[4]]
                return countryCityRegion

For every unique geopoint I am executing the above process. I then store the result so I don't have to search for the same geopoint more than once.
The issue is that it's incredibly slow. I am unfamiliar with how this process works so I'm completely lost on how to make it faster, if even possible. The shapefile is around 600 mb.


texts: Export DXF from QGIS; points' metadata as text element

texts:
I have a layer (Esri shapefile) which has points in it. Each point has 'property identifier' as metadata. (Boundaries and boundary markers are in different layers / files.)

In QGIS I changed the appearance of layer so that 'property identifier' shows as a text, for example "148-405-1-101". This way: Layer Properties > Style: Font marker > data override with Field type: string, and picked up the right metadata field. This works fine.

Now I'm trying to export/save that layer as DXF, but in DXF points will be points, not text elements. Is there a way to convert points' data into text elements in DXF?

I tried via Project > DXF Export... and Layer > Save as...

I also found this: https://www.qgis.ch/en/projects/dxf-export.
Could it be just encoding problem? Of which I don't know anything...

I'm not a GIS expert so I'm sorry if there's any clumsiness in my explanation.


texts: PostGIS point sample raster doesn't work from raster table but works when view is created

texts:
I have a table of about 90k points that represent yield values in two farm fields, and I am sampling values from a table of rasters. The raster table is a collection of many ndvi rasters, each clipped to a farm field boundary. There is a column 'year' and for each farm field there is only one ndvi raster per year.

When I use this query to sample the raster, I get about 14k points that are not assigned a value:

UPDATE pnt_table
SET ndvi_2014 = ST_Value(rast, geom)
FROM ndvi_table
WHERE ndvi_table.year = 2014 AND ST_Intersects(rast, geom);

Then running the following query confirms that over 14k points assigned 'NaN':

SELECT COUNT (*) FROM pnt_table
WHERE ndvi_2014 = 'nan';

Viewing this raster and the points in QGIS confirms this is not correct, as there are many points that lie over pixels that have values. It is interesting that all the points assigned 'nan' are clustered in one part of the field.

However, by creating a view and sub-setting the rasters according to the two fields I am working in, I can get the expected results:

CREATE OR REPLACE VIEW raster_subset_view AS
SELECT * FROM ndvi_table
WHERE field_name = 'field1' OR field_name = 'field2';

Now re-running the UPDATE line:

UPDATE pnt_table
SET ndvi_2014 = ST_Value(rast, geom)
FROM raster_subset_view
WHERE raster_subset_view.year = 2014 AND ST_Intersects(rast, geom);

Then:

SELECT COUNT (*) FROM pnt_table
WHERE ndvi_2014 = 'nan';

Now all records have been assigned the expected value and there are no 'Nan' values.

This doesn't make any sense to me. It's not a huge deal to just create a view each time I need to do this, but I don't see why I need to. Does anybody have an idea of what might be happening?


texts: How can I verify the z values of my countour line file in ArcGIS/QGIS?

texts:
I have a shapefile with countour lines  and I want to verify if the z value for each countourline  is correct when compared with their neighbors.

How can I do that in ArcGIS or QGIS?


texts:
I have a map a at the link below:

http://gis.xyz/openlayers/openlayer.html

I added one of my GeoServer WMS layers to the map, now I try to add the attributes window pop up on click. I found to examples on the OpenLayers 3 homepage, but I cannot put together a working code...

http://openlayers.org/en/v3.1.1/examples/popup.html?q=popup

http://openlayers.org/en/v3.1.1/examples/getfeatureinfo-tile.html?q=popup

Is there any chance to have a working sample somewhere?


texts: Copying only shapefiles which have specific field attribute in another folder using ArcPy?

texts:
I have 100 shapefiles that all have 2 fields in common called Type2017 and Type2018. The whole of the fields have only one number (1 or 2 or 3). fo example, Type2017  = 1 or Type2018 = 2

I would like to copy those shapefiles which one of these fields have been filled with number 1.

Usually I use this code for copying shapefiles to another folder:

NewFolder = WORK_DIRECTORY + r'/01_SHP'
copyFeatureClass = os.path.join(NewFolder, os.path.basename(Shapefile))
arcpy.CopyFeatures_management(Shapefile, CopyFeatureClass)

But I don’t know how to input the field value that I want .


texts:
I am beginner for widget customization and I am trying to have a custom widget that has fade in and fade out functionality using a slider for a tiled map. I need to have the same functionality as the swipe widget, except it is slider unlike the swipe widget.
Can I identify the basemap and the maplayer (tiled cached layer) so as to apply the fade in out function using the opacity property with dojo/javascript?


texts:
Can I extrapolate a raster with the extent of North America (current geography) to the shape of the continent during the Last Glacial Maximum?

I would want to extrapolate pH-values to the additional coastline areas present in the LGM that currently are under water.


texts: Shortest Path Analysis with Turns in FME?

texts:
Using FME and a road networks dataset, I am trying to find the shortest path from A to B, while taking into consideration of turning width (say for a truck). I have created a point feature class which contains all the turns/junctions, each of these junctions will have a maximum allowed turning width.

How would I incorporate this into my shortest path finder analysis?

Example Scenario: A truck needs to travel on certain roads with certain restrictions, and can only turn at intersections wide enough to accommodate the truck.

I have used FME to select all the roads that satisfy the criteria (i.e. clearance height, maximum weight, shortest length traveled), but how would I incorporate the turns restrictions (turn width) from another feature class?

Below is a screenshot of my FME workbench.


texts:
I am not able to do LOS analysis in QGIS. I have point Features and Raster SRTM Data. I have installed the LOS Plugin in QGIS but not getting any results.

Kindly suggest?


texts: google earth engine: access ImageCollection features list?

texts:
the console gives me a view into the 'properties' of the ImageCollection object:

ImageCollection LANDSAT/LC08/C01/T1_SR (107 elements)
type: ImageCollection
id: LANDSAT/LC08/C01/T1_SR
version: 1524816484446982
bands: []
features: List (107 elements)
properties: Object (23 properties)

and the docs describe methods to access the properties of the properties object.

problem: how do I read the features list? I want the 'id' property of each of the images (which the console can show me). I have tried:

var imageList = myImageCollection.toList(100).features;

and the error is: Cannot read property '0' of undefined

hoping for an undocumented method, I tried:

var imageList = myImageCollection.features();

and the error is: myImageCollection.features is not a function


texts: QField Sync original project path missing

texts: I have a project in QGis, configured for QField. I moved that original project and all the files to a different folder. Now when i want to sync data back from my tablet using the QField-sync plugin, it says the Original path was not found - which is true, the path does not exist anymore. I do not understand why it matters though, i have set the project to relative paths so the data on the tablet should be synced to the new project folder without problems! How can i change the "original path" within Qfield, QGis, QField Sync etc?


texts: I have 2 layers (.shp) point and polygon type. I want to create using QGIS a polygon layer or a raster around the layers' features in a fixed distance (500 meters). I found out how to make a raster around the point type layer but could not for the polygon type layer.


texts: Conversion of elevation from Wusong Datum to Chinese 1985 National Datum to WGS84

texts:
I  want to combine/compare raster data from a digital elevation model (EarthEnv-DEM90) with surface water level data from China. Unfortunately the water level data is given in a vertical Datum with the name Wusong Datum (or Wu Song) while the DEM is given in WGS84.

I may be able to convert the Wusong Datum to the 1985 National Datum (Yellow Sea Datum or EPSG 5736) by subtracting 1.84m. Although I'm not entirely sure that value is correct (because Ive found different statements about this), it should be good enough for the purpose I've in mind.

So my actual Question is:

How do I convert the water level values from EPSG 5736 to water level values that fit into my WGS84 DEM?


texts: Georeferencing CAD Files and Creating Geodatabase

texts:
I'm hoping this is an easy question to answer. I have georefrenced a CAD data set to an orthophoto that was projected in NAD 83, and turned it into a GDB.

Even though the GDB is not projected, is this still considered to be Georefrenced since its based off of the CAD files?

When I bring in a file from the GDB I get a warning saying "The following data sources you added are missing spatial reference information. This data can be drawn in ArcMap, but cannot be projected"

In order to put it into a GDB, I used the select by attributes to only take the certain features I was working with, and exported them into their own Polygon or Polyline file. I then created a GDB and imported the files into it.


texts:
I'm having trouble selecting a feature inside a KML, and wondering if this is even supported.

as a pretty basic example i have:

var vectorSource = new ol.source.Vector({
    //url: 'https://developers.google.com/kml/documentation/KML_Samples.kml',
    url:'https://openlayers.org/en/latest/examples/data/kml/2012-02-10.kml',
      format: new ol.format.KML()
  });

  var vectorLayer = new ol.layer.Vector({
    source: vectorSource
  });

  var map = new ol.Map({
    layers: [
      new ol.layer.Tile({
        source: new ol.source.OSM()
      }),
      vectorLayer
    ],
    target: 'map',
    controls: ol.control.defaults({
      attributionOptions: {
        collapsible: false
      }
    }),
    view: new ol.View({
      center: [876970.8463461736, 5859807.853963373],
      projection: 'EPSG:3857',
      zoom: 10
    })
  });
  var selectClick = new ol.interaction.Select({
    condition: ol.events.click
  });
  map.addInteraction(selectClick);

which you can find running here: https://codepen.io/vishnu4/pen/oQzVVz
you'll notice the select action does nothing, where i have a similar geojson example that works fine : https://codepen.io/vishnu4/pen/xQEMjy .

Is this as expected?  Or is there a workaround to get this to work so i can select the specific features (placemarks?) inside the KML object?

Just as an aside, I know that getting a feature using forEachFeatureAtPixel works ( https://codepen.io/vishnu4/pen/YRGbyR ) , but I'm looking to use the Select method.


texts:
How can i symbolize a map dynamically by data stored in mysql database in many tables (these tables have relationship through Primary Key and FK). For instance "Disease Name"(chosen from drop down then "Affected Number" (from another drop down) "Sex" (to chosen from another drop down) then "Year" from another drop down box, after all an ok button to be hit to produce intended map. I am using geoserver 2.2  postgresql( postgis). I am try to outline the process but ended in smoke (my first design and so so), would any one ,as samaritan ,help me, plz. Thanks in advance.

my project is as in https://docs.google.com/open?id=0B1NvUoakjNxwZ0VGSm15bE9fNkE


texts:
The introduction of Bing Imagery with the release of ArcGIS 10.0 provided a fantastic means to very quickly produce Site Location Maps.  Prior to its' introduction, I was using a GIS server connection to www.geographynetwork.ca_80 for site location maps.  I highly preferred using Bing, mostly for aesthetic reasons.

Now that Bing no longer available (without buying a licence), I wondered what others are using as a replacement.  I find the options available under the "add basemap" dropdown in ArcMap to be far from appealing.  Often, I find the placement of the citation information to be rather inconvenient.  The resolution is often very poor as well (such as when using the "Open Street Map" or "Imagery options").  The labels are often fuzzy and hard to read as well.

What are those of you who don't have a Bing subscription using?

Below, Is a picture to demonstrate what I meant by "Site Location Map".  In many cases - not always - it is faster and simpler to just use an air photo or satellite image to show where something is.  Point to the study area, label a few roads .... done!!


texts: Only selecting areas not overlapped by another layer in ArcMap?

texts: How do I create a layer that displays only the bright red areas on the map (The areas of the red layer that are not overlapped by the purple or yellow layer)?


texts: I have a *.dgb file on landcover (vector format). I am however battling with trying to edit this file. None of the edit tools are active. I have changed the file to not read only.


texts:
I have 5 different layers. I would like to change the order they are drawn. I read that (and I also used to) simply select and drag the layer from in the position I want. Nevertheless, now I am trying to do it again and when I try it seems impossible to me to drop the layer anywhere. My layers come from the same geodatabase, as shown in

Why is it not working?


texts: Iterating over column from CSV using pandas?

texts:
The problem that I am experiencing is as following:

I have a csv with the following columns: 'time' (with date and time), 'id', 'lat', and 'long'. The data shows movements and id represents a mobile device.
I need to calculate the distance and the velocity of an activity. Activity is movement that lasted above 5 minutes (or so, a conditional that I added). Some IDs have multiple movements in a single day so I need to separate first by ID and than by time.
So, far, I have managed to apply several things and I am not asking for someone to do this for me. I am a newbie and I want to learn on examples. So I am trying to do this on my own. However, developing an algorithm is really hard. I know what I want but don't know how to put it.

import geopy.distance
import pandas as pd
import math

#import csv and define dtypes and headers
headers = ['time', 'id', 'lat', 'lon']
dtypes = {'time': 'str', 'id': 'str', 'lat': 'float', 'lon': 'float'}
parse_dates = ['time']
df = pd.read_csv('C:\Posao\Python\People_Movement\python_velocity.csv', header=None, names=headers, dtype=dtypes, parse_dates=parse_dates)

df = df.sort_values(by=['id', 'time'])

# Group the sorted dataframe by ID, and grab the initial value for lat, lon, and time.
df['lat0'] = df.groupby('id')['lat'].transform(lambda x: x.iat[0])
df['lon0'] = df.groupby('id')['lon'].transform(lambda x: x.iat[0])
df['t0'] = df.groupby('id')['time'].transform(lambda x: x.iat[0])

#Calculating the distance between two points
def haversine(lon1, lat1, lon2, lat2):
# convert decimal degrees to radians
  lon1, lat1, lon2, lat2 = map(math.radians, [lon1, lat1, lon2, lat2])

# haversine formula
dlon = lon2 - lon1
dlat = lat2 - lat1
a = math.sin(dlat/2)**2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon/2)**2
c = 2 * math.asin(math.sqrt(a))
r = 6371 # Radius of earth in kilometers. Use 3956 for miles
return c * r

lat_1 = df.lat0[0]
lon_1 = df.lon0[0]
lat_2 = df.lat0[1]
lon_2 = df.lat0[1]
#while df.id[0] == df.id[::1] and (df.to[1] df.to[0] )
    #distance + distance
    #return final_distance
    #lat_1+=1
    #lon_1+=1
####this is still in work

###while lat i lon exist:
lat_1 +=1
lon_1 +=1
return haversine

#Calculating the time different between two points
t = df.t0[0]
def calc_minutes(time_end, time_start):
    duration = (time_end - time_start).seconds
    minutes = duration/60
    return minutes
time start = df.t0[0]
time_end =

###while (lat0==(lat0+=1)) and (t0==(t0+=1)):
                           calculate haversine
                           #calculate time_difference
                           #calculate velocity
                           df.id+=1
                           df.t0+=1
                           df.lat+=1
                           df.lon+=1
                           df.lat0+=1'''

### Wherever I added this is an idea of mine. It still makes no sense.

How can I iterate over a column? For example, as long as id is equal to the id in the row after it, keep calculating and adding distance between lat and long in the corresponding rows. And I need to add another conditional besides the id and that is as long as id is equal to the id in the following row and the time difference between the two rows is less than 5 minutes, keep adding.


texts: Warning 6 while loading data to MSSQLSpatial using ogr2ogr

texts:
I am loading layers from a .gdb to a database in SQL Server using the following command:

ogr2ogr -f MSSQLSpatial "MSSQL:server=server_name;database=database_name;trusted_connection=yes" "gdbsourcelocation" -a_srs "EPSG:4326" -lco "GEOM_TYPE=geography" -lco "GEOM_NAME=geog4326" -progress

This fails with the following warning:

Warning 6: Linestring has no points. Removing the geometry from the output.

Even when I use -skipfailures it fails.

Any help will be appreciated.


texts:
I am new to GIS and trying to understand how to edit an OpenDRIVE map.
How do I open and edit an OpenDRIVE map in ArcGIS Pro?


texts:
I have registered my layer as "Registered as  versioned" with the "Register the selected objects with the option to move edits to base".

I can publish the the data the feature server with "Sync" features on.

Data can be used online, but when trying to "download the map" for "offline" use.

The application logs says :

Error Domain=NSCocoaErrorDomain Code=400 "Failed to create replica." UserInfo={NSURL=http://52.28.226.146/arcgis/rest/services/yyy/yyyy/FeatureServer/jobs/j57ac7b5825c446e0a9e29de2faee7700?f=json&token=QE3SFRWbxkrxEbwfgxz3w-54a8amIH1Rv7rcIAnrfGv1mM5C8pbi1QVr0rB_KjFC, responseJSON={
    error =     {
        code = 400;
        description = "Failed to create replica.";
    };
    lastUpdatedTime = 1458648838000;
    replicaName = "";
    resultUrl = "";
    status = Failed;
    submissionTime = 1458648835000;
}, NSLocalizedFailureReason=Failed to create replica.}

And when checked "Arcgis server Manager" logs 2 error occurs :

Error 1

  Error executing tool.:
  ErrorMsg@SyncGPService:{"code":400,"description":"Failed to create
  replica."} Failed to execute (Create Feature Service Replica).

Error 2

  Error: An error occurred while adding an entry to the
  FEATUREDATASET_CLASSES or the FEATURECLASSES table. [Read only data or
  unversioned data cannot be replicated. Data versioned with the option
  to move edits to base cannot be replicated. GlobalIDs are required for
  two way and one way replica data. Creating one way replica using
  archiving requires archiving the replicated parent data.].

The achivement goal is :

Publish a feature server to use on mobile for online/offline editing with versioned data and all changes and updates shown on service immediately.

What i have already tried :

Read the prepare offline map and feature documentation and made sure
that i had Global ID
When i do not use "Register the selected objects
with the option to move edits to base" for data versioning option.
Data is not displaying on the server when sync in offline mods even if i stop the service and tried to compress the data to DB.

Why we need versioning :

To be able to change and edit on Arcmap and webservices for multiple users.

Using : Arcgis 10.3


texts:
Can you edit an ArcSDE Oracle spatial view?

By 'edit' I mean edit the geometry in the feature class, not any attributes in the joined table. Unversioned.

I really had it in my head that you can, as long as the objectid column is included in the view along with the shape column, some half-remembered reference. So I was forging ahead and decided to double-check if there was fine print, and now it seems I was whackadoodle in the first place.

I've seen nothing to suggest it is supported, and seen in Administering ArcSDE for Oracle, page 225: "Spatial
views cannot be edited through an ArcSDE client."

But then I went ahead and tried it today, and I could edit a spatial view in ArcMap, but given the above I question how safe it is.


texts: How to determine whether a map click event has been set?

texts:
The Evented section of the Leaflet documentation shows how to set an event listener for the map object:

  Leaflet deals with event listeners by reference, so if you want to add
  a listener and then remove it, define it as a function:

function onClick(e) { ... }
map.on('click', onClick);
map.off('click', onClick);

In my application, there are various functions which may run on a map click event, depending on which mode the application is in (along the lines of having an "active tool" mode).

How can I determine which function is currently set as the map's on click event?

For example, if the line map.on('click', xyz) has been run, I would like to know that the current on click function is called xyz. And if map.off('click', xyz) has been run and there is no current click event I would like to know this.


texts: Given distances from two locations how to accurately find geographic coordinates

texts:
I am working on a research problem which requires me to calculate a geographic location accurately (with error distance of few meters), given distances from two geographic locations. For instance:

Since the points are static I know their location already. Flipping of the triangle is avoided by adding one more known location. I have calculated the distances using great circle distance formulae.

double theta = lon1 - lon2;
double dist = Math.sin(deg2rad(lat1)) * Math.sin(deg2rad(lat2))
            + Math.cos(deg2rad(lat1)) * Math.cos(deg2rad(lat2)) * Math.cos(deg2rad(theta));
dist = Math.acos(dist);
dist = rad2deg(dist);
dist = dist * 60 * 1.1515; // in miles

How to find the coordinates of C?
The correct answer is 39.98, -83.03, but I am getting 8 miles off the track.
Am I doing anything wrong?


texts: QGIS 2.12 attribute table update question

texts:
I'm trying to fill in a whole lot of cells in the same column of an attribute table all at once using the field calculator. Am I on the right track with the following expression?

"Author"  || example text

I have tried " ", ' ' and ( ), none of which seem to get the result I'm hoping for.


texts: Changing thematic symbology based on an atlas feature reference

texts:
I have created a layout atlas based on a geographic layer containing defined store 'catchments' (ie as polygons).

Each atlas 'page' shows just the relevant catchment region.

For each atlas feature, I  also want to customise a thematic representation of each store's market share by postcode. This runs from a separate layer which is based on postcode regions and has a field for each store showing the market share within each postcode.

I am trying to reference the field name (in the picture below the field name is 'Albury_mkt') dynamically as a concatenation of the atlas page name ('Albury') and the text string "_mkt".

If successful, each atlas page will show the thematic only for the relevant region.

I have tried to use the function  concatenate(@atlas_pagename,'_mkt') but this doesn't seem to parse properly.


texts: calculating HRMS using NMEA sentences including GPGGA , GPGSV, GLGSV, QZGSV,QZGSA, GPGSA?

texts:
I have NMEA sentences including  GPGGA , GPGSV, GLGSV, QZGSV,QZGSA, GPGSA and GNGSA. my goal is to achieve HRMS accuracy meter. produced android devices for special intent gives GP/GN GST nmea message which includes latitude and longitude error meters but android device that I have is not produced for special aim and dosent give GP/GN GST nmea message.  I dont know How I am able to calcute HRMS value as meter. On the other hand, when Location is changed Android locationlibrary provide accuracy as meters. But I would like to calculate Horizontal accuracy using nmea messages. is there any way to calculate HRMS value using messages that I mentioned above of this writing.

nmea messages
$GPGSA
GPS DOP and active satellites

eg1. $GPGSA,A,3,,,,,,16,18,,22,24,,,3.6,2.1,2.2*3C
eg2. $GPGSA,A,3,19,28,14,18,27,22,31,39,,,,,1.7,1.0,1.3*35

1    = Mode:
       M=Manual, forced to operate in 2D or 3D
       A=Automatic, 3D/2D
2    = Mode:
       1=Fix not available
       2=2D
       3=3D
3-14 = IDs of SVs used in position fix (null for unused fields)
15   = PDOP
16   = HDOP
17   = VDOP

http://aprs.gids.nl/nmea/


texts: How to build Geometric network in MapInfo?

texts: I have sample dataset of Water network. I wanted to make it geometrically connected network. I was trying to trace the network end to end, I have used "FnodeTnode tracing" plugin, but I was not lucky enough to achieve tracing. Do I need to change/ Add the Columns in the TAB files to make that plugin work? Is there any Extension that can help me in this? Is there any Mapbasic script that I can use to create the Geometric network connection and helps e to run the Upstream and downstream tracing?


texts: How to compute a line intersection point in a street layer in PostGIS 1.4?

texts:
I have a street layer downloaded from OSM. I need to calculate intersection point as we do in QGIS from vector menu -> Analysis tool -> Line intersections.

I am new to PostGIS/PostgreSQL, I managed to download and install it spending whole day on Ubuntu machine. Using shp2psql-gui I imported data with pgAdmin3.

I know how to do in QGIS, but tell me how to do the same in PostGIS.  I have heard some st_intersects etc. but it can't help. I want to return something like this.

Below are snapshot for more clearer picture of my issue. Blue indicates street layer.

name_1 | name_2 | laitude | longitude
xyz st | abc st | #value| #value
lmn st | pqr st | #value| #value
xyz st | pqr st | #value| #value
pqr st | abc st | #value| #value
lmn st | xyz st | #value| #value
lmn st | abc st | #value| #value


texts: Can´t categorize field for unknown reason

texts:
I am trying to categorize a field without any success

The values should be symbolized from 1-5 and yet, I just get that blank symbol.
at first I thought that the problem is that the column is defined as string, so I created a new number column and after updating to the original, I got an NULL column.
Iv'e created another string column and again I got a blank column .
A weird thing is that in the field calculators output preview, the original value appears as it should

I'm using QGIS 3.4.7 and I have the same issue in other versions as well

here's the way it  looks in the field calculator