texts: Getting count of like features types based on field attributes using ArcPy?

texts:
I'm trying to create an update cursor which will take the count of the number of like features in a field and then use that number in a random number generator. The field would look something like that shown below:

**Field**
X
X
X
X
X
Y
Y
Y
Y
Y
Z
Z

As had been suggested by someone on here previously, I tried using collections however, I find that the result is not always accurate with the count total. For example I may have values greater than the total number of like features. I managed to get something working however, now I need something which constrains the random number without replacement and I feel the only way to do this is using collections. My code is shown below:

def sumdict(listvals):
    ddict = collections.defaultdict(int)
    for val in listvals:
       ddict[val]+=1

    return ddict

dic = sumdict(row[0] for row in arcpy.da.SearchCursor(MyFC, (Field1)))

with arcpy.da.UpdateCursor(MyFC, (Field1, Field2)) as cursor:
    for row in cursor:
        row[1] = sample(1, dic[row[0]])
        cursor.updateRow(row)

Below is an example of what the output should look like, naturally the random numbers will vary.

Field1  Field2
X       2
X       3
X       1
X       5
X       4
Y       3
Y       5
Y       4
Y       2
Y       1
Z       2
Z       1


texts: How to aquire a topographical map of Britain?

texts:
First of all, let me say that I know next to nothing about GISs.

I'm trying to find / make / aquire a topographical map of britain. I want to use it as a height map in Dwarf Fortress so it needs to be gradient grayscale, without labels and bmp format.

I've done some googling and so far haven't found any maps I could use, free or otherwise. Which seems absurd, so I must be looking for the wrong thing?

It seems possible to retrieve the height from Google Maps API. My first thought is to make a program that requests the height for every coordinate from Google. Leave the program running for a week and then construct an image from the data.

Though I thought I should ask here to see if people look at me like I'm crazy!?


texts: Opening GeoJSON file with mixed geometries in QGIS?

texts:
If I have a GEOJson file, for example a file named test.js with inside the following:

{ "type": "FeatureCollection",
"features": [
{ "type": "Feature",
  "geometry": {"type": "Point", "coordinates": [102.0, 0.5]},
  "properties": {"prop0": "value0"}
  },
{ "type": "Feature",
  "geometry": {
    "type": "LineString",
    "coordinates": [
      [102.0, 0.0], [103.0, 1.0], [104.0, 0.0], [105.0, 1.0]
      ]
    },
  "properties": {
    "prop0": "value0",
    "prop1": 0.0
    }
  },
{ "type": "Feature",
   "geometry": {
     "type": "Polygon",
     "coordinates": [
       [ [100.0, 0.0], [101.0, 0.0], [101.0, 1.0],
         [100.0, 1.0], [100.0, 0.0] ]
       ]
   },
   "properties": {
     "prop0": "value0",
     "prop1": {"this": "that"}
     }
   }
 ]
}

how to open it with QGIS?

I open QGIS, press "add vector layer", choose "File" but I am not able to open it.

If I choose GEOJson format, my file is not showed in browser window.

If I do not specify input type, I get the error: "your file is not valid or recognized data source".

How have I to format this file?

It seems there is a lack of documentation about this topic.


texts: Joining MS Access (accdb) table to Feature Layer in ArcGIS for Desktop?

texts:
I would like to join a table from a microsoft access database (.accdb) to a featurelayer in ArcGIS 10.1.

How do I associate tables from that database with existing layers that have a matching attribute?


texts: Extracting coordinates from shapefile using QGIS?

texts:
I have a shapefile and spatial coordinate file. I wish to extract information for each coordinate from shapefile.
Is it possible in QGIS?

I have 19 coordinate points and I'm trying to extract soil type data for each coordinate. Shapefile is polygon and coordinate file is .txt. I'm trying to extract information precisely for each coordinate.


texts: Produce background map images for plotting data with Google Earth Engine

texts:
I am submitting a manuscript for journals that do not accept background tile images from ArcGIS or similar services.

I thought it would be easy to produce your own images with GEE. I wrote this code combining examples I found on the web. I export this image and end up using QGIS to select bands, edit and produce a final pdf, svg or png.

However, I would like to:

export a thumbnail, pdf, png image that is ready to be used as background. Ideally an image that is no more than 20Mb, given any polygon.
spice-up the background image so it looks more similar to the tiles of any GIS service.

My ideal result is an image that looks pretty much like a screenshot of Google Earth (or maybe even one in black and white, greyscale, so that data plotted on top would visualize better, straight out of GEE and without having to use QGIS.

This is the code:

var geometry =
    /* color: #d63000 */
    /* displayProperties: [
      {
        "type": "rectangle"
      }
    ] */
    ee.Geometry.Polygon(
        [[[-81.29700581101451, 14.41901072292579],
          [-81.29700581101451, -7.683476214952404],
          [-63.718880811014515, -7.683476214952404],
          [-63.718880811014515, 14.41901072292579]]], null, false);

// Function to cloud mask from the pixel_qa band of Landsat 8 SR data.
function maskL8sr(image) {
  // Bits 3 and 5 are cloud shadow and cloud, respectively.
  var cloudShadowBitMask = ee.Number(2).pow(3).int();
  var cloudsBitMask = ee.Number(2).pow(5).int();

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

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

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

// Map the function over one year of data.
var collection = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR')
    .filterDate('2016-04-01', '2018-7-31')
    .map(maskL8sr)

var composite = collection.median();

var VIS = {bands: ['B3', 'B5', 'B7'], min: 0, max: 0.747};

// Display the results.
Map.addLayer(composite,VIS ,"compt2");

//export image

Export.image.toDrive({
  image: composite,
  description: 'background-image',
  scale: 10000,
  region: geometry,
  maxPixels: 1000000000
});


texts: Heat Map using a Date Range in QGIS 3.16.3

texts: I'm new to QGIS and I'm trying to make a Heat map with roughly 2000 data points using a date range. What I'm wanting is for the newer date to be the "hot" portion of my map. I've tried to figure this out on my own and I've given up, is this even possible?


texts:
I've installed two OSM maps, but I cannot use them in Garmin. See the window of Utilities / Manage Map Products:

The maps are not lockable, but I cannot use them in Garmin... What shall I do to be able to upload them to my GPS? Thanks a lot!


texts:
I have a river geometric network built that is fully operational with flow direction working.
In ArcMAP when I trace upstream I get a selected feature that I can convert into a feature class (Data > Export > Selected Features).

How do I do this in a model as I have 50,000 Flags?

Here is the model.

The output appears to be a "Ouput Group Layer Name". How do I get this to a feature class in the model above?

In the help menu it states that I should use "Make Feature Layer" and "Select Data" but I have been unable to make this work.


texts:
I would like to include some state agency data in an  application I am creating.  The agency will not share the data so I have to find another approach.  What I am interested in are the points on the webmap below along with the attributes of each point (onClick).  I know I can easily digitize the points but would prefer not to take this approach

http://fwp.mt.gov/gis/maps/huntPlanner/?bma=true&region=7

Is there an easier way instead of digitizing?  Is there a way to access and export the dataset from this webmap instead of digitizing?

thanks,

-j


texts: Data Interoperability Extension Live connection to WFS from arc10

texts:
I am trying to make use of a live non-cached Data Interoperability (DI) connection to consume a WFS services served by geoserver. The data displays in arcmap fine, and using fiddler to monitor the HTTP traffic, I can see the live http requests and responses as I pan and zoom the map. However, viewing the attribute table, or other operations such as select by rectangle are not working. Looking at the requests and responses with fiddler, the problem seems to be with the format of the filter xml embedded by arc10/DI within the request.

Arc/DI is embedding something like this:

<ogc:Filter>
<ogc:PropertyIsEqualTo>
<ogc:PropertyName>SomeProperty</PropertyName>
<ogc:Literal>100</Literal>
</ogc:PropertyIsEqualTo>
<ogc:/Filter>

The server is failing to parse this filter xml.

When I manually edit the xml in the requests within fiddler to this (ie drop the ogc namespace)

<Filter>
<PropertyIsEqualTo>
<PropertyName>SomeProperty</PropertyName>
<Literal>100</Literal>
</PropertyIsEqualTo>
</Filter>

the server understands the request and send the feature gml back.

Is this is a bug on the Arc/DI side?


texts: Error while loading shared libraries libqgis_core.2.9.0

texts:
I have complied QGIS version 2.8.1 from source, according to the instructions from https://github.com/qgis/QGIS/blob/master/INSTALL

Now after successfully running ccmake .., make and make install commands and I move the install directory, I get the following error in terminal, when executing the generated binaries

./qgis

./qgis: error while loading shared libraries: libqgis_core.so.2.9.0: cannot open shared object file: No such file or directory

Any thoughts?

Edit
What worked:

After some search on the internet I found the following procedure that works permanently

~$ cd {path}/lib/     # change directory to qgis installed app lib directory {}
lib\$ pwd > qqis.2.8.conf  # get current working directory and paste to a file
lib\$ sudo mv qqis.2.8.conf /etc/ld.so.conf.d
lib\$ sudo ldconfig

Now compiled qgis installed at a custom location, will load all the time ./qgis is executed


texts: USA Historical Satellite Imagery for an Online Application

texts:
I am looking for different options to add historical satellite imagery to our online-based application.  We currently incorporate Bing Maps into our application, but cannot get historical imagery.  I am essentially looking for something like Google Earth with its historical imagery slide bar.  I am currently looking into using services from Google Earth and Digital Globe.

I have also looked into NAIP imagery.  I tried downloading the .sid files for one state and one year, but am having trouble.  I downloaded the data from Geospatial Data Gateway in county files.  I was trying to clip the county .sid files because of the noise around each .sid file.  The noise that I am referring to is not the black border, but pixels that received none zero values in the border.  In the end I need to have a seamless image.  I am not sure where I would go from this spot of having clean county imagery.  I need an opensource way to create a mosaic dataset.  If I stick with this option, I think I would be looking at Terabytes of imagery.

What are my options?  What else is out there?  What could I be doing differently?


texts: Convert TIF to PNG without normalizing alpha

texts:
I'm trying to make some height maps with elevation data from https://earthexplorer.usgs.gov/. Since they give you assets in separate TIF tiles, I need to convert them to PNG and stitch them together to the area of land that I'm trying to extrude in Blender. The only problem is, once I run

gdal_translate -of PNG -ot Byte -scale input.tif output.png

with QGIS's terminal, the output png has all it's values 'normalized'. It end up looking like this in photoshop when I try to put the tiles together.

I partially can get around that by tediously adjusting contrast for each tile so they look normal to the eye, but it's always got artifacts in it. Is there a way to avoid this when converting with gdal?


texts:
I have a new task where the client is asking for a .dwg file of layer labels. My first thought would be to create an annotation of the labels and then try and use the Export to CAD tool in ArcMap. Does this sound like an appropriate workflow?

I am using ArcMap 10.1

Ultimately I would like to create a model and tool to do this.


texts: Unable to delineate inundation map from HEC-RAS to Hec-GeoRAS?

texts:
I am trying to delineate the flood map generated from HEC-RAS by importing them into the HEC-GeoRAS. However, I got this error:

How will I resolve this?


texts: Is the Qspatialite 6.0.8 font for SQL free text supposed to be huge?

texts:
A week ago I upgraded my Qspatialite plugin in QGIS 1.8, which I frequently use in parallel with Spatialite_GUI, to get the best of both Sqlite worlds.  I immediately noticed that my SQL free text entry in Qspatialite is huge, as I type.

I think I prefer larger fonts than many folks, especially when typing SQL by hand.  This font size happens to be more than large enough to suit my tastes.

Is the Qspatialite 6.0.8 tab for typing in SQL by hand supposed to be huge?


texts: How to transfer attribute of nearest Point to Line feature?

texts:
Feature class 1: Addresses (Point)
Feature class 2: Road centerlines (Polyline)

Which is the most appropriate tool for transferring the Address attribute (e.g. Road No. or Bldg) to the nearest Road line segment?

Tools reviewed and rejected

Attribute transfer ( Spatial adjustment)
Near/Near table


texts: I want to perform a time series analysis of Landsat NDVI. I could not find a prepared product online (except in GE Engine but it is no possible to download). Thus I would like to calculate the NDVI in batch, possibly applying a rule to keep the filename (I need the reference of the date).


texts: Alternative to Delete Polygons And Extend Lines of Production Mapping?

texts: Maybe someone knows an alternative for the Delete Polygons And Extend Lines (Production Mapping) tool or a simple way to achieve the same result but without requiring a Production Mapping license?


texts: Is there a way to use formula in filtering a dataset in Mapbox Studio?

texts: I am very new to Mapbox and was wondering if anyone knew how to filter a layer based on a formula defined for a field in that layer. For example, I have a 'depth' field in the data and I want to display all the lines that have depth <=500 and is completely divisible by 100 (with no remainder). So the result should look something like this: show only the lines that have depth 100, 200, 300, 400, and 500. There is something like this for styling the layer as a text symbol but for symbolizing the layer as a line, the only option I get is a filter where I need to choose a specific value and not use a formula.


texts: How do I turn a GeoTIFF file into a bmp file?

texts:
I need to get GeoTIFF files transformed into bmp files so that I can display them as apart of some software that I am working on. Unfortunately, I am a total newbie with GIS and really don't know much about geography. I am working in C++ and need to have this included in what I am doing. A dll file would be great, although it is my understanding that we might need to include the source code from which the dll is created. I am wide open to any other technique that could be used as well. Can anyone help me with this?

I should also note that I am working on a Windows machine.


texts: How to get Euclidean Distance to include all input features?

texts: I've run Euclidean distance on a water source dataset. Water is represented as polygons and I run the Euclidean Distance tool in ArcMap. It appears to work great but when I zoom in I notice that some of the smaller water source polygons have not been included as results show the area as far from a source. Is there anyway I can make the tool include every polygon, not matter it's size in the analysis?


texts: Changing SVG marker color using Expression String Builder

texts:
What is the proper way to use the Expression String Builder to set a specific color for an SVG symbol?
I found some geology svg symbols through this StackExchange answer that I had been looking for https://gis.stackexchange.com/a/269648/65188. The documentation for the USGS symbol library states, "In addition, modifying the size, color, and (or) lineweight of an existing symbol to suit the needs of a particular map or output device also is permitted", however, I am not able to change the color of my svg symbol. The Fill Color and the Stroke Color options in the Layer Styling panel are not automatically available, so I believe using the Expression String Builder is the way to go. Below are some of the unsuccessful expressions that I tried.

color_rgb( 179,14,238)
@symbol_color =color_rgb( 179,14,238)
'#b709f6'
'[179,14,238,0]'

Using QGIS v. 3.14 on Windows 10


texts: Coordinate Reference System issues between SHP and Raster file in QGIS?

texts:
I am having an issue with the CRS of a given Aster DEM.
I have downloaded a Aster DEM from www.oziexplorer3.com which comes with the EPSG:4326 - WGS 84 as the CRS.

This is a DEM from Gabon, and all the other data I have for the project I am working on has the CRS defined as EPSG:32632 - WGS 84 / UTM zone 32N.

I have to create some contours from the DEM in order to use them in another software for the elevation definition of a water network. When I open the DEM and SHPs I have in QGis they perfectly overlap (although they have different CRS) and so do the contours I extract from the DEM.

But then, when I import the contours into the other software, because they have different CRS, they do not overlap and so I can't use the newly created contours in order to define the water network's elevations... and if I change the CRS in QGis (for the DEM or for the contours) they do not overlap even in QGIS...

Is there a solution for this problem?

PS - please know that I am quite inexperienced with GIS softwares, but you might already noticed that...


texts:
I am trying to find a way to get certain service properties from multiple published map services. One thing I am looking for is all services with "Feature Access" capabilities enabled:

I've managed to get the errors and warnings with a Python script by creating SDDrafts as it is described here.
In code example 7 I found this passage:
# Turn on feature access capabilities
configProps = doc.getElementsByTagName('Info')[0]
propArray = configProps.firstChild
propSets = propArray.childNodes
for propSet in propSets:
    keyValues = propSet.childNodes
    for keyValue in keyValues:
        if keyValue.tagName == 'Key':
            if keyValue.firstChild.data == "WebCapabilities":
                keyValue.nextSibling.firstChild.data = "Query,Create,Update,Delete,Uploads,Editing"

So I searched the SDDraft files for "Query,Create,Update,Delete,Uploads,Editing". But this is not only included in files with Feature Access enabled but also in all other files I've checked.
How can I get the right information? Would GetCapabilities be an option?


texts: OpenLayers 3.15.1 - styleFunction is executed every time I move/change the map

texts:
I'm testing the new OpenLayers v3.15.1 and I've a simple question about styleFunction used to plot GeoJSON.
My code contains:

var vectorSource = new ol.source.Vector({
    url: <my geoJson URL>,
    format: new ol.format.GeoJSON({
        defaultDataProjection :'EPSG:4326',
        projection: 'EPSG:3857'
    })
});

var styleFunction = function(feature) {
    console.log('styleFunction');
    var stroke = new ol.style.Stroke({color: 'black', width: 1});
    var fill = new ol.style.Fill({color: 'red'});
    var geom_name = feature.getGeometry().getType();
    var styles = {
            'Point': new ol.style.Style({
                image: new ol.style.RegularShape({
                    fill: fill,
                    stroke: stroke,
                    points: 3,
                    radius: 7,
                })
            })
        }
    return styles[geom_name];
};

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

My GeoJSON contains more than 500 Point than running the map the first time the styleFunction is executed 500 times; moving the map, it is executed again 500 times and so on.

Looking the JavaScript console, I noticed that every time I move the map or I change the map zoom, the styleFunction is executed. Is this normal? Why the style should be rendered every time I move the map?


texts:
I have a dataset where multiple attributes need to be shown with the same symbol. I can achieve this using multiple 'PropertyIsEqualTo' operators nested within 'Or' logical operators, see below example;

    <ogc:Or>
    <ogc:PropertyIsEqualTo>
        <ogc:PropertyName>CURRENT_ST</ogc:PropertyName>
         <ogc:Literal>AB-LOC</ogc:Literal>
    </ogc:PropertyIsEqualTo>
    <ogc:PropertyIsEqualTo>
        <ogc:PropertyName>CURRENT_ST</ogc:PropertyName>
         <ogc:Literal>ABD</ogc:Literal>
    </ogc:PropertyIsEqualTo>
    </ogc:Or>

This seems a little inefficient, especially considering I have about 10 different symbol types each with about 10 contributing attributes.
Ideally you would just clean the attribute data up, or add a symbol field that grouped them together but that is not an option in this case.

Is there a more efficient way of doing this?

For example being able to list the attributes something like the below would be great;

    <ogc:PropertyIsEqualTo>
        <ogc:PropertyName>CURRENT_ST</ogc:PropertyName>
         <ogc:Literal>'AB-LOC','ABD','ABC','DEF'</ogc:Literal>
    </ogc:PropertyIsEqualTo>


texts: How to find the two furthest points between a group of points?

texts:
I have a table with 20000 records and with those columns Shop, longitude, latitude, geom

I would like to find a way to find the pair of ths farthest shops between a selection of shops

For example we have outletnames : A,B,C,D,E,F,G,H,I,J,K,L,M,N,O...

I would like to know what is the bigest distance between all those shops A B C D and E
A-B = 4 km
A-C- = 2 km
A-D = 1 km
A-E = 1,3 km
B-C = 20 km
C-D = 50 km ....

I would like to have an optimised query that would give me as a result 50 km

Thank you


texts: How to (automatically) extract a substring from attribute table

texts:
I have an attribute table that contains in a field/slot called "Description" data like this:

xyz: 11.02.2016 04:05 - 14.02.2016 04:03

0 km/h, 0 hauk/min , ODO 0m

xyz: 01.03.2016 04:01

0 km/h, 0 hauk/min , ODO 0m

0 km/h, 0 hauk/min , ODO 0m

In two of these three samples is a date and a time. Actually the first one has a time period.

I want to automatically extract the substrings that contain this data and if this data is not included like in the last sample then nothing should be extracted.

I was thinking that everything between > (from ) and < (from ) could be extracted.

How can I achieve that?

Please note that I have no coding experience with Qgis Python. I recently started to learn Python 3 online, but in Qgis Python is different. I use QGIS 2.12.


texts: Iterative point feature removal based on Thiessen polygon values

texts:
The background is: I have point feature class and Thiessen polygons created based on it. I wanted to know for each point Thiessen polygon area value, so I have calculated Thiessen polygon area value and added this information to point class. I have also sorted points based on Thiessen polygon areas in descending order, so the first point in a table has the biggest Thiesssen polygon area.

The task is to remove first point with the biggest Thiessen area, then rebuild Thiessen polygons (without deleted point), update Thiessen areas for points, sort them again, check first on the list, remove it and so on until the number of points that left is equal to some fixed number.

I am working in ArcGIS 10.1 and I am using ModelBuilder. Do you think I can achieve something like this in ModelBuilder?  How for example take the first point on the list?Or the point with max attribute value?And which tools should I use for iteration/loop?


texts: Calculating XY from Depth, Azimuth, Inclination using Python

texts:
I have been using the answer from here and trying to see what I am doing wrong, but just don't see anything.

dMD = D2 - D1
B = acos(cos(I2 - I1) - (sin(I1)*sin(I2)*(1-cos(A2-A1))))
RF = 2 / B * tan(B / 2)
dX = dMD/2 * (sin(I1)*sin(A1) + sin(I2)*sin(A2))*RF
dY = dMD/2 * (sin(I1)*cos(A1) + sin(I2)*cos(A2))*RF

X2 = X1 + dX
Y2 = Y1 + dX

The blue are the correct points, but the red is what i end up with.

I've tried to convert to radians as well, with similar results, offset from actual location. I've also tried several different projections

Is there something that i'm missing?

Edit.
So the mistake relating to dX and dY, which would explain the offset and although my inclination appeared to be in radians, it seems it wasn't.


texts: Adding new geoJSON polyline layer to mapbox styles using Android Studio

texts:
I am trying to add new geoJson polyline layer to mapbox in Android Studio. It works fine if I am adding it using the default style of mapbox, but if I set my own styles, the line is not drawn anymore in my app. Could anyone tell me what could be wrong?

mapView.getMapAsync(new OnMapReadyCallback() {
        @Override
        public void onMapReady(MapboxMap mapboxMap) {

            // Create the LineString from the list of coordinates and then make a GeoJSON
            // FeatureCollection so we can add the line to our map as a layer.

           if ("city".equals(spinner2)) {
                mapboxMap.setStyleUrl("mapbox://styles/johan88/cj1fdzcoo00012spgu1lldduk");
            }
            else if ("beach".equals(spinner2)) {
               mapboxMap.setStyleUrl("mapbox://styles/johan88/cj3cuj7a800032rmqsrhir1av");
            }
            else if ("adventure".equals(spinner2)) {
                mapboxMap.setStyleUrl("mapbox://styles/johan88/cj37g9zzz00092rqp7c55ufxt");
                // Customize map with markers, polylines, etc.
            }

            LineString lineString = LineString.fromCoordinates(routeCoordinates);

            FeatureCollection featureCollection =
                    FeatureCollection.fromFeatures(new Feature[]{Feature.fromGeometry(lineString)});

            Source geoJsonSource = new GeoJsonSource("line-source", featureCollection);

            mapboxMap.addSource(geoJsonSource);

            LineLayer lineLayer = new LineLayer("linelayer", "line-source");

            // The layer properties for our line. This is where we make the line dotted, set the
            // color, etc.
            lineLayer.setProperties(
                    PropertyFactory.lineDasharray(new Float[]{0.01f, 2f}),
                    PropertyFactory.lineCap(Property.LINE_CAP_ROUND),
                    PropertyFactory.lineJoin(Property.LINE_JOIN_ROUND),
                    PropertyFactory.lineWidth(5f),
                    PropertyFactory.lineColor(Color.parseColor("#e55e5e"))
            );
            mapboxMap.addLayerBelow(lineLayer, "Water labels");
    }
});


texts: aggregating using 'gdalwarp -average' resulting in incorrect values

texts:
I'm trying to aggregate data up from 0.05 to a coarser 0.25 resolution and am having trouble getting gdalwarp to work correctly. From comparing the original input data with the aggregated data and manually calculating the averages, I can see that the result I'm getting is not what I would expect. Many (not all) aggregated pixel values are just wrong (FYI, I have specifically been looking where the water meets land on the horn of Africa). Also, it's my understanding that gdal is supposed to ignore NoData pixels in the average calculation, but in some areas, the 25 pixels underneath the aggregated pixels were all NoData, yet the aggregated pixel reported a number, when it should have been NoData.

The command I used below:

gdalwarp -tr .25 .25 -r average -srcnodata -15000 -dstnodata -15000 /path/to/temp_CHIRPS_201512.tif /path/to/CHIRPS_201512_ave.tif

I can post the entire script if need be.

I'm using gdal 2.0.1 on linux.

Here is a link to the input file I used, in addition to the result I got by calling the above command: https://www.dropbox.com/sh/dk0c41638jg8fir/AAD9z9IXeEfJ2PyEUVnjG4aya?dl=0

I tried running the same command on a different dataset (I am doing the same process for several datasets) and it appears gdalwarp did in fact work.

Can anyone help me understand why it is working for one dataset and not the other?


texts: Advice on normalizing census data in choropleth maps

texts:
I'm trying to figure out the best way to normalize census data in two choropleth maps showing the population of people with hearing loss. Their geometries:

The states of the United States.
The census tracts of one state.

For the first map (the states of the United States), I have the following data fields:

Total population
Total population with hearing loss

Ways I can think of to normalize the data:

Divide the total population with hearing loss within a state by the state's total population
Divide the total population with hearing loss within a state by the total population with hearing loss in the United States

For the second map (the census tracts of one state), I have similar data fields:

Total population
Total population with hearing loss

Ways I can think of to normalize the data:

Divide the total population with hearing loss within a census tract by the census tract's total population
Divide the total population with hearing loss within a census tract by the total population with hearing loss in the state

Essentially the same question for both maps. I see mixed advice when I Google this. What is the significance of each approach, and is there a clear better choice?


texts: Is it possible to implement Geofencing in Openlayers?

texts:
We were planning to develop an application using Openlayers, which can display the incoming GPS signals as markers and create alerts , if the markers( positions) are outside a particular zone.

Is this possible through Openlayers code ? or do we have to use something JSTS Topology suite with Openlayers ?


texts:
When I try to launch my model (graphical modeler) in QGIS I get this message (in french) :

Cet algorithme ne peut pas être lancé :-(  Le modèle que vous essayez
de lancer contient un algorithme qui n'est pas disponible
:native:countpointsinpolygon

It translates roughly in the fact that the  :native:countpointsinpolygon  algorithm can't be used because it's not available.
I guess it is because I changed QGIS version (currently 3.10.13).
How can I convert it/export ? It was quite a big model, and I this feature (counting points in polygons) didn't disappear from QGIS.


texts: Are FME2013 and 64bit Background Geoprocessing of ArcGIS 10.1 SP1 compatible?

texts:
We are running ArcGIS 10.1 sp1 and have installed 64-bit background geoprocessing.  We are also upgrading to FME 2013.  I know in the past we have had to install the 32-bit version of FME due to compatibility issues - I think regarding writing to File Geodatabases.

Now that ArcGIS has 64-bit background processing available, am I able to install the 64-bit version of FME?  Or do I still need to install the 32-bit version in order to have full compatibility?  Does the 64-bit background geoprocessing even have anything to do with the compatibility?


texts: Is there a way to convert GeoTIFF to DEM?

texts:
Simple enough question. The NASA Aster data comes in GeoTIFF format, but the tool I'm using (Sketchup) only accepts DEM. Is there a way to convert?

I did a lot of searching, and only found the program GEOTIFF4, but it's obsolete, it's Windows-only, and source code is not available.


texts: Using a dictionary along with geoprocessing tools?

texts:
I am new to python and trying to create a script that will select specific parcels from a large parcel layer. Is it possible to create a dictionary with the names of the parcels as the key and the query used to find the parcels as the value and input those into a arcpy.Select_analysis with a for loop to loop through the dictionary? What would be the best way to set this up?

example

Import arcpy
import os  #importing operating system functions
arcpy.env.workspace = "c:\gis\python\other\parcels.gdb"
arcpy.env.overwriteOutput = True

properties{'johns pizza':"MUN in ('1323') and BLOCK in ('168') and LOT in ('20')",'chinese palace':"MUN in ('1326') and BLOCK in ('12') and LOT in ('100')"}`

for key in properties:
    arcpy.Select_analysis("parcels", key, value)


texts: how can I set the StartBounds, units, and tileSize correctly in OpenLayers?

texts:
how can I set the StartBounds, units, and tileSize correctly in OpenLaters?I can not make sure this parameters,for I use jsonp,and it will dynamticly read layerinfo from esri rest resources.

 var maxExtent = new OpenLayers.Bounds(102.9702000357959, 26.573810737410064, 110.21822155809812, 29.11729567697843);
    baseLayer = new OpenLayers.Layer.ArcGISCache("baseMap", layerURL, {
        isBaseLayer: true,
        //numZoomLevels: 11,
        StartBounds: maxExtent,
        zoom: 5,
        layerInfo: layerInfo,
        projection: 'EPSG:' + layerInfo.spatialReference.wkt
    });
    console.log(layerInfo.spatialReference.wkt);
    //console.log(baseLayer.spatialReference);
    map = new OpenLayers.Map('map', {
        StartBounds: maxExtent,
        initialExtent: baseLayer.initialExtent,
        maxExtent: maxExtent,
        units: (baseLayer.units == "esriFeet") ? 'ft' : 'm',
        resolutions: baseLayer.resolutions,
        tileSize: baseLayer.tileSize,
        eventListeners: { //set event handler
            "moveend": moveEvent},
        //displayProjection: baseLayer.displayProjection
        projection: 'EPSG:' + layerInfo.spatialReference.wkt
    });


texts: Trying to plot coordinates through xy data, but ArcMap plots the points in the Atlantic ocean near 0 degrees by 0 degrees every time

texts: No matter what I do when trying to plot these XY data points it is always plotted in the Atlantic ocean, basically the south of Ghana every time. I've attempted exporting the .csv file, and .xlsx file, but no luck. The points are all in the same county at roughly 38, -95. I'm using the World Topographic base map along with a shape file of county boundaries in the US that align perfectly.


texts: Seeking GIS Solutions And ROI For Local Government?

texts:
In your experience, what GIS solutions provide the greatest ROI for local government's responsibilities?

Examples: Parcel management, managing CRM cases from the public, emergency response, etc.


texts: Private Surface Water Mains In Between Public Mains Analysis

texts:
I am working in a geometric network of a surface water (stormwater) system. We have both public and private assets (pipes/mains, manholes, inlets/CB's, etc.) in the network. I am doing a bit of QC on the data trying to identify any private assets in between public assets. I've worked through this a couple times not really accomplishing what I'm trying to do although my first method found 3 items that fall in this category. I am looking for alternatives to do this analysis and hopefully by consulting the community here I can figure out a different or new approach. The method that worked but only found a few issues (maybe our data is this good but I'd like to verify this with another method) was buffering all manholes by 10 feet and then intersecting private and public mains to find the locations that meet these criteria. I also ran this same analysis method using 30 ft and 50 ft finding no new incidences in the data.

I have searched google, geonet as well as gis.stackexchange for further help with no real ideas of how to improve this analysis. I have come up with a list of over a thousand points where public changes to private but I'd have to go through the entire list manually to see if they fit the criteria and that just doesn't seem like a wise use of time especially if I need to do this again in the future say on the sanitary system side or as changes are made to the system over time.


texts: Performing Tool Validation in Python Toolbox?

texts:
I'm trying to make sure there is no '-' in the output filename, and if there is, to set an error (circle with the red 'X') so that the tool can't be run, and have a tooltip that explains the issue.  I've tried a bunch of things, but from what I've read, it seems this should work:

if parameters[7].value:
            if parameters[7].value.value.count("-") > 0:
                parameters[7].setErrorMessage("Output cannot have a '-'." + \
                    "Please change the dash and rename the output after the tool has run")
            else:
                parameters[7].clearMessage()
        return

I'm not super clear on what is returned from the list of parameters.  The paramater type is "DEShapefile".  There is no error message being thrown when I use "--" as the filename.

Can someone see why it's not working for me?


texts: Conflate (merge) private shapefile data with OSM data

texts:
Background

Downloaded Alberta OSM data from Geofabrik and have it running on a private Linux server using PostgreSQL 9.1, PostGIS 2.0, Mapnik 2.1.0, osm2pgsql, Apache 2, mod_tile, renderd, and OpenLayers:

The data was imported using osm2pgsql as follows:

osm2pgsql -W -K -S /usr/local/share/osm2pgsql/default.style -d osm alberta.osm.bz2

Problem

The OSM data for Alberta is incomplete. I was given a set of shapefiles that improves upon the OSM data:

City.dbf, City.prj, City.sbn, City.sbx, City.shp, City.shp.xml, City.shx

Plus additional shapefiles for villages, urban areas, municipal district boundaries, and so forth. I have successfully imported the shapefiles into PostgreSQL using a pgAdmin plugin. The City.prj file describes its projection as follows:

GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]]

The shapefile data and OpenStreetMap (OSM) data each have their own database, but I believe I want to import the shapefile data into the OSM database. (This is a private server and a local copy of the OSM data; the shapefile data cannot legally be shared.)

Update #1

To be clear, the shapefiles do not contain roads: only municipality boundaries (that are definitely not part of the OSM data), cities (some of which are part of the OSM data), and city boundaries for larger cities (some of which might conflict with the OSM data).

Question

How do I merge the shapefile data with OSM data so that the new cities appear on the map?

Note: My main concern is resolving duplicate data (e.g., Edmonton is listed in both OSM and the procured shapefiles).

Related Links

http://wiki.openstreetmap.org/wiki/Shapefiles
http://wiki.openstreetmap.org/wiki/Import/Shapefile
http://wiki.openstreetmap.org/wiki/Convert_shp_to_osm_using_grass_and_gpsbabel
https://github.com/wavded/js-shapefile-to-geojson
http://redmine.yellowbkpk.com/projects/list_files/geo
http://www.gdal.org/ogr2ogr.html
http://wiki.openstreetmap.org/wiki/Ogr2osm
https://github.com/andrewguertin/ogr2osm
http://svn.openstreetmap.org/applications/utils/import/shp2osm/
http://crschmidt.net/blog/archives/354/polyshp2osm/
http://www.mail-archive.com/talk@openstreetmap.org/msg12487.html

Thank you!


texts: arcpy.Intersect_analysis() returns empty result, but Intersect tool in ArcMap is succesful

texts:
I've got two large feature classes that I'm trying to perform an intersect on. Using the following python code, I get an empty output feature class:

arcpy.Intersect_analysis(["Schools","Zips"], "Intersect_output","ALL")

Yet when I use the Intersect_analysis tool from the arcMap toolbox with these same two feature classes, it performs the intersect as expected. I've used the above code as part of a script that intersects a whole bunch of feature class pairs, and it does fine with most of them, but not these two feature classes in particular.

Can anyone think of a reason that these feature classes in particular might fail where others succeeded? I can't post the data, unfortunately, as I am not authorized to distribute it.

UPDATES:

Some more information:
*The failure now appears to be intermittent. After running code above exactly as shown in the arcmap python interpreter, it failed. I entered the exact same code (changing only the output feature class name, and this time it generated the correct intersect. I didn't touch the data, nor did I change the code, nor did I even close the python window between runs. it seems completely random. Has anyone ever encountered this behavior before? Might this be an arcpy bug?

*Both feature classes are using the exact same coordinate system.

*I have a separate script that calls the line of code above, but I'm testing it in the arcmap python interpreter, and I has failed there as well.


texts: Creating a multispectral image from scratch

texts: I want to make a multispectral image from cero to do some tests on it. Something really simple like 5 completely uniform bands with salt and pepper noise on them or a square of different values at the center. Clearly this would just be a stack of matrices, a multidimensional array, which is pretty straight forward to generate. I want to achieve this using python and gdal but gdal is being pretty hermetic, I don't get the hang of it at all. I ideally would want to create a geotiff file. Could anyone help me on this? some pointers or some gdal tutorial which is very gentle? Thank you all.


texts: Batch polygon clip doesn't work [ArcGIS]

texts:
I have an OSM layers of a region, but want to extract data only for city that I am interested in for further export in Illustrator:

AYCS, I've created a polygon which includes the territory that I am interested in. But when I try to create a batch clip, I input all layers in a window and ArcGIS gives me the message "Overlay not Polygon".

What solution could it be? Really want to work only with features of selected area, because ArcGIS loads all the other features from whole region.

UPD: That's how this batch clip looks like


texts: Is there a good library for WMS "Tile Servers"?

texts:
I am looking for WMS Tile Servers that I can consume with "OpenLayers".

One: Are there any good library / resources online that collate WMS tile servers or where is a good place to start looking? (Street maps & others etc)

Two: Are there any free tile servers available?


texts: Adding simple table of text (Excel) into layout and preserving formatting and clarity?

texts:
I am trying to put a simple table of text (excel) into my layout and preserve cell formatting including alignment, borders and text sharpness.  I've tried many approaches but nothing so far is satisfactory.

I am using ArcGIS 10.1 for Desktop.


texts: Automated correction of parcel geometry based on area measure?

texts: i have a set of parcels in DWG file and i have the areas that should be in these parcels in an excel file, so i have to modify each parcel to match the area specified in the excel file with a moving the minimum of points of each parcel (polygone), i'm searching any tool to automate this process in GIS (autocad , arcgis , arcgis engine , postgis ..)


texts:
I am working on a heat map that will show by region the number of issues found in which area. The code are not displaying any error on console, my functions which parse my data are ok (I have just debbuged all of them), but for some reason I don't know, the features are not displayed.

This is my code

var listaBairros = document.getElementById("lista-bairros").value;
var lotes = [];

lotes = lotes_parser();

var bairros = assignFeaturesToPolygon();

var styles = [
    new ol.style.Style({
        stroke: new ol.style.Stroke({
            color: 'blue',
            width: 3
        }),
        fill: new ol.style.Fill({
            color: 'rgba(0, 0, 255, 0.1)'
        })
    }),
    new ol.style.Style({
        image: new ol.style.Circle({
            radius: 5,
            fill: new ol.style.Fill({
                color: 'orange'
            })
        }),
        geometry: function (feature) {
            // return the coordinates of the first ring of the polygon
            var coordinates = feature.getGeometry().getCoordinates()[0];
            return new ol.geom.MultiPoint(coordinates);
        }
    })
];

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

var source = new ol.source.Vector({
    features: bairros
});

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

var map = new ol.Map({
    layers: [osm, layer],
    target: 'map',
    view: new ol.View({
        center: ol.proj.fromLonLat([-46.539028170768447, -23.952416758520052], 'EPSG:3857'),
        zoom: 11
    })
});

function polygonToArray(entrada) {
    var removePolygon = entrada.replace("POLYGON ((", "").replace("))", "");
    var points = removePolygon.split(", ");
    var finalPoints = [];
    points.forEach(function (point) {
        var temporaryPoint = point.split(" ");
        temporaryPoint[0] = Number(temporaryPoint[0]);
        temporaryPoint[1] = parseFloat(temporaryPoint[1]);
        finalPoints.push(temporaryPoint);
    })
    return finalPoints;
}

function assignFeaturesToPolygon() {
    var lot = [];
    lotes.forEach(function (value, index) {
        var new_feature = criarFeature(value);
        lot.push(new_feature);
    });
    return lot;
};

function  criarFeature (geometria) {
    var feature = new ol.Feature({ geometry: new ol.geom.Polygon([geometria]) });
    feature.getGeometry().transform('EPSG:4326', 'EPSG:3857');
    return feature;
}

function lotes_parser() {
    var bairros_separados = [];
    var bairros_parser = listaBairros.split(" _ ");
    bairros_parser.forEach(function (lote, index) {
        bairros_separados.push(polygonToArray(lote));
    });
    return bairros_separados;
}

What am I doing wrong?

Before you ask about my projection, that is everything ok with them.

This is the input field that my first line is reading (I used pastebin because is too many characters to write here)

https://pastebin.com/ra2iWdaH


texts:
I have three shapefiles (shp); A, B, and C, for example.

Each shp contains their data. I want to merge A and B into C. I did merge A into C, and it worked--let's call it D. But when I merged B and D (become E), data on B didn't replace data on D (unsubstituted).

What do I have to do?


texts: Is it possible to create floating/"flying" 3D objects in QGIS?

texts:
Is it possible with QGIS 2.14 and/or the Qgis2threejs-plugin to create 3d objects that don't "start" on the ground level (z=0), but are 3d objects "flying/floating" in the air?

For example wind turbines (their blades), bridges or certain towers or structures with overlapping elements.

If this is possible with Qgis, how can it be done?


texts:
Trying to understand how version releases work on GitHub. In particular the DotSpatial library https://github.com/DotSpatial/DotSpatial.

If I download the master branch and compile it, the results are the DLLs are all at version 2.0.0.0.

However if I just want the compiled DLLs that are ready for use I can click on the link https://ci.appveyor.com/api/projects/mogikanin/dotspatial/artifacts/Source/bin/Release.zip?branch=master
This contains the DLLs with version 2.0.0.402.

How do the versions work on GitHub?


texts: zoom set to 4 but on get zoom shows zero ,OpenLayers

texts:
var map1 = new OpenLayers.Map('firstmap');
var wms1;
var state1="<%=state%>";
var cql_filter1="st_nm='"+state1+ "'";
wms1 = new OpenLayers.Layer.WMS( "OpenLayers WMS",
       "http://localhost:8083/geoserver/wms", {layers: 'india:india',
       'CQL_Filter': cql_filter1,center: [22.2587, 71.1924],
        zoom: 4
        } );

map1.addLayer(wms1);
map1.zoomToMaxExtent();
console.log(map1.getZoom());

The output of map1.getZoom() is 0 and map is also not zoomed so how should I set zoom?


texts: Not sure why my scale area has turned grey and is not allowing me to zoom in or out?


texts:
In PostGIS 2.0 I try:

select ST_GeomFromGML(
  '<gml:Point srsName="urn:ogc:def:crs:EPSG::28992">
       <gml:pos>275466.0 565559.0 0.0</gml:pos>
   </gml:Point>');

This gives an error:

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

ERROR: invalid GML representation
SQL state: XX000
Context: SQL function "st_geomfromgml" statement 1

What is wrong with my point?


texts: Ogr2ogr is rounding numeric field from PostgreSQL to 6 decimals

texts:
I am using ogr2ogr to convert the data from PostgreSQL to GeoJSON. PostgreSQL contains field named shape.area, type numeric, which stores values up to 17 decimal places, like 0.15627328781960101.

The conversion to the GeoJSON is successful, but the decimal digits for shape.area field are rounded to 6 decimal places producing the value 0.156273.

The coordinate values are displayed with 15 decimal places and that is good.

I am using this command for the conversion:

ogr2ogr -f "GeoJSON" "D:\test.geojson" PG:"host=localhost user=postgres dbname=xx active_schema=xx  password=xx port=5432" "table"

I did try this conversion with GDAL 2.0.2 on CentOS and 1.11.3 on Win and the rounding is the same. The conversion to shapefile preserves the precision so this must have to do something specifically with the conversion to GeoJSON.

Am I missing some configuration switch? Is there a way to preserve the precision while exporting numeric values to GeoJSON format?


texts: Why does IFeatureClass.Search fails with -2147215922 (FDO_E_SE_INVALID_NUM_PARTS)?

texts:
in a class extension, I'm trying to do search features inside a point feature class(part of a geometric network on an sql server sde database) and I get a -2147215922 (FDO_E_SE_INVALID_NUM_PARTS) error.

The geometry I'm using to create my search filter is actually a Polygon created with the ISegmentCollection.SetCircle function. Basically, I'm trying to find some features inside a circle defined by a point and a radius

void CTools::FillArrayWithFeatures(IFeatureClassPtr ipFeatureClass,
              IPointPtr ipPoint, double fRadius, CAtlArray<IFeature*>* pArrFeature)
{
    //create the geometry
    IGeometryPtr ipGeometry;
    IPolygonPtr ipPolygon(CLSID_Polygon);
    ISegmentCollectionPtr ipSegmentCollection(ipPolygon);
    HRESULT hr = ipSegmentCollection->SetCircle(ipPoint, fRadius);
    ipGeometry = ipSegmentCollection;
    ISpatialReferencePtr ipSpatialRef;
    ipPoint->get_SpatialReference(&ipSpatialRef);
    ipGeometry->putref_SpatialReference(ipSpatialRef);

    //Create the spatial filter
    ISpatialFilterPtr ipSpatialFilter(CLSID_SpatialFilter);
    ipSpatialFilter->putref_Geometry(ipGeometry);
    BSTR strName;
    ipFeatureClass->get_ShapeFieldName(&strName);
    ipSpatialFilter->put_GeometryField(strName);
    ipSpatialFilter->put_SpatialRel(esriSpatialRelIntersects);

     // Execute the search
    IFeatureCursorPtr ipFeatureCursor;
    HRESULT hr = ipFeatureClass->Search(ipSpatialFilter, VARIANT_FALSE, &ipFeatureCursor);

    //Fill the array with the results
    if(ipFeatureCursor)
    {
        IFeature* pFoundFeature;
        hr = ipFeatureCursor->NextFeature(&pFoundFeature);
        while(pFoundFeature)
        {
            pArrFeature->Add(pFoundFeature);
            ipFeatureCursor->NextFeature(&pFoundFeature);
        }
    }
}

Am I doing something wrong? The same code works on a personnal geodatabase, but not on an sde sql server

Thanks


texts: How to draw features on OpenLayers geoserver raster source?

texts:
I can successfully draw a point on an OpenLayers map with source new ol.source.Vector({features: features});. However, I am trying to show these points on the raster source that I am pulling through geoserver using ol.source.TileWMS

Everytime I try to add a point on the raster source it does not appear. I can select the point by hovering over it on the raster source. However, it does not appear on the raster source. I need to show something to indicate the clicked point. What can I do ?

The following are the relevant parts of my code :

           draw = new ol.interaction.Draw({
                  source: vectorsource,
                  type: ("Point"),
                  wrapX: true
              });

           draw.on('drawend', function (event) {
            var style= new ol.style.Style({
                          fill: new ol.style.Fill({
                                color: 'rgba(255, 255, 255, 0.2)'
                              }),
                          stroke: new ol.style.Stroke({
                                color: '#ffcc33',
                                width: 2
                              }),
                          image: new ol.style.Circle({
                                radius: 7,
                                fill: new ol.style.Fill({
                                  color: '#ffcc33'
                                })
                              }),
                          text: new ol.style.Text({
                                text: ++counter+"",
                                fill: new ol.style.Fill({
                                  color: '#fff'
                                })
                              })
                                 });

var vectorsource = new ol.source.Vector({features: features});

var rasterSource = new ol.source.TileWMS({
              url: 'http://localhost:8080/geoserver/system/wms',
              params: {'LAYERS': 'system:ia-soil', 'TILED': true},
              serverType: 'geoserver',
              crossOrigin: 'anonymous'
                });

The images below shows the effect in greater detail:

Notice How point 2 appears on the map , but point 1 is overridden by the raster source. The point 1 is there you can select it, by hovering over it. I want to show point one location on the raster source. How can I do that?

The code that shows loading layers to the map.

map = new ol.Map({
        controls: ol.control.defaults().extend
          ([  new ol.control.LayerPopup(),new ol.control.FullScreen()]),
        interactions: ol.interaction.defaults().extend([
          new ol.interaction.DragRotateAndZoom()
        ]),
        layers: [raster,vectorsource,rastersource],// raster is just OSM layer.
        target: 'map',
        view: view
      });

This is how I am loading my map.


texts: How to control line width in legend in QGIS?

texts:
I have a map in which I am defining the sizes of various features in map units.  When I come to generate a pdf of the map, using the Print Composer, the map itself is fine, but when I include a legend then the symbol elements of the lines are represented at a totally inappropriate size (e.g. a 'Main road' symbol swamps the whole legend).

I think it is the use of map units that is the cause of the problem as switching the line width to absolute units (e.g. mm) on the original map produces an appropriately sized line symbol in the legend.  It seems like the legend that appears in the Print Composer is directly reflecting the way it appears in the Layer TOC - where, again, using map units produces oversized symbols relative to the map.  However this is only true for lines.  Point symbols rendered in map units appear exactly as I would expect in the Print Composer legend.

I would be very grateful if anyone can enlighten me as to whether this is a known problem, my misunderstanding of something, or there is a solution/setting/trick that I have missed.

I'm using QGIS 1.7.4 on Mac OSX 10.6.

Thanks for any help.

Phil


texts: Importing raster layer under Ubuntu 16.04?

texts:
I am using QGIS 2.18 with Ubuntu 16.04 and with Windows 7.

I can import a raster layer under Windows 7.
It does not work under Ubuntu 16.04.

I got a message telling me that the file cannot be recognized.
Why ? Shall I import a special library under Ubuntu?


texts: replace records in layer with records form another layer by id and column name QGIS

texts: I have a two polyline layers - one with 2800 rows, the other with 150 rows.   Each row has a  unique id column. I need to replace the records in one column of the the 2800-row layer with the records of the same column in the  150-row layer by id.  What would be the fastest method to do this?


texts: Creating an Action using the Value column in Value Relation widget

texts:
I want to create an Action that would use the attributes (such as report_id and county) of a site_investigation_report feature to generate a file path to a PDF report.

The county attribute uses the Value Relation widget.

Unfortunately, the Open type action, "[% "basepath %]/[% "county" %]/[% "shz_report_id" %]", evaluates to (using the Key column):

  "/Reports/ORA/SHZ_S_KI_12475"

instead of using the Values column

  "/Reports/Orange/SHZ_S_KI_12475".

So, how can I create an such an Action?


texts: Identify the latitudinal and longitudinal point directly beneath an object in space

texts:
Given the azimuth and altitude / latitude and longitude from the observer to the satellite, how can I accurately find the latitude and longitude of the point of intersection directly under the satellite?

This problem has been plaguing me for quite some time now, thank you in advance for your help!

Edit:
I also have the L1/2 p and L1/2 c/a values.


texts: How to install statist plugin in QGIS / OS X

texts:
I am stuck with installing a plugin for QGIS namend "Statist".

I get the following message when starting QGIS: "Missing matplotlib Python module".
So basically this question relates to all plugins for QGIS using the matplotlib python module.

I have tried to stick to every advice on stackexchange in order to install the whole QGIS package on OS X correctly. However, I don't have any clue why QGIS cannot find the properly installed matplotlib /Library/Python/2.7/site-packages directory.

I installed matplotlib with sudo pip install matplotlib.

I know that some questions already relate to his problem, but none helped me with solving this problem.

I use

OS X 10.7.4 QGIS 1.8.0-Lisboa
I've installed all necessary QGIS-related stuff from kyngchaos in the latest version.
QGIS had some problems with recognized rgdal first. I used the following
advice – then it worked properly.

This is what i get when i check for all versions installed (while "USER" is my userdirectory):

      import sys;sys.version;sys.path
      '2.7.1 (r271:86832, Aug  5 2011, 03:30:24) \n[GCC 4.2.1 (Based on Apple Inc. build 5658) (LLVM build 2335.15.00)]'
      ['/Users/USER/.qgis/python/plugins/GeoCoding', '/Users/USER/.qgis//python/plugins/sextante', '/Users/USER/.qgis/python/plugins/GeoCoding', '/Applications/QGIS.app/Contents/MacOS/../Resources/python', '/Users/USER/.qgis//python', '/Users/USER/.qgis//python/plugins', '/Applications/QGIS.app/Contents/MacOS/../Resources/python/plugins', '/Library/Frameworks/SQLite3.framework/Versions/3/Python/2.7', '/Library/Python/2.7/site-packages/numpy-override', '/Library/Frameworks/GDAL.framework/Versions/1.9/Python/2.7/site-packages', '/Library/Python/2.7/site-packages/Scrapy-0.14.4-py2.7.egg', '/Library/Python/2.7/site-packages/lxml-2.3.4-py2.7-macosx-10.7-intel.egg', '/Library/Python/2.7/site-packages/w3lib-1.2-py2.7.egg', '/Library/Python/2.7/site-packages/readline-6.2.4.1-py2.7-macosx-10.7-intel.egg', '/Library/Python/2.7/site-packages/nose-1.2.1-py2.7.egg', '/Library/Python/2.7/site-packages/six-1.2.0-py2.7.egg', '/Library/Python/2.7/site-packages/pytz-2012j-py2.7.egg', '/Library/Python/2.7/site-packages/tornado-2.4.1-py2.7.egg', '/Library/Python/2.7/site-packages/pyzmq-2.2.0.1-py2.7-macosx-10.7-intel.egg', '/Library/Python/2.7/site-packages/pika-0.9.9-py2.7.egg', '/Library/Python/2.7/site-packages/Jinja2-2.6-py2.7.egg', '/Library/Python/2.7/site-packages/patsy-0.1.0-py2.7.egg', '/System/Library/Frameworks/Python.framework/Versions/2.7/lib/python27.zip', '/System/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7', '/System/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/plat-darwin', '/System/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/plat-mac', '/System/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/plat-mac/lib-scriptpackages', '/System/Library/Frameworks/Python.framework/Versions/2.7/Extras/lib/python', '/System/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/lib-tk', '/System/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/lib-old', '/System/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/lib-dynload', '/System/Library/Frameworks/Python.framework/Versions/2.7/Extras/lib/python/PyObjC', '/Library/Python/2.7/site-packages', '/Library/Python/2.7/site-packages/PIL', '/Users/USER/.qgis/python/plugins/mmqgis/forms', '/Users/USER/.qgis/python/plugins/QuickMultiAttributeEdit/forms', '/Applications/QGIS.app/Contents/Resources/python/plugins/fTools/tools']


texts: Regarding Animation of Vector data in Google Earth Engine

texts:
I want to animate the vector file of Wind Direction in Google Earth Engine.
https://code.earthengine.google.com/9bd384d4041e4426e677281704dcd407


texts: Multiple CON Expression in ArcGIS Raster Calculator

texts:
How can I write the following conditional statement with CON?

If ((VALUE = 20 OR IF VALUE=24) then assign value 0, else 1)

Con(("Raster1", 0, 1, "VALUE =20")|("Raster1", 0, 1, "VALUE =24"),0,1)

I just can´t get this to work


texts: What is "Simplify tolerance" in QGIS's Simplify geometries tool?

texts:
I am trying to use QGIS 2.8 to simplify vector geometries.  But I do not understand the significance of the Simplify tolerance value in the Simplify geometries tool.

The documentation does not explain it.


texts: jquerymobile-TypeError: Object doesn't support this action on windows phone

texts:
Obtaining

  TypeError: Object doesn't support this action

in windows phone 8.1, IE v11. The application is plotting route and drawing circles in Android & iOS devices.

Is there any compatibility issue for ArcGIS for Javascript with windows phone? I am using dojo script to load the map..

I have a requirement for finding route between to points using RouteTask on my custom map layer. A TypeError is obtained on calling

var routeTask = new RouteTask(myRouteServerURL);

The same errors are obtained on initializing Circle and graphicsUtils.graphicsExtent(features) statements...
This is the require module used in my application

require([
    "esri/map",
    "esri/symbols/PictureMarkerSymbol",
    "esri/graphic",
    "esri/SpatialReference",
    "esri/geometry/Point",
    "esri/geometry/webMercatorUtils",
    "esri/InfoTemplate",
    "esri/lang",
    "esri/layers/ArcGISTiledMapServiceLayer",
    "esri/tasks/FeatureSet",
    "esri/graphicsUtils",
    "esri/tasks/RouteTask",
    "esri/tasks/RouteParameters",
    "esri/symbols/SimpleLineSymbol",
    "esri/Color",
    "esri/geometry/Circle",
    "esri/symbols/SimpleFillSymbol",
    "esri/layers/GraphicsLayer",
    "dojo/_base/connect",
    "dojo/on",
    "dojo/domReady!"
    ],
    function(Map, PictureMarkerSymbol, Graphic, SpatialReference, Point, webMercatorUtils,
      InfoTemplate, lang, ArcGISTiledMapServiceLayer, FeatureSet, graphicsUtils, RouteTask,
      RouteParameters, SimpleLineSymbol, Color, Circle, SimpleFillSymbol, GraphicsLayer, connect, on) {
            .......
            routeTask = new RouteTask(myRouteServerURL);
            .......
            var circle = new Circle({
              center: point,
              radius: 50
            });
            .......
});

I used the code the same code provided in Find a Route.


texts:
I am looking for a free, easy to use android app or mobile friendly website to plan a route and save it in GPX format.
So far I have seen https://mapstogpx.com/ which unfortunately does not work with google maps app url and https://www.gpxgenerator.com/ which works but in a rather convoluted way for my purposes, needing to copy the resulting gpx, paste it in a new file and save.


texts: Creating loop to reproject multiple rasters in ArcPy gives ERROR 999999?

texts:
I am trying to batch reproject multiple raster files using arcpy.

Here is the code I have so far:

# Import system modules
import arcpy
import os

# Set environment settings
arcpy.env.workspace = "H:\GIS Project\Bulk Density\In Rasters" #The file where your files are
arcpy.env.overwriteOutput = True
outWorkspace = "H:\GIS Project\Bulk Density\Reprojected BDW" #where you want the new rasters to go
arcpy.env.parallelProcessingFactor = "50"
rasters = arcpy.ListRasters()

#Begin loop

for inrl in rasters:

        # Determine if the input has a defined coordinate system, can't project it if it does not
    dsc = arcpy.Describe(inrl)

    if dsc.spatialReference.Name == "Unknown":
        print ('skipped this raster file due to undefined coordinate system: ' + inrl)
    else:
        # Determine the new output feature class path and name
        outrl = os.path.join(outWorkspace, inrl)

        # Set output coordinate system
        outCS = arcpy.SpatialReference('GDA 1994')

        # run project tool
        arcpy.AddMessage("Reprojecting" + inrl)
        arcpy.Project_management(inrl, outrl, outCS)

        # check messages
        print(arcpy.GetMessages())

print(arcpy.AddMessage("Reprojection complete"))

print "Script run ok"

However, when I attempt to debug this code I get an error in mixins.py with this part of the code highlighted:

self._arc_object.createFromFile(item)

and the following generic error:

  RuntimeError: ERROR 999999: Error executing function

Now I think it has something to do with the arcpy.ListRasters() function as if I run it with arcpy.ListFeatureClasses() the code works, but the created list is empty. But that is as far as I got with troubleshooting this code.

Here is the traceback message

Traceback (most recent call last):
  File "C:\Python27\ArcGISx6410.3\lib\bdb.py", line 400, in run
    exec cmd in globals, locals
  File "H:\GIS Project\Python Script\Reprojection.py", line 17, in <module>
    import arcpy
  File "C:\Program Files (x86)\ArcGIS\Desktop10.3\ArcPy\arcpy\arcobjects\mixins.py", line 949, in __init__
    self._arc_object.createFromFile(item)
RuntimeError: ERROR 999999: Error executing function.


texts: Extract by mask with multiple shapefiles and multiple rasters

texts:
I have two folders, one which is the shapefiles which define the polygons for my extract by mask (there are 33 of them) and one which has my rasters (10 of them).  For each of the 33 shapefiles I want to do an extract by mask with all 10 of the rasters.  So my output will be 330 new rasters.

I have tried two different ways to do this each with the same result.  The first is to use a zipped list:

arcpy.CheckOutExtension('Spatial')
arcpy.env.workspace=r'G:\raster_pathway'
rasters=arcpy.ListRasters('*.tif')
arcpy.env.workspace=r'G:\shapefile_pathway'
shapefiles=arcpy.ListFeatureClasses('*.shp')
out=r'G:\output_pathway'
for ras1, shape1 in zip(rasters, shapefiles):
    outpath=os.path.join(out,shape1[:-4] + '_' + ras2[:-4] + '.tif')
    arcpy.gp.ExtractByMask_sa(ras1,shape1, outpath)

the second way is to use a nested loop:

  arcpy.CheckOutExtension('Spatial')
  arcpy.env.workspace=r'G:\raster_pathway'
  rasters=arcpy.ListRasters('*.tif')
  arcpy.env.workspace=r'G:\shapefile_pathway'
  shapefiles=arcpy.ListFeatureClasses('*.shp')
  out=r'G:\output_pathway'
  for shapefile in shapefiles:
    for raster in rasters:
        outpath=os.path.join(out, shapefile[:-4] + '_' + raster[:-4] +'.tif')
        arcpy.gp.ExtractByMask_sa(raster,shapefile, outpath)

in both cases the following error is returned:

Traceback (most recent call last):

  File "<ipython-input-11-5034b4c5c647>", line 1, in <module>
    runfile('G:/Stefano/python codes/ecoregions_extract_mask.py', wdir='G:/Stefano/python codes')

  File "C:\Users\spotter\AppData\Local\Continuum\Anaconda\lib\site-packages\spyderlib\widgets\externalshell\sitecustomize.py", line 685, in runfile
    execfile(filename, namespace)

  File "C:\Users\spotter\AppData\Local\Continuum\Anaconda\lib\site-packages\spyderlib\widgets\externalshell\sitecustomize.py", line 71, in execfile
    exec(compile(scripttext, filename, 'exec'), glob, loc)

  File "G:/Stefano/python codes/ecoregions_extract_mask.py", line 22, in <module>
    arcpy.gp.ExtractByMask_sa(raster,shapefile, outpath)

  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\geoprocessing\_base.py", line 498, in <lambda>
    return lambda *args: val(*gp_fixargs(args, True))

ExecuteError: Failed to execute. Parameters are not valid.
ERROR 000865: Input raster: CDL_2006_Group_lookup_new_prj30m.tif does not exist.
Failed to execute (ExtractByMask).

CDL_2006_Group_lookup_new_prj30m.tif is the first raster within rasters

It definitely exists though.


texts:
I am quite new to GIS, and was wondering if you might be able to help me out?

I have a raster layer for gradient as well as a line created from GPS tracking points(using the point connector tool). Is there a way to sample the raster layer according to the distance covered along each raster value by the line?


texts: Dissolving "neighbour"/bordering polygon features in QGIS

texts:
I have a polygon layer with buildings for the area of the whole city.

My question is how to dissolve/aggregate/merge only those buildings, that are "neighbours", in a sense that they share at least one mutual wall or the distance between walls is miserable, e.g. roughly saying less than 0.5 m.
My desired result should look like

I have in mind several approaches which might be applied.

In layer create a column with a unique number and then apply

Vector > Geoprocessing Tools > Dissolve (by this unique number)

With much of manual work

Select Feature by Freehand > Merge Selected Features

Unfortunately, in both cases, it requires tremendous manual work, therefore I would love to optimize the process.
For large-scale maps, it is definitely an advantage to have each building separately, nevertheless, when working with bigger scales, the appearance of too many buildings on a small piece of a map, makes it overcrowded, causing the "ruffling effect" and cartographic generalization is a strong necessity.
The following examples do not really suit my request:

Dissolve polygons by specific field- but only neighbours polygons
Dissolve polygons - group to others with touching corners
Merging overlapping polygons with Dissolve?
Dissolve polygons in QGIS (GEOS) with tolerance


texts: Geoserver Cascading WMS, Wrong Projection

texts:
I configured a geoserver store with an external WMS. I used the following WMS url

http://www2.dmsolutions.ca/cgi-bin/mswms_gmap?VERSION=1.1.0&REQUEST=GetCapabilities

Once I saved the store, Geoserver listed the layers on the screen. I clicked the "Publish" link of the "DEMO" layer and I saved it with the default options.

When I display the layer preview, the layer has a wrong projection and is repeated many times.
Here you can see the image:
https://stackoverflow.com/questions/14336018/geoserver-cascading-wms-wrong-projection

The projection of the layer seems to be wrong. I tried to force to another projection or to keep the native one but results are not better.

Have you any ideas about the problem?


texts: Satellites needed to receive proper 1pps signal

texts:
So I have read that you need 4 satellites to be able to get a proper position from your GPS. But how many do you need if you are only interested in the 1pps signal from them? And what effects play a role for the precision of this signal?

(I am using a 1pps input from a GPS receiver for an atomic rubidium clock)


texts:
I have an Excel file with three columns, IncidentID, Latitude, and Longitude.  Each row in the file corresponds to an incident and the location of the incident in units of decimal degrees.

I also have a road shapefile of the area in which these incidents occurred.  The road shapefile has the WGS_1984_UTM_Zone_42N Projected Coordinate System and GCS_WGS_1984 Geographic Coordinate System and whose whose units (as far as I can tell) are in meters.

I can use 'Add Data' to add the excel file, and then use 'Display XY Data' to import the data, but the import interprets the Latitude and Longitude values as meter values.  How do I tell ArcGIS to interpret the Latitude and Longitude values as decimal degree values?


texts: Python: best practise to ceate a large shapefile

texts: What is the best approach to create a create a very large shapefile (close to the max of 2GB) using python? Best approach means: robust and fast. Above all fast.


texts: GeoPandas won't read PostGIS Polygon Table

texts:
I want to read my data on postgis using geopandas but it won't read. It works on Linestring data, but not on Polygon. When I open the postgis data using QGIS and it works fine, the polygon reads good. Here's my code:
import psycopg2
import geopandas as gpd

con = psycopg2.connect(user = "postgres",
                      password = "admin",
                      host = "127.0.0.1",
                      port = "5432",
                      database = "spatial")
cursor = con.cursor()
sql = "SELECT * FROM public.polygonData"
df = gpd.GeoDataFrame.from_postgis(sql, con)

if(con):
            cursor.close()
            con.close()
df

and this is the error. I wonder why it won't read the polygon as enclosed? How do I fix this?
IllegalArgumentException: Points of LinearRing do not form a closed linestring
---------------------------------------------------------------------------
WKBReadingError                           Traceback (most recent call last)
<ipython-input-10-b3bfe339a21d> in <module>
      6 cursor = con.cursor()
      7 sql = "SELECT * FROM public.polygonData"
----> 8 df = gpd.GeoDataFrame.from_postgis(sql, con)
      9
     10 if(con):

~\anaconda3\lib\site-packages\geopandas\geodataframe.py in from_postgis(cls, sql, con, geom_col, crs, index_col, coerce_float, parse_dates, params)
    362             coerce_float=coerce_float,
    363             parse_dates=parse_dates,
--> 364             params=params,
    365         )
    366

~\anaconda3\lib\site-packages\geopandas\io\sql.py in read_postgis(sql, con, geom_col, crs, index_col, coerce_float, parse_dates, params)
     89             load_geom = load_geom_text
     90
---> 91         df[geom_col] = geoms = geoms.apply(load_geom)
     92         if crs is None:
     93             srid = shapely.geos.lgeos.GEOSGetSRID(geoms.iat[0]._geom)

~\anaconda3\lib\site-packages\pandas\core\series.py in apply(self, func, convert_dtype, args, **kwds)
   3846             else:
   3847                 values = self.astype(object).values
-> 3848                 mapped = lib.map_infer(values, f, convert=convert_dtype)
   3849
   3850         if len(mapped) and isinstance(mapped[0], Series):

pandas\_libs\lib.pyx in pandas._libs.lib.map_infer()

~\anaconda3\lib\site-packages\geopandas\io\sql.py in load_geom_text(x)
     77         def load_geom_text(x):
     78             """Load from binary encoded as text."""
---> 79             return shapely.wkb.loads(str(x), hex=True)
     80
     81         if sys.version_info.major < 3:

~\anaconda3\lib\site-packages\shapely\wkb.py in loads(data, hex)
     13     reader = WKBReader(lgeos)
     14     if hex:
---> 15         return reader.read_hex(data)
     16     else:
     17         return reader.read(data)

~\anaconda3\lib\site-packages\shapely\geos.py in read_hex(self, data)
    439         if not geom:
    440             raise WKBReadingError(
--> 441                 "Could not create geometry because of errors "
    442                 "while reading input.")
    443         # avoid circular import dependency

WKBReadingError: Could not create geometry because of errors while reading input.

my geopandas version is 0.7.0.
my osgeo version is 3.0.4.
My Psycopg2 version is '2.8.6 (dt dec pq3 ext lo64)'.


texts:
I have an ESRI Addin in built in Visual Studio 2013 using .NET for ArcMap 10.3
I have a custom toolbar with a few buttons. One of the buttons is not ready to be shared with the rest of the group. Can I comment out the line in the Config.esriaddinx file? For example, I want to disable the "Edit" button on my toolbar so I can distribute it without that functionality, but not remove it entirely.

<ESRI.Configuration xmlns="http://schemas.esri.com/Desktop/AddIns" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance">
  <Name>Awesome Toolbar</Name>
  <AddInID></AddInID>
  <Description>My toolbar draws</Description>
  <Version>66</Version>
  <Author>Matthew<Author>
  <Date>2016-10-12</Date>
  <Targets>
    <Target name="Desktop" version="10.3" />
  </Targets>
  <AddIn language="CLR4.0" library="blabla.dll" namespace="TOOLBAR_version_6">
    <ArcMap>
      <Commands>
        <Button id="Draw" class="btn_Draw" onDemand="false" message="Draw" caption="Draw" tip="Draw" category="My toolbar" image="Images\btn_Draw.bmp" />
        <Button id="Edit" class="btn_Edit" onDemand="false" message="Edit" caption="Edit" tip="Edit" category="My toolbar" image="Images\btn_Edit.png" />
      </Commands>
      <Toolbars>
        <Toolbar id="My toolbar" caption="version 6" showInitially="false">
          <Items>
            <Button refID="esriArcMapUI.File_PageSetup" />
            <Item refID="btn_Draw" />
            <Item refID="btn_Edit" />
          </Items>
        </Toolbar>
      </Toolbars>
    </ArcMap>
  </AddIn>
</ESRI.Configuration>


texts: Intersection of buffer with polygon layer

texts:
I have a polygon layer representing land cover, and a layer of buffers around specific points.

I want to intersect the layers to later measure the proportion of each land-use in each buffer zone. However, I get an error "'NoneType' object has no attribute 'difference'"

So I assumed I had errors in my land-use layer. I used the 'check geometries' tool and created a new "test" layer covering all my buffers to see where the error could be. After a few corrections, I can now intersect my "test" layer and my "buffer" layer; my "land-use" layer and my "test" layer, but still not my "land-use" and "buffer" layer - the error is still the same.

I can't find what could go wrong.

You can find my layers there


texts:
To clarify, I'm not trying to define the coordinate system for a raster; I know how to do that already.  However, what I'm having trouble with is actually getting the raster image to be the right scale.  For example, I know the coordinates of the corners of the image in decimal degrees, I know the width and height of the image in degrees, but when I define the coordinate system, it's blown way out of proportion and is much larger than it's supposed to be.  Is there a way to define the the scale and individual coordinates of a raster?  Or is Georeferencing my best bet?
I've got a google map static map API image download that I'm trying to overlay on a basemap.  What's the best strategy for that?  ArcPy or ArcGis desktop is helpful.


texts: Problem - ST_AsGeoJSON return text with accolade


texts: How do I know topo layer contains DEM-altitude? ArcMap 10.3

texts: How do I know if a topographic layer contains DEM-altitude?  ArcMap 10.3


texts:
I tried to import CSV files to QGIS using the "Add delimited text layer" feature. For some reason, some of the numeric variables became QString.

How do I make them become integer so I can perform heatmap calculations?


texts: In ArcObjects, the Class 'BaseCommand' has 'OnCreate' method(to handle oncreate event). Similarly is there any equivalent method available for Menu ('BaseMenu' Class)?


texts: ST_ClusterDBSCAN in Bigquery, how to get a list of points for each cluster?

texts:
The algorithm DBSCAN should be able to group points together in clusters based on their proximity.
However, the example in GCB Bigquery has only the number of clusters.
WITH Geos as
  (SELECT 1 as row_id, st_geogfromtext('point empty') as geo UNION ALL
    SELECT 2, st_geogfromtext('multipoint(1 1, 2 2, 4 4, 5 2)') UNION ALL
    SELECT 3, st_geogfromtext('point(14 15)') UNION ALL
    SELECT 4, st_geogfromtext('linestring(40 1, 42 34, 44 39)') UNION ALL
    SELECT 5, st_geogfromtext('polygon((40 2, 40 1, 41 2, 40 2))'))
SELECT row_id, geo, ST_CLUSTERDBSCAN(geo, 1e5, 1) OVER () AS cluster_num FROM
Geos ORDER BY row_id
+--------+-----------------------------------+-------------+
| row_id |                geo                | cluster_num |
+--------+-----------------------------------+-------------+
|      1 |          GEOMETRYCOLLECTION EMPTY |        NULL |
|      2 |    MULTIPOINT(1 1, 2 2, 5 2, 4 4) |           0 |
|      3 |                      POINT(14 15) |           1 |
|      4 |    LINESTRING(40 1, 42 34, 44 39) |           2 |
|      5 | POLYGON((40 2, 40 1, 41 2, 40 2)) |           2 |
+--------+-----------------------------------+-------------+

Is there any way to get a list of points for each individual cluster?
I'm not that familiar with the OVER() notation.


texts: Contour lines/ basement contour extraction

texts: My purpose: getting out area between elevation contours (like b/w 600m and 700m). The contours could either be elevation ones or the basement contours. I got elevation contours using spatial analyst tool on DEM in ArcMap but they are in the form of polyline and highly dense. I  need it for entire Maharashtra, India. Also "seismotectonic atlas of India" don't have basement contours of Deccan region.


texts: Change color vector in Google Earth Engine

texts:
I want to change color vector in Google Earth Engine.

Fill color: nothing(No color)
out line color: red

This command did not answer:

Map.addLayer(table, {color: '9bbc57'}, 'colored');


texts: Get median distance from a bunch of points from a assigned location and draw buffer accordingly

texts:
I have 2 point layers, say layer A and Layer B

Layer A with a bunch of points,
Layer B with a single point.

I first need to measure all the distances between every single point in layer A and the point in layer B,
Then draw a buffer of median distance, 20 percentile distance, 30 percentile distance etc.

How do I achieve this in bulk in QGIS?

The output will be like:


texts: Select By Attribute and Export with Arcpy

texts:
I have 4 folders with shapefiles in them, and for each folder I am selecting by attribute and then exporting the new shapefile.  My code runs with no error but the output folders have no files in them and I am not sure if the problem is the way I am selecting by attribute or the way I am exporting.  The code I am using is this:

import arcpy, os
arcpy.CheckOutExtension('Spatial')
arcpy.env.overwriteOutput = True
arcpy.env.overwriteOutput = 1
arcpy.env.qualifiedFieldNames = False

folders=r'E:\Sheyenne\Landsat_8\spatial'
out=r'E:\Sheyenne\Landsat_8\just_field_cells'
basename=os.listdir(folders)
for i in basename:
    if not os.path.isdir(os.path.join(out,i)):
            os.mkdir(os.path.join(out,i))
for folder in os.listdir(folders):
    arcpy.env.workspace=os.path.join(folders,folder)
    shapefiles=arcpy.ListFeatureClasses('*.shp')
    for shape in shapefiles:
        name=shape[:-4]
        arcpy.MakeFeatureLayer_management(shape, "lyr")
        arcpy.SelectLayerByAttribute_management ("lyr", "NEW_Selection", "FID_ = 233" or "FID_ = 604" or
        "FID_ = 887" or "FID_ = 959" or "FID_ = 1686" or "FID_ = 1731" or "FID_ = 1753" or "FID_ = 1798"
        or "FID_ = 2058")
        arcpy.FeatureClassToFeatureClass_conversion("lyr", out, name)
print "Done"


texts: How to attach latitude and longitude values to pixels values in ArcGIS?

texts: I am trying to create space time cubes of LANDSAT 8 images in ArcGIS. The space time cubes can be created using only feature classes, so I need to convert from TIFF to pixel values with Dates and Latitude, Longitudes attached.The TIFF files obtained from USGS show only values and counts in attribute table. The properties of TIFF file show a spatial reference but I can't find a way to add it to attributes table.


texts: r.surf.area command of GRASS GIS is giving a very small area

texts:
I am a beginner to GIS systems and am trying to calculate the surface area of a dem using r.surf.area command of GRASS GIS. The surface area it is giving me is unacceptable as it is very small. It is around 0.04 kmsq. What I did was: I created a track and converted it to the .gpx file. I updated elevation values in .gpx file using an online tool and took the file to QGIS. In QGIS I converted it to shapefile. Then used SAGA-> Natural Neighbour to create a DEM. I took the DEM to GRASS GIS in a new location and selected to automatically read from georeferenced data and executed r.surf.area on it and not getting required area.
I even tried to execute r.surf.area on a DEM downloaded from the internet but still, I am getting a very small area.
–

What is wrong?


texts: how to convert X, Y coordinates pixel of the screenshot image.png from the map to really Long, lat

texts:
I need to convert the x,y coordinate of the PNG image from the map into Latitude and Longitude. I took an image from google-map that map image is cut from the part of the city here is a sample image which
I have drawn two points on it. I used "pyproj" to convert x/y coordinates to lat/long. But it doesn't give me correct values, for instance, my original long and lat is(59.944321, 30.179938) but when I convert X,Y coordinate to long and lat, I received something like this value (-13.487453788374928, 0.000983113712533505), how can I convert it to the real value. The function works correctly I checked it out with the original Long and lat. I want to know that do I missed something? do I need to get extra information from Google Maps when I take a screenshot from the interest region? I am completely new to this topic and I don't know what information I need to deal with this problem. .
Here I added my script code, I will be thankful if anybody can help me.
import cv2
from pyproj import Proj, transform
import geopy.distance
import pyproj

def lonlat_to_xy(lon, lat):
    proj_latlon = pyproj.Proj(proj='latlong', datum='WGS84')
    proj_xy = pyproj.Proj(proj="utm", zone=29, datum='WGS84')
    xy = pyproj.transform(proj_latlon, proj_xy, lon, lat)

    return xy[0], xy[1]

def xy_to_lonlat(x, y):
    proj_latlon = pyproj.Proj(proj='latlong', datum='WGS84')
    proj_xy = pyproj.Proj(proj="utm", zone=29, datum='WGS84')
    lonlat = pyproj.transform(proj_xy, proj_latlon, x, y)

    return lonlat[0], lonlat[1]

def find_distance(coordinates_from, coordinates_to):

    return geopy.distance.geodesic(coordinates_from, coordinates_to).km

def draw_points(img):

    height, width = img.shape[:2]
    x, y = width // 2, height // 2
    # cv2.circle(img, (x, y), radius=5, color=(0, 255, 0), thickness=3)  # center of image

    p1 = (int(x - 80), int(y + 75))
    p2 = (int(x - 108), int(y + 165))

    cv2.circle(img, (p1[0], p1[1]), radius=5, color=(25, 255, 255), thickness=3)
    cv2.circle(img, (p2[0], p2[1]), radius=5, color=(25, 0, 255), thickness=3)

    return img, p1, p2

path1 = r'./dist.png'
img = cv2.imread(path1)

# draw points on the image
img, p1, p2 = draw_points(img)

# extract longitude and longitude from (X,Y) coordinates
lon1, lat1 = xy_to_lonlat(p1[0], p1[1])  # x-axis is longitude and y-axis lat
lon2, lat2 = xy_to_lonlat(p2[0], p2[1])

# find Distance between two points
distance1 = find_distance((lon1, lat1), (lon2, lat2))
# real lat and lon from googl map the same area.
orig_lon1, orig_lat1 = 59.944321, 30.179938  # p1
orig_lon2,orig_lat2  = 59.95343, 30.18449   # p2
original_distance_on_map = find_distance((orig_lon1, orig_lat1), (orig_lon2, orig_lat2))

print('distance1', distance1)
print('Original_distance_on_map', original_distance_on_map)

cv2.imshow('img', img)
cv2.waitKey(0)
cv2.destroyAllWindows()


texts: Getting Unicode/Codepage problem - Convert UTF-8 to ANSI with File Geodatabase

texts:
I have a Codepage issue, and for now I am a bit stumpped

My original source data contains multiple country data delivered in SHP (UTF-8). The data is then merged into a File GeoDB as the data is too large for the 2gb SHP/DBF limit. Previously I was able to use the *.CPG method with a smaller SHP dataset to correctly display the data that contains caracters that matches ANSI (ê, ë, ô, etc.). Now that the data is in a File GeoDB I am not able to use this method.

Strangely enough, I changed the codepage of a smaller SHP-UTF8 data set to SHP-ANSI using a FME process, then imported the data into a File GeoDB and the data then displayed correctly in ArcMap.

Now transforming all the source SHP-UTF8 data to SHP-ANSI is quite a timely process (many files that makes up a big dataset), so my question is if it is possible to set or convert the CodePage of the data in a File GeoDB so that the merged dataset would display correct?

The UTF8 based data looks like this : KoÃ¶perasie Street

Where it should look like (ANSI): Koöperasie Street


texts:
I am using Leaflet with Custom GeoJSON Layers to display my map. I noticed that the map off-screen is only loaded after releasing the mouse button during a drag.

Hopefully, this demo makes it more clear what I mean:

Is there a way to make the map draw before or during a pan movement?

I found the Leaflet.EdgeBuffer plugin which would do what I want but it appears to only work for tile layers.

As you can probably tell, the map I am using is fairly simple so memory or performance shouldn't be impacted much by even buffering the whole map.


texts: I am a newbie in R. I want to classify images using random forest in R. I am trying to learn the codes. But in the mean time, is there a push button (like a gui) in R that can be used for classification. Looks model map package can do it, but i cannot figure out how?


texts: Hillshade result after mosaicing multiple rasters appears irregular

texts:
I am trying to mosaic 19 rasters which i already reprojected in WGS 1984 World Mercator and then hillshading them. These are the steps i follow please correct me if i am doing something wrong.

I am am using the mosaic to new raster tool , i have been referring to this question  Mosaic to New Raster tool fail so i am not creating the raster dataset in the same location where the rasters are stored (used Desktop in this case). I made sure that the pixel type and the number of bands are the same.

The result is acceptable.

Here comes the problem. Now i hillshade this raster. At first look appears to be ok but when i zoom in as you can see in the image below the hillshade comes with some grid like lines. I don´t understand what am i doing wrong.


texts: Problems with adding CSV layers ETRS89 (no latitude/longitude but X and Y coordinates)

texts:
I am dealing with some data I need to get on QGIS.  It's CRS should be ETRS89, but instead of
longitude/latitude I only have X and Y coordinates, and as soon as I add the data as CSV layer it ends up at a place where it should not belong.

I use this CSV here: https://data.gov.uk/dataset/b5d2f859-9a06-47bc-8f09-89171cfe97ba/shropshire-brownfield-site-register
It includes the X and Y coordinates and states it is ETRS89.

Within the CSV I have the coordinates.
CordX: 3853505 and CordY -205876

and it should be ETRS89
but if I look at the area where the points should end up, as ETRS89 the coordinates should be more like.

X: -3.059671 and Y: 52.862264

Is there a way to import the CSV with these coordinates instead? (in QGIS)
Or do I have to convert them prior manually? How do I do that in the best way? (it's a file with a total of ~1500 records).


texts: Total distance to drive all roads in network (with some driven multiple times)

texts:
We need to drive all the main roads in Iceland. We are wondering how many KM we need to drive. We see that some roads need to be driven twice to be able to reach all roads.
We have ArcMap basic - Spatial analyst - ArcScan - Network analyst
I've tried OD and Route - it fails and I don't know what I'm doing. The internet information is just costing days of trying no luck.
Does anybody know this?


texts: Getting name of subtype field in feature class in ArcPy?

texts:
I have a geodatabase with many feature classes. Some of the feature classes have subtype fields. Their names vary feature class by feature class. I now want to find the names of those subtype fields.
What I have so far is:

for fc in arcpy.ListFeatureClasses("*","",dataset):
    subtypefields = arcpy.da.ListSubtypes(fc)
    for subtype in subtypefields:

Now, I don't know how to continue from here.


texts: Error when trying to edit the header of an envi file?

texts:
I have a file with the following Metadata in Envi 5.3

when I click on Edit Metadata or when I try to access the option through Raster Management > Edit ENVI Header, I get the error:

Does anyone know what's the cause of the error?
Is it a bug in 5.3 version?
Or normally I can't access header of files that are imported from PolSARpro?
Do you know any other way to edit the header (within ENVI or outside)? especially map_info part?


texts: Getting UK local authority boundaries in OpenStreetMap through R

texts:
I'm trying to plot some infrastructure data that comes from the OpenStreetMap using R. I'm using the osmdata package with ggplot2. I'm trying to pull the local authority boundaries of a few local authorities in Scotland, mainly islands in the west and north. I can't locate the local authority boundaries at an administrative level ("admin_level"). According to the OSM wiki, there appears to be a bespoke tag for local authorities. My question is: how can I use R to identify local authority boundaries for the UK?
My code, which returns zero features is:
### cycling through the admin levels 1:n doesn't provide the correct boundary

outer_heb_1 <-
      getbb("Outer Hebrides") %>% ## bounding box for Outer Hebrides
       opq()   %>% # build the Overpass query
    add_osm_feature(key = "admin_level", # Admin, level 4
                    value = "4") %>%
  osmdata_sf() # return a simple features

### trying boundary/local authority in add_osm_feature
outer_heb_1 <-
      getbb("Outer Hebrides") %>% ## bounding box for Outer Hebrides
       opq()   %>% # build the Overpass query
    add_osm_feature(key = 'boundary', # Admin, level 4
                    value = "local_authority") %>%
  osmdata_sf() # return a simple features


texts: Extracting GeoTIFF data with coordinates using R?

texts:
I have downloaded the following data as a GeoTIFF and I'm interested in determining the value from this raster at specific locations.

raster = readGDAL("filename.TIFF")
summary(raster)
Object of class SpatialGridDataFrame
Coordinates:
   min max
x -180 180
y  -90  90
Is projected: FALSE
proj4string :
[+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0]
Grid attributes:
  cellcentre.offset cellsize cells.dim
x            -179.5        1       360
y             -89.5        1       180
Data attributes:
     band1
 Min.   :  1.0
 1st Qu.:255.0
 Median :255.0
 Mean   :226.9
 3rd Qu.:255.0
 Max.   :255.0
 NA's   :245

head(data[7:8]) #Xp = longitude, Yp = latitude
              Xp        Yp
    7  -50.41777   1.33138
    8  -59.25075 -19.16780
    9  -67.25043 -27.83412
    10 -54.91759 -15.00130
    12 -67.58375 -28.00078
    13 -65.08385 -34.33386

However, when I run the following, I get an error message:

sp <- SpatialPoints(data[7:8])
summary(sp)
Object of class SpatialPoints
Coordinates:
         min       max
Xp -76.08341 -35.25171
Yp -38.66702   1.33138
Is projected: NA
proj4string : [NA]
Number of points: 184

data$ndvi<-extract(raster, sp)
    Error in .subset(x, j) : invalid subscript type 'S4'

If I use simulated data for the raster, I don't seem to have the same issue:

simr <- raster(ncol=36, nrow=18)
simr[] <- 1:ncell(simr)
simndvi<-extract(raster, sp)
head(simndvi)
[1] 301 373 408 373 408 444

How can I figure out why this isn't working for the GeoTIFF?


texts: I want to assign bus routes with bus stations where this routes passes. How could I do this knowing that in one bus stations passes several bus routes? I would like to handle in QGIS Environment, could I?


texts:
I have a question regarding converting a vector forest cover layer to a raster layer.

The basis of what I am trying to do is determine what percent of forest cover is within a set of 1x1 km grid cells in arcmap. I originally attempted to do this using fishnet, tabulate area tool, and a vector layer for forest cover, but I got the error message "maximum record length has been exceeded" (if anyone knows a way to get around this without converting the vector layer to a raster, please feel free to also let me know). So, I decided to convert the vector forest cover layer to a raster and try to use the zonal statistics tool. Here, I run into problems where it seems that the numbers do not add up correctly.

First off, I originally merged two vector files for secondary and primary forest into one vector file I called "forest"; using the merge tool; is it possible I executed this wrong and lost important information in this process?

Next, should I be using the feature to raster tool or the polygon to raster tool? I am wondering when I used the feature to raster tool if I am not inputting the correct information to provide the proper data related to the forest cover pixels. For example, what should be in the VALUE field, can I just put FID? Should it be something for specific? Do I need to be specifying a cell size, and should it be 30? (Since the map is resolution is 30x30 m pixels).