texts: automation help within an attribute table ArcMap 10.3

texts:
I have 100+ catchments within a watershed.
There are 6 different types of land use classifications within in the watershed.
There are 5 different types of soil types within the watershed.

My example shows Catchment ID 1 is 216 acres total with 4 land use classifications and 4 types of soils. In my example I would start editing, merge the 5 highlighted, then do the same thing for the rest in catchment 1 and so on for the remaining catchments.

Question: Without individually merging same land use types with same soil types in each catchment is there a way I can automate this process.


texts: Sum vector data correlating with a second vector layer

texts:
There are a lot of questions very similar to mine, but I couldn't find one that accurately described what I want to do.  I am working in QGIS.

I have one vector layer which, among other attributes, has a production volume for oil wells.   I have a second layer which represents townships.  Some townships have five wells, some have none etc., all with different levels of production.

I want to create a new layer which displays all the townships with a summation of the total production volume within the township.  Can anyone think of an efficient way to do this? EDIT: I would also like to be able to the sum of the production volumes to the township layer.

Thanks!


texts: Implementing "conformal mapping" in OpenLayers?

texts:
I want to preserve shape of objects, so how do I implement "conformal mapping" in openlayers 4.6.5?

My current objects:
I have several geoJSON objects (WGS84) around the globe in long/lat with a given shape. These objects are converted to 3857 this way:

ol.proj.transform([long, lat], 'EPSG:4326', 'EPSG:3857').

My current map: Using this: map.getView().getProjection();
I get the following output:

xb {wb: "EPSG:3857", a: "m", i: Array(4), oe: Array(4), b: "enu", …}

My previous post : Why do I get different resolution between horizontal/vertical directions, from lon/lat to openlayers map?

also explains what i want to achieve, however how do i do it in practice?

Post and pages already read:

geoJSON: https://tools.ietf.org/html/rfc7946
8 years old: In what coordinate systems do true circles appear round?
Proj4js http://proj4js.org/


texts:
When I capture the screen and save the file it also generates a .pwg file, I do not need it, only the .png

self.iface.mapCanvas().saveAsImage("doc.png")


texts: How to make length difference between two network analysis?

texts:
I have one quite challenging problem and can't find solution for it.

I am doing closest facility analysis via network analysis. I have bus stations as facilities and houses as incidents. I have run analysis for this scenario and saved it. Then I did the same again but I included also flood as polygon barrier into analysis. Now I want to find out which incidents (houses) have longer distance when there is flooding compare to scenario when there is no flooding. Note, that during this analysis there were some houses which were excluded as errors as they are sooo flooded that they can't access road.

I was thinking to export both tables in Excel and simply minus them from each other and then import them back to include difference. Problem is that numbers of "Shape_Length" are not in same row so in the end there are less rows for scenario with flooding. Nevertheless ID of incidents is matching for both scenarios so I can find difference for each house manually. Problem is that I have approx. 9 000 of these houses... so I need to find for some smarter solution.

Do you know please about some good solution how to make difference between these two scenarios so I could see what difference of length is there?


texts:
I have a CSV data sheet that I would like to import in QGIS as a set of polygons. The CSV file contains a column of geometries with long lines of data (see example below). Can somebody tell me how I can import this into QGIS in order to see the polygons?

  0106000020E61000000100000001030000000100000007000000465C28D96628E43F86DF2C8898F8494006B64134B928E43F5A21943A98F849401B5F5A8CF028E43F7E4EAEA299F849407AE3D3EE2D29E43F261E4D329BF849405F597722D828E43F41BF64829BF849408283070E9B28E43F2CF0BEF299F84940465C28D96628E43F86DF2C8898F84940


texts: Labeling polygons along their boundaries using QGIS?

texts:
When labeling polygons (not polylines!) in QGIS, is it possible to label them along their borders rather than their area?

I know I could convert polygons to polylines and then label them, but I am looking for a solution, where I would only label a specific selection of polygons from a polygon shapefile, therefore I would like to keep them in the same shapefile (i.e. the same shape).

Below is the example of a desired result:


texts:
I am drawing a bearing line using a line with max vertices of 2, I have the calculations for figuring out the actual bearing from the lat lon of the original point and the mouse latlon

I want to place the label on the actual mouse before the second point is placed and the line feature is drawn.

Is this actually possible?

getMouseBearing : function(firstLatLon)
{
    this.map.addControl(new OpenLayers.Control.MousePosition());

    this.map.events.register("mousemove", map, function(e) {
        var mousePosition = this.getLonLatFromPixel(e.xy);
        var bearing = getLineData(false, mousePosition, firstLatLon);
        //Need something here to add label to mouse cursor
    });
},


texts: Reload the dimension values of a layer with OpenLayers 4

texts:
I'm looking to reload the dimension values of a layer with OpenLayers 4.

My problem is when my App calls the GetCapabilities of the layer, the dimensions received are for a certain datetime, for example:

<Dimension name="time" units="ISO8601" default="2018-07-27T14:20:00Z" nearestValue="0">2018-07-27T11:20:00Z/2018-07-27T14:20:00Z/PT10M</Dimension>

but after a couple of hours, those datetime values are expired, so when I call the GetMap with these date parameters, the WMS tell me:

Invalid model run or date or timestamp value out of time domain : 2018-07-27T11:30:00Z/2018-07-27T14:30:00Z/PT10M

I try:

layer.getSources().refresh()

but it doesn't work. Is it something possible to do without delete/add the layer?


texts: What am I not getting about JPEG 2000 coordinate systems?

texts:
I need some help with understanding JPEG 2000 imagery and coordinate systems.

I've downloaded NAIP files from USGS in jp2 format and am attempting to place them over a vrt height map I've created using QGIS 2.6.  The .vrt shows the proper lat/lon decimal degrees such as -77.15, 41.91but when I add the .jp2 raster, it is completely off the map with coordinates like -8,575,110, 5,164,330.  The jpeg 2000s are delivered in WGS 84 projection.

The jpeg 2000 Projection was "WGS_1984_Web_Mercator_Auxliary_Sphere".  I used gdal warp to change it to a WGS 84 projection to the jpeg 2000 which coverted it to a GEOTIFF.

Comparing the information on the vrt and new GTIFF, the two files are nearly identical, the GEOGCS, DATUM, SPHEROID, AUTHORITIES, PRIMEM are identical.  However, the Origins are different with vrt containing the proper decimal degree origin and the jpeg 2000 having the wacky -8.5 million, 5.1 million as the origin.

How do I align these files?


texts:
I have made some web search but coultn't find any answer regarding how to perform viewshed analysis in R without the use of external tools (e.g., GRASS). I do know that a number of interesting packages are available in R to perform viartually all types of GIS-related operations, from visualization to spatial analysis with all which lies in-between, but I seem to not have found anything regarding viewshed.

I need to have a package which does not require the use of external tools, such as GRASS or ArcGIS for instance.


texts: Convert Oracle schema to GeoDatabase using Direct Connect only

texts:
Either this is a gap in the documentation or I'm really missing something!

I know that I can direct connect to an SDE GeoDatabase without needing ArcSDE to be installed on the database server, but I can't figure out how I'm supposed to create a GeoDatabase or convert an existing schema to a GeoDatabase without going through the ArcSDE installation process.

I'm using ArcGIS Server 10 which is installed on a web server. I don't have sufficient licensing to install ArcSDE on a separate server but the license does allow me to direct connect to a GeoDatabase on a separate machine.

Is the conversion accomplished through the ArcGIS Server post-install? I started testing but I don't currently have the existing SOM and SOC service passwords so I don't want to change them in the post-install and break something. I can't get to where I think the database configuration stuff is withouth setting these accounts.

Is this something I would do through ArcCatalog? I can't see the relevant option.


texts: Having download/latency issues with ESRI JS-API recently?

texts: Is anyone using the ESRI JS-API having issues over the last week or two with either latency or performance in the client download and operation? Since about 7/12/2011 we have seen across numerous ISPs and Browser/OS platforms in WA state issues with users getting the application to load and operate normally.  The only workaround I have had is to switch my app to using the compact api, to reduce the download footprint.


texts: What is most suitable DEM to delineate watershed (from spot height, ASTER and SRTM)?

texts:
I have three DEMs of our country which are spot height DEM, SRTM DEM and ASTER DEM.

SRTM has 90m resolution and ASTER has 30m resolution.

My supervisor told me to delineate watershed using these three DEMs separately, and also using merged DEMs.

What is the main differences between the three types of DEM that I mentioned above?
How do I merge these three DEMs?
What is the most suitable DEM to delineate watershed?


texts:
I have a bit of SQL code that look similar to this:

SELECT
sum(case when latitude = '0' then 1 else 0 end) as count_zero,
sum(case when latitude is NULL then 1 else 0 end) as count_null,
sum((case when latitude = '0' then 1 else 0 end) + (case when latitude is NULL then 1 else 0 end)) as total_zero,
count(latitude) as count_not_nulls,
count(*) as total
FROM sites_database;

Is there a "cleaner" way to write this same query. I have tried using the "sum" expression using the column alias, something like:

Sum(count_zero + count_null) as total_null

But this doesn't seem to work for some reason


texts:
I've been working with the append tool inside arcpy and discovered that, even though I have domains constraining the values of some of my columns, I never receive an error when the domains are violated.

I've read through the documentation for the append tool, and it doesn't discuss this in any way. Other google searches have failed me in finding an article about this issue.

Have others had this problem?

Does anyone have references to articles that have discussed this problem?

Are there any strategies to ensure that domains are enforced properly?

UPDATE

Thanks for the great answers. It unfortunately seems that there isn't an easy way to validate features using arcpy. I can't make the switch to 10.2 because this script has to run on a variety of versions that aren't 10.2 (with no plans to upgrade soon).

Thanks to everyone's help, I was able to find this article about the sdelayer -o feature_info command. It seems like this command will at least give me a list of invalid features, but I haven't had the chance to try it yet. I'll update when I've played around with it.


texts: Snapping doesn't work when many layers loaded!

texts: I am using Arcmap 10, when I have only limited number of layers loaded, the snapping experience is perfect but when I have many layers loaded specially polygons, snapping process is painful. My cursor turns to "Busy" shape (the sand clock) and can't snap to the feature. I have no idea why it happens. Any solution?


texts: QGIS open tif file larger than 10 GB, '*.tif is not a supported raster data source'

texts: When trying to use QGIS to open a tif file larger than 10GB, I encountered an error says "*.tif is not a supported raster data source." Is there any way or other software to open large tif files?


texts: Cannot use GDAL algorithms in processing modeler (QGIS 3.0)

texts:
I am trying to run a fairly complex model in QGIS 3.0, but when I try to use a GDAL algorithm which uses the output of another algorithm as input, I get error messages about an unexpected '(' character in the ogr2ogr command. The command it tries to run is:

ogr2ogr "/tmp/processing_8e14b7ffb2794d6985ee9b732c45d59b/a1896a49f4a648ef80fef64fd3679b70/OUTPUT.shp" memory?geometry=LineString&crs=EPSG:4647&field=Klasse:string(1,0)&field=Breite:double(12,11)&field=Bemerkung:string(254,0)&field=offset:double(10,3)&field=layer:string(100,0)&field=path:string(200,0) -dialect sqlite -sql "SELECT ST_Line_Interpolate_Point(geometry, 0.5) AS geometry, Klasse,Breite,Bemerkung,offset,layer,path FROM 'None'" -f "ESRI Shapefile"

so it seems to be the part where memory comes into play in the command. I don't know if this is a bug in the processing modeler or if it is something that I am doing wrong. I have also tried with a simple input to GDAL "Points along lines" and then output, and this throws errors as well.

My question is - am I doing something wrong with trying to use GDAL algorithms in the Modeler, or is there some way of getting it to accept the input from one of the other algorithms that are being run?


texts: Generate a GPX file that matches a road network from an imported drawing

texts:
I'm looking for a tool that'll help create a GPX file which will help me run a track and create "GPS Art" or "GPS Drawing".

I would like to draw simple shapes, like animals or texts by following a GPX path on my runs. I need a tool that'll take a lines, and given a specific area, will match these lines to a real world road, sidewalk or path network. The tool needs to best match the drawing to the road network as possible to maintain the idea of the shape.

Does such a tool already exists? If not, is there an existing API or algorithm that would help me create the tool?

I'd prefer to write this in Python. I'm not asking for code, but help finding the right existing tool, or algorithm.


texts: Upright/Horizontal labels when labeling circular polygons at perimeter?

texts:
I have several polygons (circles) I want to label with their ID number. The label should be outside the circles so I use the positioning options "Use perimeter" and the tickbox "Next to the line (or similar, my qgis talk spanish...).

QGIS now automatically aligns the labels with the curvature of the line. Is there an option or a way to make them oriented horizontally?


texts:
The below code is working perfectly on click event to retrieve Lat/Long and Current Zoom level. How do I get the same when I hover mouse/drag on map Window (not on submit button)?

    <!DOCTYPE html>

<input type= "text" ID=latlng>
<input type="text" ID=Zoom>

<button onclick="copyText()">Submit</button>


texts: Create point in intersections of two layers

texts:
I have two layers. One of them is 'Road_line' and another one is 'forest_shape'.
I want to create a point layer that has points in intersection of that layers. I couldn't find any tools for it.


texts:
I used the following link to create labels with leaders, but since I have a polyline layer I used a different expression for the Geometry Generator:

make_line(
  centroid($geometry),
  make_point("auxiliary_storage_labeling_positionx", "auxiliary_storage_labeling_positiony")
)

It appears to work, but looks quite ugly, where some of the leaders go past the geometry centroid as follows:

And others fail to reach it:

Why does this occur and how can I fix it?


texts: "Attribute only valid on v:image" error on popup

texts:
I'm implementing a script based on the popup on vector layer example into a cognos report. It all works fine in firefox and i dont get no errors. But if i open my script with internet explorer, as soon as I click on a feature, I can't close the popup anymore and internet explorer crashes with the failure notice: "Attribute only valid on v:image". This happens in both IE 8 & 7.

When i open the regular example i can open/close the popups without an error.

If not for my work i wouldnt care at all that the script can't be used by IE users :)

I hope someone can help.

EDIT: Yes, i've googled before i posted here.


texts: Query of line direction histogram plugin

texts:
I was wondering if anyone could explain how this plugin works I have had a play around with it and I do like how it visualizes how vectors are orientated. However I ran one and opened the CSV file showing the total vector length and the number of vectors within a segment (in total 16 segments/bins).

My question is how does this plugin identify which vectors go into each segment and is there a way of identifying which vectors(in particular the vector id) fall into particular histogram segments as I would like to then use this information for further analysis these vectors and accompany the results with the line direction histogram graph


texts: How to find nearby objects using Easting and Northing range

texts:
I have a bunch of objects which have associated Eastings and Northings in format OSGB36

for example:

265397  869562
547297  179002

If I have a starting Easting and Northing, and wanted to find related objects within X km is it possible to compute a Start/End Northing and Start/End Easting range that I can search for programatically?

I am using C#.

I.e. searching for objects using pseudo-code :

Easting >= 265390 && Easting <= 265400 &&
Northing >= 869560 && Northing <= 869564

would this find nearby objects? What resolution do these coordinates operate to?
I can operate in miles or km, just need to know the scale - if such a thing is possible.

Sorry this is new to me.


texts: Clipping a raster with shapefile and export the raster with shapefile extent

texts:
I am clipping a raster file with a shapefile grid and I want to save each raster subset individually with the polygon spatial extent.

To select the extent of each polygon, first I select the polygon with arcpy.SelectLayerByAttribute_management and then I get the polygon extent and I use this extent in arcpy.Clip_management function.

But when I run the code each raster subset don't have the extent that I want.

The code is:

def execute(self, parameters, messages):
        """The source code of the tool."""

a= parameters[0].valueAsText
b= parameters[1].valueAsText`
#-----------------------------

with arcpy.da.SearchCursor(a, ["FID_malla"]) as cursor:

    for row in cursor:

        #select a polygon
        arcpy.SelectLayerByAttribute_management(a, "NEW_SELECTION", '"FID_malla" ='+str(row[0]))

        #extent
        shp_ex = arcpy.Describe(a)
        exshp = shp_ex.extent

        xmin=exshp.XMin
        ymin=exshp.YMin
        xmax=exshp.XMax
        ymax=exshp.YMax

        ExtStr = "{} {} {} {}".format(xmin, ymin, xmax, ymax)

        #Execute clip
        arcpy.Clip_management(b, ExtStr, "C:\Users\mvinuales\Desktop\Image_Chips\imagen\\" + str(row[0]) +".tif", a, "256", "ClippingGeometry","NO_MAINTAIN_EXTENT")

    arcpy.SelectLayerByAttribute_management(a,"CLEAR_SELECTION")


texts: PostGIS intersected polygon fails on contains

texts:
I have two functions, first I check if a polygon is within a defined polygon. If it's not, but on the border, I cut it with ST_Intersection, and get a new polygon. If I check the result, I get that it is still intersecting, so it gets cut again and again.

BEGIN
   RETURN QUERY

SELECT
CASE
WHEN (SELECT ST_Contains(a, b)
    FROM (SELECT (SELECT geom from big_poly) As a, ST_SetSrid($1,4326)  As b)
    As tbl) = false
THEN (select * from geometry_check($1))
WHEN (SELECT ST_Contains(a, b)
    FROM (SELECT (SELECT geom from big_poly) As a, ST_SetSrid(ST_GeomFromText($1),4326)  As b)
    As tbl) = true
THEN null

END;
END

BEGIN
   RETURN QUERY
SELECT ST_Intersection(a, b)
    FROM (SELECT (SELECT big_poly.geom from big_poly) As a, ST_SetSrid($1,4326)  As b) as foo;

END


texts: On python, how to generate an array of coordinates between two points

texts:
Perhaps I'm looking to a solved question, but while searching I didn't find any solution. I need to build a matrix of lat, lon values, based on two points, and a sampling rate (distance between points).

A possibility is implement a solution using numpy, with all the rules. However some exception treatment is needed, specially if my region is on an pole or mid-pacific zones, where a WGS-84 coordinate based matrix would change from 179.9 to -179.9, for example.

I would like to know if already exists a library that is able to generate interactively an array within those two points, and able to control those issues.


texts: CalculateField_management - Updating a field with value from another field?

texts:
I have automatically generated X/Y coordinates using arcpy.AddGeometryAttributes_management but want to rename these fields. I've planned to do this by adding a new field for MidX and MidY then using CalculateField_management to update this new field with the value from the old one ("MID_X" and "MID_Y")

This isn't working:

Error message:


texts:
I have several DEM rasters in different projections and at different resolutions.  The rasters are not rectangular, but they have appropriate NoData values.  They overlap enough that I should be able to make a seamless dataset across the entire area, but is there a way I can use the GDAL tools to warp/mosaic them to get that seamless dataset?

I'm currently using the following command line to warp and mosaic the datasets:

gdalwarp -t_srs epsg:3857 -tr <calculated-x> <calculated-y> -r cubicspline output.tif input1.tif input2.tif ...

The problem is that the edges of the raster datasets don't always match each other, so I get effective cliffs at the edges.  Is there a way to average them together for a few pixels near the edges of data for the dataset on top?

I'd really prefer to be able to do this in batch via command line or Python bindings.

Edit:

Here's an (extreme) example of what I'm trying to avoid:

I'm not hugely fond of the ridges at the edge there.

I'm using the USGS NED 1/3 data as a seamless starting point, so I'm guaranteed to have something at every pixel (and the higher-resolution datasets that don't cover the entire area will always have something to be blended with).


texts: Calculating distance from a point to all nearby polygons?

texts:
I'm using ArcMap 10.3 on Windows 10 and I have the 3D analyst, ArcScan, Geo-statistical analyst, network analyst, publisher, schematics, spatial analyst and tracking analyst extensions.

How do I calculate the distances from my points (n = 157) to the polygons (n = 1891) that are within a 2km buffer around the points with the intent to use the distances in SPSS analysis and to generate a habitat suitability model using ArcMap (which will involve rasters)?

This is what the points and polygons within a 2km buffer look like:

The Generate Near Table tool produces the following:

This isn't very useful, as the object_ID field is blank so there is no way to know which polygon the distance refers to. I need to know this as the polygons represent habitat types and want to find out the distance from the point to all the habitats within the buffer. Is there a way to do this? I've tried inverting the calculation (i.e. select the polygons that are closer by and run the near distance from the polygon to the point rather tan from the point to the polygon) but that didn't fix the problem of not knowing what habitats the polygons are.

I have tried selecting by attribute to select pertinent habitats (e.g. woodlands) but doing this then causes a generic 999999 error when using the Near or Near Table tool. I've also tried selecting by attribute and exporting the table, but didn't realise the resulting table can't be used in Near.

I used the "planar" method when using the Near tool and my data are in the Irish grid coordinate system.

I've tried the euclidean distance tool, but I've no idea what to do with the output (i.e. the "heat map" type output and the lack of an attribute table).


texts:
How to achieve the statusbar as the same effect with the arctoolbox‘s tool called in ArcMap add-in?
i do not know how to use the progressbar as a result of my project involves in a large mount of data calculate and it is very complicated.So I choose the startbar the same effect
with the arctoolbox‘s tool(when the tool was called ,the stausbar show ,while the work is finish ,the statusbar is hidden).
I want to achieve the same effect as below (picture link)：


texts: In-line variable %Name% problem in Model Builder

texts:
Having a problem with the inline variable structure of %Name%.  Attached is my model builder program. I want to be able to take the Name of the file (which is produced as a variable in the iterate raster loop) and use it to name the final output of the model.  The main idea is that I want to retain the final name.  The naming structure within my work space is the date on which the raster was made (satellite imagery).  I wrote a few batch files to eliminate erroneous characters in the file name leaving a YYYYMMDD schema.  I want to be able to insert the final work space of rasters into a raster catalog in order to show time (using the Name as the temporal discriminator).  I just need to name the final output like the input file.

When I use the output name E:\ArcGIS\GIS_Data\Kuwait_Clipped_Normalized.gdb\%Name% I get an error.

Not sure why I can not reference %Name% as the file name as the Reclass by Table output.

I have even tried a few other options, but they did not work.  I used the Parse Path on the Raster output and referenced %Value% as the file name (instead of %Name%).  I also calculated a field as well, with no success.

Or, Is there a way I can carry over the file name attribute with the raster?  I just need to retain the file name variable to preserve the time series.


texts: Load a GeoJSON layer with mixed geometries via QGIS's Python module

texts:
I have been able to connect my collection in MongoDB , write and store the collection to as .json file and  load the GeoJSON file as a layer by QGIS's Python module unless there are mixed geometries.
However, I have a problem when the collection in MongoDB has more than geometry type. It loads only one geometry to QGIS canvas. But when I drag and drop the file, I can load both geometries. So .json file, that has been created, is fine.

client = MongoClient()
client = MongoClient('localhost', 27019)
db = client['pb']
cursor=db.section.find()
feature = QgsFeature()# create layer

os.chdir( "C:\Users\Public" )
feature =[]

for data in cursor:
    geom=(data['geometry'])
    id  = (data['properties'])

    feature.append({"geometry":geom, 'properties':id})

json_str = json.dumps(feature)
final={"type": "FeatureCollection","features":feature}

with open('data_3.json', 'w') as f:
     json.dump(final, f)

vlayer = QgsVectorLayer('C:\Users\Public\data_3.json',"mygeojson","ogr")
QgsMapLayerRegistry.instance().addMapLayer(vlayer)

Any advice on what I could do?


texts: Error while placing labels within polygons using ArcGIS JS API and a external Geometry Service

texts:
I have published in my ArcGIS server a map service with a table containing the polygons of all countries on Earth. So this table contains more than 200 rows and the geometry field is of type  type: esriFieldTypeGeometry , alias: Shape.

What I would like to do is creating a label inside this polygon, so what I did was modifying this ESRI example (Create points for labeling) that, instead of allowing the user drawing a polygon to be sent to the external Geometry Service, expects a click from user's side, retrieves clicked-country geometry with a Query operations and draws a predefined label on the country. Of course, last step is not done yet, but this would be the idea.

The problem is that after building a new Graphic object with the new geometry retrieved from the layer and call ESRI Geometry Service as in the example, I got an error on Firebug like this:

TypeError: _c94 is undefined
...eateSegments:function(_764){_764.shape.path=_764.vmlPath;_764.segmented=false;_7...

And well, I don't know why is this happening as error is not verbose at all. I noticed that the error is caused by:

http://serverapi.arcgisonline.com/jsapi/arcgis/?v=3.2 (On Line 34)

So these are my doubts:

The ESRI example comes, out of the box, with the 1.6 version of ArcGIS JS API. And I changed it to 3.2, because if not, I got this error just when starting the application because I imported the Feature Layer with dojo.require:

"NetworkError: 404 Not Found - http://serverapi.arcgisonline.com/jsapi/arcgis/1.6/js/esri/layers/FeatureLayer.xd.js"
uncaught exception: Could not load cross-domain resources: esri.layers.FeatureLayer
Is it possible that I need a proxy to pass the country geometry to Geometry Server? I mean, if you check the following image, you can see that there are about 2800 rings. I am not sure if these data are sent to the geometry service, but if so, I guess that the GET request will be huge and therefore, a proxy should be involved, isn't it?

Please, can anyone provide some help?

Thank you very much!


texts:
I want to run python script. I don't know the name of toolbar in python script language : menu toolbar, browser panel toolbar, and label toolbar.
The screenshots of the toolbars for which I don't know the python names of are:

To test the name of toolbar is true, I try to type in python script like this :

iface.menuToolBar().setVisible(False)

iface.browserpanelToolBar().setVisible(False)

iface.labelToolBar().setVisible(False)

But, the toolbar is still be visible. So the name of these toolbar is false. So, what the name of these toolbar?


texts: Wrong direction in r.watershed slope length calculation?

texts:
Or have I completely misunderstood something now? As a programmer I took it for granted that the row number would change if I traversed the map up and down in the program along the North-South axis while the columns would change when I traversed the map sideways from west to east. In the raster, the (0,0) coordinate would represent the north-western-most cell, with values increasing towards the south-east.

And this fact seems supported by the following line taken from init_vars.c:

max_length = 10 * nrows * window.ns_res + 10 * ncols * window.ew_res;

But then I found this in the code for calculating the length of the slope for a given cell in slope_len.c:

if (r == dr)
    res = window.ns_res;
else if (c == dc)
    res = window.ew_res;
else
    res = diag;

if (top_ls == half_res) // half_res = initialization value
    top_ls = res;
else
    top_ls += res;

Here we are dealing with a cell given by coordinates (r,c) which flows into another cell (dr,dc). So the code reads to me that if the two cells are on the same row, add to the slope the map's North-South resolution. Should it not be  West-East since the slope is horizontal rather than vertical and vice versa for the columns?

Either I misunderstand how the length is calculated or this seems like a bug in the program for me. Also on quadratic cell sizes this will not be noticeable but it is important that I get it right for other cases.


texts: Calculating total time, last value minus first value in ArcPy?

texts:
For the script tool I've written (below), one of the things I want to do is after user input has pulled tabular data (into a feature class) I want to take the last value of the 'logdt' and subtract it by the first value of the 'logdt'. 'logdt' being a datetime value in this format '2/1/2014 10:00:11 AM'. The reason being is that I want to determine the total length of time that has passed since the first and last pings for a specific vehicle in a given day. Once I have made this minor calculation I will then use it to determine how many, if any, lost pings occurred that day.

I'm sure it's very simple but I am newer to arcpy/python and cannot think of a way to do the logdt.lastValue - logdt.firstValue subtraction.

*EDIT: Previous code has been re-written and only pertinent code below is shown.

# Calculates the total number of rows for newly added column "Total_Pings"
numRows = int(arcpy.GetCount_management(queryLayerName).getOutput(0))

# UpdateCursor that will write the value of numRows to the "Total_Pings" column
cursor = arcpy.da.UpdateCursor(queryLayerName, "Total_Pings")
for row in cursor:
    row[0] = numRows
    cursor.updateRow(row)

# SearchCursor that will read the values of "LOGDT" and return the first value (earliest time)
firstLogdt = [row[0] for row in arcpy.da.SearchCursor(queryLayerName, "LOGDT")][0]

# UpdateCursor that will write the first (earliest time) "LOGDT" value to the field "First_Time"
cursor = arcpy.da.UpdateCursor(queryLayerName, "First_Time")
for row in cursor:
    row[0] = firstLogdt
    cursor.updateRow(row)

# SearchCursor that will read the values of "LOGDT" and return the last value (latest time)
lastLogdt = [row[0] for row in arcpy.da.SearchCursor(queryLayerName, "LOGDT")][-1]

# UpdateCursor that will write the last (latest time) "LOGDT" value to the field "Last_Time"
cursor = arcpy.da.UpdateCursor(queryLayerName, "Last_Time")
for row in cursor:
    row[0] = lastLogdt
    cursor.updateRow(row)

I am receiving the following error in Python:


texts: QGIS radius buffer question, can it be more circular?

texts:
When I request a radius buffer in QGIS (Dufor). my buffer comes out as a 20-sided polygon rather than a true "circle".  This has lead to some missed blocks when clipping underlying vector layers for population density analyses (e.g., they would be within a circle, but this 20-sided polygon's straight sides preclude fringe blocks from being included in the radius).

is there a way I can go into the radius buffer command lines and either increase the number if sides to make the area more circular, or make it truly a "radius" buffer?


texts:
I had symbolized different layers as you can see here:

if (vector_layer.name == "Buffer"){
vector_layer.style = {fillColor: "blue", pointRadius: 7, fillOpacity: 0.3,
strokeColor:"#F2F2F2"};

vector_layer.addFeatures(geojson_format.read(result.getData()));
legend.loadLayer(vector_layer);

It works fine, but the problem is when I try to change the color with Ext JS Library 3.4.0. The fact is with this type of layers, they are results from turf operation and ext js don't allow me to change the color. I try to choose a different color one but I doesn't work.

Code:

var cp = new Ext.ColorPalette({value:'993300', "renderTo" : cpid});
                    cp.on("select", function(palette, selColor){
                        console.log("selected: "+selColor);
                        colorPicker.hide();
                        var newStyle = new OpenLayers.Style(
                            {"fillColor" : "#"+selColor,
                            "strokeColor" :"#"+selColor,
                            "fillOpacity" : 0.4,
                            "strokeOpacity" : 1,
                            "strokeWidth" : 1,
                            "strokeLinecap" : "round",
                            "pointRadius" : 6});
                        layer.styleMap.styles["default"] = newStyle;
                        layer.redraw();

Any idea how to solve it?


texts: Using Range in Python expression of ArcGIS Field Calculator?

texts:
what is the python expression for defining a range of numbers?

I have tried numerous combinations but cannot get a simple range of angles (100 -200) to compute correctly.

Screen shot :


texts: Algorithm to find the boundary grid points

texts:
Let us consider a set S consisting of a collection of grid locations (i,j),i,j=1,\ldots,n,  in 2-dimensional gridded plane. Here the point (i,j) corresponds to a rectangular grid in the X-Y plane.

This set has the property that one can reach any element (i,j)\in S from any other element ({i'},{j'})\in S through a series of neighbouring grid points which are also in S.

I am interested about the shape of the set S in the plane, for this reason I need the boundary grid locations of the set S.

Note: R-code

data=read.csv(file="trial.csv")

require('ggplot2')

pred.points <- ggplot(data = data,
                      aes(x = lat,
                          y = lon,
                          colour = variable)) +
  geom_point()

print(pred.points)

I need an algorithm to find the boundary grid locations for the set?

Here I don't know whether the set is concave or convex.


texts: Rotating point feature on line feature in ArcMap

texts:
I have a line and a point feature class in ArcMap. I want to rotate the point feature on the line feature. I will use the Symbol Rotate Option with a defined field. I have a script for Field Calculator:
def dreh(az,ab,az1,ab1,typ):
    if typ == "tunnel":
        return 0
    elif ab == 0:
        return 90+az
    elif ab1 == 0:
        return 90+az1
    else:
        return 0`

dreh(!azimut!, !Abstand!, !azimut_1!, !Abstand_1!, !typ!)

How should I define azimut, azimut_1, Abstand, Abstand_1?
The left image is not turned (current state), the right image is turned (desired state).


texts: Standardized Precipitation Index (SPI) equation using ArcGIS Desktop

texts:
I usually calculate monthly or 3-monthly Standardized Precipitation Index (SPI) using spreadsheet, based on daily rainfall from weather station. Attached is the spreadsheet example link.

And now I would like to use satellite precipitation estimate (raster data) to calculate the SPI. The problem is in the spreadsheet I use some function like: FREQUENCY, GAMA.DIST, NORMSINV.

Is there any similar above function in ArcGIS Desktop?


texts:
I am new to QGIS and have a problem that some might be able to help me with.

Background: Windows 10 and all OS is up to date,
QGIS standalone version 2.18.13 (64 bit)

Learning tutorial http://www.qgistutorials.com/en/docs/making_a_map.html

When following this tutorial I can get up to step 21 but I can’t get the frame style to show up on the map?

My original post did not mention anything in regards to adding a grid and zebra border I can get that.
I don't know how to upload a screen caption but if you would have looked at the link I provided at step 21 you will see a pink rectangular box around the area of Tokyo which is a frame style.
This is what and can get.


texts: Snapping Elements in Layout to Neatlines using ArcGIS Desktop?

texts:
I was wondering if there is a way to snap the center of an element to a neatline using ArcGIS 10?

For example I have a neatline running down the center of the map in layout view and I want the center of the title element (text box) to snap onto that neatline. Then I can have the Title and subtitle aligned even though they are different sizes.

GIS for the government entity I worked for last had this function. For the life of me I cannot figure out how to make it work on this place’s GIS.


texts: Web service instance configuration (concurrent users) and its performance

texts:
I was recently reading about Server Software Performance and the different settings and recommendations presented by ESRI over their wiki for System Design and architecture strategies. I came across a graph that depicts the web service instance performance based on the multiple incoming requests.
The graph below shows the Display Response Time as Performance in seconds on the left, and the Displays Per Minute (DPM) can be seen on the right side of the graph. One thing it does not show is the Host platform utilization percentage. How can we calculate it?

The authors mention the percentage usage in the discussion but the graph doesn't seem to show it


texts: Creating a 3D TIN or similar 3D Delaunay triangulation in R

texts:
Background: I have topographical datasets with height information as a starting point. The endgame is to create a dataset of triangulations which I can use to import into Parasolid to be used further downstream. The takeaway point is that I need triangles.

More about the data: The data is of the type Large SpatialLinesDataFrame which represents topographical maps. I have some experience with converting this to spatialpoints/interpolating/rasterising etc., so this is not the crux of the problem. I understand the problems that arise (or don't) whether the data is regularly or irregularly spaced and am confident I would be able to coax the data between the two as required, so this is also not the crux of the problem.

The problem: Regardless of what the data format, I can't seem to find a library in R that would be able to calculate Delaunay (or any suitable) triangulations. Yes, there is deldir and RTriangle but both these are limited to datasets in the same plane i.e 2D. Mountains tend to be 3D, so therein lies my problem.


texts: Build Raster Attribute Table with GDAL Utilities

texts: Is there a way to build a raster attribute table using one of the GDAL utilities?  I know you can view it using gdalinfo and I have successfully created one using Python, but it would be so much easier to just run a command line utility to build a raster attribute table on a raster.


texts: Changing NULL values of raster to 0 in QGIS

texts:
I'm having problems with QGIS. I have some rasters that appear to have NULL values, even though I reclassified them to obtain a Boolean raster.
How can I change the NULL values to 0?


texts: I'm trying to install QWAT but when I get to the step where I have to create the pg_service.conf file I get confused I don't where to put it.


texts: Removing property page from Layer Property Sheet using ArcObjects?

texts:
I need to remove definition query tab of the Feature Layer or Table property page. Is it possible to do this with arcobjects at all?

Here is what I am trying to do:

UID uid_comCategory = new ESRI.ArcGIS.esriSystem.UID();
//Component Category of Esri Layer Property Pages
uid_comCategory.Value = "{1476C782-6F57-11D2-A2C6-080009B6F22B}";

UID uid = new ESRI.ArcGIS.esriSystem.UID();
// GUID of esriCartoUI.LayerDefinitionQueryPropertyPage
uid.Value = "{1476C788-6F57-11D2-A2C6-080009B6F22B}";

var ccm = new ComponentCategoryManager();
ccm.SetupObject("C:\\Program Files (x86)\\ArcGIS\\Desktop10.2\\bin\\FeatureLayerUI.dll", uid, uid_comCategory, false);

I get this error at the last line:

  Value does not fall within the expected range.

Furthurmore I tried to remove it manually using Component Category Manager utility, but it doesn't work there, too.


texts: Matching polygons to intersecting polygons using ArcMap

texts:
I am trying to match a set of polygons (A) to another set of polygons (B) using ArcMap. Both sets overlap themselves, and each other. Here's an example of what it looks like in my data:

The green polygons are the (A) set (about 12000 polygons), while the blue polygons are the (B) set (about 30000 polygons). Ultimately, the dataset that I want to create should have, for every polygon in (A), the set of polygons in (B) that it overlaps, with the data from (B)'s attributes table and the area of the overlap.
I've tried using the intersect tool, but it has crashed on me twice now after several hours.
Is there a better way to do this?
I haven't used ArcGIS Pro before, but have access to it and would be happy to use some tool from there.


texts: Getting gdal_translate 0ERROR 1, but the raster has a data table

texts:
I'm trying to get a large number of GeoTiff files into a geopackage for a project we'll be using QField for. Ideally, I'd like to batch the process, but in testing small samples out I'm running into problems with single images from the original data. Here's an example of the problem:

The command:

    C:\OSGeo4W64\bin>gdal_translate --config OGR_SQLITE_SYNCHRONOUS OFF
-co  APPEND_SUBDATASET=YES -co TILE_FORMAT=PNG_JPEG -a_srs EPSG:26917
-of GPKG "C:[Path]\ORTH175670482202017FBS.tif" "C:[Path]\ORTH175670482202017FBS.gpkg"

The result:

Input file size is 10020, 10020

0ERROR 1: RASTER_TABLE must be set to a non empty value

I am fairly certain that the rasters have data tables - the output for gdalinfo on the file looks legitimate (see below), but I really wouldn't know any different since this is a new level of GIS data manipulation for me.

What is this problem, and is there any way to fix it?

I'm running in a Windows environment (Windows 7), using OSGeo4W/QGIS 2.18.13.

gdalinfo output:

C:\OSGeo4W64\bin>gdalinfo "C:[Path]\ORTH175670482202017FBS.tif
"
Driver: GTiff/GeoTIFF
Files: C:\Users\James\Milton_Tiff\ORTH175670482202017FBS.tif
Size is 10020, 10020
Coordinate System is:
PROJCS["NAD83 / UTM zone 17N",
    GEOGCS["NAD83",
        DATUM["North_American_Datum_1983",
            SPHEROID["GRS 1980",6378137,298.257222101,
                AUTHORITY["EPSG","7019"]],
            TOWGS84[0,0,0,0,0,0,0],
            AUTHORITY["EPSG","6269"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4269"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",-81],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH],
    AUTHORITY["EPSG","26917"]]
Origin = (566999.000000000000000,4823001.000000000000000)
Pixel Size = (0.100000000000000,-0.100000000000000)
Metadata:
  AREA_OR_POINT=Area
  TIFFTAG_RESOLUTIONUNIT=1 (unitless)
  TIFFTAG_XRESOLUTION=0
  TIFFTAG_YRESOLUTION=0
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  566999.000, 4823001.000) ( 80d10'13.73"W, 43d33'25.04"N)
Lower Left  (  566999.000, 4821999.000) ( 80d10'14.18"W, 43d32'52.56"N)
Upper Right (  568001.000, 4823001.000) ( 80d 9'29.07"W, 43d33'24.71"N)
Lower Right (  568001.000, 4821999.000) ( 80d 9'29.53"W, 43d32'52.24"N)
Center      (  567500.000, 4822500.000) ( 80d 9'51.63"W, 43d33' 8.64"N)
Band 1 Block=10020x1 Type=Byte, ColorInterp=Red
Band 2 Block=10020x1 Type=Byte, ColorInterp=Green
Band 3 Block=10020x1 Type=Byte, ColorInterp=Blue
Band 4 Block=10020x1 Type=Byte, ColorInterp=Undefined


texts: SQL Geometry Query similar to Select by Centroid

texts:
Does anyone know if it's possible to query SQL geometry by something like 'centroid of polygon'...similar to ArcMap's select by location?

I'm using SQL Server 2008 R2 as the DBMS.

I have two polygon layers that are not exactly coincident and if I use STIntersects on the target polygon I get two different values where the polygons are not coincident.  I was wondering if there is a query that can be used instead of STIntersects that would be similar to selecting the centroid of the polygon...?

Thanks in advance!


texts: I have created a grid using the point density function in VM. I would like to crop this image to a region created from buffering a polyline feature set. Is this possible? Alternatively I have used the set clip region tool, however this only hides the rest of the image where I would like to save the cropped region as a new layer/image.


texts: ArcMap running scripts slower if you import Arcpy

texts:
I am finding ArcMap runs slower, if you have written a script, which imports arcpy.

You may wish to write a script that can be operated in ArcMap, through the python window, or as a standalone python script.  If so, you'll have to import arcpy.  However, if you do this and run it in the python window, it will run slower.


texts: Subtracting one polygon from another in MapInfo 15

texts:
I want to make a hole in a polygon.
In older versions of MapInfo I did it with ease, but how do I do it in the current version which is 15?


texts: Leaflet: Bounds from getBounds are changed after rendering

texts:
Leaflet Version: 1.4.0. At the time of render, we are saving bounds by calling getBounds() on the layer.
After the map is rendered, zoomend event is called and the bounds captured are different from previous getBounds(). This looks weird.

We want to persist the bounds and apply them after reopening. So we are capturing bounds on zoomend event. As these bounds are different, we see variation of the extents.
Code sample:

polygonLayer.getBounds(); //This will give bounds of the layer
-----
mymap.on('zoomend', function(event) {
  var bounds = event.target.getBounds();
});

These two are different as shown below. If I need to save my extent and apply again what should be done.


texts: Seeking PostGIS query to find LineStrings (Route) entering or crossing circle area?

texts:
I need to know Postgis query to find LineStrings (Route) entering or crossing circle area. I tried using ST_Crosses function but it is not resulting with correct output.

Any suggestions?


texts: How to accomplish conditional labeling based on map scale and attribute in QGIS 1.7

texts:
How can I label with qgis 1.7, based on attribute and scale of map?

For example scale>1:1000 then label all those features whose area>100000 at font size 12 and those features whose area<99999 at font size 5.


texts: Downloading subset of MODIS data through HTTPS

texts:
I need to download 'MYD04_3K - MODIS/Aqua Aerosol 5-Min L2 Swath 3km' data for India. I found two ways to download it:

Go to the website and select product details etc. I have been able to automate this process using Python. But, one of the drawbacks is that there is no way to subset the data. There are no location specific ids (e.g. tile code etc.) in the file names. So, one will have to download the files for the entire world, process the data, and then subset the data.
Go to the website and select product details etc. Select the product, time, location and process the order. Here, I am basically ordering the files using website's web-service. It does not allow you to download more than files 2000 files at a time. The links/details are then sent to my email and I can use the links to download the data.

Is there any way I can use an API or write a script to download data for a particular location without availing website's order/processing feature? I want to automate the download rather than going to website and order the data every time I download.

Edit: This question is not a duplicate of the question suggested. As I mentioned in the original question, I want to avoid subsetting the files manually on the website, and also, I have not found any location specific details in the filename. Here is a sample filename: MYD04_3K.A2004001.0015.061.2018012003848.hdf

Some datasets from MODIS have tile information, which can be used to download only subset of the data. In my case, I couldn't find any such thing. So, I asked if anyone knows how to download the subset of data automatically, without manual selection.


texts: UTM to WGS84 in QGIS 3.18 circle not round

texts:
I am using QGIS 3.18.
I have a GeoJSON file with circles in UTM CRS. I exported it to WGS84 and my circles became not round, more like an ellipse.
What is the problem?


texts:
I am new to OSM. I wonder what a "relation" is in an OSM file. Is it a road segment? In addition, I want to convert an OSM file to the following form:

s_1, t_1, lat_s_1, lng_s_1, lat_t_1, lng_t_1
s_2, t_2, lat_s_2, lng_s_2, lat_t_2, lng_t_2
......

where each line represents a road segment, s_i and t_i are the endpoints of the segment, and lat_s_i and lng_s_i are the latitude and longitude of s_i. What part of an OSM file should I look at to get this information?


texts: How do you use OSGeo-Live with Windows and Mac Installers?

texts: I'm familiar with installing from bootable USB, and have done so with Ubuntu, Linux Mint, and OSGeo-Live. But there's a version of OSGeo-Live which includes Windows and Mac installers. I would guess that this uses Wubi, but I don't see that explicitly spelled out anywhere, and Wubi is an EXE whereas the download is an ISO. So, what do the Windows and Mac installers actually do, and how do I use them?


texts: How to create an "arcgis-style" hillshade?

texts:
I came across some high res DEM and would like to create a quality hill shade.
Target quality should be ESRI's world_hillshade as shown here:
https://services.arcgisonline.com/arcgis/rest/services/Elevation/World_Hillshade/MapServer
that is not the normal hill shade I would obtain using the relevant function in my GIS (I use QGIS).

I am sure it is a composite hill shade made with multiple different illumination images blended together, still I cannot obtain a similar result. I see mountain ridges should appear white, so apart from the composite blending I suspect there might be some slope layer involved. Has anyone developed a procedure for obtaining similar "artistic" hill shades out from a raster DEM?


texts: What is the expression for getting a subset of characters in QGIS Field Calculator?

texts: I have a column of 10 character names e.g. L3NAM1988K. How can I delete the first 5 characters and the last 1 character using an expression in the QGis Field Calculator?


texts: Convert polygon shapefile geometry to text file in a specific layout, including number of nodes per polygon

texts:
I'm trying to convert a large shapefile (more than 20,000 polygons) to a text file to use in a modelling software. It needs to be in a specific format which is number of nodes, followed by all x coordinates followed by all y coordinates with each polygon on a new line (Pictured below) (first line is total number of polygons).

Currently I am using QGIS 3.6.3 and saving the shapefile as a CSV using WKT geometry. This however does not include the number of nodes and coordinates are formatted as xy xy. Please note that I am not constrained to QGIS.


texts: View attachments in popup when using the editor widget

texts:
I am using ArcGIS javascript API to construct an editor app for a layer. I am using the default out of the box editor widget to edit and display data for an app. The layers I am editing contain attachments, which show up as a link within the editor's attribute viewer. Attachments are by default opened in a new window. I would like a dialog box with the attachment ( which is images) to load via a ajax request once the link is clicked.

Has anybody been able to get attachments, such an image to open in a popup of some type, while using the default editor widget?

My logic was to trigger an event handler that changes the link to the attachment, after the attribute inspector loads but I can't seem to find an event that is triggered during that time.

Am I going at this from the right approach, or would there be a better way to accomplish this with the default editor? I know you could probably write a custom widget to to do this, but I would prefer to stick with the functionality that the editor widget provides.


texts: Compare two layers and find the percentage of the shared areas that are different between them

texts:
I have two shapefiles that identify flood zones for the same general area but in different years. The maps, however, do not cover exactly the same geographical areas (e.g., shapefile 1 covers counties A, B, C, and D while shapefile 2 covers counties B, C, D, and E) and the maps change from one year to the next as some areas are reclassified.
I need to do two things:

Restrict the exercise to the areas that are covered by both shapefiles (e.g., counties B, C, and D). The shapefiles do not, however, have a county identifier. They only have polygons and the flood zone designation (e.g., A*, X, VE).
For the area in which they do overlap, I need to compute the percentage of areas designated as, say, "A" in shapefile 1 that are still "A" in shapefile 2, and the percentages that were reclassified to other designations (e.g., X, VE). And I need to do this for all possible classifications that exist in shapefile 1.

The expected outcome has two parts. First, the percentage of the overlap area that has the same designation in both shapefiles. Second, the percentage of the overlap area in each category in shapefile 1 that is reclassified into different categories in shapefile 2 (as in a transition matrix).
I'm using QGIS for this (3.10 or 3.16) on Linux.
Two examples follow:

This shows shapefile 1 with the polygons identified by black lines and shapefile 2 with the polygons in red. It shows areas that are only covered by shapefile 1 (white area and black polygons), only by shapefile 2 (white area, red polygons), and the overlap area.

Adding labels (black for shapefile 1 and red for shapefile 2) for a zoomed in area to show reclassification.


texts: Why would Dissolve stop working in ArcGIS for Desktop?

texts:
I  have an ArcMap Add-In that completes various geoprocessing steps, including Dissolve (using ESRI.ArcGIS.DataManagementTools.Dissolve).  This tool worked at 10.0, but I just upgraded to 10.1 SP1 and it no longer works.  I have run the tool on the same feature class in 10.0 and 10.1 and it worked in 10.0 but does not work in 10.1 SP1.

I'm attempting to dissolve on 2 fields, single part, in case that's relevant.

I can still complete the dissolve manually inside of ArcMap, but every time I try to run it programmatically, it crashes with the ever-helpful "Error HRESULT E_FAIL has been returned from a call to a COM component."

Upgrading to 10.2 is not an option because the rest of my organization won't be upgrading to 10.2 for at least another year.

Has anyone else encountered this?  Any suggestions on how to solve this problem?  Knowledge of something that changed at 10.1 that would require some extra conditions in order for dissolve to run correctly...?


texts:
I have 2 projects, with data from different regions and with different EPSG, 1 with EPSG 31983 and the other with EPSG 31984, according to the folder structure shown below, in Windows explorer.

In both projects I have the same folder structure and file names. I created a QGZ project in QGIS 3.12, using project data 01 (EPSG 31983), included the data (vector and raster), defined the project's symbology and visual coordinate system as being EPSG 4674 (figure below).

I saved the project with the name “Project.qgz”). I went back to Windows explorer and copied the QGZ to project 2 folder (EPSG 31984) Crtl + c / Crtl + v.

I opened QGIS again, now choosing the QGZ project that I previously copied into the folders of project 2 (EPSG 31984) The project opened normally, however when reading the two themes QGIS kept the EPSG of the previous themes (project 1), or that is, he fetched the data correctly, applied the symbology but kept EPSG 31983 for project 1 not updating to EPSG 31984 for project 2.

I need to apply this example to several projects, just by copying and pasting the “.QGZ” from one folder to another.

How do I get QGIS to automatically update the EPSG of the data being presented, in the same way as it updates the symbology for the new data being loaded?

The shapefile file inside the folder 'Project EPSG 31984' is set to EPSG: 31984. I checked it. I am referring to CRS. The problem is that QGIS is memorizing the EPSG of the shapefile within the QGZ project, and then, when I make a copy of the QGZ project (copy / paste) to the other project, QGIS does not change the EPSG when reading the new shapefile.
I am trying to show that when copying the .QGZ project from one folder to another, QGIS updates all data to contain the shapefiles in the new folder, however, it does not update the coordinate systems (CRS), causing problems.
There any suggestion to solve this?


texts: How to convert .ps to .png in a ps.map script?

texts: I already managed to create a .ps file. How do I convert it to .png using a GRASS script?


texts: calculate bearing from one MapPoint to another MapPoint

texts:
I am using Esri.ArcGISRuntime package inside a C# program we are writing.

I need to create specific shapes from a reference point and project them in different directions.  This is easy to create a new point 5000m away on a 45 degree bearing.  ie

var startPoint = new MapPoint(125.4, -33.2, SpatialReferences.Wgs84);
var endPoint = GeometryEngine.GeodesicMove(startPoint, 5000, LinearUnits.Meters, 45);

Later on when I need to reverse engineer the points, I can easily recalculate the distance between these points using

var distance = GeometryEngine.GeodesicDistance(startPoint, endPoint, LinearUnits.Meters);

However, I cannot find anyway to determine what the original angle or bearing was between these 2 points.

I have tried using Arctan2, but the MapPoint coordinates are skewed by earth curvature (I think) and therefore straight forward calculation is incorrect.

var degrees = Math.Atan2(endPoint.Y - startPoint.Y, endPoint.X - startPoint.X) / Math.Pi * 180;

The above returns a figure of roughly 39 degrees (PS I tried changing X & Y around as well, no better).


texts: How to set QGIS mapCanvas to exact values

texts:
My QGIS 2.8.2 standalone app's mapCanvas fills a central widget with a fixed quadratic size x=y=1000 pixels. When I set the mapCanvas extent() to an area defined quadratic as well given by a QgsRectangle in UTM coordinates, it happens that the extent() gets wider then appropriate. The y-axis fits well. The mapScale() has a significant to great value (about 20% to both sides). The only feature I want to show at this time is a polygon reflecting the extent (kind of background), so this mis-fit is very obvious.

After the first scaleChanged() event, eg. zooming in and out again, the extent() reflects the correct value with width equal to height, and the scale() is correct as well.

Has this something to do with the 5% extra on the extent or with the initial view (without any features)? How can I set the map extent such that it fits well into the widget without any extra white spaces?


texts: Creating quadtree index on latitude, longitude data in PostGIS

texts: I have a table of locations in which one of the columns represents the lat-long for the location.  These coordinates have been extracted from the OpenStreetMap file. How can I create a quadtree index on these coordinates using PostGIS?


texts: Is there a standard st_relation_readable() function for translate DE-9IM strings?

texts:
Spatial databases, such as PostGIS, implement the DE-9IM model with function  st_relate(),
which returns a string with domain of {0,1,2,F}, meaning 0=point, 1=line, 2=area, and F="empty set"...
See DE-9IM at Wikipedia.

The DE-9IM codes are not understandable for programmers and mortal users,
but each code  have the semantic of a list of predicates, then,
it can be translated into a human-redable list of predicates.

Example: when st_relate(a,b)='FF1F0F1F2' or st_relate(a,b)='FF10F0102' or st_relate(a,b)='FF1F0F1F2' the DE-9IM code can be understood as "Intersects & Touches" or "Intersects & Touches with a point".

Is there a standard function or library that implements this translation, from a DE-9IM code to a human-readable string?

REPORT EXAMPLES: instead Report-1 I need Report-2.

Report-1 (have a standard string, but not human readable)

 0F1F00102 |       2 | group 1
 0F1FF0102 |     616 |
      total      618

 1FF00F102 |       1 | group 2
 1FF0FF102 |       1 |
      total        2

 1FFF0FFF2 |      23 | group 3
      total        23

 F01FF0102 |      31 | group 4
 F01FFF102 |       1 |
 FF10F0102 |      42 |
 FF1F00102 |   19673 |
 FF1F0F1F2 |      25 |
     total     19772

Report-2 (human readable... There are a SQL function to do this?)

  618 | point-Intersects & point-Crosses & line-Overlaps
    2 | point-Intersects & line-Crosses & line-Within & point-CoveredBy
   23 | line-Equals
19772 | point-Intersects & point-Touches


texts: How do I retain all attribute data when clipping two polygons in R?

texts:
I would first like to say that my question is very much similar to this answered one, but the proposed solutions don't seem to be working for me. I think that my question is slightly different, as I'm performing a clip/intersection rather than a difference. Onto it!

My input can be found here at the FEWSNET website, and an example country outline can be downloaded here at the GADM site. I've been using the Ethiopia ESRI shapefile.

What I'm trying to do is quite simple. I first clip the country shapefile to the FEWSNET shapefile. Once the clip is performed, I need to symbolize the ML1 variables from the latter shapefile. In QGIS, I perform the vector clip, then with the newly clipped shapefile, I go to Properties > Style, change it from single to categorized, then choose my symbology for the ML1 category.

The code that I have so far has also been quite simple:

library(rgdal)
library(rgeos)

#read in ML and Ethiopia shapefiles
ea.ml1 <- readOGR(dsn = "D:/Mapping-R/Ethiopia", layer = "EA201507_ML1")
eth <- readOGR(dsn = "D:/Mapping-R/Ethiopia", layer = "ETH_adm0")

#perform clip
eth.clip <- gIntersection(ea.ml1, eth, byid = TRUE)

The problem is that the data that I need to symbolize (ML1) exists under ea.ml1@data$ML1. When I perform the clip, it appears that ea.ml1@data is lost altogether. How do I retain that data with my clip so that I can properly symbolize it?


texts: How do I perform a spatial query on this area?

texts: I am in the process of trying to perform spatial analysis bike on bicycle trip data publicly found on Citibike's website. I want to extract data for each square in the rectangle (the boundary-lines shapefile that encompasses the bike stations. In the box group, I have a lines layer that defines the grid of the rectangle. What I wish to do is extract data of each square in this rectangle, ex. find the middle point, find the number of stations in each square and more. How would one attempt do this? Do I have to make use of raster data instead?


texts:
I am  doing a number of calculations on a timeseries of data and want to create a 2d netcdf file with the results of those calculations.
However, when saving the calculated values in the new array I get the following message:

IndexError: too many indices for array.

I am using python 3.6 with the NetCDF4 and numpy libraries.
val_list=[]
av_array=[]
#create netcdf
out_file=netCDF4.Dataset("test.nc","w", format="NETCDF4")
tempgrp=out_file.createGroup(var)
tempgrp.createDimension('lon', 575)
tempgrp.createDimension('lat', 541)
tempgrp.createDimension('time', None)

longitude = tempgrp.createVariable('Longitude', 'f4', 'lon')
latitude = tempgrp.createVariable('Latitude', 'f4', 'lat')
temp = tempgrp.createVariable(var, 'f4', ('time', 'lon', 'lat'))
time = tempgrp.createVariable('Time', 'i4', 'time')

longitude[:] = 575
latitude[:] = 541

for i in range(540):
    for j in range(574):
                del val_list
                val_list=[]
                for t in range(30):
                    #create list of i,j coordinate
                    val_list.append(data[t,0,i,j])
                    #calculations
                av=sum(val_list)/float(len(val_list))#example
                temp[0,j,i] = av #this is where the error occurs

This is the output of print(temp):

<class 'netCDF4._netCDF4.Variable'>
float32 Chlx(time, lon, lat)
path = /Chlx
unlimited dimensions: time
current shape = (0, 575, 541)
filling on, default _FillValue of 9.969209968386869e+36 used

The idea is to in the end create a 3d array (that is why the time dimension is there).
Can anyone give me a hint to why I get this index error?


texts: Sorting numbers into categories in QGIS?

texts:
I need to sort a column of population figures into categories in a new column so that I can split them out and create layers for each one. Something like 0-99, 100-249, 250-500 etc.

Struggling to figure out a way to do this? See screen grab to show what I mean.


texts: How convert OSM relation to GeoJSON -based all kind of relation member roles

texts:
I just want to convert OSM data to GeoJSON where it is pretty simple with node and way tags of OSM data, but relation which is giving complex can't find the logic to convert OSM relation to GeoJSON format, especially how to form coordinates parameter of GeoJSON based on their relation member roles (ex. Inner, outer,admin_centre,riverbank,side_stream,main_stream,from,via,to,forward,backward,stop_entry_only,stop_exit_only and platform_exit_only).

I have tried many documentation but no logics or code found in Java and Python.

In Python https://pypi.org/project/osm2geojson/ osm Official https://wiki.openstreetmap.org geojson Official https://geojson.org/ sample data: https://download.geofabrik.de/asia/maldives-latest.osm.bz2 https://download.geofabrik.de/asia/bhutan-latest.osm.pbf

As simple, how to form relation member's latitude and longitude in coordinates parameter of GeoJSON based on their roles?

EXample code

public void processRelations() {
    try {
        LineIterable lineIterable = LineIterable.openGzipFile(dir + OsmJoin.REL_ID_COMPLETE_JSON);
        try (JsonWriter writer = createJsonWriter(OsmType.RELATION)) {
            Processor<String, JsonObject> p = compose(entryParsingProcessor, jsonParsingProcessor, new Processor<JsonObject, JsonObject>() {
                @Override
                public JsonObject process(JsonObject input) {
                    // FIXME see if we can extract some useful things from relations
                    return null;
                }
            });
            try (ConcurrentProcessingIterable<String, JsonObject> concIt = processConcurrently(lineIterable, p, 10, 9, 100)) {
                for (JsonObject o : concIt) {
                    if (o != null) {
                        writer.add(o);
                    }
                }
            } catch (IOException e) {
                throw new IllegalStateException(e);
            }
        }
    } catch (IOException e) {
        throw new IllegalStateException(e);
    }
}

I have tried in Java and Python but reason is I could not find the logic behind the conversion of the relation tags from OSM to GeoJSON

I am expecting the logic behind the conversion of the relation tags from OSM to GeoJSON and sample code for that.


texts:
How is R used with ArcGIS for Desktop?

I have seen posts talking about combining rasters while using R.

I have not had much experience with R so not 100% certain of its capabilities.

I had thought that R was a scripting language for data analysis with tabulated data and data bases, but its seems much more than that.

Are there particular R libraries to incorporate it with ArcGIS for Desktop?


texts: Convert dwg (CAD) to shp in any open source GIS software?

texts:
I would like to do what's needed in a GIS (GRASS, QGIS, etc) alone, no other converters-only.

The process is to read the dwg files that would be located in a folder and read them as shp and perform calculations like snapping later or using a geodatabase to store the files and later export them as shp.


texts: Transfering from traditional GIS software developing to internet area

texts:
I am a GIS developer, formerly doing some secondary development with GIS software. I now want to obtain some internet-based GIS area in my career.
Somebody advised me to first read some networking/linux books. When prepared, then transfer to some company to provide experience and so on.
Is there any opensource project belonging to this type: server side provide some service, for example, provide some navigation service(or other GIS area service).  I want to learn something about these software platforms and how to realize the functionality provided by server to client, Scheme about the whole thing. Ideall6 such type of thing is using in c++/c. Easier for me to learn.
My former area and expertise has been using desktop software.
Does anyone have any advice?


texts: Buffering selected objects in the map in Graphical modeler

texts:
Process: Buffer selected point to 3000m then unselect that point and grab all records from a particular layer that intersect with the buffer.
I can't find a way of specifying that the buffer must use only the selected object in the map. Any ideas?


texts: Previewing layer on GeoServer instead of downloading WMS file?

texts:
I am new to Geoserver 2.8.0. I just installed it and passed some layers from PostGIS.

I used the existing, default cite workspace, then added a myProject Store, with type PostGIS. I loaded the layers and set the Bounding Boxes, I enabled Response Cache Headers, set some default styles and set the min/max zoom levels from 0 to6.

In this new Geoserver interface, how can I set or at least check what type of layer is this? Is it a WMS, a WFS? My PostGIS layers are some OSM data, that I converted and saved on PostGIS

My main question is that I can preview the planet_osm_line, planet_osm_road and planet_osm_point but not the planet_osm_polygon. When I hit OpenLayers to preview it in Geoserver, instead of opening a new tab with OpenLayers it automatically downloads a wms  file.

How do I fix this? Is there something wrong with my settings?


texts: Solving the Point in Polygon Problem Using Google Maps and Fusion Tables

texts:
The Point-In-Polygon problem is fairly well documented on this forum. My question is if there is a way to use Google's free tools, i.e. Google Maps, Fusion Tables, and possibly Google Refine, to solve this problem.

I am ultimately attempting to create a heat map of neighborhood boundaries. In basic language here are the steps that I wish to perform:
1. Geocode a list of home addresses using Fusion Table's built-in feature. This yields a set of lat/long coordinates for each address.
2. Import a Shapefile of polygons into Google Fusion to store the boundaries (polygon coordinates).
*I am able to get this far.. The next steps are where I am looking for some guidance.
3. Apply a point-in-polygon algorithm to identify which addresses fall within each boundary. What is the simplest way to do this considering that I am ultimately attempting to color-code each polygon depending on the number of points that lie within its boundary lines.

Thanks,


texts:
I have a DTM file (created from OS terrain data) which I would like to split into 4 smaller files, based on the features in a polygon shapefile, in order to perform further processing in ArcScene and ArcMap.

I am using Data Management Tools --> Raster --> Raster Processing --> Clip to clip the .tif file.

However, the output of the smaller DTM does not exactly match up with the input polygon feature, as shown here:

I have tried other methods (i.e. split) and also tried in QGIS but I still get the same problem.

Can anyone suggest why this is happening and help me navigate around it?


texts: The coordinates measure are out of bounds in Android route task

texts:
I am a beginner in ArcGIS and I was following this tutorial. I am trying to find a route between two points. I created a sample application, in which first, i get the two points on map single tap as follows

@Override
        public boolean onSingleTap(final MotionEvent e)
        {
            Point point = map.toMapPoint(new Point(e.getX(), e.getY()));
            points.add(point);
            if (points.size() == 2)
            {
                new FindRouteTask().execute();
            }
            return true;
        }

and my FindRouteTask is as follows

private class FindRouteTask extends AsyncTask<Void, Void, RoutingResult>
    {

        ProgressDialog  dialog  = new ProgressDialog(FindeRouteActivity.this);
        RoutingResult   mResults;

        @Override
        protected void onPreExecute()
        {
            dialog.setMessage("Finding routes...");
            dialog.show();
        }

        @Override
        protected RoutingResult doInBackground(Void... params)
        {
            RoutingParameters routeParams = new RoutingParameters();
            routeParams.setImpedanceAttributeName("Length");

            String routeTaskURL = "http://tasks.arcgisonline.com/ArcGIS/rest/services/NetworkAnalysis/ESRI_Route_NA/NAServer/Route";
            RoutingTask routeTask = new RoutingTask(routeTaskURL, null);

            NAFeaturesAsFeature naFeatures = new NAFeaturesAsFeature();
            // Create the stop points from point geometry
            StopGraphic startPnt = new StopGraphic(points.get(0));
            StopGraphic endPnt = new StopGraphic(points.get(1));
            // set features on routing feature class
            naFeatures.setFeatures(new Graphic[] { startPnt, endPnt });
            // set stops on routing feature class
            routeParams.setStops(naFeatures);
            try
            {
                mResults = routeTask.solve(routeParams);
            } catch (Exception e)
            {
                e.printStackTrace();
                return null;
            }
            return mResults;
        }

        @Override
        protected void onPostExecute(RoutingResult result)
        {
            dialog.dismiss();
            points.clear();
            if (result == null || !result.getRoutes().isEmpty())
            {
                Toast.makeText(FindeRouteActivity.this, "There was some problem in fetching the routes. Please try again!", Toast.LENGTH_SHORT).show();
            } else
            {
                findRoute(result.getRoutes().get(0));
            }
        }

    }

But the problem is that routeTask.solve(routeParams) always throws a exception

Failed to convert the input propertyset into a recordset to load into the NAClass "Stops".  The coordinates or measures are out of bounds.

Could anyone please tell me what am I doing wrong and what the problem is??
-Thanks in advance


texts: QGIS Atlas pages with different raster base layers?

texts:
Is it possible to enable raster layers only for specific Atlas pages, without using multiple map Windows?

Unfortunately, there is no way to set transparency of raster layers by an expression function.


texts: Identifying flow of a Wastewater Network

texts:
I am trying to identify the flow path of a wastewater system.

I am using Arcgis 10.3 and I have set up a geometric network system and identified the sources and the sink (sources = blue dots, sink = green dot).  The network is definitely connected as the below picture showing the result of a "find Connected" Analysis shows.

I have used the Set Flow Direction on the Utility Network Analyst toolbar but it is showing the flow of every pipe as indeterminate flow.  I have also tried the ArcToolbox geoprocessing tool called Set Flow Direction but this hasn't changed the results.  Here is what I would expect the flow to look like (sorry about my poor arrow drawing abilities).

Is there a step that I have missed or should I be creating this network a different way?

I understand that the tool with have issues determining flow within loops.  Is there some way that I could integrate a contour layer to assist in flow direction?

This is just a small portion of several networks I have, hence I need an automated solution.


texts: What does the srid code mean in an sde.st_geometry?

texts:
In case of an sql insert query like :

INSERT INTO OUTBREAKS (OID, address, city_code, loc)
   VALUES (1,
   '1420 kirchestrasse',
   43,
   sde.st_geometry (0.00003, 0.00051, null, null, 4326)
   );

Does the srid (4326) mean that the decimal coordinates used for x and y are of the spacial reference specified from the srid?
Or  it says that the geometry will keep its coordinates in that spatial reference?
In other words can I store a geometry with decimal coordinates (in WGS_84) for latitude/longitude ( like 50.91212 , -0.1221) and use an srid of a projected system which uses meters? Will they be automatically convented from wgs_84 to my new spatial reference and be kept in that way in the database ?


texts:
I'm trying to find the erase geoprocessing method as I am familiar with in ARC, but I can't seem to find it in OGR using Python. I have looked through the module OGR documentation but perhaps I missed it.

An example of the function would be even better.


texts: Converting protected areas polygons into raster

texts:
I am working with protected areas for vulnerability analysis, I have values of 1 for protected areas being vulnerable, and 2 to non-protected areas being less vulnerable. I have transformed this into a raster file.
My question is: for standardizing this data should I run the following equation: ( x - min(x) ) / ( max(x) - min(x) ) in order to aggregate this layer with other data?


texts:
I was georeferencing a map. (it’s now a working map in the GeoTIFF format) Then I used
gdal2tiles -z5-11 myFile.tif myFolder

but the resulting tiles were incompatible with Mapbox and compared to working ones it seems as if it would be using some different projection.
I used gdalwarp to convert the GeoTIFF to the correct ESPG:3857 format but the tiles are in the same format. I also tried to add the gdal2tiles -p mercator option but it didn’t change it the expected way.
How can I generate ESPG:3857 tiles?


texts: "Failed to load system tools" error from ArcGIS Pro Geoprocessing

texts:
I updated my ArcGIS Pro to the latest version, and after opening a project and clicking on the Geoprocessing tab I get an error that starts

Failed to load system tools

which is then followed by a long traceback along the lines of

\Toolboxes\Intelligence Tools.pyt Traceback (most  recent call last):
File "<string>", line 28, in <module> File "c:\program
files\arcgis\pro \ResourcesArcToolbox\Scripts
\GenerateBlindSpotAreasToolClass.py", line 24, in...

and so on (the error text was long and not copy-able so I've just typed in the first few lines from my screenshot).  After I click OK I'm presented with a blank Geoprocessing Pane where there would normally be a list of the Toolboxes available in ArcGIS Pro

I found a Technical Support page for the error on the Esri website - Error: Failed to load system tools - which suggests the following differing solutions or workarounds:

Modify the proenv.txt file to match the environment specified in
ArcGIS Pro Clear cache
Reinstall ArcGIS Pro

Each step had instructions which I followed, but even after trying all three steps (including resinstalling ArcGIS Pro 4 times!) the error persisted and I was still getting the error.
What can be done to get my ArcGIS Pro Geoprocessing Pane to work correctly again?


texts: QGIS Load WFS everytime different Output

texts:
I want to load a WFS in QGIS and execute some queries. I tried different versions of QGIS.

QGIS 2.18 number of objects: 160036 (first query)
QGIS 2.18 number of objects: 80018 (second query)
QGIS 3.2.0 number of objects: 81103
QGIS 3.4.5 number of objects: 80018

To check the correct number of objects I load the WFS in FME and get 80018.

I tried other WFS Clients and there I get in every version the correct number of objects.

Unfortunately I can´t give you the WFS because the access is restricted.

Does anyone have an idea why each version spends a different number of objects and why reloading in the same version will dump different objects?


texts: How to set the layer attribute for IDW interpolation tool within a QGIS Python script

texts:
I have a custom Python script in QGIS 3.10 that integrates some QGIS tools. One of them is the IDW interpolation tool. My issue is I don't know how to properly set the layer attribute to be used for the interpolation.

# file connections
infile = os.path.join(scratch_folder, "lakes.shp")
lakes_meta = QFileInfo(infile)
lakes_title = lakes_meta.baseName()
lakes_layer = QgsVectorLayer(infile, lakes_title, 'ogr')

# run IDW interpolation tool
processing.run('qgis:idwinterpolation',
                {'INTERPOLATION_DATA': lakes_layer,   ### < --- how to specify layer attribute???
                 'DISTANCE_COEFFICIENT': 2,
                 'EXTENT': lakes_layer,
                 'PIXEL_SIZE':1,
                 'OUTPUT': os.path.join(scratch_folder, 'IDW_output.tif')})

The QGIS documentation alludes to the use of either '::|::' or '::~::' to assign the correct column for interpolation, but I've never used anything formatted like that before and am not sure how to proceed.


texts: GeoPandas: Calculate coordinates of a point within a linestring

texts:
I have two datasets, one is a GeoJSON of a highway with around 1509 linestrings.
The other is a list of incidents that have occurred along said highway, but, instead of coordinates, the file specifies the kilometer where the incident happended.
I.e

incident_id
km
incident_type

1
10
vehicle theft

2
152
dui

What I need to do is get the coordinates of each incident.
I know that you can calculate the distance between two coordinates using geopy.distance, but is there any way to calculate the coordinates of a point which is x kilometers away from a given set of coordinates?


texts:
I've been trying to deploy GeoServer in a WAMP/Apache environment via a .war file. (GeoServer being downloaded as a web archive)

When I try to enter the address in the web browser, I get this:

Has somebody already encountered the same issue ? Do you know the reason for such a failure ?


texts: Getting error while importing GeoPandas in Python console

texts:
I have installed GeoPandas and its dependencies using pip install package_name. Here is the pip list output:

Package         Version

attrs           18.2.0
Click           7.0
click-plugins   1.0.4
cligj           0.5.0
cycler          0.10.0
descartes       1.1.0
Fiona           1.8.4
geographiclib   1.49
geopandas       0.4.0
geopy           1.18.1
kiwisolver      1.0.1
matplotlib      2.2.3
munch           2.3.2
numpy           1.15.4
pandas          0.22.0
pip             18.1
psycopg2        2.7.6.1
pyparsing       2.3.0
pyproj          1.9.5.1
PySAL           1.14.4.post2
python-dateutil 2.7.5
pytz            2018.7
Rtree           0.8.3
scipy           1.2.0
setuptools      40.6.3
Shapely         1.6.4.post2
six             1.12.0
wheel           0.32.3

But as soon as i am importing it in Python console as

  import geopandas

I am getting an error as

  Traceback (most recent call last):   File "", line 1, in
     File
  "/home/rajesh/panda_projects/env/lib/python3.4/site-packages/geopandas/init.py",
  line 1, in
      from geopandas.geoseries import GeoSeries   File "/home/rajesh/panda_projects/env/lib/python3.4/site-packages/geopandas/geoseries.py",
  line 5, in
      from pandas import Series   File "/home/rajesh/panda_projects/env/lib/python3.4/site-packages/pandas/init.py",
  line 26, in
      from pandas._libs import (hashtable as _hashtable,   File "/home/rajesh/panda_projects/env/lib/python3.4/site-packages/pandas/_libs/init.py",
  line 4, in
      from .tslib import iNaT, NaT, Timestamp, Timedelta, OutOfBoundsDatetime   File "pandas/_libs/src/numpy.pxd", line 865, in
  init pandas._libs.tslib ValueError: numpy.ufunc has the wrong size,
  try recompiling. Expected 192, got 216

What is the solution for this?


texts: Printing attributes of selection in ArcPy?

texts:
I'm trying to use ArcPy to print a list of addresses of land parcels crossed by a power line. I'm using the following code to accomplish this:

# Loop through features in parcels.shp
for feature in parcel:
    # Search cursor for address field
    cursor = arcpy.da.SearchCursor(parcel, 'SITUSADDR')

    # Loop through cursor
    for row in cursor:
        # Select features crossed by powerline.shp
        matchcount=arcpy.SelectLayerByLocation_management(parcel, 'CROSSED_BY_THE_OUTLINE_OF', power)

        # Loop through matchcount, print
        for i in matchcount:
            print((row[0]))

    # Give feedback if there are no matches
    if matchcount == 0:
        print('No features matched spatial criteria.')

The code is able to print addresses, but isn't filtering them the way I'm expecting it too. Rather than printing several addresses from 'matchcount' it prints the addresses of all parcels and seems to repeat them three times.

How can I edit this code to print only the selected features?


texts: How to polygonise a binary raster with python?

texts:
I already asked this in StackOverflow, but maybe it's better here. I have a ndvi raster that I binarized (I put my threshold at 0.35) with a script python and I wanted to polygonize this raster, in order to get ndvi into polygones, to store them after in a database. I found and used this code :

import rasterio
import rasterio.env
from rasterio.features import shapes
import fiona
mask = None
with rasterio.Env():
    with rasterio.open('binarized_ndvi.tif') as src:
        image = src.read()

    results = (
    {'properties': {'raster_val': v}, 'geometry': s}
    for i, (s, v)
    in enumerate(
    shapes(image, mask=mask, transform=src.transform)))

    with fiona.open(
    'ndvi_shapefile.shp', 'w',
    driver='Shapefile',
    crs=src.crs,
    schema={'properties': [('raster_val', 'int')],
    'geometry': 'Polygon'}) as dst:
        dst.writerecords(results)

But all I get is a shapefile, with only one polygon whose shape is the same as my raster, like if only the first elment of the raster was taken into account ...

I can't send you my raster, it's too big.


texts:
I have a Dijkstra’s algorithm. When I run this portion of my script:

graph = {}

graph["begin"] = {}

graph["begin"]["a"] = 6

graph["begin"]["b"] = 2

graph["begin"]["c"] = 3

graph["begin"]["d"] = 1

graph["a"] = {}

graph["a"]["c"] = 3

graph["b"] = {}

graph["b"]["a"] = 3

graph["b"]["d"] = 5

graph["d"] = {}

graph["d"]["fin"] = 1

graph["c"] = {}

graph["c"]["d"] = 1

graph["c"]["fin"] = 3

graph["fin"] = {}

infinity = float("inf")

costs = {}

costs["a"] = 6

costs["b"] = 2

costs["c"] = infinity

costs["d"] = infinity

costs["fin"] = infinity

parents = {}

parents["a"] = ["begin"]

parents["b"] = ["begin"]

parents["c"] = ["begin"]

parents["d"] = ["begin"]

parents["fin"] = None

processed = []

def find_lowest_cost_node(costs):

  lowest_cost = float("inf")

  lowest_cost_node = None

  for node in costs:

    cost = costs[node]

    if cost < lowest_cost and node not in processed:

      lowest_cost = cost

      lowest_cost_node = node

  return lowest_cost_node

def run_dijkstra():

  node = find_lowest_cost_node(costs)

  while node is not None:

    neighbors = graph[node]

    for neighbor in neighbors.keys():

      newcost = neighbors[neighbor] + costs[node]

      if newcost < costs[neighbor]:

        costs[neighbor] = newcost

        parents[neighbor] = node

    processed.append(node)

    node = find_lowest_cost_node(costs)

  return

  def print_path(begin, fin):
print ["cost to destination", str(costs["fin"])]

     print ["fin"]

  nextNode = parents[end]

  while nextNode != begin:

    print ["nextNode"]

    nextNode = parents[nextNode]

  print ["begin"]

run_dijkstra()

print_path("begin, end")

I get this error message:

  Traceback (most recent call last):   File "python", line 70, in
   NameError: name 'print_path' is not defined

Does anyone know a decision of this problem?)


texts: Leaflet only allows polygon layer from GeoJSON to be displayed one time

texts:
I have used an AJAX call to load a Black Sea GeoJSON that I want to display as a polygon on the map.
var blackSea = [];
var blackSeaData = new L.GeoJSON.AJAX("./blacksea.json");
blackSea.push(blackSeaData);

blackSea[0].addTo(map);

The first time I use blackSea[0].addTo(map) it works, the polygon is displayed, and I can even change its style, but if I clear the map with clearMap() and try it again it does not display. If I clearMap() before the initial call it still displays, but after I clearMap() again it does not.
What could be the problem here?
edit: the clearMap() function:
    function clearMap() {
    for(i in map._layers) {
        if(map._layers[i]._path != undefined) {
            try {
                map.removeLayer(map._layers[i]);
            }
            catch(e) {
                console.log("problem with " + e + map._layers[i]);
            }
        }
    }
}


texts: From Structure From Motion (.ply) to GIS (i.e. grass gis)

texts:
i'm not an expert about 3D reconstruction but i wanted to ask, what is the best method to import a 3D model (gained by using Structure from Motion tecniques).
I used VisualSFM to generate the 3d reconstruction, and I got a .ply file as output.
I now want to import it and analyse its elevation model, by using GRASS Gis, but I get this error:

(Fri Aug 19 11:38:27 2016)
g.extension.py extension=v.in.ply svnurl=http://svn.osgeo.org/grass/grass-addons/grass6
Fetching <v.in.ply> from GRASS-Addons SVN (be patient)...
Compiling...
Makefile:5: /usr/lib/grass64/include/Make/Script.make: No
such file or directory
make: *** No rule to make target
`/usr/lib/grass64/include/Make/Script.make'.  Stop.
ERROR: Compilation failed, sorry. Please check above error messages.
(Fri Aug 19 11:38:31 2016) Command finished (4 sec)


texts: QGIS 2.14 Processing plugin crashes - TypeError: must be string, not None

texts:
I've recently installed 2.14, and it worked fine for the first couple of days (Appart from Grass tools), including several program and computer reboots.

Today I had this error show up

Couldn't load plugin processing due to an error when calling its initGui() method

TypeError: must be string, not None

The full Pyhton error report:

The Processing plugin version is 2.12.99. I've tried reinstalling and downgrading the plugin, that didn't work.


texts: How to measure areas using an official CRS which is not available in QGIS?

texts:
I am very new to QGIS and doing a project on Glacial Lake Outburst & Flooding (GLOF).

I am trying to find the area of a lake using images from Google Earth 7.0.3.8542. I am using the WGS84 CRS from my Google Earth image for the area calculation.

I am getting some messy area values like 0.000103145e-5. I can't even guess in what unit it's being displayed.

I have an official Local Coordinate System for our country in .prg format but I can't find it in the list of available CRSs in QGIS 2.0 or how to load it as a custom CRS.

How can I load the local CRS in my computer into QGIS?.

The contents of the .prg file are:

  PROJCS["DRUKREF 03 local grid
  GASA",GEOGCS["GCS_DRUKREF03",DATUM["D_DRUKREF03",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",250000.0],PARAMETER["False_Northing",-2500000.0],PARAMETER["Central_Meridian",90.0333333333333],PARAMETER["Scale_Factor",1.0],PARAMETER["Latitude_Of_Origin",0.0],UNIT["Meter",1.0]]


texts: Disabling service recycling in ArcGIS Server

texts: Is it possible to disable recycling for a service? I am using Server 10.6.1. I am thinking there is some hack to this because its clearly not an option using the Manager interface.


texts: Which polygons intersect rectangle in PHP/Linux GDAL?

texts:
I'm developing a website where I want to test which polygons from a shapefile intersect a given rectangle. Using PHP, I thought about using exec() function, calling a GDAL program, like ogrinfo. However, the following command line

ogrinfo -spat x0 y0 x1 y1 my_shape_file.shp

is not working how I thought it should. It always returns the name of the layer (after "Open of ... using driver 'ESRI Shapefile' successful."), despite the polygon is inside the given rectangle (x0,y0,x1,y1) or not. Anyway, I want the polygons intersecting the rectangle, not the whole shapefile. What is wrong in my usage? The instructions say nothing! (Of course I used actual coordinates instead of x0,y0...)

Other question: even when I get this to work, maybe using one of PHP GDAL extensions (version 5 or version 6) is faster, is it? And is it safer?

EDIT: I wasn't using the layer name (as returned from the previous command) after the shapefile name:

ogrinfo -spat x0 y0 x1 y1 my_shape_file.shp layer_name

Now I only need to know if PHP extensions are preferred over external commands.


texts: Unexpected result converting a vector to a raster in GRASS 7.0.0

texts:
My purpose is to create a DEM raster in GRASS 7.0.0 starting from a shapefile containing contour lines.

The attributes table has a field named QUOTE_INT containing the elevation as integer.

I started importing the shapefile with the command v.in.ogr.

Here I'm stuck because when I run the command v.to.rast I obtain a raster filled monochromatically instead of a raster containing the contour lines coloured based on their elevation.

What am I doing wrong?


texts:
I am using File Geodatabase API into my own MFC C++ application (which only reads data from .gdb folders). I have include all the headers and libraries required for the api in the project Properties section. (except FileGDBAPI.dll and FileGDBAPID.dll files because i don't know the the including directory for them). there is no compile time error but when i create any object from FileGDBAPI (e.g GeoDatabase geodatabase;) then i am getting an error "the appliaction was unable to start correctly (0xc000007b). Click OK to close the application"

Can someone please tell me what i am ignoring or anything that can make this error go away.
How can i add *.dll files related to File Geodatabase API?


texts: How can QGIS Field Calculator sum values from a spatially joined table?

texts:
Can the Field Calculater be used to sum values from a related table via a Spatial Join?

I have two vector polygon layers A and B, each with a numerical column VALUES.  I want to update VALUES in layer A with the sum of VALUES from layer B where geometry from B is contained within geometry A.

So far I have the expression geomcontains('Layer B','VALUES') and what I want is to Sum this.

I have seen other ways to do this via plugins like SpatialJoin etc, but they all involve modifying either saving new layers or modifying the structure of Layer A, which I think is fairly cumbersome - is there a simple solution? (In MapInfo this is extremely simple, using Update Column.)


texts: Auth Error Cascading MapBox WMTS with GeoServer

texts:
I'm trying to cascade MapBox WMTS with Geoserver using the following GetCapabilities url:

https://api.mapbox.com/styles/v1/mapbox/dark-v10/wmts?access_token=pk.

But Geoserver is getting a 401 {"message":"Not Authorized - No Token"} response from MapBox:

WMTS Connection test failed: Server returned HTTP response code: 401 for URL: https://api.mapbox.com/styles/v1/mapbox/dark-v10/wmts?REQUEST=GetCapabilities&VERSION=1.0.0&ACCESS_TOKEN=pk.

It seems that access_token is case sensitive as I get a 200 response when using curl:

curl "https://api.mapbox.com/styles/v1/mapbox/dark-v10/wmts?REQUEST=GetCapabilities&VERSION=1.0.0&access_token=pk"

StatusCode: 200

Is there a workaround to prevent Geoserver from capitalizing the request parameters, or perhaps manually creating the store?

I see the 'Key Authentication' module does something similar for incoming requests.