texts: ArcGIS map service not resizing for print-friendly page

texts:
( I asked this same question on ArcGIS forums, with no response...I hope someone here can help!)

I am using ArcView, ArcServer 10.0 and the JS API v3.1 (so no print function) ... I have an app where the user chooses parameters from a few dropdown menus, and a print-friendly map page opens in a new tab. From this page, the user can click a 'print' button, which uses window.print() to open the browser's print dialog. I have a print.css page, which is working well with my 'normal' divs, text and images. I have this tag in my HTML:

<link rel="stylesheet" href="print.css" type="text/css" media="print" />

My problem is that I cannot seem to reset the size of the Dojo map object. The map in my printed page (or PDF, if I choose CutePDF) is the same size as whatever the size was in the print-friendly page. I would like my print.css to hard define the printed version / PDF as an exact height and width, regardless of the window size in the print-friendly page. I cannot define the size of page to be printed in the print-friendly page, b/c for one, my client would like to pan and zoom the map after the print-friendly page is generated. I have no problem setting the size of the other elements in my print.css rules (text box, hard coded legend with images)...but any rules I try to set for the map itself in my print.css page are ignored.

I tried using 'exportMapImage', thinking that maybe I could make an image of the map and its contents, but this loses the basemap which I need, and also I do not want the user to open yet another page or tab. (I am having the same problem with my scalebar, but I this isn't so important for this project).

My goal is to create a PDF or printed page that is the same size no matter what the screen size or current window size., and I have it working except for the map itself.

Any ideas? I hope I explained my problem clearly. Thanks , J

Also,  I am hoping to find a client-side solution; I don't have much experience with server-side programming.  If this can't happen, please let me know any server-side ideas, thanks!)


texts: Create a report with filter from different layer using Atlas

texts:
I have the following scenario:

Layer ‘grid’, with a field called ‘number’
Layer ‘log’ with field called ‘number’ in common with ‘grid’ layer and another field called ‘date’
Layer ‘problem’ with field called ‘number’ in common with    ‘log’ and ‘grid’ number and another field called ‘type’

Now I’d like to make layout with Atlas enabled; in every page I need to print some fields from layer ‘log’ filtered by the common field ‘number’ from ‘grid’ layer.

On the same page I need also some field from layer ‘problem’ filtered by the common field ‘number’ from layer ‘grid’.

Unfortunately I cannot use $atlasfeatureid as filter since  I deleted some grid on the ‘grid’ layer and there’s no more correspondence.


texts: How do i add online GeoServer layers on android mobile app

texts:
I am trying to develop an android app, and in this app I want to display a Google map and on this map I also want to add some layers on the Map which I have created from GeoServer.

Basically my requirement is that I want a Google map with some layers of GeoServer in my android app. which layers are online on my server not on my localhost

I am new for the android app as well as GeoServer.


texts:
I've used the following code. Output can be seen in map.
How do I add a title to my legend?
my.palette <- brewer.pal(n = 7, name = "RdGn")
spplot(AGCmap1, "AGCmap_1", main = "Carbon density per ha", zlim=c(0,80),
       col.regions = my.palette, cuts = 6, col = "transparent")```
[![Carbon Map][1]][1]

  [1]: https://i.stack.imgur.com/0Ubn0.png


texts: Is tileJSON supported by ArcGIS for Server?

texts:
Does ArcGIS 10.2 for Server support tileJSON? I mean, can I retrieve in some way the tileJSON from and ArcgisTiledLayer ?

I need this so as to use it with wax library and have a grid along with the tiled image which will contain the encoded location of the features that the image shows.


texts:
About a year ago I installed the newest version of QGIS, but it kept freezing after  a few minutes of use, so I uninstalled it and installed 2.18 because it was supposed to be more stable. It seemed to work at the time. Not until the last few weeks have I really started to use QGIS intensively, and now it has started freezing again. I haven't been doing anything complicated, just creating polygons and editing their attributes. It doesn't seem to matter exactly what i'm trying to do, after a certain time or number of clicks the whole thing just freezes and won't function no matter how long I wait.

I've seen similar questions, but without any definitive answers. I'm running Windows 10 on an MSI 72. The computer should definitely be fast enough, I never had a problem with ArcGIS or AutoCAD. I haven't installed any plugins, and I don't see how it could be a problem with the data as I'm creating it in QGIS and it doesn't matter if I'm editing the data or just trying to zoom in/out.


texts: I'm very new to QGIS and am having some issues using the SAGA IDW tool. I'm attempting to perform an interpolation using the 10 closest sample points and an inverse distance power of 2. I would like to specify a cell size of 5000 meters; however, the cell size appears to be much larger than 5000 meters on the output grid. What am I doing incorrectly?


texts: Which versions of python3 can we use with current QGIS3 version

texts:
Mac OS 10.04, QGIS 3.10.

I have been a way from QGIS for a while and I am now setting up a "new" Mac and having my usual battle with python.  Great to see that we now have signed package for QGIS but are we still constrained to using python 3.6.5?

Mac Os now installs 3.7 by default and I would rather not have to learn about virtual environments if I don't need to!  (But I fear I will ;)


texts: Conditionally omitting field values while calculating new field in ArcGIS Desktop?

texts:
I have six fields in a table and I want to calculate a new field based on a formula using these six fields. Suppose the formula is simply:

newField = field1 + field2 + field3 + field4 + field5 + field6

how can I omit a field from this formula if the field value is 0? For example if for a row field1 = 0 then I want the new field to be calculated as:

newField = field2 + field3 + field4 + field5 + field6

Any field value can be 0 in any row. That field needs to be ommitted during calculation for that row. Multiple fields can also be 0, then all those fields are omitted. The fields cannot contain null values and are numeric. The formula needs to change conditionally based on the presence of 0 values.
Is there anyway to do this with python or normally?

PS: Okay it looks like nothing needs to be changed in the above formula since adding 0 doesnt change anything. But its just an example. The actual formula is complex and if 0 value is used, the result will become false.


texts: Show point, based on coordinates, on a map in QGIS

texts:
I would like to enter coordinates (e.g. UTM longitude 323519 & UTM latitude 3777727) by hand into QGIS and have it identify that point (either by centering the map on that point or, better yet, by putting a marker at the corresponding location on the map).

I know how to generate a CSV file with coordinates but don't want to go through the process of generating a CSV file each time I want to find a location on a map.

Seems like a simple question but I have been struggling to find an answer.


texts: How to quickly remove all WMS-Services in QGIS

texts:
I uploaded a lot of WMS-Services in QGIS via a .XML file.

When I modify the .XML file, for instance delete a couple of Services or add new ones, my changes are not changed in QGIS.

As an alternative I would like to just delete all of them and upload my new .XML file but I don't find the option to delete all at once. I don't just want to delete one by one.

Where are my services saved and how can I delete all of them at once?

I'm using QGIS 3.4


texts: Total of several values in field for points within polygons

texts:
I have ~4000 insurance claims for a city. I have made a fishnet and would like to count the number of claims within each polygon. It is easy to count the total number of claims, but I would like to find out how many claims there are within each cell for every claim date.

This is how the data set looks like:

    FID Shape       Name    Date
83  Point ZM    A2528   20030718
84  Point ZM    A2529   20100814
85  Point ZM    A2530   20140831
86  Point ZM    A2531   20140831
87  Point ZM    A2532   20140831
88  Point ZM    A2533   20110607
89  Point ZM    A2534   20120116
90  Point ZM    A2535   20100915
91  Point ZM    A2536   20100813
92  Point ZM    A2537   20140831

As you can see, the flood claims are from different dates (~50). Does anyone of you have an idea how to solve this problem? I have read some similar Q&A's (Compute total area of several polygons within a polygon; Total of attribute field of point layer per polygon layer; How to write FIDs of points within each polygon to polygon field?), but could not find any answer that works in this case. The resulting table will consist of ~4000 rows (one for each cell) and ~50 columns (one for each claim date). The value in each field will tell how many flood claims there are in each cell for a certain date.

I work in ArcMap 10.1.


texts:
after using zonal statistics by hand which works perfectly fine.
I want to use it in Batch-mode to extract Mean/SD/Median for different rasters in the same polygon. I can not select a rasterband, there is no option.

Unfortunately the output is empty, see below.

If anyone can help me would be great.


texts:
I am trying to create a feature in a polygon layer in a file geodatabase but I get this error:

Create: An unexpected failure occurred.
The spatial index grid size is invalid.[test]
The spatial index grid size is invalid

I can't figure out what my mistake?

update:

I can't delete this field in Arcatalog -- I  get this error:

The ObjectID index cannot be deleted, it is a required index


texts: Multi-band raster color model in QGIS. Why RGB (668, 983, 1007) in Uint16 appears green?

texts:
I created an RGB Uint16 composite raster in QGIS .
These are band designations:

Band 1 - blue
Band 2 - green
Band 3 - red

When i query value of a pixel that appears green - it shows (B3: 668, B2: 983, B1: 1007). Since the pixel is green, I expected DN in the green band to be the highest, although this is not the case.

Why is that?

The raster is stretched and clipped using 2 - 98% cumulative count, although the values queried are raw DNs.


texts: Who can teach me how to reach the StackProfile method in ArcEngine10.2?

texts:
I try my best to write code to get it.

The codes list below, but I think I don't understand these parameters so that I can't run the gp Tool successfully. So I want to consult someone who can tell me how to write correctly.

 private void AddSection(IFeatureLayer linFeatureLayer, IRasterLayer pRasterLayer,string name)
    {
        Geoprocessor gp = new Geoprocessor();
        gp.OverwriteOutput = true;
        ESRI.ArcGIS.Analyst3DTools.StackProfile sfsf = new ESRI.ArcGIS.Analyst3DTools.StackProfile();

        sfsf.in_line_features = linFeatureLayer;
        sfsf.profile_targets = pRasterLayer;
        sfsf.out_table = @"C:\Users\Haker\Desktop\section\StackProfile5.dbf";
        sfsf.out_graph = @"C:\Users\Haker\Desktop\section\newgraph";
         try
        {
            gp.Execute(sfsf,null);
        }
        catch (COMException e)
        {
            string message = "";
            for (int i = 0; i < gp.MessageCount; i++)
            {
                message += gp.GetMessage(i) + "\r\n";
            }
            MessageBox.Show(message + e.ToString());
        }}


texts: How to create a 3D Map with multiple layers?

texts:
This is a Generic question, if it is not applicable please delete.

I was wondering if it was possible to make a 3D map (such as seen below) on QGIS version 3.10.0

The image is basically showing sediment samples and the estimated thickness of oily sediment at each sample.


texts: Entering data down column using tab or return key in QGIS for Mac?

texts:
I have many shapefiles that require I manually populate data fields in an attribute table.  I am using QGIS in a Mac system.

For speed and efficiency, I prefer to use the tab or return carriage to move from cell to cell within a column.

Currently the tab key takes me from cell to cell moving across a row and it seems the return key gives me the option to delete the cell.

Is there a setting I can change somewhere to force the tab function to move from cell to cell within a column?


texts: Atmospheric Corrections for Landsat 8 OLI Using Orfeo Toolbox

texts: I am trying to correct some Landsat 8 OLI imagery using the Orfeo Toolbox Optical Calibration tool via QGIS.  I've retrieved data from Aeronet as well as relative spectral response data for the OLI sensor from here.  However, I can't seem to get the tool to run.  I'm having difficulty finding any documentation beyond that I need the RSR file and the Aeronet file.  What I'm unsure of is A: what format the relative spectral response data needs to be in, and B: what file counts as the Aeronet file?  Aeronet gives me several datasets to choose from here.  Alternatively, if anyone has any suggestions for alternative methods or tools, I'll gladly hear them.


texts: North arrow not syncing with map QGIS 2.18

texts:
I read about clicking on the 'sync with map' option after adding a north arrow (Layout -add image)(in Print composer), and I've done that, but when I rotate the map, the north arrow does not rotate with it. I've tried refreshing.
What have I missed?

edit _ this is still an issue.

edit - solved it.
There are TWO places where you have to tell it which map the N arrow should rotate with.
1 is in the layout tab for the image, under general settings, make sure reference map has the right map selected. THEN go to item properties, expand image rotation and select the sync with map option and tell it which map.


texts: Why does ArcPy Excel To Table give AttributeError: 'module' object has no attribute?

texts:
I'm trying to import an excel workbook using arcpy. ExcelToTable, but I get the following traceback:

ERROR

updateParameters Execution Error: Runtime error Traceback (most recent call last):
File "c:\program files (x86)\arcgis\desktop10.2\ArcToolbox\Toolboxes\Conversion
Tools.tbx#ExcelToTable_conversion.UpdateParameters.py", line 4, in File
"c:\program files (x86)\arcgis\desktop10.2\ArcToolbox\Toolboxes\Conversion Tools.tbx
#ExcelToTable_conversion.UpdateParameters.py",
line 40, in updateParameters AttributeError: 'module' object has no attribute
'get_sheet_names'

I have imported arcpy, have I missed something else?

EDIT: Here's my code:

from win32com.client import Dispatch
import arcpy, os
from arcpy import env

arcpy.env.workspace = "C:\\Users\\Ant\\Documents\\ArcGIS\\ITN.gdb"
env.overwriteOutput = True

x = 1
with arcpy.da.UpdateCursor("Sites",["SizeScore","SizeWt","UPRNScore","UPRNWt"]) as cursor:
for row in cursor:
    SQL = "FID = {0}".format(x)
    arcpy.MakeFeatureLayer_management("Sites","Siteslyr",SQL)
    arcpy.TableToExcel_conversion("Siteslyr","C:\\Sites.xls")
    def RunExcelMacro(name):
        myExcel = Dispatch('Excel.Application')
        myExcel.Visible = 0
        myExcel.Workbooks.Add('c:\\Data\\UpdateScores.xlsm')
        myExcel.Run(name)
        myExcel.DisplayAlerts = 0
        myExcel.Quit()
    if __name__ == "__main__":
        RunExcelMacro('MacroEnableDBF')
    def RunExcelMacro(name):
        myExcel = Dispatch('Excel.Application')
        myExcel.Visible = 0
        myExcel.Workbooks.Add('c:\\Data\\UpdateScores.xlsm')
        myExcel.Run(name)
        myExcel.DisplayAlerts = 0
        myExcel.Quit()
    if __name__ == "__main__":
        RunExcelMacro('AppendScores2')
    arcpy.ExcelToTable_conversion("C:\\Sitescores.xls", "Carescores")
    with arcpy.da.SearchCursor("Carescores", fields) as curs:
        for rows in curs:
            row[0] = rows[0]
            row[1] = rows[1]
            row[2] = rows[2]
            row[3] = rows[3]
            cursor.updateRow(row)
    x = x + 1


texts: Calculate a field in a GetFeatureInfo request with MapServer 7.0.4

texts:
I would like to be able to do some simple calculations on fields returned by a GetFeatureInfo with MapServer 7.0.4.

My case is that I am querying a layer (OGR connection to an ESRI File Geodatabase) and one of the fields returned is the area of a polygon. The area is in square meters but I would like to do an on-the-fly calculation to show the value in acres or hectare.

I didn't find anything in the official documentation indicating that this is possible.

Is it possible, without having to modify the source layer?

Here is the template I use to output the result of the GetFeatureInfo:

<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN"
       "http://www.w3.org/TR/html4/transitional.dtd">

    <title>MapServer Template Sample</title>

                    Attribut
                    Valeur

                    Numéro de lot
                    [no_lot_renove]

                    Prénom
                    [prenom_proprio]

                    Nom
                    [nom_proprio]

                    Superficie
                    [shape_area] m2


texts:
I've styled some census tract tilesets and set filters in Mapbox Studio, where everything looks fantastic, but I want the tracts to appear shaded irrespective of the zoom level, and Studio won't let me set a minzoom lower than 8.  I've tried changing the minzoom after load in javascript:

map.setLayerZoomRange ('ba_percent', 4, 22)

And I see the change reflected if I log

map.getLayer('ba_percent')

But that doesn't affect what appears in the browser.  There must be a way to show census tract data at whatever zoom level I choose, right?  Just for reference, this is how the map looks when tracts are shaded.


texts: Feature Service related table not showing up in Survey 123 form

texts:
Using Arcmap 10.3.1 I have published a featured service using My Hosted Services with a related table. I have confirmed the relation was working before publishing the service and can see the related table attributes using the web viewer of the feature service. I also confirmed I chose feature access instead of tiled mapping when publishing from service editor.

Now, I am trying to use Survey123 connect to build a survey off of this feature service (with related table). The survey creates questions from the standard database just fine, but nothing from the related table. I am expecting a "begin repeat" to show under the "type" xlsx form where the related table questions begin.

Is there a way to build the survey off of the feature service WITH the related table included in the Survey123 form?


texts: Creating column in attribute tables which sequentially numbers features based upon a direction using QGIS?

texts:
I am using QGIS 3.6.

I'm trying to create a map which labels each parcel according to it's position along the red line, from West to East. I've tried creating a column which uses each features position within the attribute table to do this, however as you can see in the picture, they are not in any sensible order. I normally do this manually, however I have over 3000+ features to edit in that case, it would be good to find an automated solution.


texts: Creating other shared objects like map (CMV)

texts: We need to create an object which is accessible from the controller and widgets. Rather like how we can currently access the map object. Could anybody shed light on how this can be done?


texts: Topo To Raster executing but not creating the raster

texts:
I am using Topo to Raster by File (ArcGIS version 10.1) via arcpy on a large stream and point dataset. Arcpy runs fine and executes without issuing any error message, however it does not create the raster as it normally does. When I look in the project workspace I can see that the S File that Topo to Raster creates temporarely (before deleting it when execution is terminated) have been created and its size is slightly over 2GB.
Not sure what is the S File used for as it is a binary and can't inspect it.
Is there a limit on the size of the S File created by Topo to Raster?
What can I do to get around this issue in addition to reduce the size of the point dataset?

Note that:
There is enough space in the hard disk thereore the issue is not space; and
The tool runs ok when I reduce the size of the dataset (option which I would like to avoid).


texts: jpg to geotiff transformation with gdal_translate and gcps fails

texts:
I am trying to create a geotiff from a jpeg by inputting gcps with gdal_translate:

gdal_translate -of GTiff -a_srs EPSG:4326 -gcp [1616.0 0 -121.459869532 38.5822533549] -gcp [0 0 -121.460081357 38.5831052365] -gcp [0 1080.0 -121.460807872 38.5829948432] -gcp [1616.0 1080.0 -121.460596039 38.582142963] pict20140910_131103_0.JPG tmp.tif

followed by gdalwarp to set the geotiff projection:

gdalwarp -s_srs EPSG:4326 -t_srs EPSG:4326 tmp.tif geo_pict20140910_131103_0.tif

however the gdalwarp process fails with:

  ERROR 1: Failed to compute GCP transform: Transform is not solvable.

The gcp coordinates are the image footprint vertices and were calculated using the image location, altitude, orientation and camera attributes.  The gcps are in lat/lon (wgs84 epsg:4326), and the image size is 1616x1080.

Here is the gdalinfo for the file after gcps have been added, but before gdalwarp:

Driver: GTiff/GeoTIFF
Files: tmp.tif
Size is 1616, 1080
Coordinate System is `'
GCP Projection =
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
GCP[  0]: Id=1, Info=
          (0,0) -> (-121.459869532,38.5822533549,0)
GCP[  1]: Id=2, Info=
          (0,0) -> (-121.460081357,38.5831052365,0)
GCP[  2]: Id=3, Info=
          (0,1079) -> (-121.460807872,38.5829948432,0)
GCP[  3]: Id=4, Info=
          (0,1079) -> (-121.460596039,38.582142963,0)
Metadata:
  AREA_OR_POINT=Area
  EXIF_BrightnessValue=(4.275)
  EXIF_ColorSpace=1
  EXIF_ComponentsConfiguration=0x1 0x2 0x3 00
  EXIF_CompressedBitsPerPixel=(3)
  EXIF_Contrast=0
  EXIF_CustomRendered=0
  EXIF_DateTime=2014:09:10 13:11:02
  EXIF_DateTimeDigitized=2014:09:10 13:11:02
  EXIF_DateTimeOriginal=2014:09:10 13:11:02
  EXIF_DigitalZoomRatio=(1)
  EXIF_ExifVersion=0230
  EXIF_ExposureBiasValue=(0)
  EXIF_ExposureMode=0
  EXIF_ExposureProgram=4
  EXIF_ExposureTime=(0.001)
  EXIF_FileSource=0x3
  EXIF_Flash=16
  EXIF_FlashpixVersion=0100
  EXIF_FNumber=(4)
  EXIF_FocalLength=(16)
  EXIF_FocalLengthIn35mmFilm=24
  EXIF_GPSAltitude=(66)
  EXIF_GPSAltitudeRef=00
  EXIF_GPSImgDirection=(79)
  EXIF_GPSImgDirectionRef=T
  EXIF_GPSLatitude=(38) (34) (57.45)
  EXIF_GPSLatitudeRef=N
  EXIF_GPSLongitude=(121) (27) (37.22)
  EXIF_GPSLongitudeRef=W
  EXIF_GPSVersionID=0x2 0x3 00 00
  EXIF_Interoperability_Index=R98
  EXIF_Interoperability_Version=0x30 0x31 0x30 0x30
  EXIF_ISOSpeedRatings=2000
  EXIF_LightSource=0
  EXIF_Make=SONY
  EXIF_MaxApertureValue=(2.96875)
  EXIF_MeteringMode=5
  EXIF_Model=NEX-5T
  EXIF_Orientation=1
  EXIF_PixelXDimension=1616
  EXIF_PixelYDimension=1080
  EXIF_ResolutionUnit=2
  EXIF_Saturation=0
  EXIF_SceneCaptureType=0
  EXIF_SceneType=0x1
  EXIF_Sharpness=0
  EXIF_Software=NEX-5T v1.01
  EXIF_WhiteBalance=0
  EXIF_XResolution=(350)
  EXIF_YCbCrPositioning=2
  EXIF_YResolution=(350)
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 1080.0)
Upper Right ( 1616.0,    0.0)
Lower Right ( 1616.0, 1080.0)
Center      (  808.0,  540.0)
Band 1 Block=1616x1 Type=Byte, ColorInterp=Red
Band 2 Block=1616x1 Type=Byte, ColorInterp=Green
Band 3 Block=1616x1 Type=Byte, ColorInterp=Blue


texts: How to import a 1x1 degree text file into QGIS?

texts:
I have a text file with 360 columns and 180 lines, each number indicating a value for a one-degree cell.

1.32    1.43    1.67    1.61    1.57    1.39    1.22    ...
1.19    1.28    1.22    1.07    0.94    0.76    0.69    ...
0.83    0.83    0.77    0.67    0.58    0.55    0.62    ...
...

I haven't found any indication on how such a file can be imported into QGIS.
How can I do this?


texts: How to split OSM roads into individual segments at intersections?

texts:
I want to create a road-network for use with pgRouting using OpenStreetMap data. I loaded a shapefile from GeoFabrik into a Postgres table (with PostGIS enabled). However, one problem I had was that the roads do not always end at intersections, so I decided to split them all at every intersection or crossing.

To identify all the interesections where roads crossed or intersected I used the following SQL (similar to a previous question):

CREATE TABLE split_points as
SELECT DISTINCT
   ST_GeometryN(ST_Intersection(a.geom, b.geom),1) as geom
FROM
   roads as a,
   roads as b
WHERE
    ST_Touches(a.geom, b.geom)
OR
    ST_Crosses(a.geom, b.geom)
    AND a.gid != b.gid
GROUP BY
   ST_Intersection(a.geom, b.geom);

I now want to split the roads using these points. I used the following approach:

CREATE TABLE split_roads as
SELECT
    ST_GeomFromEWKB((ST_Dump(ST_Split(g.geom, blade.geom))).geom) As geom,
    generate_series(1,ST_NumGeometries((ST_Split(g.geom, blade.geom)))) as gid
FROM
    split_points as blade,
    roads as g
WHERE
    ST_Intersects(g.geom, blade.geom);

The problem with this split approach is that the full road length remains in addition to all of the split pieces. To remove these un-split road geometries that were included I used the ST_Equals() function to identify them, and to delete them:

DELETE FROM split_roads USING roads
WHERE ST_Equals(split_roads.geom, roads.geom)

However, this approach does not remove all of the original unsplit geometries (although it does remove some of them). Is there a better approach for deletion (or overall) so that I have only the split geometries in a table?


texts: Unable to use a Toolbox with a input from DataInteroperability

texts:
I preface all of this by explaining that I'm not a Python developer (I usually use PHP and Java) nor do I use the ArcGIS product.  I've been attempting to use Python to assist a colleague but I'll try to provide all the information I can.

Our aim is to use a Python script that utilises toolboxes to update data.  There are two sources of data for the model: one from DataInteroperability and the other files on a mapped drive.  My current code is below:

import arcpy

if arcpy.CheckExtension("DataInteroperability") == "Available":
   arcpy.CheckOutExtension("DataInteroperability")
   print "Checked out \"DataInteroperability\" Extension"
else:
   print "Data Interoperability license is unavailable"

arcpy.ImportToolbox("G:/data/Toolbox/ScheduledLandChargesConversion.tbx","Land")
arcpy.Land.ScheduledLandCharges("Interoperability Connections\\ScheduledConversionLandCharges.fdl\\Assets_of_Community_Value Polygon","Q:\\Shp_Files\\Land_Charges")

When this is run I receive the following error:

  Traceback (most recent call last):
    File "Jonathan2.py", line 10, in
      arcpy.Land.ScheduledLandCharges("Interoperability Connections/ScheduledConversionLandCharges.fdl/Assets_of_Community_Value Polygon'","Q:\Shp_Files\Land_Charges")

  File "G:\data\Toolbox\ScheduledLandChargesConversion.tbx", line 25, in ScheduledLandCharges

  arcgisscripting.ExecuteError: Failed to execute. Parameters are not valid.

  ERROR 000732: Assets_of_Community_Value Polygon: Dataset Interoperability Connections/ScheduledConversionLandCharges.fdl/Assets_of_Community_Value Polygon' does not exist or is not supported

  The value does not exist.

  ERROR 000732: Input Features: Dataset 'Interoperability Connections/ScheduledConversionLandCharges.fdl/Assets_of_Community_Value Polygon'' does not exist or is not supported

  Failed to execute (ScheduledLandCharges).

I think this is because I'm passing the model the wrong kind of input, as it's after a feature layer.  I've copied the text for the DataInteroperability "path" from the toolbox's Python export, so that may be what's going wrong.  The model in question (ScheduledLandCharges) requires the following:

Help on function ScheduledLandCharges:

ScheduledLandCharges(Assets_of_Community_Value_Polygon=None, Land_Charges=None)
    ScheduledLandCharges_Land(Assets_of_Community_Value_Polygon, Land_Charges)

       (null)

    INPUTS:
     Assets_of_Community_Value_Polygon (Feature Layer)
     Land_Charges (Folder)


texts: Error importing shp into postgis using ogr2ogr?

texts:
I run an automated script that copies a number shapefiles into PostGIS. I have a problem with one dataset which fails during the export. As the roads.shp dataset is maintained by someone else, I automate making a copy of their shapefile into my PostGIS database. My error is;

ERROR 1: COPY STATEMENT FAILED
ERROR: numeric field overflow
DETAIL: A field with precision 13, scale 11 must round to an absolute value less than 10^2.
CONTEXT: COPY ccbc_roads, line 2, column length :"222.05900000000"

My code is as follows;

ogr2ogr.exe -f PostgreSQL PG:"dbname=postgis host=localhost user=postgres password=mypassword active_schema=ccbc_data" "roads.shp" -overwrite -a_srs "EPSG:27700" -lco GEOMETRY_NAME=geom -nln "ccbc_roads" -nlt PROMOTE_TO_MULTI --config PG_USE_COPY YES

The field that is causing the issue is the road length. Is there a way to change the scale from 11 to 2 in my code so that the script will work?


texts: How to sum up values of overlapping polygons in QGIS?

texts:
I have large single layers with many features overlapping.

Simplified example with 50% transparency:

Where the features overlap I want the attributes added together. I manually created the output I want below.


texts: LiDAR intensity image from multiple flight lines?

texts:
I'm working with multiple raw las tiles, and I've used Fusion to create an intensity image, but it shows the flight line overlap (see below). Any idea how can I smooth the image so that the overlap is blended into the rest of intensity values. Do I need to do any calibration or correction?

Can you please advise what tool I can use?


texts:
I am trying to merge 3 satellite images of 3 times. The images are median images. I have used the following code but there is an error.

var maskL8 = function(image_m) {
  var qa = image_m.select('BQA');
  var mask = qa.bitwiseAnd(1 << 4).eq(0);
  return image_m.updateMask(mask);
};

var collection_1 = ee.Image(ee.ImageCollection('LANDSAT/LC08/C01/T1_TOA')
.filterBounds(roi)
.filterDate('2018-06-16', '2018-10-15')
.map(maskL8)
.median());
Map.addLayer(collection_1, {bands: ['B5', 'B4', 'B3'], max: 0.3}, 'collection_1');

var maskL8 = function(image_s) {
  var qa = image_s.select('BQA');
  var mask = qa.bitwiseAnd(1 << 4).eq(0);
  return image_s.updateMask(mask);
};

var collection_2 = ee.Image(ee.ImageCollection('LANDSAT/LC08/C01/T1_TOA')
.filterBounds(roi)
.filterDate('2018-03-16', '2018-05-15')
.map(maskL8)
.median());
Map.addLayer(collection_2, {bands: ['B5', 'B4', 'B3'], max: 0.3}, 'collection_2');

var maskL8 = function(image_w) {
  var qa = image_w.select('BQA');
  var mask = qa.bitwiseAnd(1 << 4).eq(0);
  return image_w.updateMask(mask);
};

var collection_3 = ee.Image(ee.ImageCollection('LANDSAT/LC08/C01/T1_TOA')
.filterBounds(roi)
.filterDate('2017-11-16', '2018-02-15')
.map(maskL8)
.median());
Map.addLayer(collection_3, {bands: ['B5', 'B4', 'B3'], max: 0.3}, 'collection_3');

var collectionFromConstructor = ee.ImageCollection([collection_1, collection_2, collection_3]);
print('collectionFromConstructor: ', collectionFromConstructor);

var collectionFromImages = ee.ImageCollection.fromImages(
  [ee.Image(3), ee.Image(4), ee.Image(5)]);
print('collectionFromImages: ', collectionFromImages);

var mergedCollection = collectionFromConstructor.merge(collectionFromImages);
print('mergedCollection: ', mergedCollection);
Map.addLayer(mergedCollection, {bands: ['B5', 'B4', 'B3'], max: 0.3}, 'mergedCollection');

Error says:

mergedCollection: Tile error: Expected a homogeneous image collection, but an image with incompatible bands was encountered:
  First image type: 12 bands ([B1, B2, B3, B4, B5, B6, B7, B8, B9, B10, B11, BQA]).
Current image type: 1 bands ([constant]).
          Image ID: 0
Some bands might require explicit casts.

Error in using collection_1.merge(collection_2):

"collection_1.merge is not a function."


texts: Geoserver WMS tile rendering is too slow

texts:
I have following geoserver status. I have installed JAI but memory used by jai is 0 and image rendering for wms is very slow.

Is it issue with memory allocation? I am using Sql server 2008.


texts: USGS/FGDC SVG geoscience point symbols for use in QGIS

texts:
I have over 700 symbols created and tested for use on geological maps, using the USGS symbol documents as a source (https://ngmdb.usgs.gov/fgdc_gds/), converting the original PDF using Inkscape.

I am using only a few of these on my own geological maps, so wondered if they would be of more general interest or use in the QGIS community. They are at http://members.iinet.net.au/~richard.langford/USGS.zip, and if they are fit for purpose I would like to know.

Has anyone already done this for USGS symbols as SVGs, or is this new?


texts:
I have an ESRI ASCII file with usual 6 rows of headers followed by elevation values in metres. The NO_DATA value in the header is set to -9999.

I am trying to convert this ASCII to XYZ using gdal_translate and do not want to include nodata values of -9999 in the output XYZ File. I am using the a_nodata command to ignore -a_nodata values in the output XYZ file however have so far been unsuccessful in getting this to work. I'm hoping someone can help me out as to where I am going wrong with the -a_nodata command. I've included my current code below which creates the XYZ file but does not ignore no data values.

gdal_translate -of XYZ D:\WORK\COASTAL_PRODUCTS\P_10853\ASCII_GRIDS\DSM\tq6575_20161005lidaru.asc D:\WORK\COASTAL_PRODUCTS\P_10853\ASCII_XYZ\DSM\tq6575_20161005lidaru.txt -a_nodata -9999


texts: How to get nanmean of each pixel in multiple raster image?

texts:
I wrote a code based from rasterio and numpy that gets the nanmean of every pixel from multiple of raster images. I first read it using rasterio, then I convert the array to float32 so that I can set the NaN (-32768) values to np.nan. I then add each image to a list.
When the list is filled. I then compute it using np.nanmean from numpy. Real enough that NaN values were ignored in computing the output image. However, I do not understand the values from the resulting image. The input raster images have a value from 0 to 5. My output has a value 0 - 20000.
Note: The input images are really Float32 in QGIS. I dont know why it was converted to Int16 in rasterio.

Why was it converted to Int16?
How do I correct this? I want an output with nanmean from the inputs. The real values should be between 0 - 5 also.

import rasterio
from glob import glob
import numpy as np

path = '/my/directory/*.nc'
output = 'output/AOT_L2_Mean_MAM.tif'
datasets = glob(path)

rasters = []
for data in datasets:
    ds = rasterio.open(f'netcdf:{data}:AOT_L2_Mean')
    band = ds.read(1)
    #print(band.dtype)
    band = band.astype('float32')
    band[band==-32768] = np.nan
    rasters.append(band)

out = np.nanmean(rasters, axis=0)

with rasterio.open(
    output,
    'w',
    driver='GTiff',
    height=out.shape[0],
    width=out.shape[1],
    count=1,
    dtype=out.dtype,
    crs='EPSG:4326', # +proj=latlong
) as dst:
    dst.write(out, 1)


texts: Count the number of buffers intersecting with zip code polygons

texts:
I have a short R script in which I read two sets of data:

a shapefile showing the zip codes of a European country
an Excel showing multiple points coordinates throughout the country

I used the points to define buffers around them. These buffers naturally intersect with the zip codes polygons and one polygon can intersect multiple buffers.
Here is the code so far:

F_A <- read.xls("F_A.xlsx", sheet = "FA")
F_A$lon <- as.numeric(levels(F_A$lon))[F_A$lon]
F_A$lat <- as.numeric(levels(F_A$lat))[F_A$lat]
F_A <- na.omit(F_A)

coords <- cbind(F_A$lon, F_A$lat)
pts <- SpatialPointsDataFrame(coords = coords, data = F_A, proj4string = CRS("+proj=longlat +datum=WGS84"))
pts_proj <- spTransform(pts, CRS("+init=epsg:3347"))

Buffers <- gBuffer(pts_proj, width = 7000, byid = TRUE)
Buffers_xy <- spTransform(Buffers, CRS("+proj=longlat +datum=WGS84"))

zips <- readOGR(dsn = "Z", layer = "PLZ")
zips.wgs84 <- spTransform(zips, CRS("+proj=longlat +datum=WGS84"))
zips.wgs84.df <- fortify(zips.wgs84)

Polygons <- readOGR(dsn = "Poly", layer = "POLYGON")
Polygons.wgs84 <- spTransform(Polygons, CRS("+proj=longlat +datum=WGS84"))
Polygons.wgs84.df <- fortify(Polygons.wgs84)

I would like to add a new attribute to the zip code object counting the number of buffers intersecting with each zip code but fail to find the right command.
How would one go about doing this?

To be complete, that is a follow up question on a somewhat similar issue I had already asked about (QGIS creating a heatmap) and managed to solve using QGIS.


texts:
I am  having trouble writing to an excel file from ArcGIS Pro Sdk add-in that I created, I have a basic understanding of multithread, I am using task but the main thread does not wait for the excel to finish before trying to move ahead.

here is the helper I created.

public class ReportHelper
{
    public ReportHelper()
    {

    }

    public async Task WriteReport(string path)
    {
        await Task.Run(() =>
        {
            Microsoft.Office.Interop.Excel.Application excel = new Microsoft.Office.Interop.Excel.Application();
            excel.Visible = false;
            excel.DisplayAlerts = false;
            Workbook wrkbook = excel.Workbooks.Add(Type.Missing);
            string filename = Path.GetFileName(path).Replace(".gdb", "_importreport.xlsx");
            string sheetname = new string(filename.Take(11).ToArray());
            string directoryname = Path.GetDirectoryName(path);
            string savename = Path.Combine(filename, directoryname);
            Worksheet wrksheet = (Worksheet)wrkbook.Sheets.Add(After: wrkbook.Sheets[wrkbook.Sheets.Count]);
            //MessageBox.Show(sheetName);
            wrksheet.Name = sheetname;
            int x = 1;
            List<DataRecord> dRecords = Model.FeatureRecord;
            if (dRecords == null || dRecords.Count() < 1) return;
            foreach(DataRecord drecord in dRecords)
            {
                wrksheet.Cells[x, 1] = drecord.FeatureType;
                wrksheet.Cells[x, 2] = drecord.FeatureCount;
                wrksheet.Cells[x, 3] = drecord.ImportCount;
            }
            wrkbook.SaveAs(savename, Microsoft.Office.Interop.Excel.XlFileFormat.xlWorkbookDefault, Type.Missing, Type.Missing, false, false, Microsoft.Office.Interop.Excel.XlSaveAsAccessMode.xlNoChange, Type.Missing, Type.Missing, Type.Missing, Type.Missing, Type.Missing);
            wrkbook.Close();
            excel.Quit();
        });

    }
}

I have found async streamwriter methods but I cannot find any awaitable methods for excel.

I am just trying to call these from button clock, everything works fine until I try to write the excel.

private  void button1_Click(object sender, EventArgs e)
{
    panel1.Enabled = false;

    Loader loader = new Loader();
    ReportHelper rpthelper = new ReportHelper();

    QueuedTask.Run(async () =>
    {
        await QueuedTask.Run( async () =>
        {
            foreach (string path in Model.DailyGdbs)
            {

                await loader.ReturnMasterPointList(Model.MasterGdb);

                await loader.Load(path);

                await rpthelper.WriteReport(path);
            }
        });
        await AsyncMessageBox();

    });
    //await AsyncMessageBox();
    //MessageBox.Show("The data has been loaded");
}


texts:
I am using the code below to find a country (and sometimes state) for millions of GPS points. The code currently takes about one second per point, which is incredibly slow. The shapefile is 6 MB.

I read that geopandas uses rtrees for spatial joins, making them incredibly efficient, but this does not seem to work here. What am I doing wrong? I was hoping for a thousand points per second or so.

The shapefile and csv can be downloaded here (5MB): https://www.dropbox.com/s/gdkxtpqupj0sidm/SpatialJoin.zip?dl=0

import pandas as pd
import geopandas as gpd
from geopandas import GeoDataFrame, read_file
from geopandas.tools import sjoin
from shapely.geometry import Point, mapping,shape
import time

#parameters
shapefile="K:/.../Shapefiles/Used/World.shp"
df=pd.read_csv("K:/.../output2.csv",index_col=None,nrows=20)# Limit to 20 rows for testing

if __name__=="__main__":
    start=time.time()
    df['geometry'] = df.apply(lambda z: Point(z.Longitude, z.Latitude), axis=1)
    PointsGeodataframe = gpd.GeoDataFrame(df)
    PolygonsGeodataframe = gpd.GeoDataFrame.from_file(shapefile)
    PointsGeodataframe.crs = PolygonsGeodataframe.crs
    print time.time()-start
    merged=sjoin(PointsGeodataframe, PolygonsGeodataframe, how='left')
    print time.time()-start
    merged.to_csv("K:/01. Personal/04. Models/10. Location/output.csv",index=None)
    print time.time()-start


texts: How do I ensure polygons are properly displayed when using 'on the fly' CRS transformations in QGIS 2.6?

texts:
QuantumGIS 2.6.0 (Windows 7 64-bit)

I have created a PostgreSQL POLYGON as WGS84 (CRS EPSG:4326). The polygon is a series of points around (but not actually on) the edge of the map. This is a rectangle (see the first image below).

POLYGON((-179.9 -89.999,-179.8 -89.999,-179.7 -89.999...), 4326

I have changed the CRS 'on the fly' to World Robinson (CRS ESRI:54030). In the previous version of QGIS this would now display my polygon as the racetrack shape of the Robinson projection. In QGIS 2.6 the correct shape is only displayed when I toggle 'Toggle Editing'.

I believe that only the extreme points are being displayed, until I toggle editing, when it is using all of the points. I have also observed this behaviour in linestrings.

Why is this, and how can I get my desired behaviour? Any help gratefully received!

Update
The polygon starts to render correctly as you zoom in.


texts: Converting Coordinates from wgs84 to UTM?

texts:
I have a geographical (longitude/Latitude) coordinates and
i want to Convert it to UTM.

11° 13' 47.034 N , 125° 1' 30.492" E

Can anyone help me?


texts: Why are there holes in these radar blobs?

texts: I was recently looking at a live weather radar online and I noticed several "holes" in the radar blobs usually near the center of the blob. I also noticed that in a few cases that the holes stayed in roughly the same place even when the blob moved around it. What is the cause of these holes? Is it some sort of radar error (ie. if the radar is situated inside the holes and can't detect anything that is too close to it? If not, then what is the cause?


texts:
I have found a great source for historical hurricane data, but I can't seem to find a way to download the track information.  It appears that you can only view the data on the website.  I was wondering if anyone else has had this issue, and if anyone knows a way to get access to the data (either as a .shp or .kml).

Here's the website: http://www.csc.noaa.gov/hurricanes/#

Yes I did see this, I was hoping there was some data already made.  Wishful thinking.  Guess I will make the files.


texts:
I have two raster layers that I would like to add together (raster1+raster2) but although they have the same resolution(15x15 nautical miles), they are not perfectly overlapping because they have been created using two different grids and cropped at the end using the same polygon. One is the result of a Utilisation Distribution evaluation via the adehabitatHR package for R and the other is the result of a longer procedure which includes R for data gathering and analysis, ArcMap 10.4 for spline interpolation with barriers and QGIS for further editing/geoprocessing.

I have looked at different proposed solutions via ArcMap (e.g. https://community.esri.com/thread/22270), R (e.g. Sum rasters with different extent in R) but none of them seems to work for me. Two of my rasters can be downloaded from https://www.dropbox.com/sh/aq746s58hbgerru/AADci-IGR1bgqIFWKQkOY7ata?dl=0. I have several rasters, about 120, that cover 10 different areas. And Ideally I would like an automated way to “snap” one to the other by area.

EDIT:
I have actually found a solution very easily by using the basic raster calculator in QGIS. For some reasons the SAGA and GDAL raster calculators were not producing an output while the basic QGIS one works fine. Furthermore I was using an older version of QGIS that I still prefer to the newest 3.x and this prevented me from using the basic raster calculator in a model. So, in other words, I just opened my project in QGIS 3+ and put the raster calculator in a model. This works fine with me.


texts:
I have a problem importing Terra MOPITT data into QGIS (2.18) correctly.
The data is in hdf5(HDFEOS) format. QGIS/GDAL imports the data, but the image is rotated by 90° and shifted.

QGIS error: CRS was undefined: defaulting to CRS EPSG:4326 - WGS 84

Do you know which CRS I have to use?


texts: Overlaying Raster image with Google Maps in OpenLayers3

texts: I am trying to overlay Google Maps with a raster image layer in Openlayers3 and am new to this. Can anyone help by providing the steps involved in performing this task?


texts:
I have a large amount of point data that I have in CSV form. I am able to insert it all into a PostGIS table no problem. What I want to do though, is split the data up on insert into multiple separate tables. The data is statewide and I would like to have separate tables for each county within the state.

I know I could insert all of the points into a temp table and then INSERT INTO each county based on an ST_Intersect() but is there a better way of doing this? Some sort of trigger function that will run faster? Or is there a way to use a CONSTRAINT on each table where anything within/not within a certain polygon is rejected?

Edit: Sorry for the change but we ended up changing how we wanted to store the data: One table for our county of interest and another for the rest of the data because we occasionally need it. The other answers were not quite what we were looking for, we really wanted to avoid temp tables and adding fields and I did come up with a solution for it.


texts: How do you create a rainbow polygon hover effect in Mapbox?

texts:
We're working on a fun website called VicePrice. It's a happy hour website using GIS instead of a list directory to show happy hour info. We've run into a snag and perhaps someone could give us their thoughts.

You can play with our current design here. In Washington DC, Logan Circle has played a large role in the gay community. It's served as a place of pride for decades and I'd like our map to tip our hat to that.

We're thinking that when a user hovers over Logan Circle, instead of a gold overlay, a rainbow overlay would show instead. Like with the rainbow flag here.

So far the only possibly solution we've come up with has been to create a horizontal polygon for each color in the rainbow, and when a user hovers over one, it triggers all polygons. But that's been a huge challenge applying the same rule to those polygons.

I was wondering if someone else has dealt with this issue and how they would suggest tackling this problem.


texts: Python GDAL showing incorrect metadata for GRIB file

texts:
I'm trying to read a GRIB file containing ERA-5 single level data that I downloaded from the Climate Data Store (CDS, https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-single-levels?tab=form).

Opening the grib-file in the NOAA Weather and Climate Toolkit shows me that all the correct data is inside the file, i.e., the variable names match the data that I selected in the CDS and I can also check the dates/times of the variables, they are correct. See this screenshot:

However, when I open the same file using python/GDAL, i see metadata that is completely different. I'm running the following little script to see which variables are inside the file:

import gdal
import numpy as np

filepath = r"F:\era5\sfc_2008_01.grib"

dataset = gdal.Open(filepath)

var_names_long = np.array([])

for layer in range(1, dataset.RasterCount):

    meta = dataset.GetRasterBand(layer).GetMetadata()

    var_names_long = np.append(var_names_long, meta['GRIB_COMMENT'])

print(np.unique(var_names_long))

This gives as output:

['Convective cloud cover [%]' 'Density [kg/m^3]'
 'Evaporation (of water) [m]' 'Ice cover (ice=1;no ice=0) [proportion]'
 'Precipitable water content [kg/m^2]' 'Salinity [kg/kg]'
 'Surface pressure [Pa]' 'Total cloud cover [%]'
 'Total column water [kg/m^2]' 'Total precipitation [m]'
 'Water runoff [kg/m^2]']

These are completely different variables, that I've never downloaded/selected?! Any ideas what might be causing this? It seems to me like it's a problem inside GDAL, since the NOAA toolkit gives the correct names... I'm using gdal 2.4.0 (which is the most recent version).


texts: How to incrementally add points to a polygon and eventually increasing the size of the polygon in postgis?

texts:
I want to implement a clustering algorithm in postgis in which new points are incrementally added to an existent polygon and the polygon subsequently increases it's size to a the minimum possible value so as to cover all the points.

One solution that I came up with is:-

1) When there is first point then nothing to do.

2) When the second point comes we simply make a linestring out of it.

3) When the third point come we make a closed linestring (by making any one point the starting point as well as the last point using ST_AddPoint) using the given three points and ultimately make a polygon (using ST_MakePolygon).

4) Now when the third point comes we use the following algorithm:-
If the new point is outside the existing polygon and within a threshold distance from it, we find two of closest points of the polygon from the new point. Now we connect the new point to these two points and eventually form the required polygon.
Now, how to find the two closest points? We can find distances of all the points constituting the existing polygon from the new point, order them in decreasing order and limit by 2.


texts: How to open a layer in QGIS from another application?

texts:
Let's say I have some metadata about a layer stored in some custom application, and I would like to open or view the layer in QGIS.

What, if any, are my options?

Edit: formats include shp, wms, wfs, geotiff.


texts: Render Mapbox tiles with Mapzen's Tangram

texts: We've built a few projects using Mapzen's vector tiles and rendering them through Tangram. However, now that Mapzen is shutting down, we need to use other tile servers (we would like to maintain the rendering through Tangram, since we've put some work into getting the styles right).Is there any way of rendering the tiles provided by Mapbox http://a.tiles.mapbox.com/v4/mapbox.mapbox-streets-v7/{z}/{x}/{y}.mvt using Tangram and their .yaml files?


texts: How to add a new permanent column that contains result from ST_AREA select queries

texts:
I have a bunch of polygons in my CartoDB map, and I want to create a new column called AREA in the table MY_TABLE that contains values from my area calculation query:

select ST_Area(the_geom::geography) as area_sqm FROM My_TABLE

this statement works, but the column area_sqm that contains the area value clears every time I clear the query, how can I have this permanently as part of my table so I can display the area value in the info table of my map?


texts: Pycharm cannot find reference 'SearchCursor' in 'arcpy.da'

texts:
Whenever I try to use a search cursor from arcpy in pycharm, I get this unresolved reference error. I should note that it was working fine earlier in the night, I have no idea what happened.

I've tried everything in this Q&A already:

https://stackoverflow.com/questions/11725519/pycharm-shows-unresolved-references-error-for-valid-code

Other arcpy functions seem to work fine. I've also tried restarting computer and reinstalling pycharm. I can use the search cursor when I use IDLE.


texts:
I am trying to build QGIS 3.10 on Windows x64 and follow INSTALL file.

I am getting following error

2>------ Build started: Project: version, Configuration: RelWithDebInfo x64 ------
2>------ Build started: Project: qgis_core_autogen, Configuration: RelWithDebInfo x64 ------
2>  Automatic MOC for target qgis_core
3>------ Build started: Project: crssync, Configuration: RelWithDebInfo x64 ------
3>  crssync.vcxproj -> C:\project\QGIS-release-3_12_vs14\ms-windows\osgeo4w\build-qgis-test-x86_64\output\bin\RelWithDebInfo\crssync.exe
3>  crssync.vcxproj -> C:/project/QGIS-release-3_12_vs14/ms-windows/osgeo4w/build-qgis-test-x86_64/apps/qgis-test/pdb/crssync.pdb (Full PDB)
4>------ Build started: Project: synccrsdb, Configuration: RelWithDebInfo x64 ------
4>  Running crssync
4>C:\Program Files (x86)\MSBuild\Microsoft.Cpp\v4.0\V140\Microsoft.CppCommon.targets(171,5): error MSB6006: "cmd.exe" exited with code -1073741515.
========== Build: 3 succeeded, 1 failed, 2 up-to-date, 0 skipped ==========

Earlier solutions as described in http://osgeo-org.1560.x6.nabble.com/Build-qgis-errors-td5219922.html unable to solve the problem as it is based on Qt4.

I had also tried Debugging QGIS from Visual Studio and https://www.shaeffer.co/compiling-qgis-on-windows    but no solution.


texts: Creating auto generated number field as like ObjectID in ArcGIS Desktop?

texts: I have a FGDB with many feature classes. I need a "FeatureID" field it should be auto generated with a sequence number while adding a new feature in Edit session. It is like the existing field "OBJECTID". It is like auto numbering while adding new feature


texts: Erroneous location when converting from a GCS_WGS_1984 to UTM

texts:
I have a shapefile map of Belgium divided by cantons that I want to re-project in order to have units in meters and to be able to calculate distances. I know how to do it with ArcGIS Projection tool, but the result always gets wrong map location.

I will detail my steps. First of all, my original map had unknown projection and coordinate system. However, it is pretty clear that it is in GCS_WGS_1984, with decimal degrees, since the Left/Top/Right/Bottom limits are precisely equal to the real ones of Belgium, when a map of this country is projected that way. So far, so good. I defined the GCS_WGS_1984 projection with DD units in ArcGIS, using the Define Projection tool.

The problem is in the next step. I then want to convert from that to another projection. One from which I can have units in meters. I tried a bunch of UTM options, using the Project tool after finding out that Belgium is around UTM 31N. However, no matter which UTM projection I choose, the location of the country after re-projected becomes wrong.

For instance, re-projecting from thw GCS_WGS_1984 to WGS_1984_UTM_Zone_31N gives me the following extent limits: Left:468861,233828 m; Top:5707355,258593 m; Right:742547,667605 m; Bottom:5485665,513781 m. That's very off the target.

Could anyone point me what am I doing wrong or what is missing?

PS: Another detail is the following. I tried opening the original map among other country maps with the GCS_WGS_1984 projection in decimal degrees. All countries are positioned correctly. Then, I do the same operation described above to convert some countries besides Belgium to check what happens. All of them get located wrong - what means, my problem doest not seem to be in my Belgium files, but rather in me doing something wrong.


texts: I am working in ArcMap modelbuilder, doing a cost path anaylsis. In my analysis I am running a pipeline. it will cost 100 dollars to cross or parallel a road. I have buffered my roads in modelbuilder (which created about 14 polygon records in the buffered roads attribute table). I have added a new field to the attribute table called "cost" in modelbuilder. Now I want to make the value of all 14 records in the "cost" field equal 100 --- before I convert to raster. How do I make the value of all 14 records in the "cost" field 100?


texts: Layer disappears when changing from CRS in degrees to CRS in meters in QGIS

texts:
The problem

I have tens of thousands of small polygons for which I want to create buffer/doughnuts for, but I don’t get the units right. I want meters and not degrees, and changing to (what I think is) the right CRS my layer doesn’t longer show. The polygons are small farming plots. I work in QGIS 3.4.15.

Let me take you through it

I start by adding a vector layer. I opt for ‘Directory’ as a source type (Encoding: UTF-8 by default) and under Source (Type: UK. NTF2 by default) I point to the folder in which I store the shapefile.

The default CRS is set to EPSG:4326 – WGS 84, for both layer and the project (checked in bottom right), and unit is thus in degrees. I add the layer world map (EPSG:4326 – WGS 84 by default) by taping “world” in the coordinate box in the middle bottom, I also add Open street map (EPSG:3857 – WGS 84 / Pseudo-Mercator by default) from XYZ Tiles. The polygons appear nicely on the map considering some are just adjacent to the nearby lake, many nicely adjacent to roads, etc.

When attempting to do buffers the distance is obviously set to degrees. I have tried to change the CRS, both for the layer(s) and the projects, but without luck. After researching (e.g. here: https://epsg.io/?q=Kenya%20kind%3APROJCRS) I believe that the correct CSR would be Arc 1960 / UTM zone 36S EPSG:21036 (it is my understanding that either 36N or 36S would work as the polygons appear both south and north of, but close enough to the equator). When setting the mentioned CRS for my polygon layer, the polygons disappear from the map, and remain gone even after setting the project CSR (by bottom right) to the same. It seems that whenever I try with different CRS in UTM, the layer doesn’t show. Does this mean that the polygons are somehow not compatible with units in meter? Is there an issue with importing the .prj file?

Other things I've tried:

I have made a good attempt searching for solutions. I have ensured that the polygon layer is on top of the layers list. I have ensured the file names are identical and tried all different suitable CRS I can think of. I have seen similar questions been solved by saving the .shp 'as', but I don't think that's my issue as it's already a .shp file. Can it be the case that I need to perform a datum transformation? I gave it a shot using a transformation code (1122) provided by above link, but then the polygons were projected on the Antarctica on the world map, so obviously something went quite wrong...

Note that when I have changed to what I believe is the right CRS I can execute buffers in units of meters, however, they do not appear on the window…


texts:
I'm currently having a problem with Geotools.

I try to generate a buffer (in meter using Lambert 72) on a map that has WGS84 as it's coordinatesystem.
(problem is a bit like this one: http://www.mail-archive.com/geotools-gt2-users@lists.sourceforge.net/msg08046.html, however I can't manage to solve my problem with it).

I have a couple of points in a map (in WGS84), and I want to add buffers around it using meters. As said here ( http://docs.geotools.org/latest/userguide/library/jts/operation.html ) you should transform the Geometry into the other CRS, buffer and then go back.

So what I do:

CoordinateReferenceSystem WGS = DefaultGeographicCRS.WGS84;
CoordinateReferenceSystem lambert = CRS.decode("EPSG:31300");
MathTransform convertToMeter = CRS.findMathTransform(WGS, lambert,false);
MathTransform convertFromMeter = CRS.findMathTransform(lambert, WGS,false);

//Start
targetGeometry=(Point) simpleFeature.getDefaultGeometry();
//transform to meter
targetGeometry = JTS.transform(targetGeometry, convertToMeter);
//add a buffer of 1000meters
targetGeometry= targetGeometry.buffer(1000d);
//go back to WGS 84
targetGeometry = JTS.transform(targetGeometry, convertFromMeter);

So now I expect to get a buffer of 1000meters, however what I get is a buffer of around 750meter, and I don't have any clue why.

Anyone could give me some advice (or could tell  me the stupid error I'm making)?


texts: Google Earth Engine - getting download URL of reprojected image not working?

texts:
I am using the Google Earth Engine Python API.

I'm having a strange issue where when I try to generate a download URL for an image, it works fine.

But when I reproject the image before generating the URL it does not work (gives me a "Site can't be reached" error. Also, when I pass in the desired CRS to the actual function that generates the URL it results in an error and doesn't even run. Here is my code so far:

This works:

test_url = test_img.getDownloadURL({
    'image': test_img.serialize(),
    'region': polygon,
    'bands': out_bands,
    'filePerBand': 'false',
    'name': img_name,
})

and outputs a URL that works fine.

However when I try this:

test_img = test_img.reproject('EPSG:32628')
test_url = test_img.getDownloadURL({
    'image': test_img.serialize(),
    'region': polygon,
    'bands': out_bands,
    'filePerBand': 'false',
    'name': img_name,
})

I get a "Site can't be reached" error when I click on the url, it still generates a url but it doesn't work.

When I try this:

test_url = test_img.getDownloadURL({
    'image': test_img.serialize(),
    'region': polygon,
    'bands': out_bands,
    'filePerBand': 'false',
    'name': img_name,
    'crs': 'EPSG:32628'
})

I get the following error:

ee.ee_exception.EEException: The requested image band is too large: was 14801736 bytes per scanline; must be no more than 262144 bytes per scanline.

Oddly enough, after doing the reprojection, I can verify manually that the reprojection took place, as all the bands in the image now have the CRS I specified.

Weird...not sure what I'm doing wrong as it seems pretty straightforward.

Any suggestions anyone?


texts: Configuring qgis map server to apache web server platform?

texts:
I installed Apache Web Server from OSGeo4W. here is the image shown below.
 But when i put URL
http://localhost/cgibin/qgis_mapserv.fcgiSERVICE=WMS&VERSION=1.3.0&REQUEST=GetCapabilities
i get an error which is shown below

Can you help me how to work in web platform using qgis software?. i am new user of QGIS currently having a project on it.


texts: Why not measure elevation from the center of the earth?

texts:
In trying to understand Geoid vs Ellipsoid height, I ran across this page.

One of the claims is:

  Zero elevation as defined by Spain is not the same zero elevation
  defined by Canada

Basically, they are indicating that MSL is not constant and that makes comparing things harder.

So my question is: why not measure elevation as a radius from the center of the earth?

I guess I'm proposing {lat, lon, radius} to define a point on or near the earth. This could eliminate many problems in my thinking. (Let's use inertial center of the earth as our reference point...)

Who cares that the topography of the earth is not constant? It wasn't constant in the first place, and to get height above ground, you just need some terrain data and a simple subtraction. I guess I'm proposing using one common reference point that doesn't change.


texts: Extracting NDVI for large extensions using ArcGIS for Desktop?

texts:
I'm trying to calculate NDVI in ArcGIS for Desktop using MODIS to entire neotropical region but I can't find MODIS values for all neotropical.

Instead, I foud only for small regions.

Does someone know how to extract NDVI for large extensions at the same time?


texts: exporting table in to a drive from Google earth engine returns blank rows

texts:
I am trying to find the mean precipitation of a watershed. When I print the feature collection BasinMean, I can see the mean precipitation in the console. A lot of them are zeros but the first two rows have non-zero precipitation.

When I export the table in to a drive, I expect a CSV file which is a timeseries with 3 columns - Station, date(timestamp) and mean precipitation. However, with the below code, the generated CSV has empty cells for mean precipitation. Any idea why the export function does not work?

Link to the workspace : https://code.earthengine.google.co.in/7f21f81224e88bed21e2db64fb7fda30

// Filter the image colelction for a certain timeframe
var filtered = Chirps_daily.filterDate('1981-01-01', '2016-12-31');

// Display the basin
Map.addLayer(Basin,{color:"#008080"},"Basin");
print(Basin,"Basin");

// Load the first image of the collection
var Chirps_first = ee.Image(filtered.first());
print(Chirps_first,'Chirps_first');

// Load the scale of the chirps image
var Chirp_scale = Chirps_first.projection().nominalScale();

var BasinMean = Basin.map(function(f) {
    return filtered.map(function(i) {
       var mean = i.reduceRegion({
       geometry: f.geometry(),
       reducer: ee.Reducer.mean(),
    });
    return f.setMulti(mean).set({date: i.date()})
  })
})

// flatten
BasinMean = BasinMean.flatten()
print(BasinMean.limit(20))

Export.table.toDrive({
  collection:BasinMean.limit(20),
  folder: "Google EE results",
  selectors:(["Station","date","precipitation"]),
  });


texts: MapBasic IDE Switch from Code view to Design view

texts:
I am brand new to mapbasic and I am working on a mapbasic project already developed. The previous developer has left and is not available for questions. I have searched mapbasic documentation and cannot find an answer to my question.

I am using MapBasic IDE beta 1.7

I have an mb file and want to confirm it relates to the windows dialog form I see in MapInfo. I can view the code for the mb file but I want to view the form itself to confirm the which function the text bar is referring to when the user types in information. Is there a method to switch from code view to design view? Or should there be a different file that contains the form design?


texts: ogr2ogr losing all ExtendedData attributes of KML

texts: I'm trying to convert a KML to a shapefile using ogr2ogr, but whenever I do, I end up with only blank Name and Description fields on the elements, losing the ExtendedData that is attached in the KML. A number of old answers here suggest that this should've been solved with GDAL 1.9.2; I'm on GDAL 1.11.2, so that shouldn't be a problem. How do I get the ExtendedData information into my shapefile?


texts: Calculating the distance between polygons and polygons (or polygons and lines) in QGIS

texts:
I'm new to working with GIS data and QGIS. I'm trying to find an efficient way to measure/calculate the width of streets in a city.
I have two datasets: buildings (polygons) and streets (lines).
I'm interested in the average distance from the facade to facade for each side of a block (from 1 crossing to the next).
Is there a way to calculate the distances between the polygons (distance A in the picture) or between the polygons and the lines (distance B)?

The nodes of the buildings are manly at the corner points (with some exceptions, as seen in the picture). The lines of the streets end at every crossroad. Very few of the lines have additional nodes (e. g. at the bottom left corner of the picture).


texts: Discrepancy between CAD and GIS route lengths

texts:
I have imported a route from AutoCAD which is marked every kilometre and converted it into a shapefile. I put markers every 100 m using X-Tools. There is a difference between the two sets of markers, very small at first but growing substantially over the 380 Km route. I'm thinking that the CAD route will also be offset in the same way against the underlying DTM meaning that I can't be sure of the exact distance of any given feature along the route line.

Why is this happening and how do I correct it?


texts: Set custom extent in ArcGIS JSAPI 3.4/3.5

texts:
I am using ESRI suit 10.1 and ArcGIS JSAPI 3.4/3.5.

I want to set custom extent (for India) for this map and trying to assign the coordinate in center but its not getting proper map.

In previous version (3.0/3.1) we used dojo.toJson(map.extent.toJson()); with the help of firebug to identify the extent of the specific area(x min,y min, x max,y max) then we used that extent for base map.

What I tried :

-I checked the API reference but I did not get any information.

-Is there any Dojo problem (might be few references gets updated?)

-Tried to use lat,lon values at Center property but not getting perfect map for India

In new version I am not getting this so is there any way to identify the custom extent using previous tool and use those in center  property ?

Any help will be great !!!


texts: Definition Query to exclude certain values and return remaining values

texts:
I have been trying to create a Definition query to symbolize a set of water features. I cannot edit the .shp file and I want to exclude a multiple values. I have tried using NOT but have had no luck.

This is what works for one value:

NOT "OGF_ID" = 214620160

I don't know how to add more values to exclude.

What am I missing?


texts: Using geopandas GeoDataFrame in gdal.Grid for spatial interpolation viz IDW, nearest neighbour etc

texts:
I want to interpolate the point data using IDW with GDAL osgeo.gdal.Grid.
IDW_gdal = gdal.Grid(Output, input, zfield="rs",
               algorithm = "invdist",
               outputBounds = [ulx,
                               uly,
                               lrx,
                               lry],
              width = width, height = height)

It works when I am passing the shapefile in 2nd argument
but it does not work when I pass the geopandas GeoDataFrame
and
it throws the error
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-23-a7776d5dffb0> in <module>
     18                                lower_right_x_coordinate,
     19                                lower_right_y_coordinate],
---> 20               width = raster_width, height = raster_height)

C:\ProgramData\Anaconda3\lib\site-packages\osgeo\gdal.py in Grid(destName, srcDS, **kwargs)
    986         srcDS = OpenEx(srcDS, OF_VECTOR)
    987
--> 988     return GridInternal(destName, srcDS, opts, callback, callback_data)
    989
    990 def RasterizeOptions(options = [], format = None,

C:\ProgramData\Anaconda3\lib\site-packages\osgeo\gdal.py in GridInternal(*args)
   3295 def GridInternal(*args):
   3296     """GridInternal(char const * dest, Dataset dataset, GDALGridOptions options, GDALProgressFunc callback=0, void * callback_data=None) -> Dataset"""
-> 3297     return _gdal.GridInternal(*args)
   3298 class GDALRasterizeOptions(_object):
   3299     """Proxy of C++ GDALRasterizeOptions class."""

TypeError: in method 'GridInternal', argument 2 of type 'GDALDatasetShadow *'

Alternately, I have passed the GeoDataFrame into ogr as GeoJson as given here
df_to_pass_in_ogr = ogr.Open(df_new.to_json())
print(df_to_pass_in_ogr)

output is
<osgeo.ogr.DataSource; proxy of <Swig Object of type 'OGRDataSourceShadow *' at 0x0000000BCDA48A80> >

but it did not work too and throws the same error.
Is there any way to use geopandas GeoDataFrame in the 2nd argument where you have passed the point feature shapefile like input_shape_file_name here?


texts: ArcGIS export layout with inserted EMF graphics

texts:
using ArcGIS 10 SP3 for some simples mapping tasks. But::: I have a layout with some inserted EMF graphics. These graphics are in fact formatted Excel-tables. The first rows in the tables are shaded with a light gray background.

The EMFs are drawn correctly in the ArcMap layout, but dependent on the position on the page the shades disappear in the exported maps: if the graphic is positioned at the right page border the shades disappear, at the bottom or left side, the output will be correct.

I have copied and pasted the same graphics a dozen times and placed it side by side. It seems, that all graphics which lay mainly outside the DIN A4 portrait format are not correctly exported, independent of the page settings (in my case landscape). In the meanwhile I have deleted ALL OTHER layout elements and the effect still is the same.

Has anybody any idea? Is this one more of ArcGIS strange features?


texts: How to use ArcPy to iterate through each grid cell and make it the current extent?

texts:
I have a topology that covers a big area. I cannot validate all at once 'cause I get errors. So I've created a grid that covers the area.
Through a script runned from my .mxd I intend to:
1. itterate through each grid cell
2. set it as the current extent of the map,
3. use the tool "Validate Topology", having checked the option "VISIBLE_EXTENT".

The workflow sounds correct to me, but I'm blocked at the first two steps. Here's my code:

import arcpy
mxd = arcpy.mapping.MapDocument("CURRENT")
df = arcpy.mapping.ListDataFrames(mxd, "Layers")[0]
gridLayer = arcpy.mapping.ListLayers(mxd, "Caroiaj_RO", df)[0]

# Create SearchCursor
with arcpy.da.SearchCursor(gridLayer, ("GridName2")) as cursor:

# Start a loop that iterates through all rows
    for row in cursor:

    # Clause to select only the current record
        whereClause = ' "GridName2" = ' + str(row.GridName2)

        # Select Layer by attribute
        arcpy.SelectLayerByAttribute_management(gridLayer, "NEW_SELECTION", whereClause)

        df.extent = gridLayer.getSelectedExtent(True)

        arcpy.RefreshActiveView()

del cursor, row
del mxd

Thanx a lot!


texts: Setting a default value in a QGIS creates empty field

texts:
I use Qgis 2.18.3 and tried to set default values in new column in a existing shapefile via the layer properties as described here: Setting a default value in a QGIS field
or here
How to set a default value for attributes in shapefiles?.

But fields are not populated, tried several field combinations: text, real, integer but always no result. Tried some combinations with editing layer/not editing layer.
Is there kowen issue or bug with this function? Do I miss an crucial step?


texts: GeoJSON vector layer created using OpenLayers is broken

texts:
What i want to do here is to create a vector layer (using openlayers) from a GeoServer JSON Response. Actually i'm doing this:

var url = "http://"+ip+":"+puerto+"/geoserver/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=tiger:poi&outputFormat=text/javascript&format_options=callback:getJson";

jQuery.ajax({
                    jsonp: false,
                    jsonpCallback: 'getJson',
                    type: 'GET',
                    url: url,
                    async: false,
                    dataType: 'jsonp',
                    success: function(jsonObject) {
                        alert(jsonObject);
                        //handler(jsonObject);
                        create_GeoJSON(jsonObject);
                    }
});

With the code above (after a long time of headaches with cross domain issues) we can retrieve data from GeoServer (through WFS) using JSONP.

The function below handles the response, tries to create the vector layer and put that vector over the map, using GeoJSON.

function create_GeoJSON(json){

    var in_options = {
            'internalProjection': window.map.projection,
            'externalProjection': window.map.projection
             };

    var jsonString = JSON.stringify(json);

    var geojson_format = new OpenLayers.Format.GeoJSON(in_options);
    var vector_layer = new OpenLayers.Layer.Vector("Consulados");
    vector_layer.addFeatures(geojson_format.read(jsonString, "FeatureCollection"));
    window.map.addLayer(vector_layer); }

Now, the problem is that the map isn't showing the layer components properly. As you can see i only see a point (or maybe a lot of superposed) drawed over the (0,0). This seems to be an issue related with the coordinates, projections or something like that.

By the way, here is the map init:

var options = {
   projection: "EPSG:900913",
   displayProjection: "EPSG:4326",
   controls: [],
};

window.map = new OpenLayers.Map('map_canvas', options);

What you think about that? I "waste" a lot of time over this right now and still don't have a solution.

Thanks.

PS: Feel free to edit the post, my english is a little poor.


texts: In ArcGIS Pro, how can I transfer the color scheme to other maps so that the same color is maintained for each class?  For example 2,1-4,0 degrees Celsius "deep blue". Class 4,1-6,0 "light blue".


texts: QGIS world imagery does not repeat itself

texts: Is it possible to make raster imagery from a service like OSM repeat itself like google maps does? I want to create a map of the world that shows Australia at the center, but the imagery stops at new zealand? It doesnt repeat the imagery from the other side of the dateline. I imagine this is just a tick box somewhere, any ideas. (i know i could do this manually by exporting and joining, but i am more curious if this is possible in qgis in a more dynamic way)


texts: Locating where different polygons (buffer rings) from different layers are intersected based on same attribute value (QGIS)

texts:
I need to locate where different polygons (buffer rings) from different layers are intersected based on a same attribute value.
As shown in the figure, there are different buffer rings from different layers (blue, red and black) at different distances with a different attribute values in every ring, let's say for example 'temperature' values.

Is it possible to make a spacial query to select where same temperatures values from the different buffer rings are intersenting? I mean, for example, where the 3°C buffer ring of every layer are intersecting.

Spacial query just interacts with whole layers not specific parts of them.


texts: arcpy.mapping.Layer error when working with files over multiple dataframes

texts:
I can get the following code to run in an mxd with a single dataframe. I need it, however, to run and loop through an mxd that has multiple dataframes with identically named files and data sources to those in my lyrlist.

EDIT Additional Info - The intent of the last if statement was to provide a workaround for the inability of arcpy to set the Display Expression in the layer properties dialog by re-adding a file with a name field. So, the crux of the question is how best to update a layers Display Expression in an automated a way as possible (preferably with python as I have no familiarity with ArcObjects).

import arcpy
lyrlist = ['Inset Bridges','Bridges', 'Transportation Points', 'Transportation Lines']
mxd = arcpy.mapping.MapDocument("Current")
dflist = arcpy.mapping.ListDataFrames(mxd)
for df in dflist:
   if df.scale <= 126720:
      for lyr in (arcpy.mapping.ListLayers(mxd, "", df)):
         if lyr.supports("DEFINITIONQUERY"):
            if lyr.name in lyrlist:
                NewFile = arcpy.mapping.Layer(lyr.dataSource)
                OldFN = lyr.name
                NewFile.name ="NewFi"
                lyr.name = "OldFi"
                arcpy.mapping.InsertLayer(df, lyr, NewFile, "AFTER")
                #the following line results in the ValueError when run in an
                #mxd that has multiple dataframes with filenames from lyrlist
                updateLyr = arcpy.mapping.Layer("NewFi")
                sourceLyr = arcpy.mapping.Layer("OldFi")
                arcpy.mapping.UpdateLayer(df, updateLyr, sourceLyr, True)
                updateLyr.name = OldFN
                arcpy.mapping.RemoveLayer(df, lyr)
arcpy.RefreshTOC()
arcpy.RefreshActiveView
print "script complete"
del mxd

The code works well until it gets to the where I assign the updateLyr variable.  That gives me the following error:

  Runtime error Traceback (most recent call last):

  File "", line 44, in

  File "c:\program files
  (x86)\arcgis\desktop10.1\arcpy\arcpy\arcobjects\mixins.py", line 389,
  in init
      super(LayerMixin, self).init(lyrfile)

  File "c:\program files
  (x86)\arcgis\desktop10.1\arcpy\arcpy\arcobjects_base.py", line 47, in
  init
      for arg in args))

  ValueError: Object: CreateObject Layer invalid data source


texts:
geo2=# create extension postgis;
ERROR:  PostGIS is already installed in schema 'public', uninstall it first
geo2=# create extension postgis_topology;
ERROR:  required extension "postgis" is not installed

It seems I'm getting an error that postgis is installed, but when I try to install postgis_topology it says that it isn't. When I try to drop the extension I get:

geo2=# drop extension postgis;
ERROR:  extension "postgis" does not exist

How can I drop postgis from the public schema so I can add it to my database?


texts: Modify blob representation field data with ArcPy

texts:
I want to modify the data stored in the "Override" blob type field in ArcGIS 10.7 with ArcPy.
When I use the search cursor, it finds nothing:
data = row.getValue("Override")
if data:
   print("Gotcha!")

else:
   print("Nothing...")

I somewhere read that search cursors doesn't recognise these fields.
Is there a way to get the stored data, perhaps by modify the blob "geometry"?


texts: arcpy Polyline.within(Polygon) returns False for all lines between pairs of Polygon vertices?

texts:
I'm trying to develop a function which will find the longest straight line length that can be drawn within any of the polygons in a feature class.  (To keep things simple, I'm only including lines that can be drawn between any two vertices of a polygon.)

My logic so far, is (for each polygon):

create Polyline geometries between all possible pairs of vertices of the polygon
for each Polyline created, test if it is fully within the polygon (eg, does not cross the boundary of a concave part of the polygon)
if it is fully within the polgyon and it's length is the largest found so far, update the maximum length variable.

Unfortunately, my arcpy line.within(poly) always returns false, even though the lines are created from pairs of points of that polygon.

I've commented out the test of line length to ignore that for now, but the "WITHIN" message (below) is never printed and the function always returns 0 (zero).

What I'm doing wrong in this code, below?

I can confirm that the function does create over 72,000 polylines for my test polygon, and that they all have realistic looking (greater than zero) lengths.  (Confirmed by using arcpy.AddMessage() for these values.)

def PolyLength(polyFC):
    maxLength = 0
    with arcpy.da.SearchCursor(polyFC, ["SHAPE@"]) as cur:
        for row in cur:
            poly = row[0]
            featureMaxLength = 0
            points = []
            lines = []
            for part in poly:
                points.extend(part)
                arcpy.AddMessage("Points count:  {}".format(len(points)))
            for a in points:
                for b in points:
                    if not a.equals(b):
                        lines.append(arcpy.Polyline(arcpy.Array([a,b])))
            arcpy.AddMessage("Lines count:  {}".format(len(lines)))
            for line in lines:
                #if line.length > featureMaxLength and line.within(poly):
                if line.within(poly):
                    featureMaxLength = line.length
                    arcpy.AddMessage("WITHIN and greater length:  {}".format(line.length))
                else:
                    arcpy.AddMessage("NOT WITHIN or NOT greater length:  {}".format(line.length))
            if featureMaxLength > maxLength:
                maxLength = featureMaxLength
    return maxLength


texts: select 8 rows in the table and calculate statistics for them, then select next 8 in the same table and calculate statistics for them

texts:
I have a large table which structure looks like that:

This is just to show you what I want to achieve, because my original table has nearly 1 million records.

Now imagine that I have temperature values calculated for every region which is represented by column geoid.

In this example I have 4 regions. For each region I have temperature values calculated 2 x daily for every day in a year.

What I need to get, is to calculate max temperature value for each day, thus in this example I need something which let me:

select first 8 rows -> calculate max temp. value for each geoid -> then select next 8 rows -> calculate max temp. value for new selection.

I need this to repeat for every day in a whole year.

How can I do this?


texts: qgis refresh button for attributeactions

texts:
QgsVectorLayer* vLyr=static_cast<QgsVectorLayer*>(aktLayer->qgsLayer);
vLyr->actions()->addAction(QgsAction::ActionType::OpenUrl, "Link im Webbrowser öffnen", "[% "+sAttribute+" %]" , ":/images/images/hotlink.png", false);
vLyr->actions()->setDefaultAction(vLyr->actions()->size()-1);

In my plugin after loading a layer to mapregistry i add with the code above an attribute-action to the layer. But i can't use it, because the button in attributes-toolbar is disabled. Looking at the layerproperties-dialog all seems to be ok and only after clicking 'OK' (not 'Apply') in the dialog the button in toolbar gets enabled and all works fine.
What function do I have to call to reload/refresh the attributes-toolbar in QGIS? Does anybody know this?

EDIT:

found the method "void QgisApp::refreshFeatureActions()" which maybe could solve my problem but how can i access this class from my cpp-plugin?

EDIT:

The final solution for me is the following:

QgsVectorLayer* vLyr=static_cast<QgsVectorLayer*>(aktLayer->qgsLayer);

QToolBar* attributeToolbar= mQGisIface->mainWindow()->findChild<QToolBar*>("mAttributesToolBar");

vLyr->actions()->addAction(QgsAction::ActionType::OpenUrl, "Link im Webbrowser öffnen", "[% "+sAttribute+" %]" , ":/images/images/hotlink.png", false);
vLyr->actions()->setDefaultAction(vLyr->actions()->size()-1);

QgsMapLayerAction * mapAction=new QgsMapLayerAction("Link im Webbrowser öffnen", attributeToolbar, aktLayer->qgsLayer, QgsMapLayerAction::Target::AllActions);//, new QIcon(":/images/images/hotlink.png"));
QgsMapLayerActionRegistry::instance()->addMapLayerAction(mapAction);
QgsMapLayerActionRegistry::instance()->removeMapLayerAction(mapAction);

I have to remove the QgsMapLayerAction because it doesn't know that it is an OpenURL an so it is empty. But only adding an entry to the QgsMapLayerActionRegistry emits a signal to check the toolbar availibility.
I don't know why I have no access to QgisApp-class, can't include qgisapp.h.


texts: How to smooth "internal" polygon borders

texts:
This is my first foray into GIS and GIS software. I am using ArcMap currently.

I have a shapefile of the world divided into different ecological regions. I would like to reduce the file size by generalizing/smoothing the edges of polygons inside the continent land masses. I want to perform these operations on edges that are adjacent to polygons for other ecological regions, but not on edges adjacent to oceans or waterways (no data).


texts: Move features from one layer to another existing layer

texts:
How can I move Selected features from one layer to another existing layer?
I have tried Copy (Data Management), Copy Features (Data Management), But I got the error:

  ExecuteError: ERROR 000258: Output (Layer name) already exists


texts: ogr2ogr getting ESRI File Geodatabase (FileGDB) work

texts:
I'm using ogr2ogr to load ESRI .gdb into PostGIS DB. I'd also like to export to .gdb the other way round.
But the included OpenFileGDB driver only has read access.

Now there are hints to use ESRI's FileGDB API SDK but I cannot get it work. What do I have to copy and where or what are the steps to get it work?

EDIT: It's MS 10 and GDAL 2.1.3. Installed it via OSGeo4W64.


texts: Get timeline for multiple nodes in Overpass

texts:
I am trying to use the timeline function to get information about the history of the nodes of an area, instead of the history of a single node as shown in the documentation.

This code works, and I am able to refer to the node id in the query filter condition by using u(id()):

[bbox:45.4,9.2,45.5,9.25];
node["amenity"="vending_machine"]({{bbox}});
foreach (
  if (u(id()) == 4612313390) {
    out;
  }
);

I try to pass the id to timeline, but I get a syntax error and an additional error saying that id should be a positive integer.

[bbox:45.4,9.2,45.5,9.25];
node["amenity"="vending_machine"]({{bbox}});
foreach (
  timeline(node, id());
  out;
);

Error:

Error: line 4: parse error: ',' or ')' expected - '(' found.
Error: line 4: parse error: ';' expected - ')' found.
Error: line 4: static error: For the attribute "ref" of the element "timeline" the only allowed values are positive integers.

What is the correct way to pass an id as argument to a function like timeline?


texts:
How can I get a 3x3 or 4x4 transformation matrix between two CRS, using the proj4 library in C/C++?

It's OK if it makes some approximations.


texts: PostGIS ST_Transform results shifted by about 25 metres

texts:
I have a road network downloaded from OSM and imported into PostGIS with pgrouting. I convert it with a specific local EPSG ST_Transform(geom, 26717), but the result is shifted about 20 metres to the east and 5 to the north. The new coordinates are instead equivalent to a different local EPSG 2029. The only real difference between these two seems to be NAD27 vs NSD27(76).

Example:
https://epsg.io/transform#s_srs=4326&t_srs=26717&x=-79.3545864&y=43.7094892 returns x=632549.24, y=4840701.20

but

SELECT ST_AsEWKT(st_transform(st_setsrid(st_makepoint(43.7094892, -79.3545864), 4326), 26717))

returns

SRID=26717;POINT(632566.193912585 4840706.93374429)

I'm at a loss what could be causing it or how to resolve...

Using PostGIS version
POSTGIS="2.5.1 r17027" [EXTENSION] PGSQL="110" GEOS="3.7.0-CAPI-1.11.0 673b9939" PROJ="Rel. 5.2.0, September 15th, 2018" GDAL="GDAL 2.3.1, released 2018/06/22" LIBXML="2.9.1" LIBJSON="0.12" LIBPROTOBUF="1.3.0" RASTER
on Amazon RDS Postgres 11.4-R1


texts: Downloading map from website such as metrolinx

texts:
https://rem.info/en#map
I want to download this map and bring it into ArcGIS Pro or ArcGIS Earth.
How can I download it?


texts:
Why is longitude also written as B and latitude as k?

I see that abbreviation when working with Google maps API.


texts: Create an 'heatmap' based off of a numeric attribute in a point feature class

texts: I have a shapefile with several points, each of which has a numeric attribute 'price'. I am trying to find areas that have the most potential revenue concentrated based off of the price in the feature class. How might I do this?


texts:
I work with ArcGIS 10.3.1. I my map i have 200 layers. some of them repeat several times. Is it possible with python code to remove all duplicated layers from table of content, that have the same group layer name?

all the layers are in different group layers- each group layer has different name except the duplicated. I know the "Delete Identical" tool- i have only basic license.

My question differ from Remove Duplicate items in legend because i seek for duplicated group layers (with the same name).


texts: Should I use ISO 3166-1 Numeric codes as Primary Keys for Countries in a database table?

texts: My main concern is, will one country ever adopt a numeric code that has been used for another in the past, such as can happen with alphanumeric 3166-1 codes? If not, then it seems like a perfect identifier for a primary key in a database...


texts:
I try to add the gxp styler plugin, to my suit sdk, but I cannot make it work perfectly.

I had to set up a digest authentication in geoserver (and the user has to log in to see any layer), and I have two users (admin, client). The styler plugin icon is inactive after log in and the browser says:

PUT http://91.109.10.190:9080/geoserver/rest/styles 401 (Unauthorized)

It will change to active only one case, when I log in to geoserver and to the sdk viewer on the same browser.

Can somebody provide me a solution what sould I change in the code or in geoserver? (probably something in the geoserver)


texts: Copy an attribute from one polygon shapefile to another with ArcGIS

texts:
I'm very new to GIS, and have been toying with something that intuitively seems simple, but I can't make it work. I have an ArcView license (small municipal government).

I have a parcel polygon shapefile and a building polygon shapefile. I'm trying to copy parcel ID numbers to the building shapefile. The buildings are all wholly contained within the parcels. I cannot seem to find a way to make this work.

Once I know how to do this (which I'm assuming will be easy and obvious once I see it), I will try to do the same for the driveway shapefile, which aren't always completely contained within the parcels.

My end goal is to come up with the area that is neither driveway nor building for about 8000 parcels.


texts: How to get specific area boundaries on a map?

texts:
I want to get specific area boundaries highlighted on the map or by coordinates. I tried to search on google maps with this structure (district, city, country), sometimes it works and give me the area boundaries like this :

and sometimes it shows the area but without highlighted boundaries. How can I do it ??


texts:
Is it possible to create a MapBasic SQL command to subtract one date column from another and get the number of days between the two (positive or negative values)?

I have found the following process in FME and ArcGIS but need a MapBasic solution

https://knowledge.safe.com/questions/4349/subtract-dates.html and
Add or subtract days to date field using the field calculator


texts: Is it possible to render a large number of shapefiles to a pdf document using Python?

texts:
I'm extremely new to GIS, and recently made a map by loading a bunch of shapefiles into QGIS using a python script. However, it took 15 mins to scroll in or out on the image, and the pdf export crashed the program (I attribute both to the number of shapefiles, it was several thousand).

Is there a way to take all the shapefiles in python and output them to a pdf document? I only want black lines and I only want one line style so the only complicated part is dealing with the sheer number of them.

Alternatively, is there a trick I missed with exporting massive pdfs in QGIS?

====================

Edit for details:

The dataset I'm trying to render is the TIGER/Line USGIS data, specifically the All Roads set. There are 3,108 shapefiles, each in a subdirectory.

I have successfully pulled all this into QGIS using a Python script, set the line styles using another little tool I found, and then exported it as a png file at 1200 dpi (it refused to go higher). However, while the png file does look cool, I am unsatisfied with the quality of the print I can create from it. When the printer prints a raster graphic, the subpixel rendering makes result look washed out and blurry.

My goal is to send a vector document to the printer so that I can print this map without the greyed out look caused by subpixel rendering in a raster graphic.

So, to phrase my question better, what is the best way to create and print a ridiculously complicated document made of a huge pile of shapefiles? I use PDF just because I'm used to abusing it with LaTeX to make cool stuff, but is there a better way? I'm happy to go as deep into the code as possible, GUI just gets in the way when there's this much data and the result is visually so simple.

I was able to create a workflow where I exported the file from QGIS and then turned it into a vector pdf with a potrace script and then printed it with an lp script (I'm on Mac) which looks amazing for my state but does not work at the scale of the whole country.


texts: Geopandas intersects/contains fails between GeoDataFrames

texts:
I am struggling to make the GeoDataFrame functions: intersects/contains work between GeoDataFrames instances.

Let t be a GeoDataFrame with contents:

    id                                          geom
0  230  POINT (154860.3254511474 -468176.1141802836)
1  231  POINT (149247.4170036312 -436634.9002944119)
2  232  POINT (192530.4468407042 -465274.3331522844)
3  233  POINT (138125.4083732415 -419237.1333624315)
4  234  POINT (195314.0952219047 -468429.8881146177)
5  235  POINT (196320.0493430012 -479448.5877129403)
6  236  POINT (199780.1489409541 -446899.1169124786)
7  237  POINT (161437.1203666388 -440485.6615666826)
8  238  POINT (192413.1253873391 -465423.3888242305)
9  239  POINT (196039.9555616449 -473696.9492107867)

I extract the 5 last points in another GeoDataFrame named prt:

    id                                          geom
5  235  POINT (196320.0493430012 -479448.5877129403)
6  236  POINT (199780.1489409541 -446899.1169124786)
7  237  POINT (161437.1203666388 -440485.6615666826)
8  238  POINT (192413.1253873391 -465423.3888242305)
9  239  POINT (196039.9555616449 -473696.9492107867)

executing the command:

t.intersects(prt)
Out:
0    False
1    False
2    False
3    False
4    False
5     True
6     True
7     True
8     True
9     True

However when I create a polygon from the points of prt and repeat the process it fails:

import shapely.geometry as sh
polygon = sh.polygon.Polygon([[p.x, p.y] for p in prt.geometry.tolist()])

import geopandas as gpd
pds = gpd.GeoDataFrame([polygon], columns=['geom'], geometry = 'geom', crs = {'init':'epsg:3309'})

pds.geom_equals(polygon)
Out[873]:
0    True
dtype: bool

t.intersects(polygon)
Out[875]:
0    False
1    False
2     True
3    False
4     True
5     True
6     True
7     True
8     True
9     True
dtype: bool

t.intersects(pds)
Out[874]:
0    False
1    False
2    False
3    False
4    False
5    False
6    False
7    False
8    False
9    False
dtype: bool

Can you please provide me some help/guidelines? I am pretty new in Python, so any advice is welcome.


texts: Areas where DGPS with post processing (not RTK) is still required

texts: I am wondering in which domains or countries where Differential GPS (DGPS) with post processing (not Real Time Kinematic [RTK]) is still required? It's known that RTK has been rapidly adopted by cadastral services in several countries, but I wonder if DGPS with post processing is still required? And why we still need it when we can accomplish the same task using RTK?