texts: Procedure for converting NAD83 to NAD83 1992 epoch (1991.35)?

texts:
I am attempting to convert a shapefile from NAD83 HARN / egm96 vertical to NAD83 1992 epoch (1991.35) / NAVD 88 vertical. I have attempted to use NADCON to no avail. I know I can use the Define Projection tool in ArcMap to convert to NAVD 88, but am stumped on transforming to NAD83 (92).

Is anyone aware of a utility that can facilitate the conversion to this specific, rarely used datum?


texts:
The project I am working on currently contains 10+ boolean features and I would like to create an interactive checkbox SQL filter where only the features checked are applied as AND filters. For example on the dashboard,
[✓] buildings built in 1901 and 1930
[✓] residential use buildings
[ ] buildings with more than 3 units
Would yield a filter for buildings built in 1901 and 1930 AND residential use buildings (but not buildings with more than 3 units). Currently my code yields both filters even with one checkbox marked.


texts: Can I find the area of use for a given projection in GDAL?

texts:
We make an application (using Python and open source libraries like GDAL) that stores raster data from around the world. We wrote some code that acts as a WMS server, that can create images from the raster data in various projections.

Now in the GetCapabilities response we want to list some projections that can be used with a given raster. We only want to show projections for which the data falls within the right area of use, e.g. we only want to show EPSG:28992 if the raster falls entirely within the Netherlands. But we don't want to hardcode the area of use for each projection we may want to show.

Is there a way to find out what the official area of use is for any projection that GDAL can handle? E.g. to find out what the bounding box of the area of use of EPSG:28992 is, preferably in terms of its own coordinates.


texts: Digitizing polygons by cut or split feature tool in QGIS

texts:
For saving digitizing time i created one polygon feature cover all satellite imagery (study area) and styled it outline, tried many tools inside qgis 2.18 to split / cut polygon features (i don't mean split polygon start from edge polygon feature) but digitize building inside polygon by split it.

I did this task by ArcGIS for Desktop!.


texts: Google Earth Engine: Multiply an image collection by a different number depending on the month

texts:
I am trying to multiply a monthly imageCollection by a different number depending on the month. What I want is to multiply all the pixels of the first 12 months by 0.85, the following 12 months by 0.95, the following 12 by 0.65. All this for a total of 468 months from 1981 to 2019 (the example contains only 2 years, 24 months).

var monthly_sum = ee.List.sequence(0, 2*12).map(function(n) { // .sequence: number of years from starting year to present
  var start = ee.Date('1981-01-01').advance(n, 'month'); // Starting date
  var end = start.advance(1, 'month'); // Step by each iteration

  return ee.ImageCollection("UCSB-CHG/CHIRPS/DAILY")
        .filterDate(start, end)
        .sum()
        .set('system:time_start', start.millis());
});

print(monthly_sum)

var my_list = ee.List([]);

var otherList = ee.List([0.884748080310752,
0.884748080310752,
0.884748080310752,
0.884748080310752,
0.884748080310752,
0.884748080310752,
0.884748080310752,
0.884748080310752,
0.884748080310752,
0.884748080310752,
0.884748080310752,
0.95,
0.95,
0.95,
0.95,
0.95,
0.95,
0.95,
0.95,
0.95,
0.95,
0.95,
0.95
])

Which would be the easiest way to multiply each monthly image by the corresponding multiplicative factor creating a new imageCollection?

Thanks,


texts: How to import CSV with addresses (no coordinates) into QGIS?

texts:
I have a CSV file with many address and I need to find the coordinates.
I think to import the CSV into google-mymap hoping that it find those points, then I would like to export a KMZ to import it in QGIS?

Am I wrong in this thought? If this is the only way to solve it, how can I find info to create this CSV file?

I imported the file into mymaps and I can see the points on the map.
Then I exported it in KML/KMZ.

When I go to layer/add layer, QGIS import only a table (with the orginal info from CSV file) not a point shape file.

In the picture below you can see information from my imported addresses KMZ (this is not a geometry point shape file) and a correct test create into mymaps.

I have to add that into mymaps app the dataset looks identical.


texts: Is it possible to integrate a route planner into QGIS?

texts:
I am looking for an open source GIS and now asking myself if it is possible to integrate a route planner into QGIS.

The route planner should be able to calculate the distance from one position in the map in all directions (e.g. 30 minutes driving time) and mark it as a catchment area. It would also be helpful to export the data of the marked area to excel (e.g. name of the administrative districts and also demographic date).


texts: Is it possible to use PostGIS on gisgraphy database?

texts:
I have the entire gisgraphy database inside my server. Is it possible to use PostGIS on it to query such as "get locations within 10 miles from this coordinate"?

I tried reading about spatial queries and PostGIS but I am not sure if I could directly use the tables in gisgraphy to obtain results.


texts: Assign points within 150 meter buffer of each other the same unique identifier?

texts:
I have a bunch of points on a QGIS map. Each point is a unique sampling event. However, some of these points have been sampled multiple times but there is no common identifier in the attribute data to determine this. The only way we can see this is by looking at the map and seeing points very close to each other.

I am wondering if anyone has a method to assign points within a 150 meters of each other a common identifier?


texts: Size based scaling for symbology in StoryMaps journal

texts:
Is there a way to have size based scaling for symbology StoryMaps journal?
Scale based sizing of symbols works in ArcGIS Pro, and web map has been uploaded to my organization, however when viewing in map view or storymaps symbols stay the same size


texts: Find country from GPS coordinates, offline using python

texts:
I'm aware that there are several similar questions already but I have been searching for days and not been able to solve my problem:

All solutions I've found require GDAL, which I have not been able to install (even doing all the environment variables tweaking etc.) except by using OSGEO4W, which allows me to use osgeo, but only in its own shell, rather than in the shell I use for everything else.

I have the country shapefiles, I have the shapely package. Is there a way to do this without GDAL? (or as an alternative, an easy way to get GDAL packages in my own shell?)


texts:
For a certain project I'm required to run QGIS processes on a large data set on a Linux server as my PC does not have the sufficient memory to execute them (snapping points to lines + Steiner tree extraction using national road data).

Is it possible to run QGIS without the aid of its GUI only using the Command Line Interface?


texts: GeoExt WMS GetFeatureInfo in proxy server. Problem - popup is empty in chrome and not responding in firefox?

texts:
I have python2.7 and apache installed in my computer.
In proxy.cgi, I included in topline :

  #!/usr/bin/env python.exe -u;

I also added 'localhost:8080/geoserver' and 'localhost:80 in allowedhosts.
I moved the OpenLayers files to htdocs folder of apache.
I copied the cgi file to apche/cgi-bin.
Now, when I click on feature, popup appears but is empty.

Please help!! where my code went wrong?
Given below is my code:

    <script type="text/javascript" src="http://localhost:8080/recipes/resources/ext/ext-all.js"></script>
    <link rel="stylesheet" type="text/css" href="http://localhost:8080/recipes/resources/ext/resources/css/ext-all.css" />
    <link rel="stylesheet" type="text/css" href="http://localhost:8080/recipes/resources/ext/resources/css/xtheme-gray.css" />
    <link rel="stylesheet" href="http://localhost:8080/recipes/resources/openlayers/theme/default/style.css" type="text/css" />
    <link rel="stylesheet" href="http://localhost:8080/recipes/resources/app/theme/style.css" type="text/css" />
    <link rel="stylesheet" href="http://localhost:8080/recipes/resources/geoext/resources/css/gxtheme-gray.css" type="text/css" />
    <script src="http://localhost:8080/recipes/resources/openlayers/build/OpenLayers.js"></script>
    <script src="http://localhost:8080/recipes/resources/geoext/build/GeoExt.js"></script>
    <script>
        OpenLayers.ImgPath = "http://localhost:8080/recipes/resources/openlayers/img/";

        // the `panel` and `popup` variables are declared here for easy debugging
        var panel, popup;

        Ext.onReady(function() {
        OpenLayers.ProxyHost = "/cgi-bin/proxy.cgi?url=";
            // queryable WMS layer
            var zoning = new OpenLayers.Layer.WMS(
                "railroads",
                "http://localhost:8080/geoserver/usa/wms",
                {layers: "usa:states", format: "image/gif", transparent: "true"}
            );

            panel = new GeoExt.MapPanel({
                title: "MapPanel",
                renderTo: "map-id",
                height: 300,
                width: "100%",
                layers: [
                    new OpenLayers.Layer.WMS(
                        "Global Imagery",
                        "http://maps.opengeo.org/geowebcache/service/wms",
                        {layers: "openstreetmap", format: "image/png"}
                    ),
                    zoning
                ],
                center: [-122.87, 42.34],
                zoom: 13
            });

            // create a control to get feature info from queryable layers
            var control = new OpenLayers.Control.WMSGetFeatureInfo({
                url: "http://localhost:8080/geoserver/usa/wms"
                title: 'Identify features by clicking',
                layers: [zoning],
            });
            panel.map.addControl(control);
            control.activate();

            // register a listener for the getfeatureinfo event on the control
            control.events.on({
                getfeatureinfo: function(event) {
                    // close existing popup
                    if (popup) {
                        popup.destroy();
                    }
                    popup = new GeoExt.Popup({
                        title: "Popup",
                        map: panel.map,
                        lonlat: panel.map.getLonLatFromPixel(event.xy),
                        width: 250,
                        autoScroll: true,
                        collapsible: true,
                        bodyStyle: {padding: 5},
                        html: event.text
                    });
                    popup.show();
                }
            });

        });

    </script>

After changing the path to wamp/www and still I am getting this new following error:


texts:
I am trying to append ~10,000 tables together using an ArcPy script.

I believe the method I am using does not work as the script is attempting to append too many tables at the same time. Would a possible work-around be a script that appends using looping?

import os
import sys
import string
import arcpy
from arcpy.sa import *
arcpy.CheckOutExtension("spatial")

WD = "C:/workspace"
arcpy.env.workspace = WD

#selects the 10,000 tables and tries to append them all together at the same time

arcpy.CreateTable_management(out_path=WD, out_name="append", template="C:/Workspace/template", config_keyword="")
tables1 = [table for table in arcpy.ListTables() if not table.startswith('append')]
arcpy.Append_management(inputs=tables1, target="ndvi_append", schema_type="TEST", field_mapping="", subtype="")


texts: Geodatabase feature class not drawing in inset maps, but draws fine in main map

texts:
I am working with some census data, identifying minority populated blocks (shown in blue), in a geodatabase.  The blocks draw fine in my main map but do not draw in either inset. They can be selected (highlighted polygons in both insets) from the both inset Data frames, as shown in the screen capture below.

Neither inset has any scale restrictions set to the data frame.

The other weird thing about this situation is when I turn the basemap on or off in either inset (but not my main data frame, I get the following error:

I am sure this is related to the feature class is not rendering, but I have no idea why or how to fix it.


texts:
I have a question about the style files used for osm2pgsql.

These files have a structure as follows:

# OsmType  Tag          DataType     Flags
node,way   access       text         linear
node,way   addr:housename      text  linear
node,way   addr:housenumber    text  linear
node,way   addr:interpolation  text  linear
node,way   admin_level  text         linear
node,way   aerialway    text         linear
node,way   aeroway      text         polygon
node,way   amenity      text         polygon

My question is this, how do I know if I should use polygon, linear, hstore, etc.. for the Flags values?  The style file gives a brief explanation, but I am looking for further information.


texts: Javascript Function for transforming WGS84 to EPSG:5514

texts:
I'am searching for JavaScript function which converts WGS84 coordinates to this special coordinate format EPSG:5514 https://epsg.io/5514.

Function should looks like function wgs2sjtsk(lat,long) and output will be coordinates in  EPSG:5514 (sjtsk).


texts: Getting pixel values at single point using rasterio

texts:
To get a single pixel value at a point in a raster using rasterio, there is an example here: https://github.com/mapbox/rasterio/pull/275

However, is there a direct API within rasterio (and not the cli) which can be used to extract value at a single point in a raster?

-- EDIT

with rasterio.drivers():

    # Read raster bands directly to Numpy arrays.
    #
    with rasterio.open('C:\\Users\\rit\\38ERP.tif') as src:
        x = (src.bounds.left + src.bounds.right) / 2.0
        y = (src.bounds.bottom + src.bounds.top) / 2.0

        vals = src.sample((x, y))
        for val in vals:
            print list(val)


texts: I have to restart GeoServer when I change the mapping file?

texts: I'm trying to use GeoServer with app-schema plug-in to map a complex features. I have some trouble but my question is: I have to restart GeoServer every time that I modify the mapping file so that the wfs request work properly?


texts:
I'm new to GIS and wanted to start play around with some data about the appalachian trail (AT) and PostGIS.

However I already failed at the first and most simple task to calculate the length of the AT. I downloaded the data here. And imported it using the pgShapefileLoader. After loading I noticed that the srid of at_centerline wasn't set. Inspecting the table with select st_astext(geom) from at_centerline limit 1 showed that it was in lat/lon format. Therefore I updated the table with select updateGeometrySRID('at_centerline','geom',4326).

So, time to calculate the length

select sum( st_length( geography( geom ) ) ) / 1000 * 0.621371
from at_centerline

yielded 2122.61779398273 miles. However the AT currently is 2190 miles long.

I also tried to directly calculate the length in miles with this query

select sum( st_length( st_transform( geom, 2877 ) ) ) / 5280
from at_centerline;

which yielded 2126.50802521765 miles.

So how to explain the discrepancy? Is it simple the data provided by the ATC isn't accurate enough? I would rather doubt that. So what am I doing wrong here?


texts: I have several rasters that I want to merge but after selecting the layers, in the raster > misc > merge tool, the panel window just closes without any error and of course no new merged raster has been created. I tried with both qgis3.10 and 3.12 and the SAGA Mosaic raster layers as well. With no success, nothing happens.


texts:
I know that this question has already asked but no concrete answer has given!!
My environment is Ubuntu 12.04 and I'm using eclipse-pydev as IDE.


texts: Using script parameter validation in ArcPy

texts:
I am developing a script to Add, Delete, and Update rows within a table. The script basic function works well using InsertCursor and UpdateCursor. My questions or troubles come when I attempt to use script validation to auto populate some of the parameters.

When the users choose to Update a Row, they pick a country from a list. Then the validation prepoulates 5 parameters from the existing row in the table. I am using SearchCursor within script validation to auto populate those fields. That generally works well.

The problem is the script validation never seems to stop. So after the user picks a country code, 5 parameters are populated in the tool dialog with data from the fields in the table. Ideally I would like the user to be able to type over those entries and then resubmit to the table for changes. However every time the user changes one of the parameters the script validation overwrites it with the original value from the table. Does anyone have any suggestions on how I can make the parameters populate with the field data in the table, and still allow the user to be able to manually change those parameters for update purposes?

Here is a code snippet:

if self.params[1].value == True or self.params[2].value == True:

  self.params[4].filter.list = [str(val) for val in sorted(set(row[0] for row in arcpy.da.SearchCursor(Table,Field)))]

  wClause = "COUNTRY_CODE ='{}'".format(self.params[4].value)

  # SO FAR EVERYTHING WORKS FINE

  # HERE IS THE VALIDATION THAT USES SEARCHCURSOR AND IT WORKS

  if self.params[4].altered:

    values = next(arcpy.da.SearchCursor(Table, Fields, where_clause=wClause))

    self.params[5].value = values[0]

    self.params[6].value = values[1]

    self.params[7].value = values[2]

    self.params[8].value = values[3]

    self.params[9].value = values[4]

    self.params[10].value = values[5]

So after the SearchCursor populates the params, how can I stop the validation and let the user manually edit those params?


texts: Convert feature class to raster based on attribute field values

texts:
As part of a tool, I’m trying to convert a feature class to a raster based on attribute field values.  I have accomplished this conversion using IRasterConvertHelper, but this is done without specifying an attribute field.  So I then was able to do the conversion using the IGeoProcessor to execute ConversionTools.FeatureToRaster , but was only able to use numeric fields (i.e., tool didn’t work when a text field was selected).

Going another route, I tried dissolving the feature class based on a selected attribute field using IBasicGeoprocessor and Dissolve, but this fails with:
{System.Runtime.InteropServices.COMException (0x80004005): Error HRESULT E_FAIL has been returned from a call to a COM component.  (thanks ESRI, lots of help!)

Does anyone have any idea on how to accomplish the same thing as the ArcToolbox Conversion Tools Feature to Raster tool and convert a feature class into a raster based on an attribute field?


texts: Summing Values for Each Item in a Dictionary

texts:
# Create dictionary with segmentid as key, calculated length as value1, measured length as value2
fiberCableLengths = defaultdict(list)
with arcpy.da.SearchCursor(fiberCable, ['segmentid', 'calculatedlength', 'measuredlength']) as cursor:
    for row in cursor:
        fiberCableLengths[row[0]].append(row[1])
        fiberCableLengths[row[0]].append(row[2])

# Create slackloop dictionary with segmentidfkey as key, total up the length values of each slackloop under that segmentidfkey, and store as the dictionary value
slackloopLengths = defaultdict(list)
with arcpy.da.SearchCursor(slackLoop, ['segmentidfkey', 'length']) as cursor:
    for row in cursor:
        if row[1] is not None:
            slackloopLengths[row[0]].append(int(row[1]))

sumValueDict = {k: sum(v) for k, v in slackloopLengths.items()}

# For each fibercable segmentid in dictionary, If calculated length + total slackloop length = measured length – Add to good list, Else – Add to bad list
goodList = []
badList = []
with arcpy.da.SearchCursor(fiberCable, ['segmentid', 'calculatedlength', 'measuredlength']) as cursor:
    for row in cursor:
        for k, v in sumValueDict.items():
            if row[0] == k:
                if row[1] and row[2] is not None:
                    if int(row[1]) + int(v) == int(row[2]):
                        goodList.append(str(row[0]))
                        print "Good List: " + row[0]
                        print int(row[1]) + int(v)
                        print int(row[2])
                    else:
                        badList.append(str(row[0]))
                        print "Bad List: " + row[0]
                        print int(row[1]) + int(v)
                        print int(row[2])
print "Good List Size: " + str(len(goodList))
print "Bad List Size: " + str(len(badList))

The code above creates a dictionary where each item looks like the following:
'1f169eae-5d7c-4ebb-90cd-abbe0e1b72dd': [50.0, 50.0, 50.0, 50.0]

How do I get it to sum the values for each key? The output I'd be looking for here would be:
4ebb-90cd-abbe0e1b72dd': [200]

Resources I've looked at:

https://stackoverflow.com/questions/4880960/how-to-sum-all-the-values-in-a-dictionary/4881100
https://stackoverflow.com/questions/34849651/how-to-add-the-values-of-dictionary-with-same-keys-in-python
https://stackoverflow.com/questions/34849651/how-to-add-the-values-of-dictionary-with-same-keys-in-python


texts: Map render complete event is not triggered - Openlayers 5

texts:
I have written a function in Map render complete event to hide the progress icon if all map tiles are rendered. But the map render complete event is not fired after I added a empty vector layer on map. Below is the code snippet.

var osmLayer = new ol.layer.Tile({
    source: new ol.source.OSM({
        "url": "http://{a-c}.basemaps.cartocdn.com/light_all/{z}/{x}/{y}.png"
    })
});
osmLayer.getSource().on('tileloadstart', renderStart);

var extentTransform = ol.proj.transformExtent([minx, miny, maxx, maxy], "EPSG:7844", "EPSG:3857");
osMap = new ol.Map({
    layers: [osmLayer],
    target: 'map',
    view: new ol.View({
        zoom: 11,
        maxZoom: 20,
        minZoom: 5,
        extent: extentTransform,
        center: ol.extent.getCenter(extentTransform),
        controls: []
    })
});
osMap.on('rendercomplete', renderComplete);

var putStyle = new ol.style.Style({
    stroke: new ol.style.Stroke({
        color: '#900C3F',
        width: 3
    }),
    fill: new ol.style.Fill({
        color: 'rgba(144, 12, 63, 0.5)'
    })
});

vectorSource = new ol.source.Vector({
    features: []
});

var vectorPolygon = new ol.layer.Vector({
    source: vectorSource,
    style: putStyle
});

osMap.addLayer(vectorPolygon);

function renderComplete(event) {
    $(".MapLoadPopup").hide(); // Tile loading progress hide
}

function renderStart(evt) {
    $(".MapLoadPopup").show(); // Tile loading progress show
}

function test(){
    vectorSource.clear();
    var geojsonObject = {
        'type': 'Feature',
        'geometry': {
            "type": "Polygon",
            "coordinates": JSON.parse(polygon),
        }
    };

    var features = new ol.format.GeoJSON().readFeatures(geojsonObject);
    vectorSource.addFeatures(features);
  }

If I commented the osMap.addLayer(vectorPolygon); line, the render complete event is triggered correctly. When sourcing this issue in internet, it is came to know that the map render complete event is not fired if we added a vector layer without url.

Myself also tested this case with some url and it is working fine. If so, what could be the solution to overcome this issue.


texts: Confusion matrix from supervised learning

texts:
I am performing a supervised classification in GEE and want to verify how well my classifier has got trained.
For this, I found a method confusionMatrix() and my respective program is

var training = india_image.sampleRegions(ft,['class'],30);
var trained = ee.Classifier.cart().train(training, 'class', new_bands);

var trainAccuracy = trained.confusionMatrix();

var exportAccuracy = ee.Feature(null, {matrix: trainAccuracy.accuracy()});

// Export the FeatureCollection.
Export.table.toDrive({
  collection: ee.FeatureCollection(exportAccuracy),
  description: 'exportAccuracy',
  fileFormat: 'CSV'
});

The exported csv is here

Can you please explain it or give some suggestion about how to find presision, recall and accuracy


texts: Creating a multiring buffer in ArcObjects

texts:
I am currently writing an an addin for ArcGIS in C#. The premise of the application is, the user will select a polygon, a multi ring buffer is created from the centroid of the polygon. The buffer distances in this case being 250m and 500m.

I have found a bit of code which works, in that the multiring buffer is created

            IBufferConstructionProperties bcp = new BufferConstructionClass();
            ESRI.ArcGIS.esriSystem.IDoubleArray dis = new ESRI.ArcGIS.esriSystem.DoubleArrayClass();
            dis.Add(250);
            dis.Add(500);
            bcp.DistanceOffsets = dis;
            IBufferConstruction bc = bcp as BufferConstruction;

            IGeometryCollection originalGeometryBag;
            IEnumGeometry originalGeometryEnum;
            IGeometryCollection outBufferedGeometryCol;
            IGeometry resultGeo = null;
            originalGeometryBag = new GeometryBagClass();
            originalGeometryBag.AddGeometry(centroid, Type.Missing, Type.Missing);
            ISimpleFillSymbol simpleFillSymbol = new ESRI.ArcGIS.Display.SimpleFillSymbolClass();
            IFillShapeElement fillShapeElement;

            simpleFillSymbol.Style = esriSimpleFillStyle.esriSFSHollow;
            originalGeometryEnum = originalGeometryBag as IEnumGeometry;
            outBufferedGeometryCol = new GeometryBagClass();
            bc.ConstructBuffersByDistances2(originalGeometryEnum, dis, outBufferedGeometryCol);
            IEnumGeometry resultingGeometries;
            resultingGeometries = outBufferedGeometryCol as IEnumGeometry;
            if (resultingGeometries != null)
                resultGeo = resultingGeometries.Next();
            while (resultGeo != null)
            {
                fillShapeElement = new ESRI.ArcGIS.Carto.PolygonElementClass();
                fillShapeElement.Symbol = simpleFillSymbol;
                element = (ESRI.ArcGIS.Carto.IElement)fillShapeElement; // Explicit Cast

                element.Geometry = resultGeo;
                gc.AddElement(element, 0);
                resultGeo = resultingGeometries.Next();
            }

            mxDoc.ActiveView.PartialRefresh(esriViewDrawPhase.esriViewGraphics, null, null);

However the end result is not what I was expecting

As you can see from the hollow buffer, it has created the 500m buffer but has then added the 250m buffer onto the end, when I was expecting something similar to the coloured buffer, where the origin of the buffer is where the centroid of the point is.

Does anyone have any idea on creating a buffer that looks like the coloured one, which was done using the Create Multiple Rings Buffer tool in ArcToolbox.

EDIT**
I should have added, from my final buffer, the larger buffer needs to have a donut hole in the middle.

Thanks


texts: Adding time enabled feature to layer running newest to oldest

texts:
I am seeking information about the Time Enabled feature of my map created in ArcGIS 10.1, Basic (ArcView).

Specifically, I am wondering if the layer can be made to run from newest to oldest.  I update my layer a few times a year adding new information but I need the layer to open at the newest data and progress to the oldest data.  I have data spanning from 2007-2014, so I need the 2014 data to be displayed first then to progress backwards to 2007.  My published map is NC Shellfish Lease Siting Tool (UNCW).


texts: How do we call the plane a latitude circle is in?

texts:
I was just wondering if there is a proper term to refer to the plane that contains a circle of latitude.

I know for the special case of equator, that plane is called the equatorial plane. What about the circle of latitude for any latitude.

I googled "latitudinal plane" and "latitude plane" without finding anything


texts: Zonal statistics for several rasters in ArcGIS

texts:
I am running this Zonal Statistics as a Table code for a single raster and it works:

import arcpy, os
destination_path= r'\\Mac\Home\Desktop\arc_stuff'
shapefile= r'\\Mac\Home\Desktop\arc_stuff\subfolder1\adm_shp_all_selection.shp'
arcpy.env.workspace=r'\\Mac\Home\Desktop\arc_stuff'
raster='lights_1992.tif'
outZSaT = arcpy.gp.ZonalStatisticsAsTable_sa(shapefile, "unique_id", raster,'lights_1992.dbf',"DATA","SUM")

I'll have to process many rasters, hence I tried building the loop below, but I always get this error:

ERROR 999999: Error executing function. Create output table failed Failed to execute (ZonalStatisticsAsTable).

This is the loop:

import os
import arcpy, os
destination_path=r'\\Mac\Home\Desktop\arc_stuff\output'
shapefile= r'\\Mac\Home\Desktop\arc_stuff\subfolder1\adm_shp_all_selection.shp'
arcpy.env.workspace=r'\\Mac\Home\Desktop\arc_stuff'
rasters = arcpy.ListRasters("lights_*.tif")
for raster in rasters:
    print('Destination Path: %s\nRaster: %s\nDestination File: %s\n' % (destination_path, raster, destination_file))
    destination_file = os.path.join(destination_path, raster + ".dbf")
    outZSaT = arcpy.gp.ZonalStatisticsAsTable_sa(shapefile, "unique_id", raster, destination_file,"DATA","SUM")

I am aware of this related thread: Batch Zonal Statistics as table?, but I think the issues pointed out there do not apply to my problem.

Who can help?


texts: Seeking open source geospatial viewer of attachments stored in file geodatabase?

texts:
I have a file geodatabase that contains imagery and photo attachments linked to each record of a feature class.

I am looking for a non-proprietary solution that allows a user without Esri to view these attachments.

I do not have ArcGIS Publisher, hence the question.


texts:
i have a an array of polylines in epsg:3857 projection, and i want to calculate length of each polyline. i am new to arcgis and arcobjects, i couldn't find the functions that i need.

Should i use IGeometryServer2.GetLengths2 or something else?

note: if it matters i am developing a server object extension on C#

another note: if i need to use IGeometryServer2.GetLengths2 what should i use for IlinearUnit parameter?

Edit:

Administrator marked question as possible duplicate but this link
Programmatically calculating the length of all polyline features in a featureclass
shows calculating length in map units not in meters

this is the code i am currently using but i think it doesn't correctly calculates length

IProjectedCoordinateSystem myCoordSys = (IProjectedCoordinateSystem)myGeom.SpatialReference;

ILinearUnit myLinUnit = myCoordSys.CoordinateUnit;
IPolylineArray ploughList = new PolylineArrayClass();

ESRI.ArcGIS.Geometry.IPointCollection4 eachPloughPointColl = new ESRI.ArcGIS.Geometry.PolylineClass();

ESRI.ArcGIS.esriSystem.WKSPoint[] aWKSPointBuffer = new ESRI.ArcGIS.esriSystem.WKSPoint[2];
WKSPoint ploughStartPoint = new WKSPoint();
WKSPoint ploughEndPoint = new WKSPoint();
ESRI.ArcGIS.Geometry.IGeometryBridge2 pGeoBrg = new ESRI.ArcGIS.Geometry.GeometryEnvironment() as IGeometryBridge2;

//create lines
for (int j = 0; j < ploughCount; j++)
{
    //check point collection - clear
    if(eachPloughPointColl.PointCount > 0)
        eachPloughPointColl.RemovePoints(0, eachPloughPointColl.PointCount);

    //compute point coords for plough line
    ploughStartPoint.X = myGeom.Envelope.XMin;
    ploughStartPoint.Y = myGeom.Envelope.YMax - myFieldList[i].ploughWidth * (j + 1);

    ploughEndPoint.X = myGeom.Envelope.XMax;
    ploughEndPoint.Y = ploughStartPoint.Y;

    aWKSPointBuffer[0] = ploughStartPoint;
    aWKSPointBuffer[1] = ploughEndPoint;

    //set points into point collection
    pGeoBrg.SetWKSPoints(eachPloughPointColl, ref aWKSPointBuffer);

    //cast to polyline
    IPolyline4 bigPloughLine = (IPolyline4)eachPloughPointColl;

    ITopologicalOperator myTopo = (ITopologicalOperator)bigPloughLine;

    IPolyline4 finalPloughLine = (IPolyline4) myTopo.Intersect(myGeom, esriGeometryDimension.esriGeometry1Dimension);

    //rotate back to original position -- unnecessary remove!!! REMOVE!!!!!!
    rotationCNT = (ITransform2D)finalPloughLine;
    rotationCNT.Rotate(geomOrigin, Math.PI * 2 - myFieldList[i].ploughDirection.Angle);

    //add plough line to array
    ploughList.Add(finalPloughLine);
}
IGeometryServer2 myGeomS = new GeometryServerImplClass();

IDoubleArray myLengthList = myGeomS.GetLengths2(myGeom.SpatialReference, ploughList, myLinUnit);

I am trying to draw lines on a polygon and calculating it's lengths
like this post
Filling a polygon with lines using Postgis

But I think I can't get lengths correct.


texts: Converting raster stack or brick to ts time series object to run in BFAST


texts:
Background: I'm using GIS tutorial basic workbook (ArcGIS 10.1) with ArcGIS 10.5 software, and my MS office is 2010

I was practicing building a file geodatabase, and I have imported shapefile without issue, but I couldn't import the data table (.xlsx). The error message keeps popping up saying "Failed to connect to database. An underlying database error occurred. Class not registered".

I'm still fairly new to using the software so I'm not sure what went wrong.


texts: Reproject in QGIS doesn't work the way like in ArcGIS

texts:
Sometimes I get polygon data that digitizes buildings, but they are slightly displaced like in the following example

However, when I use the ArcGIS Tool "Project" to project the layer from a geographic to a projected CRS the displacement is corrected.

I have tried to do the same with the QGIS Tool "reproject Layer" but it doesn't work. The tool does change the unit of the extent from degrees to meters but it has no effect on the representation.

Honestly, I don't really know what is happening in the ArcGIS tool. In the description it says that "project" changes the unit from degrees to meters, thus it does the same as QGIS "reproject" but the results are not the same. Does anyone know what is happening here and can explain it to me?

EDIT

Actually, I made a mistake. I wasn't projecting from unprojected to projected but from a DHDN based to a ETRS_1989 based projection


texts: QGIS referring to a features length with python in function editor

texts:
I tried to modify the code given in https://gis.stackexchange.com/a/212762/7849 so the value from [field_name] would be divided by the length of the feature.  I tried to change the last line in

return sum( f[field_name]/f.geometry().length() for f in features if geometry.contains(f.geometry()) )

But I get an error: 'builtin_function_or_method' object has no attribute 'length'


texts: I am facing this problem how to fix it. How can I fix it? I already uninstall and reinstall numpy but that did not work.


texts: Bootstrap 4 popover is not fixed on the map

texts:
I am using bootstrap 4 to   show some popovers on the points of the map, but they are not fixed, when you scroll the map the popover changes position.

In bootstrap 3 it works fine though. I don't know if there is any solution in OpenLayers to this.

https://jsfiddle.net/70p621w9/


texts: How to Import Circular Polar Orbit (HDF5) data into ArcMap with WGS 1984 Projection

texts:
I'm currently trying to retrieve data from http://smap.jpl.nasa.gov/data/ and use it in ArcMap. The issue is that the data is curved towards the poles and so the bottom left of the extent isn't the same latitude as the bottom right. This means that even if I define the projection as WGS 1984, it doesn't look as it should look. It comes up as a rectangular raster that isn't in the right geographic location.

Some things to note...
The HDF5 contains a grid of latitude values as well as a grid of longitude values along with the soil moisture subdatasets I'm interested in working with. These are separate subdatasets but maybe I could save them as text files, import them into excel, and then merge them somehow into file readable by ArcMap. This way when combined with the soil moisture data I'll be able to see it in ArcMap in the right geographic location.


texts:
I don't have any answers (yet) from from previous question on GDAL memory raster + OpenLayers, so i'll try to ask an easier question :)

Anyone knows how to access GDAL memory with command lines when using the -of MEM output ? As ex, after this command :

gdaldem color-relief -of MEM -alpha C:\tmp\myImg.tiff C:\tmp\myLgd.txt myMemOuput


texts: I am working with Topology plugins in Qgis 3.8, but I don't see tolerance in topology. I hope to get help. I have attach images.


texts: Pyscripter and arcpy not working together

texts:
I have installed pysripter on my machine and am not able to run any arcpy scripts from it despite being able to import arcpy.  Models that work in model builder do not transition to pyscripter here is an example of an error I got

Traceback (most recent call last):
  File "S:\Shared\1 Survey Operations\Data Management\Databases\DB\Downloads\Testing\ftof.py", line 13, in <module>
    arcpy.ImportToolbox("Model Functions")
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\__init__.py", line 124, in ImportToolbox
    return import_toolbox(input_file, module_name)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\toolbox_code.py", line 438, in import_toolbox
    toolbox = gp.createObject("Toolbox", tbxfile)
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\geoprocessing\_base.py", line 379, in createObject
    self._gp.CreateObject(*gp_fixargs(args, True)))
IOError: The toolbox file Model Functions was not found.


texts: QGIS 2.18.6 How to avoid the refresh of map when you export as image

texts:
In the version 2.18.6 if I make a layout with many thematic maps coming from the same shape, how can I avoid the refresh of the maps?

I am doing a layout with the same shapefile split in four maps with different thematic design but when I export it all maps are refreshed with only design.


texts:
I was having some trouble when trying to zoom to certain point in ESRI base map. So basically I got an ESRI base map and added a map overlay on top of it:

function setMap() {
function init() {
    require(
        [
            "esri/map",
            "dojo/dom-construct",
            "esri/geometry/Point",
            "esri/symbols/SimpleMarkerSymbol", "esri/symbols/SimpleLineSymbol",
            "esri/graphic", "esri/Color",
            "dojo/domReady!"
        ],
        function
        (
            Map, domConstruct, Point,
            SimpleMarkerSymbol, SimpleLineSymbol,
            Graphic, Color
        ) {
            map = Map("map-canvas",
            {
            });
            map.setZoom(0);
            coreFunctions();
        });
}
dojo.ready(init);
gsvc = new esri.tasks.GeometryService("http://tasks.arcgisonline.com/ArcGIS/rest/services/Geometry/GeometryServer");
}

Here is the part for me to add map overlay with wkid:3414:

function addLayersToData() {
var layer = new esri.layers.ArcGISTiledMapServiceLayer("https://www.onemap.sg/ArcGIS/rest/services/BASEMAP/MapServer");
mapLayers.push(layer);

}

And then I got a function to zoom to certain point of the map:

function zoomPostal(postalCode) {
 $.getJSON("http://www.onemap.sg/API/services.svc/basicSearch?token=qo/s2TnSUmfLz+32CvLC4RMVkzEFYjxqyti1KhByvEacEdMWBpCuSSQ+IFRT84QjGPBCuz/cBom8PfSm3GjEsGc8PkdEEOEr&searchVal="
                        + postalCode
                        + "&otptFlds=SEARCHVAL,CATEGORY&returnGeom=0&rset=1", function (data) {

    var loc = new esri.geometry.Point({ "x": data.SearchResults[1].X, "y": data.SearchResults[1].Y, "spatialReference": { "wkid": 3414 } });
    console.log(data.SearchResults[1].X);
    console.log(data.SearchResults[1].Y);
    map.centerAndZoom(loc,1);
    });
}

When I print out the X and Y, I am getting: 23948.7465,47505.7986. And I am getting these error message after I tried to zoom to certain point:

Map: Geometry (wkid: 3414) cannot be converted to spatial reference of the map (wkid: PROJCS["SVY21",GEOGCS["SVY21[WGS84]",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",28001.642],PARAMETER["False_Northing",38744.572],PARAMETER["Central_Meridian",103.8333333333333],PARAMETER["Scale_Factor",1.0],PARAMETER["Latitude_Of_Origin",1.366666666666667],UNIT["Meter",1.0]])

I have no idea how to solve this. Any ideas?

Basically if I added these code to see the conversion of the coordinate by plotting a marker, it works as in the marker did shows up. Just that I cannot zoom to that point:

var mrtIcon = [];
    var symbol = new esri.symbol.PictureMarkerSymbol('img/transportation/mrt_marker.PNG', 30, 30);
    var PointGraphic = new esri.Graphic(loc, symbol);
    map.graphics.add(PointGraphic);

    var graphic = PointGraphic;
    graphic.setSymbol(symbol);
    mrtIcon.push(map.graphics.add(graphic));
    map.centerAt(loc);

EDIT

var x =  data.SearchResults[1].X;
     console.log(x);
     var y =  data.SearchResults[1].Y;
     console.log(y);
    var loc = new esri.geometry.Point({
        "x": x,
        "y": y,
        "spatialReference":'PROJCS["SVY21",GEOGCS["SVY21[WGS84]",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",28001.642],PARAMETER["False_Northing",38744.572],PARAMETER["Central_Meridian",103.8333333333333],PARAMETER["Scale_Factor",1.0],PARAMETER["Latitude_Of_Origin",1.366666666666667],UNIT["Meter",1.0]]' });

    map.centerAndZoom(loc, 5);

So I got variable x and y which hold the same value as your hardcoded x and y coordinates. It does not work if I put it this way, it only work when I hardcoded the x and y with coordinates like your solutions. Any ideas?


texts: Simple point in polygon query at PostGIS & PostgreSQL

texts:
Using PostGIS & PostgreSQL (Amazon RDS instance), I need to create a (quite large) table with points (just x, y coordinates)
Then what I need to do is very simple point-in-polygon queries. Points and polygons are in 2D space (not geography types).

I actually managed to do this with Amazon Athena with this simple query:
SELECT ST_WITHIN(ST_POINT(x, y), ST_POLYGON('polygon((1  1, 11  44, 44  44, 44 11))')) FROM "my_table"

How it would be possible to have this kind of query working at PostgreSQL?


texts: ArcPy script causes Python crash when run outside of ArcGIS Python window

texts:
I have a Python script that will successfully complete when run from the Python window in ArcMap or ArcCatalog, however, when I run as a standalone it will cause python.exe stopped working crash once it completes.  I've tried from a command prompt and a couple different IDEs with the same result.

The target files and updates complete before the crash and when I compare the results from the stand alone script to the one run within the Python window they appear identical.  I've placed print statements at various locations including on the last line to help with debugging and all of them print.  Since the crash is happening after the all lines in the script have been run, I'm not sure how to further debug this. I'm running ArcGIS 10.5.1.

If I comment out the arcpy.CreateRapidGraphic line then the crash does not occur.

I'm primarily running via a command prompt

C:\Python27\ArcGIS10.5\python.exe D:\Projects\py\createmap.py

but have also tried double-clicking and running in visual studio code.

The last lines of my code are similar to this:

arcpy.CreateRapidGraphic_defense(geodb, areaOfInterestLyr, productType, version, mapDocumentWork, rasterLyr)
print "graphic done"
arcpy.CheckInExtension('defense')
arcpy.CheckInExtension('foundation')
arcpy.CheckInExtension('spatial')
arcpy.CheckInExtension('3D')
print "checkin done"


texts: Create a blank raster from scratch using QGIS?

texts:
I'm wondering if there is any way to create a raster with no input data aside from an extent and cell size (and possibly a default value).

I'm using QGIS as a way to try my analysis and eventually script it, but just am looking for how to do it in QGIS at the moment.


texts: Software or script for Batch Shapefile to GPX

texts: I am looking for a software that does batch processing of shapefiles to GPX files.  I have seen many including Arcmap that will batch gpx to shapefiles but none vice versa.  I currently am using DNR Garmin to complete the task one file at a time but it is time consuming.


texts: How do I get mobile AR mode working in qgis2threejs?

texts:
I've created a map in qgis and exported using qgis2threejs and the mobile template. Everything renders fine in the browser (laptop and phone). Rotation, pan and zoom work as expected. When I load the page in the browser on my phone (iPhone 6s Plus, iOS 12.1.4), the map renders and I get a slider control for AR (augmented reality) in the upper right, a cross-hairs button for GPS location, and settings control. If I tap the cross-hairs, the current GPS position is displayed in a dialog, but the location is not indicated on the map, nor does the map shift position/orientation. If I enable the AR slider, I get a popup asking permission to use the camera. I click yes, and the screen goes black.

I've tried exporting in the map CRS (UTM 13N) and also enabling lat/lon conversion. Same results either way.

Any clues on what supposed to happen and how to make this work?

qgis 3.4.2, current version of qgis2threejs from plug-ins repository


texts:
I only have basic GIS skills and wanted some advice or direction.
I'm only practising with a smaller dataset at the moment though my other files are much larger 10,000+ so don't know if it's possible or if I would have to separate the data by years.
Essentially I have a CSV file with the names of ancestors name, birth year, birth location, mother, father, and the coordinates of the birth location.
I wanted to create a time series / animation of their chronological birth order. I've got them in a basic map at the moment on QGIS as that is all I have access to (I don't know if that makes a difference) and span across America and Europe.
I had a look at TIME MANAGER but that seems outdated (I couldn't install it) and read that the plugin temporal profile may be good for this, however I need to convert my CSV file to a raster I believe, but I cannot understand how the rasterize would work on my CSV file because it requires output extent (I don't understand what this is) I had a look at my file and it states 0,0,0,0.
Am I going about this the right way or am I completely off track. Also I've seen people do this on YouTube with time manager but with times, can I do this just with years?


texts: Using Microsoft Footprint GeoJSON for spatial analysis with GDAL?

texts:
I am trying to use the Microsoft Footprint data (available as GeoJSON) for spatial analysis to roughly identify buildings over a certain size within a specific range of major railroads and highways. I am confident that I know how to perform the analysis, but I can not find a way to use a GeoJSON file as large as these (e.g. Alabama approx. 500 mb of txt) for spatial analysis.

I can very slowly load the GeoJSON into QGIS, but I cannot export as a shapefile or run any search-by-locations without the whole system crashing.

I have tried to use ArcMap's JSON to Feature, but it will not read a GeoJSON of that size.

I have also tried to create individual shapefiles from each feature within the GeoJSON to avoid having a shapefile over the size limit. My intention was to then include these all in a GDB for analysis, but this process too constantly crashed my Python instance and even unmounted a disk.

from osgeo import ogr
import fiona
import json

text = open('Alabama.txt', 'r')
text = text.read()
search = re.findall('{"type":"Polygon","coordinates"[^}]+', text)

i = 0

for line in search:
    i += 1
    poly = line + '}'
    polyogr = ogr.CreateGeometryFromJson(poly)
    schema = {'geometry': 'Polygon','properties': {'fld_a': 'str:50'}}
    with fiona.open('%s.shp' % i, 'w', 'ESRI Shapefile', schema) as layer:
        layer.write({'geometry': json.loads(poly), 'properties': {'fld_a': 'test'}})

Can I use GeoJSON files this large for analyst work in GDAL?


texts:
I'm trying to use Geoserver's RESTful API within limited success. I've found the REST API's documentation, but what I seem to be missing is the documentation that outlines the XML/JSON interface. From what I gather, there is a whole XML/JSON interface that is used for GET'ing/POST'ing new constructs, such as workspaces, data stores, layers, etc... Here are plenty of examples

Is there other documentation that fully outlines this XML? The above examples are great, but I'm having a hard time extrapolating from these examples. A complete interface documentation (like the REST API) would be helpful.


texts: Contour generation using Interpolation - Correcting edge Effects

texts:
I utilised the Interpolation feature in QGIS to prepare contours from a DEM I had generated.  Sometimes when contours intersect the edge of the DEM area, they run along the edge for some distance.  This doesn't seem to occur all the time, but only at select locations.

Is there a method to avoid this?  It seems if I cut 1m of the edge of the contour they would be perfect.


texts:
I'm writing a custom GP tool for printing with 10.1. I need to pass some additional parameters into GP tool. I'm using printTask from jsapi and I couldn't add more parameters just by adding to PrintParameters.

ptemplate.preserveScale = false;
params.template = ptemplate;
params.MyNewParameter = 'testing';
printTask.execute(params, printComplete);

I already have a workaround by adding it as URL parameter gpurl?MyNewParameter='testing'.

In reality my parameters are much complex than this, and I'm seeking for better solution.

Technologies are:

Arcgis server geo processing
Esri JSAPI
Arcpy


texts: How to get Perpendicular cross section along the river using python and Arc Map 10?

texts: I need to draw the perpendicular line segment along the river at a regular interval. How can i obtain like this output using Python and Arc Map 10 ? I would be very grateful if you give some idea..


texts: Calculating Geometry in ArcGIS without a Projected Coordinate System

texts:
For my project, I want to calculate the geometry of individual polygons in my layer to find area and perimeter, however it is greyed out when I try to use this tool.

Beforehand, I set the data frame for map and display units to be in meters. I also did not put a projected coordinate system. However, I was able to find which polygons are within a certain distance (ex. 1000 m) from roads using the 'Select by Location' tool.

Why is it that ArcGIS can calculate distance but won't allow me to calculate the geometry of my individual polygons. Am I missing a step here?


texts: Receive this error log when running Join attributes by location

texts:
When I try running a join attributes by location in QGIS 2.18.13 I receive this error message which I don't understand:

2017-11-22T09:39:55    2   Uncaught error while executing algorithm
            Traceback (most recent call last):
              File "C:/PROGRA~1/QGIS2~1.18/apps/qgis/./python/plugins\processing\core\GeoAlgorithm.py", line 203, in execute
                self.processAlgorithm(progress)
              File "C:/PROGRA~1/QGIS2~1.18/apps/qgis/./python/plugins\processing\algs\qgis\SpatialJoin.py", line 147, in processAlgorithm
                index = vector.spatialindex(join)
              File "C:/PROGRA~1/QGIS2~1.18/apps/qgis/./python/plugins\processing\tools\vector.py", line 219, in spatialindex
                idx = QgsSpatialIndex(layer.getFeatures(request))
            Exception: unknown

Both of my layers are relatively large shapefiles and I have created a spatial index for each of them.


texts: Orfeo Toolbox KMeans Classification failure

texts:
I am trying to process a Landsat8 OLI image using the Orfeo Toolbox  Unsupervised KMeans Classification application.

I obtain the following error message:

Error when loading model 20140915_E1B1234_km.tif_model.txt

It works fine on Sentinel2 MSI images.

Any idea what could cause the application to fail on my Landsat image, although it works on the Sentinel image? Has anybody already encountered this issue?


texts: Looking for annotated vector map of world?

texts:
I'm looking for an annotated vectorial (PS/SVG/...) map of the world (+countries).

Wikipedia contains many vectorial maps (e.g: http://commons.wikimedia.org/wiki/Category:Blank_SVG_maps_of_the_world )but the paths are defined like that:

<path fill="#ccc" d="M98.171,199.174c-0.505-0.076-0.938-0.322-0...

... and I'd like to find an identifier so I can colorize a given country (using XSLT or another method....). Something like that:

<path id="JAPAN" fill="#ccc" d="M98.171,199.174c-0.505-0.076-0.938-0.322-0...


texts: Which is the best Android tool to visualize JPG+JGW images?

texts:
I've got some georeferenced JPG+JGW images, and would like to visualize them on my Android device. Which is, in your opinion, the best tool?

I've tried OSMdroid (after converting JPG+JGW with MapTiler), but I can't create points to mark some places. Are other alternatives?


texts: Confusing offsets in shapefile written by QGIS/OGR

texts:
I think I'm stuck here.

A shapefile has been updated in QGIS. Records have been deleted and those records have been marked with an asterix in the first byte of the record in the dbf file.
I open the shapefile in a C program to convert it to postGIS. The deleted records are skipped. Things happen to the next record after a deleted record, not after the first deleted record though. The shx file says that the offset for the next record after the deleted record is 30948 and the accumulated offset in the shp file says 28740. Neither is right, or rather when you start reading in the shp file at the two offsets you land up in the middle of a point stream and not at the start of a record header. The offset 28740 should be right because that is what you get when you add up the offset of the deleted record and the contentlength of the deleted record.

The conclusion is that the shape is corrupt. But what is driving me nuts is that you can load the shape file into QGIS and it is shown correctly.

What have I missed.How can QGIS handle this.


texts: Getting vertices of segment of polyline using ArcGIS API for JavaScript?

texts:
I have a poly line with 3 co-ordinate, say (x1,y1), (x2,y2), (x3,y3). Can I get the vertices (x1,y1), (x2,y2) if clicked on line-segment drawn by (x1,y1), (x2,y2)?

I have a poly line on map that with path [[x1,y1],[x2,y2],[x3,y3]]. And I have electrical box at the corresponding paths.

When clicking the segment drawn by the points x1,y1 and x2,y2 I want to get the point (x1,y1) and (x2,y2). Same for path (x2,y2), (x3,y3).

I'm using Javascript API.


texts: Remove saving notification when use python QGIS at startup

texts:
I try to run python at QGIS startup. The script :

from qgis.core import QgsProject

from PyQt4.QtCore import QFileInfo

project = QgsProject.instance()

project.read(QFileInfo('D:/test/qgis/1.qgs'))

I save as startup.py in C:\Users\USER.qgis2\python. The project is successfully opened. There are two layer. A vector and a raster. But, there is a notification : "Do you want to save the current project?". If I click 'Save' or 'Discard', the 1.qgs is saved , and then there is a blank project. If I click 'Cancel', the 1.qgs is still opened.
It is disturb a lay user. While I don't want to close the project. I just want to load the project.

How do I remove this notification?


texts:
Similarly to this question (Is there a source for free CSV files with geo data and other desired fields?) I am working on testing an open source set of bindings for the GEOS project and would like to boot strap the process with some existing data and expected outcomes. (this polygon intersect that polygon gives this polygon) type of thing.

Anyone know of any good reference implementations of testing for geo libraries??


texts: Creating individual routes from home to multiple destinations in ArcMap ModelBuilder with ArcGIS Network Analyst

texts:
I am trying to create individual routes from Homes to multiple Destinations. The unique ID is respondent in both the Homes and Destinations tables. I've created an iterative model in ModelBuilder where the Routes are calculated based on the RouteName in this case Respondent, the problem is that I have multiple destinations (stops) for many Respondents and when running the Route model, the model adds all the stops together to create one long route for that individual Respondent.

For example, I have Home stop, and I want to calculate the routes from Home stop, Respondent 101 to multiple Destination stops that have the same Respondent 101 ID.

Can you help/advise what to change in the model to make it do this?

See the photos attached. I've got lots of data and want my model to iterate through the list of Homes, create a route with multiple paired destinations and then output to one file. The attribute table shows the Route analysis has worked but only with paired data that has 2 stops only, it didn't work for multiple stops.


texts: Arcpy Copy_management does not copy all records

texts:
I need some help with a python script that is copying a shapefile to a new folder location. I am using ArcGIS 10.3. I have a shapefile of imagery tile outlines that contains 17425 polygons. I need to make a copy of the shapefile in a new folder location. I am using arcpy.Copy_management() to copy the shapefile into the new location.  After I run the script the copied shapefile only contains 2826 of the original records. Some how I am losing tiles during the copy. If I run the tool manually from inside a map session, all the records copy over completely. I think the problem might be my shapefile and not the script, but have a look. Here is my script:

# Import system modules
import arcpy
from arcpy import env

# Set workspace
env.workspace = "C:/Documents/Working"

# Set local variables
in_data =  "Imagery_Tiles.shp"
out_data = "C:/Documents/Tiles/Tile_Shapes.shp"
data_type = ""

# Execute Copy
arcpy.Copy_management(in_data, out_data, data_type)


texts: Invalid geometry in point layer when running 'Count points in polygon'

texts:
I have a persistent problem I can't seem to fix. I'm trying to run Count points in polygon. I keep getting this error:

  Feature (747) has invalid geometry. Please fix the geometry or change the Processing setting to the "Ignore invalid input features" option.
  Execution failed after 30.86 seconds

I've run the multipart to singlepart tool, and the Check Geometry tool now comes back with no errors. Same for the Check Validity tool, with both GEOS and QGIS options. And I tried to use the Fix Geometries tool. But I'm still getting the error when I go to do the Count Points in Polygon.

Here is a link to the point and polygon files:
https://www.dropbox.com/sh/vz0ap75tab09tht/AAAUxjhR3oWZ8tUxRLclwMbPa?dl=0


texts: Copying attribute table definition to one SHP to another using QGIS?

texts:
I have a SHP with 3D vectors and in the other hand, a SHP with attribute table definition. I need one SHP with the 3D vectors and the attribute table definition.

It's important to understand that the attribute table is empty. That's actually my job to fill it.

I tried to copy/paste the 3D vectors to a SHP to the other with attribute table definition but I noticed that Z information get lost.

So I have two solution:

find out a solution to copy/paste the vectors without losing the z coordinate
find out a solution to import the empty attribute table definition into the shp containing the vectors.


texts: Is Mars on SpatialReference.org (srid 49901) broken?

texts:
Just noticed there were seemingly duplicate SRSIDs for Mars 49900 and 49901 on spatialreference.org. I saw they were both generating the same PostGIS insert statements with only a different SRIDs changed. Wanting to know more I found Interoperability in Planetary Research for Geospatial Data Analysis, it shows

  IAU2000 49900 Mars2000, aerocentric latitudes, positive East longitudes

  IAU2000 49901 Mars2000, aerographic latitudes, positive West longitudes

That seems to indicate they're different. Shouldn't that generate a +axis=wnu in the pro4j for 49901? Currently, it only shows,

+proj=longlat +a=3396190 +b=3376200 +no_defs


texts: How do I create more makers with more unique popups?

texts:
I am attempting to create a map with Leaflet that has unique files attached to different popups. I can generate a marker with the popups I want but I can not generate a second marker with different popups.
Here is my code, how do I create more markers with more unique popups?
<!DOCTYPE html>

    <title> Bowers First Web Map</title>
        <link rel="stylesheet" href="https://unpkg.com/leaflet@1.3.1/dist/leaflet.css"
   integrity="sha512-Rksm5RenBEKSKFjgI3a41vrjkw4EVPlJ3+OiI65vTjIdo9brlAacEuKOiQ5OFh7cOI1bkDwLqdLw3Zg0cRJAAQ=="
   crossorigin=""/>


texts: Getting elevation at lat/long from raster using python?

texts:
I was wondering if anyone has some experience in getting elevation data from a raster without using ArcGIS, but rather get the information as a python list or dict?

I get my XY data as a list of tuples:

xy34 =[perp_obj[j].CalcPnts(float(i.dist), orientation) for j in range (len(perp_obj))]

I'd like to loop through the list or pass it to a function or class-method to get the corresponding elevation for the xy-pairs.

I did some research on the topic and the gdal API sounds promising. Can anyone advice me how to go about things, pitfalls, sample code?

GDAL is not an option as I can't edit the system path variable on the machine I'm working on!

Does anyone know about a different approach?


texts:
I have installed Mapnik-2.0.1rc0 in windows 7 32bit machine. I then tested the installation by import mapnik from * and it did not return any error. But when I tried to input the config version checking by mapnik-config -v, it returns a name error. When I tried to run rundemo.py, the following error occurred:

  the mapnik library and python bindings must have been compiled and installed successfully before running this script

This probably happened due to improperly set PATH and PYTHONPATH.

How to add PythonPath?

Thanks


texts: Distributing maps at same distante in one axe

texts: I have many maps to put into print composer. I would like to align them at the same distance each other throw vertical axe. Is there tool that align them at the same distance automatically?


texts: Update "value relation" field problem in QGIS multi-edit mode

texts: I tried to update a value relation field in multi-edit mode (more than one feature selected). I select a field value (from combo) for all features but it doesn´t work properly. The field is only updated for one feature. I am using QGIS 2.18.9 and two tables in SpatiaLite: one for line features and another for key-values (without geometry).


texts: Find missing number from sequence of row data in FME

texts:
I have 100 tables.

In each table, there will be 9 row of records stored with unique ID like { 1,2,3,4,5,6,7,8,9 }

Need to log if there is any missing row in those tables with Unique ID.

For eg:

Below the table, there is missing row 4 and 7. need to find 4 and 7 row.


texts:
I can't use my ArcGIS Pro software. The font color is white, when the background is white.

How can I fix this?


texts:
I'm desperately trying to swap out my basemap in cartodb.js.

Here's the snippet of code, the part I switched out successfully with some basemaps is in bold:
L.tileLayer('https://dnv9my2eseobd.cloudfront.net/v3/cartodb.map-4xtxp73f/{z}/{x}/{y}.png'

When I inspect element or go here, I get a wide variety of options, like so:
http://leaflet-extras.github.io/leaflet-providers/preview/

For example, this works:
L.tileLayer('http://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x}'

But this doesn't:

L.tileLayer('http://{s}.{base}.maps.cit.api.here.com/maptile/2.1/{type}/{mapID}/hybrid.day/{z}/{x}/{y}/{size}/{format}?app_id={app_id}&app_code={app_code}&lg={language}'

How do I swap out this code??


texts: create a weight raster for a hydrology analysis

texts: I determine affecting parameters of flood flow ( evaporation,vegetation,soil type, rainfall etc.). Then, I start to analysis using Arc Hydro. How do I relate among of affecting parameters? How can learn affecting percent of parameters? I dont know how I can create weight raster of parameters.


texts: Convert shapefile data to NMEA - simulate NMEA positions from shape data

texts:
To generate waypoints to be read by the agricultural GNSS-RTK-steering system CEREA, I am looking for a way to simulate NMEA strings of a GNSS receiver by clicking the edges of the boundaries of polygon shapefiles within QGIS.
Is there a way to do that edge by edge of the polygon? The data should be exported via TCP/IP to be read by the navigation system as if the data where actual position data of a GNSS receiver using NMEA strings $GNRMC / $GPRMC and $GNGGA / $GPGGA.
If there is no way within QGIS, a workaround could be a program, that sends NMEA data from coordinates that I copy and paste from a CSV file of the shape data.


texts: Easiest way to read roads, railroads, buildings from an OSM map

texts:
I'm prototyping a Java application, which uses an OSM map of a small city.

I need to extract from the OSM map the coordinates of

roads,
railroads,
power lines and
buildings.

At the moment I'm aware of only one way to do it - setup a PostgreSQL database with PostGIS, import the map there and then run SQL queries to get the data.

Are there any alternatives, ways to implement it (reading coordinates of objects a-d) in a more lightweight way (without the nontrivial installation of PostGIS, except Portable GIS) ?

For now, I don't care about performance, I just want to learn how to work with OSM in Java.


texts: Representing overlapping polylines in QGIS

texts:
I have several shapefiles representing bus lines with different colors for each segment, however, because the lines are exactly superposed, only the shapefile on top can be seen. Is there a way to display the lines underneath for example side-by-side on these parts of segments where several lines converge (there might be 4 or 5 of them)? The polylines are made up of series of segments.

I know there are some tools in ArcGIS, but I'm using QGIS 2.8.0


texts: What worldwide raster basemaps does MapBox provide, and how much do they cost to use?

texts:
There are a number of MapBox basemaps that show up on various websites (MapBox Streets, MapBox Outdoors), but weirdly, I can't find a list of them anywhere. I have in the past seen a big list of basemaps, including a few of Afghanistan, but they weren't necessarily all worldwide, and I can't find that list now.

There is this short list of half a dozen pretty esoteric maps: https://www.mapbox.com/design/

And various old blog posts that announce new basemaps don't point anywhere anymore: MapBox Streets, MapBox Light

Logging into MapBox I can see this set of basemaps, but is that comprehensive?

And secondly, how would you go about using a standard MapBox basemap in an application? What's the pricing? Where do you pay?

Do they just not do this anymore?


texts: Non-constant transformation of photogrammetric point clouds in R?

texts:
I have a process which attempts to register two forested point clouds using CloudCompare's Iterative Closest Point (ICP) algorithm. Often times, warping (especially in the vertical direction) occurs in the generation of UAV derived photogrammetric clouds.

Using a moving window, my process attempts to match the sample clouds found in each moving window (top of figure). ICP returns a matrix of transformation values, of which the estimated X,Y and Z transformation values are then stored and associated with the center coords of each window (bottom of figure - the coloured points).

With the resulting point dataset from the moving window ICP process, we can fit models to predict the X, Y and Z shift values as a function of 2D space. I have experimented with low order polynomial functions for example:

Z shift = b0(UTMx) + b1(UTMx)^2 + b2(UTMx)^3 + b3(UTMy) + b4(UTMy)^2 + intercept

where b0-b4 are model coefficients and UTMx/y refers to the 2D coordinates.  2 other models are fit for X and Y shift.  So far the models prove to be accurate (rsquared ~0.9) and applying the shifts to the sample subsections of the original cloud show very close alignment visually.

With these models we can go back to the original Photogrammetric Cloud and apply X,Y and Z shifts to correct for warping (at least so it matches the LiDAR cloud which we know has better georeferencing procedure than UAV).  These corrections ideally would be applied to each point of the original Photogrammetric Cloud...  I am looking for an efficient way to apply the shifts (which vary according to the models) in a computationally efficient way because in some cases we deal with up to a billion points. My first instinct is to using some sort of tiling procedure, running these in parallel... I prefer using R and the lidR package however I'm open to other suggestions/comments/questions regarding any part of this process!


texts:
I'm working on an arcpy script that removes a feature class and then recreates it.  The feature class is part of a replica, and I need to recreate the replica as part of the script, however you cannot unregister a replica using arcpy.

If I manually remove the replica using Distributed Geodatabase > Manage Replicas > Unregister, the script works fine recreating the feature class and the replica, however if I don't manually remove it first the script fails:

Message: ERROR 000775: Replica with the same name already exists.
Failed to execute (CreateReplica).

I have found my replica in SQL using:

SELECT *
FROM [GDB_ITEMS]
WHERE Type = '4ED4A58E-621F-4043-95ED-850FBA45FCBC'
AND Name = 'MyReplicaName'

Can I "safely" remove this entry in GDB_ITEMS to unregister my replica?  If so, I could probably use arcpy to send the SQL to remove it just so I can continue with my script.

By "safely" I mean it's just a single reference and I don't need to remove it from 10 other tables as well to unregister the replica.

If it's relevant, the replica is two-way full SDE to SDE database replica.


texts: Mapnik wouldn't render nodata as transparent

texts:
I'm attempting to use TileStache with a Mapnik source to output what would largely be transparent raster tiles from a GEOTIFF. I'm having difficulty getting Mapnik to produce images with the proper transparency, though. My Mapnik XML file looks like

<Map srs="+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs">
    <Style name="style">
        <Rule>
            <RasterSymbolizer />
        </Rule>
    </Style>
    <Layer name="layer" srs="+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs">
        <StyleName>style</StyleName>
        <Datasource>
            <Parameter name="type">gdal</Parameter>
            <Parameter name="file">/home/ainlaen/test.tif</Parameter>
            <Parameter name="nodata">0</Parameter>
        </Datasource>
    </Layer>
</Map>

I have tried setting the nodata values both in the XML and with GDAL. The image will appear with proper transparency when viewed in QGIS after having it set with GDAL, but will have a black background when rendered by Mapnik, regardless of nodata setting. I have tried rendering images directly using Mapnik, and the result is the same as when using TileStache, so I don't believe that's the problem. In addition, the tiles are not displayed in their proper positions, but that's probably best saved for a separate question.


texts: Authenticating WFS from GeoServer with Keycloak in QGIS

texts:
I have a Geosever with the Keycloak Community extension (https://docs.geoserver.org/latest/en/user/community/keycloak/index.html) installed.
Authentication with the Keycloak Server in the Browser works fine.
Now I need to authenticate my WFS against Keycloak out of QGIS.
Is there a way to do this?


texts: How to use Zonal Stats as Table tool with in-line variable substitution

texts:
I have created a model to iterate a folder of 5 shapefiles to apply the Zonal Stats as Table tool to a raster.  My output table is saved as table%i% so i'm expecting it produce 5 tables (table0, table1, table2, table3, table4, table5).  However it keeps saving over the same table so I just get 'table0'.  Can anyone suggest where i'm going wrong?  Thanks.

--

Import arcpy module

import arcpy

Check out any necessary licenses

arcpy.CheckOutExtension("spatial")

Load required toolboxes

arcpy.ImportToolbox("Model Functions")

Script arguments

v5_River_Shapefiles = arcpy.GetParameterAsText(0)
if v5_River_Shapefiles == '#' or not v5_River_Shapefiles:
    v5_River_Shapefiles = "M:\Dissertation\Data\5 River Shapefiles" # provide a default value if unspecified

Local variables:

FeatureClass = v5_River_Shapefiles
table_i_ = FeatureClass
Name = v5_River_Shapefiles
Standard_Average_Annual_Rainfall_SAAR = "SAAR 5km (F)\Standard Average Annual Rainfall (SAAR)"

Process: Iterate Feature Classes

arcpy.IterateFeatureClasses_mb(v5_River_Shapefiles, "", "POLYGON", "RECURSIVE")

Process: Zonal Statistics as Table

arcpy.gp.ZonalStatisticsAsTable_sa(FeatureClass, "ID", Standard_Average_Annual_Rainfall_SAAR, table_i_, "DATA", "ALL")


texts:
I have a simple layer in PostGIS created as follows:

CREATE TABLE layer (
    id serial PRIMARY KEY,
    geom geometry(LINESTRING, 31370),
    ...
)

To avoid duplicate information, I created a view with some information (eg length) derived from the geometry:

CREATE VIEW myview AS
    SELECT *,
    ST_LENGTH(geom) AS length
    FROM layer;

I can access and edit the view myview in QGIS but whenever I want to save my edits, I got this error message: "Could not commit changes to layer. PostGIS error while adding features: ERROR:  cannot insert into column "length" of view "myview". DETAIL:  View columns that are not columns of their base relation are not updatable."

Is there a way to edit and save a view of a PostGIS table in QGIS?


texts: Leaflet populate L.divIcon marker's html option from feature properties

texts:
I am using Leaflet.js to make a timeline-based playback of marker locations, using a Leaflet Playback listed on Leaflet Plugins.

I am defining an L.divIcon earlier in my code:

var busIcon = L.divIcon({className: 'bus-div-icon', html:"1"});

And later on I have to pass it to the plugin's function:

var playbackOptions = {
  //...
  marker: {
    icon : busIcon,
    getPopup: function(featureData) {
      // returns something using featureData.properties.name, so here a function works
    }
  },
  //...

I want to populate the html variable in the divIcon dynamically, similarly to how getPopup has been defined (this code I have copied and customized from the plugin's documentation example).

But when I try to pass a function to the icon:

icon: function(featureData) {
  return L.divIcon({className: 'bus-div-icon', html:
  featureData.properties.busnum });
},

It did not work; I got this error in my console:

  Uncaught TypeError: t.icon.createIcon is not a function
     at e._initIcon (leaflet.js:5)
     at e.onAdd (leaflet.js:5)
     at e._layerAdd (leaflet.js:5)
     at e.whenReady (leaflet.js:5)
     at e.addLayer (leaflet.js:5)
     at e.addTo (leaflet.js:5)
     at e.addTrack (LeafletPlayback.min.js:1)
     at e.addData (LeafletPlayback.min.js:1)
     at e.setData (LeafletPlayback.min.js:1)
     at e.initialize (LeafletPlayback.min.js:1)

I tried another strategy, and no errors this time:

icon: L.divIcon({
  className: 'bus-div-icon',
  html: function(featureData) {
    console.log("Hello!");
    return '1';
  }
}),

..  but the results were hilarious!

What happened here is that instead of executing the function, it has displayed that whole code as text on the map! See screenshot above.

How to populate html dynamically?


texts: Calculating walking distance from train stations in ArcMap?

texts:
I am trying to find a way to make a buffer(?) around the train stations which follows the roads nearby. My thought is to make a visualization of how far from the trainstation you can get in 5, 10 and 15 minutes from the stations in normal walking pace - when following the roads.

Do anyone know how or if this can be done? I have managed to create a buffer with a linear value.

I am using the dataset called "Kartdata_Norge_UTM33_N50_FGDB" from kartverket (Norwegian mapping authority) in ArcMap 10.3.

The Network Analyst extension and license are installed.


texts: How does Google Earth calculate the inverse geodesic problem, in particular the forward azimuth?

texts:
Comparing the result of Vincenti’s formulae as implemented in https://geographiclib.sourceforge.io/cgi-bin/GeodSolve with the output of Google Earth Pro I find a striking difference in results for azimuth. GeodSolve will calculate the distance between the base of the flagpole on top of the White House at N 38° 53.862 W 077° 02.193 and the base of the flagpole on top of Buckingham Palace at N 51° 30.081 W 000° 08.513 as 5,912,893.79 meters with a forward azimuth of 49.387°, whereas Google Earth Pro (using the ruler) will give 5,912,887.89 meters with a forward azimuth of 49.36°. While the difference in distance is negligible, the difference in azimuth is not.
Taking the first point (WH) as start and making a projection in Google Earth using the results of GeodSolve will yield a point in Battersea Park at approximately N 51° 28.748 W 000° 09.222 more than 2,600 meters away from Buckingham Palace. I checked on the impact of height settings but taking an absolute height of 3m as opposed to “at ground” does not make a noticeable difference on the azimuth and only very few meters on the distance. Could it be that the azimuth is calculated as in great circle spheric geometry?


texts: How to change getLegendGraphic option default?

texts: I have several geoserver layers that have only one rule.  By default, they do not get a text label in the legend.  I know if you use LEGEND_OPTIONS=forceLabels:on in the getLegendGraphic call, it will label it.  Is there any way to make that the default for my geoserver installation?


texts: Tool to convert Mike Urban network (shapefiles) to EPAnet network

texts:
This might be a too localized question, but I hope I'm not the only one here doing pipe network analyses in Mike Urban.

I have a model in Mike Urban, and the client wants it in EPAnet format (since the rest of their models are in EPAnet). The network representation is not very similar. MU uses shapefiles and EPAnet, from what I understand, uses (at least) three different tables in text file format. One for node ID's and Z coordinates, one for node ID's and coordinates and one for link ID's and their corresponding breakpoint coordinates. There's some more stuff in each of these tables, but lets not worry about that for now.

Is there a tool, of any kind, that converts shapefile networks (from MU) to EPAnet style networks?

If not, what could be an efficient way of extracting necessary information from the shapefile? I have no problem extracting xyz coordinates for points (eg. nodes), but how can I get the breakpoints for lines into a text file/other table? Since there will be a lot of pipes, Python solutions are preferred.

Thank you!


texts: How to stop and then start a map service using the model builder

texts: Is it possible to stop a map service using a geoprocessing tool in the model builder ? If not How can I stop a map service using a model builder ?


texts: Using QGIS how to extract single polygon co-ordinates in multiple polygon data?

texts:
I have polygon coordinates for complete India at district level and I want to just use 2 or 3 districts out of complete polygon.

For example : I want to just select district Pune and Nasik and need this two polygon coordinates, but there are 100 district polygon boundaries in .shp file.

How to extract only selected district polygon coordinates ?


texts: Adding two raster datasets that overlap while maintaining extent of larger raster dataset?

texts:
I am using ArcMAP 10.0

How do I add two raster datasets that overlap while maintaining the extent of the larger raster dataset.  I have thought of adding zero values to the smaller raster to then be able to add them both together in raster calculator.

Is anyone familiar with a better way to achieve this?

My data is deposit thickness data, and I want to preserve the extent of the larger raster which the smaller raster is confined within when adding their elevations together - as currently when adding the data together I just get an output raster where the data overlaps.


texts: Is there a Guide to Translate ArcPy code to GDAL equivalent?

texts:
Due to cost-cutting measures I've been asked if I can migrate ArcPy code to a GDAL equivalent.

I'm fairly new to GIS, ArcPy, and GDAL.  Is there a website or guide which caters to providing equivalencies between ArcPy and GDAL?

I am not asking for any specific script examples tailored to geoprocessing shape files.  I'm looking for an overall guide.  Or something in that direction.


texts: Integrate MarkerClusterer with Google Map V3

texts:
I'm having trouble getting the MarkerClusterer into my current Google Map (which has taken a long time to get this far!!). How can I combine the two? I'm using V3 of the api.

Here's the example I'm working on http://www.mediwales.com/mapping

Here's the MarkerClusterer code:

var center = new google.maps.LatLng(37.4419, -122.1419);
var options = {
  'zoom': 13,
  'center': center,
  'mapTypeId': google.maps.MapTypeId.ROADMAP
};

var map = new google.maps.Map(document.getElementById("map"), options);

var markers = [];
for (var i = 0; i < 100; i++) {
  var latLng = new google.maps.LatLng(data.photos[i].latitude,
      data.photos[i].longitude);
  var marker = new google.maps.Marker({'position': latLng});
  markers.push(marker);
}
var markerCluster = new MarkerClusterer(map, markers);

Update: I've attempted to add the clusterer to my current code but it doesn't seem to work. Places[i] doesn't seem to feed into the clusterer.


texts:
I have published a wfs layer in geoserver, and displayed that on the google map using openlayers. By using ol3 i am able to select the polygon, but iam not able to get the coordinates of the polygon. For selecting the polygon my code is

    var interaction;
    var select = new ol.interaction.Select({
        style: new ol.style.Style({
            stroke: new ol.style.Stroke({
                color: '#FF2828'
            })
        })
    });

              interaction = new ol.interaction.Select({
        style: new ol.style.Style({
            stroke: new ol.style.Stroke({color: '#f50057', width: 2})
            })
    });

    map.addInteraction(interaction);

How to get the cocordinates of the selected polygon in geojson format?


texts: Clipping shapefile with DEM using ArcMap

texts:
I have a layer with shapefiles that represent watersheds. I'd like to exclude watershed polygons (or sections of these polygons) that fall below 2500 meters.

Is there a way to do this by clipping portions of the shapefile with a digital elevation model?

This would be for ArcMap 10.7.1.


texts: Changing number of points extracted on a linestring using `raster::extract` in R

texts:
I am trying to use a linestring to extract values from a raster using raster::extract. In the example below, I can only extract 123 points. Is there a way to extract more or less points without changing the line. I can't seem to find an argument that does this.
library(raster)
r <- raster(system.file("external/rlogo.grd", package="raster"))
x=c(5, 95)
y=c(20, 50)
line = SpatialLines(list(Lines(Line(cbind(x,y)), ID="a")))

ext <- raster::extract(r, line)[[1]]
nrow(ext)
# [1] 123


texts:
I'm trying to publish a map created in qgis desktop via qgis-web-client.

The left hand bar loads with the layers I would expect visible but I never actually see any mapping data just a grey screen with "Loading Map....".

If I make the call to the qgis server manually via a web-browser e.g.

http://adddress/qgis/qgis_mapserv.fcgi?map=/var/www/qgis-web-client/projects/test/filename.qgs&&_dc=1394809120649&SERVICE=WMS&VERSION=1.3&REQUEST=GetCapabilities

I get what looks like valid XML returned.

In the chrome dev/javascript console I'm seeing error:

Uncaught TypeError: Cannot read property 'abstract' of undefined

From webgisinit.js line 1680

var thisAbstract = wmsLoader.projectSettings.service.abstract;

Anyone seen this before or got any clues how to fix?


texts: Reading data from raster band with Java and GDAL?

texts:
After using GDAL with Python for a year or so now, I'm getting back into Java and trying to figure out how to use GDAL with it. I am wondering if anyone has had experience with reading data from a raster band with GDAL and Java. I can read the dataset and band in as follows:

Dataset theDataset = gdal.Open("k:/dem.asc", gdalconst.GA_ReadOnly);
Band theBand = theDataset.GetRasterBand(1);

...but I am not sure how to read data from theBand. From the API documentation it seems that I have to set up an array first, which means I have to know how many elements are in the block of data that I want to extract from the band. Like this:

double[] arr = new double[3];
theBand.ReadRaster(0, 0, 3, 1, arr);

I guess it works fine if you only want to extract one row of data, as above, because you end up with a one-dimensional array. Referencing elements is pretty easy. But if I want to extract a multi-row block of data, it has to be forced into a one-dimensional array. Something like this:

double[] arr = new double[9];
theBand.ReadRaster(0, 0, 3, 3, arr);

Then if I want to reference an individual raster cell in the array I just extracted (for example, if I want to process the block cell-by-cell), I need to convert my two-dimensional grid coordinate (proxy of easting and northing) into a one-dimensional integer array index.

I guess it's not too difficult to write a little method to convert from a grid reference to essentially a vector reference, but it just seems kind of convoluted, especially after being used to using grid references directly in Python-numpy-GDAL.

Is this the only way to extract data from a band in Java/GDAL?

From what I can tell, it is, but I am hoping there is a more intuitive way of processing large blocks of data.


texts: How to improve the displaying of a multiband raster in QGIS?

texts:
I have a multiband (WV3) satellite image that I'm displaying in false color (R=NIR, G=Red, B=Green) on QGIS 2.14. The quality of the results I'm getting via simply changing the band rendering in Layer-> Properties -> Style -> Band rendering --> Render Type: Multiband color are substantially different, and lower than using the SCP-plugin:

Regular rendering via layer properties:

Band rendering via SCP:

I need to do this programmatically via PyQGIS, so I can't use the plugin. I have also tried to play around with raster enhancement buttons but no success. Does anyone know how to improve the quality of band rendering?

UPDATE As pointed out by @user30184, the issue was the fact the max values were not updated in the band rendering (pressing the 'Load' button in the band rendering window). The problem remains - how to update the rendered min/max values via PyQGIS? My current code:

lyr = iface.activeLayer() #An multiband (8 band) raster

lyr.renderer().setRedBand(7)
lyr.renderer().setGreenBand(5)
lyr.renderer().setBlueBand(3)

if hasattr(lyr, "setCacheImage"):
    lyr.setCacheImage(None)
lyr.triggerRepaint()

The 'Load' action that needs to be completed in PyQGIS:


texts:
I am trying to add many layers to a project in QGIS. After adding approximately 80 layers, QGIS returns this error: "not a valid or recognized data source". There is nothing wrong with the shapefiles themselves, they will all open if done individually. It seems to just be after a certain number that it won't allow any more. The error message seems strange to me. Is there a limit to the number of layers a single project can have?

Thank you for your help


texts: Building pie chart map based on square shape, starting with land use soil shape in QGIS?

texts:
I'm a beginner and I usually use Qgis.
I would build a map made by a pie chart pattern, as in the picture attached. The pie chart returns the values in percentage of the land use soil. I would like that the system shows the percentage of the land use coverage in each quadrant of the grid.
E.g.  Square XXX 10% urban, 21% forests etc. for each quadrant.

I decided to do it this way:
1 - I build a grid of squares;
2 - Intersect squares with the use of the soil;
3 - ????

As output result I obtain a shape that is the intersection between the grid of the squares and the geometries of the land use soil.
I don't know how go ahead working.
I need to figure out how to reassemble the square keeping the land-use values.


texts: Joined tables and now field calulator is greyed out

texts:
In ArcMap (10.2), I accidentally calculated the wrong field outside of an edit session, so there was no way to reverse my action. I restored a previous version (shadow copy) and joined the tables together using the ObjectID. However, the field calculator was greyed out when I tried to restore that field's data. Why and how do I fix it?

Fundamentals of Field Calculations have a section for troubleshooting why these buttons are greyed out, but these do not fix the problem.


texts: Is it possible to access the Bing Maps "dark" mode in the ArcGIS JS API?

texts:
Stamen Designs have made a beautiful Hurricane Tracker application which uses a "dark" mode of the Bing Maps backdrop.

This is perfect for directing people's attention to the overlay rather than the background.

Is it possible to access this "dark" version of Bing Maps within the ArcGIS Server JS API?


texts:
Is there a way to export a map from QGIS with popup images?

Within QGIS, I can create a point layer with popup images using the html maptips (using Nathan Woodrow's method https://nathanw.net/2012/08/05/html-map-tips-in-qgis/ ).

I would like to export a map to a common file type that preserves this functionality.

My end goal is to be able to share maps with a non-gis users who could view the photos.

I created something like this in a PDF file that used buttons, but it was pretty clunky, and wasn't as smooth as exporting from QGIS.  Each button had to be created individually, then the photo was pasted into it.


texts: ArcMap: How to flip digitized directions of geometry based on route measures

texts:
I have a feature class that has a schema like:

routeId, fromMeasure, toMeasure, attr1,attr2..etc

within a given routeId, the segments flip directionaility many times and are not consistent. I would like to have all the geometry directions match increasing the measure of the route.

I am pretty sure I want to do the  "Set direction as M" function here https://desktop.arcgis.com/en/arcmap/10.3/guide-books/linear-referencing/about-setting-route-measures.htm

I am not use to arcMap so it's not clear to me how I can process all 20,000 segments at once and have it flip the geometry. I believe I need to convert my segment file to a route layer, and then somehow run the route editor and then somehow go from a route layer back into a segment file, perhaps with events or something.

Any guidance would be appreciated


texts:
I'm making my first ModelBuilder and my issue is about trying to get my firstly made "dissolve.shp" into an ADD JOIN like an input layer name not like a join table. (as you can see in the picture there is only that option)

Where my normal way is: Dissolve -> Add Join -> Export data -> Save like a new .shp... I don't know how to do this process by ModelBuilder.


texts: Setting up ArcSDE for Oracle Spatial connection parameters?

texts:
Just testing out the new QGIS2 install here and trying to connect it to our huge SDE/Oracle spatial database. Im struggling with the connection parameters. Here is the connection properties to the databse from within ArcGIS10.1:

Instance = sde:oracle10g:g01l10g

I can't figure out how this maps across to the connection fields in the QGIS2 popup. Can anyone help me please?


texts: Include landmarks in routing instructions based on OpenStreetMap data in ArcGIS for Desktop?

texts:
I created an OSM network dataset in ArcGIS 10.3.1 with the tool from the OSM toolbox. Now I like to include landmarks within the routing instructions. If the network dataset is built with the help of the network analyst there is a tab where the landmark layer can be defined. But in the OSM tool there is no such tab.

Does anybody have an idea?


texts:
I'm having trouble getting my SQL query to work with three tables. I'm new to SQL to I've already spent considerable time trying to understand the documentation and making some sense of it.
I'm using Postgres 8.3 with POSTGIS.
I have three tables:

A table holding polygon geometry called subshed2
A table holding point geometry called withdrawalsites with a primary key ic_site_id
A table holding values I want to summarize on a monthly basis, both within a single specific year as well as the average over all years of each month. Table name is withdrawals.
This table has a record with the same ic_site_id as found in the second table above, but there are many records for each ic_site_id. This table has my dates in timestamp form. Here's a few records:

I need my final output to look like this:

The shed_id parameter has already been supplied from another function, so it's being supplied to this query as ${param.shed_id}.
I'll be letting my user select a specific year to summarize to get the first table (and likely suppling it using another ${param.requested_year} ), but then I also want to average the values by month over all the years in the dataset. I might have to run a separate slightly altered query to get the second table, but I was hoping I could do it in one query.

My query looks like this:

SELECT
    to_char(withdrawals.begintime, 'YYYY') AS year,
    to_char(withdrawals.begintime, 'MONTH') AS month,
    sum(withdrawals.withdrawalvalue_mgd) AS mgdsum,
    withdrawals.begintime AS date
FROM
    withdrawalsites
    INNER JOIN subshed2
    on (st_contains(subshed2.the_geom, withdrawalsites.the_geom))
    INNER JOIN withdrawals
    on withdrawals.ic_site_id = withdrawalsites.ic_site_id
where
    subshed2.id = '${param.shed_id}',
    withdrawals.intervalcodes_id = '2'
GROUP BY to_char(withdrawals.begintime, 'YYYY MONTH')
ORDER BY YYYY MONTH

When I run the query, it gives a syntax error on a comma, but it doesn't tell me which one, but I suppose there's a lot more wrong with this query. Any ideas on how to query these three tables correctly, or at least point me in a direction to investigate? It's also unclear to me if I should instead be using the date_part function.

Final Solution:
I solved this problem using some of the thoughts Mike provided. I couldn't get the first function to work, but realized I was getting years with decimals (i.e. 2008.0) and was getting the month number, not the month name. I then wanted to sort those month names correctly. Here's what my worked in the end:

SELECT
    (EXTRACT(YEAR FROM withdrawals.begintime))::int AS year,
    to_char(date_trunc('month', withdrawals.begintime), 'MON') AS month,
    sum(withdrawals.withdrawalvalue_mgd) AS mgdsum
FROM
    medford.withdrawalsites
    JOIN medford.subshed2 ON
    ST_Contains(medford.subshed2.the_geom,medford.withdrawalsites.the_geom)
    JOIN medford.withdrawals ON medford.withdrawals.ic_site_id =
        medford.withdrawalsites.ic_site_id
where
    subshed2.id = '${param.shed_id}'
AND
    withdrawals.intervalcodes_id = '2'
GROUP BY (EXTRACT(YEAR FROM withdrawals.begintime))::int, date_trunc('month',
    withdrawals.begintime)
ORDER BY (EXTRACT(YEAR FROM withdrawals.begintime))::int, date_trunc('month',
    withdrawals.begintime)


texts: Why does the qgis LTR release date lag behind the current release?

texts: This is just a simple curiosity but I haven't been able to find an answer. Qgis 2.14 will be the LTR, but when it was released it was posted as the current release with 2.8.7 still as the LTR. Will the LTR be 2.14.something in the future? What is the reason that it isn't posted as the LTR along with the current release?


texts: How to decide geohash length based on zoom level?

texts:
I want to decide geohash length based on zoom level, how can I do that because so far I am unable to find any standard way to do so?

  For example, current zoom level is 10 and geohash is "qyu7652hsze2" of
  length 12. And, if I zoom out then my zoom level becomes 7 then I
  want to trim my geohash length by triming 2 ending characters and
  becomes "qyu7652hsz" of length 10.

I want to apply best practice for this, does anybody know the same?


texts: point.in.poly (spatialEco) returns non-intersecting points

texts:
I just updated the spatialEco package for R, and now the point.in.poly function behaves differently from previous behavior. Previously, it returned only points that intersected at least one of the polygons in the polygon layer. Now, it returns all points, including those that do not intersect with any polygons (the polygon fields are NA for the non-intesecting points). It seems that an extra step is now necessary to retain only the intersecting points, which doesn't seem desirable to me. It is not too burdensome to add the extra step, but it does break old code that uses point.in.poly. Is there a way to retain only the intersecting points in the output from the point.in.poly function?

Below is some example code that illustrates what I mean.

library(rgeos)
library(sp)
library(spatialEco)
library(raster)

coords <- matrix(c(0,0,0,10,10,10,10,0,0,0),5,2,byrow=TRUE)
SA <- SpatialPolygonsDataFrame(SpatialPolygons(list(Polygons(list(Polygon(coords)),1))),data.frame(pol=1))

pts <- spsample(SA,50,"regular")
pts <- SpatialPointsDataFrame(pts,data.frame(ptid=1:length(pts)))
smallSA <- gBuffer(SA,width= - 1,byid=TRUE)
grd <- raster(extent(smallSA),resolution=2)
grd <- as(grd,"SpatialPolygons")
grd <- SpatialPolygonsDataFrame(grd,data.frame(pol=1:length(grd)))
plot(SA)
plot(grd,add=TRUE)
plot(pts,add=TRUE,col="blue",pch=20)
pts2 <- point.in.poly(pts,smallSA) # includes non-intersecting points!
plot(pts2,add=TRUE,col="red",pch=20)
length(pts)
length(pts2)

# remove non-intersecting points
pts3 <- pts2[-which(is.na(pts2$pol)==TRUE),]
length(pts3)
plot(pts3,add=TRUE,col="blue",pch=20)