Latest Contributions by KimOllivier

Since writing my solution I have found an ArcView 3.x solution that shows I have not covered all the cases. You could have a line running through a pair of vertices on a polygon with rentrant sides that extent out each way to the polygon boundary between vertices. I assumed that the longest line would always be coincident with vertices at the ends. http://www.jennessent.com/arcview/longest_lines.htm Well done to Jeff Jenness for handling this case in Avenue. ... View more

It is really annoying to have text fields 255 chars wide containing 10 chars which result from loading in spreadsheets or CSV files etc. I have written a scanner in awk that counts the maximum width of the contents of each field so you can round them up and resize the fields. The obvious way is to use FME (or the Data Interoperability Extension) to reprocess the inputs to standardise the field widths. You can also reorder the fields easily. You can set the field widths in a SCHEMA.INI (see Microsoft help) before you load text files. If the sources are shapefiles then it is possible to copy the file using Conversion/FeatureclassToFeatureclass with a Field Mapping to redefine the output field widths. This is harder to program. Note that in really wide fields you can change the height of the rows in the table and the contents will wrap! I only just discovered this after 10 years of using tables. The trouble with formatting tables is that its all manual settings and does not persist when you close the table. ... View more

To save you writing a lot of scripts have a look at the ET Tools third party scripts. They have extensive tools to clean up dangling lines. ET_GeoWizards. There are some new snapping tools at 10.2 in the Editing Tools that do not use the Coverage Tools, but you still need a Standard licence. ... View more

Move all the spare points to a new featureclass. It doesn't sound like a clean dataset to be building nodes from a cluster of points. You could use Spatial Statistics>Utilities>Collect Events tool to remove the duplicates and just keep one for the network. Maybe have a relate to a table so that you can have multiple attributes for the other points? Why do you need all the points in the same featureclass? Selecting one of many points to keep is a tricky problem not available in the standard tools, but you can make a tool to create a unique ID for a location using a Geohash. Or use ET Tools, a great set that fills in all sorts of little gaps in the Esri toolset. Once you have a unique location ID you can number the points within the group and offset them a little. Needs a bit of Python programming to extend the toolset for adventurous analysts like yourself. ... View more

Duplicate points are not an error, but duplicate vertices are. However you can only have vertices in a polyline or polygon, so there is nothing wrong with your sample. To find coincident points you need another tool. It is not a useful thing to have any featureclass using play coordinates in the range 0 - 10 to test this sort of thing because of rounding errors during storage using integers and offsets. It is recommended that you always have a coordinate system defined and use real coordinates such as UTM with lots of figures if you are going to stress test Esri algorithms for bugs. There are lots of issues with the interpretation of valid geometry between systems that check geometry does not report. See FME for a more comprehensive testing tool. Something that is valid in Esri may not be allowed in OSGeo simple features. For example a vertex may be allowed to touch the edge of the polygon or not. Why not repair geometry instead of just checking it? The main issue that repair geometry handles is polygon vertex order is back to front in KML sourced polygons that need to be reversed to be Proper Polygons in ArcGIS. ... View more

Each time you upgrade ArcGIS a new folder is created in your login directory so you lose your old personal toolboxes until you move them to the new directory. ... View more

Why can't you join the layers? It does seem a bit of a lost cause if you can't at least do that. However some joins may be very inefficient if there is no index and you do not have permission to build it. Copy the data so that you can, or extract the data to enable the comparision. Don't think 'loop through' if at all possible, use SQL type queries on whole tables for efficiency. If the data is small enough to fit into memory there are Python structures to make comparisons fast. For example you could use a cursor to read in the fields from each table into a list and a dictionary and run a set comparision to get a list of differences. If that is small compared to the total, then there is a huge benefit. ... View more

Instead of drawing a graphic to use as a clip, use a tool to create a featureset which is created interactively just like a graphic when setting up parameters to run the tool. If you have hundreds of graphics that should have been features, there is a menu item to turn graphics to features. You can then select a real feature from a featureclass to use it as a clip. Just don't use graphics as features, they are really most suitable for illustrations of non GIS data in a layout. ... View more

I was intending to explain the only obscure part of your question as best as I could guess from your description. Getting an interactive point using a standard ArcTool is easy with a custom tool that you would need to write. I assume you were intending to write a script for the tool yourself? I can only guess that you want a volume from the seabed to a water surface. Here is a sample tool, the script that reads the parameters, and images before and afterwards. You will have to create your own template point featureclass to use as a template to indicate that you want points collected, not polylines or polygons. #------------------------------------------------------------------------------- # Name:        VolumeCalculator.py # Purpose:     Calculate a volume under an arbitary rectangle # # Author:      kimo # # Created:     05/06/2014 # Copyright:   (c) kimo 2014 # Licence:     Creative Commons 3.0 NZ #------------------------------------------------------------------------------- import arcpy import math fds = arcpy.GetParameter(0) width = arcpy.GetParameter(1) length = arcpy.GetParameter(2) angle = arcpy.GetParameter(3) raster = arcpy.GetParameter(4) tin = arcpy.GetParameter(5) # the fds is an in_memory featureclass, get the first point with arcpy.da.SearchCursor(fds,["SHAPE@X","SHAPE@Y"]) as cur:     row = cur.next()     x = row[0]     y = row[1] # feedback for demo msg = "{} {} {} {} {} {} {}".format(x,y,width,length,angle,raster,tin) arcpy.AddMessage(msg) # calculate volume.... ... View more

You need to add a Node layer that is now at ArcGIS a point layer with references added to each line to link the node number. This was built into ARC/INFO so it was really easy. Now you need a third party tool such as Renode Polylines to make this easy such as  in  ETGeotools.  http://www.ian-ko.com/. Just adding points to the ends does not add a single point at the ends, you get two unrelated points. Sigh, this takes a lot of extra steps to cleanup. When you have properly related point-nodes and polylines you can select all the unsealed polylines and tag the nodes with an 'unsealed end' tag. Now you can select all sealed segments to see if both end nodes have an 'unsealed' tag. Each road segment will have a field FnodeID and TnodeID. The node points will have a NodeID. Segments that share a node will have the same id. At T junctions the logic is a bit more complex. You will have to tag those as 'unsealed' if all roads are unsealed connecting it. Another useful tag on the nodes is a valence (the count of connectors) to help identify these by running a frequency on the segment end IDs. All of these steps can be done with an SQL query or a join and calculate, don't step through each segment with a cursor! ... View more

You can control the table of contents in ArcMap using arcpy.mapping module commands. So you would need to create a tool that toggled the layer and then add the tool to ModelBuilder. You need to build the tool with input and output parameters. ... View more

You don't have to wait over the weekend to know that you have not approached the problem the best way. I have a rule of thumb called the "Cup of Coffee Rule". If an single process has not finished by the time I have finished my coffee, then I interrupt the process and find a better way. It must be something seriously wrong since computers are now so fast. Maybe a missing index, working across a network instead of local data, exceeding the size that can fit in memory. Just using a tool with inappropriate tolerances or settings can turn a simple process into something impossible. In this case I can see the problem at once - you are running geoprocessing tools inside a cursor! They are not designed that way, you are supposed to frame your problem to run a tool over the whole dataset in one pass. Never in a cursor. It is just too slow to restart a process and cleanup. Since the tools are not designed that way, it is likely that garbage collection is not being worked either so you will run out of memory long before you get any useful results. I have never run a tool inside a cursor, there is always a way of avoiding it. I know it is easy to visualise a problem by considering a feature and then applying the same thing a million times, but this is not the best way usually. Think of the way that SQL queries work. You specify what is to be done with a pair of tables and you let the database engine decide how to do it, you don't even have the luxury of a cursor in the SQL language. It would be nice if GIS tools had a similar command to SQL called Explain which shows how things are going to be solved. Instead you have to think it through yourself and recognise the patterns that will not work. You need to multiply the numbers of operations together by the expected time to see if it will finish before you retire. Any 40,000 updates should finish in a few seconds or minutes, to give you a feel for what can be achieved. I know this has not reworked the problem and supplied you with a neat alternative. While I am reviewing your code, I suggest that you don't use CalculateField in a script. It is better left for ModelBuilder. Since it is just wrapping a cursor around an expression, it would be better to just use an expression in the cursor. It is much easier to understand, can be debugged more easily and you can trap unexpected data and recover with if statements. It would also be much faster. If you are using shapefiles, consider using file geodatabases. Then the  area is automatic, you don't have to calculate it. ... View more

You can do this with the standard toolbox using the featureset input type. This will capture one or more mouse clicks to draw an interactive graphic of whatever type template you have used for the featureset. It will also have fields that you can fill in the properties to be included in the new feature. In your example you would have a point feature and some properties, maybe filled with defaults. When you execute the tool you capture the input parameter as a geometry object to either convert to a featureclass or process further. In this case you would just get out the XY of the point and calculate the boundaries of a polygon to create a polygon geometry object for output. But maybe you would just create the polygon directly with a polygon template. http://resources.arcgis.com/en/help/main/10.1/index.html#//002w00000023000000 You could make a really slick interface that picked a coordinate using python AddIns, but that does not easily allow a form to input additional information like a featureset. Addins are much harder to debug. ... View more

At ArcGIS 10.2.1 there are new python modules included called XLRD and XLWT which are excellent modules for reading spreadsheets and extracting the data. It can handle multiple sheets and can pick out any cell or range of cells. Don't bother trying to use the new Excel to Table tool unless the data is already in a suitable format with valid field names etc. Here is an example of reading the 4th sheet that was not formatted for a database, no header field names and lots of comments all over the sheet. I needed three columns in the middle of the sheet. By the way watch for new uses of Sqlite with spatial extensions as Spatialite and the OSGEO recent standard http://www.opengeospatial.org/standards/geopackage #------------------------------------------------------------------------------- # Name:        loadXLS #               read Table 4 from 68 spreadsheets to get Unoccupied Dwellings # Purpose:      read Stats spreadsheets into a Sqlite database that cannot #               be read by FME, ArcGIS or anything other than Excel #               saves writing out to CSV first # Author:      kimo # # Created:     08/05/2014 # Copyright:   (c) kimo 2014 # Licence:     Creative Commons #------------------------------------------------------------------------------- import xlrd # installed with ArcGIS 10.2 by default import sys,os import re import glob import sqlite3 import arcpy # only needed for final message pat = re.compile("^[0-9]{6}") # pattern of AU code in AU name in Column 0 ("A") dType = {0:'Empty', 1:'Text', 2:'Number', 3:'Date', 4:'Boolean', 5:'Error', 6:'Blank'} def extract_tla(filename):     '''extract Unoccupied Dwellings'''     workbook = xlrd.open_workbook(filename)     worksheet = workbook.sheet_by_name('Table 4')     num_rows = worksheet.nrows     num_cells = worksheet.ncols     rec = 0     for curr_row in range(num_rows):         curr_col = 0         cellAn = worksheet.cell_value(curr_row, curr_col)         cell_type0 = worksheet.cell_type(curr_row, curr_col)         if cell_type0 == 1 and pat.match(cellAn):             a = str(worksheet.cell_value(curr_row,curr_col))             val = [a[0:6],a[7:]]             for curr_col in range(1,4):                 count = worksheet.cell_value(curr_row,curr_col)                 if count == '-':                     n = 0                 else :                     n = int(count)                 val.append(n)             rec+=1             #print val             insert_ud_row(udTable,val)     # print "Records found",rec     # fields au2013,auname,ud2001,ud2006,ud2013     return rec def create_ud_tab(conn,udTable):     # (re-)create a table     c = conn.cursor()     c.execute('''DROP TABLE IF EXISTS '''+ udTable)     c.execute('''CREATE TABLE '''+udTable+''' (         au       text NOT NULL,         auname       text ,         ud2001       integer,         ud2006       integer,         ud2013       integer )''')     conn.commit()     c.close()     return def insert_ud_row(outTab,row):     '''insert a row     requires a sqlite database to be connected     and a cursor to be open     '''     try:         cmd = "INSERT INTO " + outTab + " VALUES (" + (",?"*(len(row)))[1:] + ")"         # implied transaction already opened         c.execute(cmd,row)     except Exception,msg:         print msg         print row         print cmd         conn.rollback()         conn.close()         sys.exit() # ------------------------------------------------------------------------- if __name__ == '__main__':     try:         source_folder = sys.arg[1]     except:         source_folder = 'C:/data/census2013/source'     os.chdir(source_folder)     conn = sqlite3.connect('C:/data/census2013/dwelling.sqlite')     udTable = 'DwellingEmptyAU'     create_ud_tab(conn,udTable)     # Open a cursor and insert records from all the spreadsheets     c = conn.cursor()     print "sqlite database successfully opened"     trec = 0     lstXLS = glob.glob('*.xls')     for x in lstXLS:         rec = extract_tla(x)         trec+=rec         print rec,x     conn.commit()     c.close()     conn.close()     msg = '{} Total records inserted into {}'.format(trec,udTable)     print msg     arcpy.AddMessage(msg) ... View more