Latest Contributions by KimOllivier

Worked fine for me, with a small dataset. Maybe your environments are different. Scratch workspaces, swap space, full disk, default workspace not a filegeodatabase... ... View more

I have played around a little and have discovered that to delete a FeatureLayer you need to name the type properly, ie "FeatureLayer", not "Layer" import arcpy arcpy.env.workspace = "d:/workspace" fc = "addmar06.shp" if not arcpy.Exists(fc):     print fc,"not found" arcpy.env.overwriteOutput = True # so I can rerun to recreate a layer # make an in-memory layer result = arcpy.management.MakeFeatureLayer(fc,"fc_layer") print "this is the string-name of the fc:",result.getOutput(0) # but it is not the object # find the datatype of the layer with the name of "fc_layer" desc = arcpy.Describe("fc_layer") print "Datatype is:",desc.datatype # now delete the featurelayer using the proper datatype try:     result2 = arcpy.management.Delete("fc_layer","FeatureLayer")     print "Has it gone?",result2.getOutput(0)     # prove it has gone     arcpy.Describe("fc_layer") except Exception,errmsg:     print "Error because:",errmsg     # but fc still intact     if arcpy.Exists(fc):         print "Yes, fc still there" >>> this is the string-name of the fc: fc_layer Datatype is: FeatureLayer Has it gone? true Error because: "fc_layer" does not exist Yes, fc still there I also found that if you used arcpy.management.Delete(fc,"ShapeFile") it deleted it for me even with a Layer defined. ... View more

When you use a cursor to extract the shape outShape = row.shape it is not a geometry object unfortunately. You have to unpack the parts and then unpack the points inside the parts. Then you can rebuild a geometry object. Here is my FillDonut script as an example. (I have since found that it does not handle all cases of donuts, but I have not fixed it yet.) # fillDonut10.py # frequent forum request to fill holes # use an Update Cursor to edit polys with donuts # 22 May 2010 # Kim Ollivier [email protected] # (c) Creative Commons 3.0 (Attribution) # Method # ------ # Polygon features are stored as: # [[pnt,pnt,pnt,pnt,,pnt,pnt,pnt,pnt][pnt,pnt,pnt,pnt]] # Parts are separate lists inside the feature object array # null point in a part indicates donut ring(s) follow # So just step through and truncate at first null pnt per part, #   re-assemble parts into a polygon # Update feature if a donut found, otherwise continue # Input is a layer so can use selection in Arcmap # could also put on menu at ArcMap 10 # can be back-ported to 9.3 # 4 June 2010 # added iter lambda function for part array because isn't interable # cleared polyOuter buffer earlier to fix logic for multiple parts # 7 Nov 2011 # handle double donut parts, if a second+ parts is enclosed # then have to drop them as well. import arcpy,sys def filldonut(inlay) :     '''Edits a layer in-place     fills all features with donuts     requires a layer to honour     selected features     '''     desc = arcpy.Describe(inlay)     shapefield = desc.ShapeFieldName     rows = arcpy.UpdateCursor(inlay)     n = 0     m = 0     polyGeo = arcpy.Array() # to hold edited shape     polyOuter = arcpy.Array() # to hold outer ring     for row in rows:         feat = row.getValue(shapefield)         qInterior = False         id = row.objectid         msg1 = "%d has %d parts" % (row.objectid,feat.partCount)         arcpy.AddMessage(msg1)         for partNum in range(feat.partCount) :             part = feat.getPart(partNum)             for pt in iter(lambda:part.next(),None) : # iter stops on null pt                 polyOuter.append(pt)             if part.count > len(polyOuter) :                 qInterior = True             g = arcpy.Polygon(polyOuter)             if partNum == 0:                 polyGeo.append(polyOuter) # reassemble first part                 gOuter = arcpy.Polygon(polyGeo) # make a geometry object             else :                 if g.disjoint(gOuter) :                     polyGeo.append(polyOuter) # reassemble separate part                 else :                     m+=1                     msg2 = "%d has island at part%d" % (row.objectid,partNum)                     arcpy.AddMessage(msg2)             polyOuter.removeAll() # ready for next part         if qInterior : # in any part of this feature, otherwise skip update             n+=1             row.setValue(shapefield,polyGeo)             rows.updateRow(row)         polyGeo.removeAll()     del rows,row     msg = "Features with interior ring filled %d and %d with island donuts" % (n,m)     arcpy.AddMessage(msg)     #---------------------------------------------------- if __name__ == '__main__' :     inlay = sys.argv[1] # arcpy.GetParameterAsText(0)     desc = arcpy.Describe(inlay)     if desc.shapetype != 'Polygon' :         arcpy.AddError(inlay + " Not a polygon class")     else :         filldonut(inlay)         arcpy.RefreshActiveView() ... View more

There is a lot of help in the Python documentation on the Decimal module. Numbers are stored as strings and converted to a Decimal datatype. Basically you avoid binary arithmetic by storing all numbers as decimal strings and using the Decimal datatype which has some basic arithmetic operations. >>> import decimal >>> a = decimal.Decimal("1234.12") >>> be = decimal.Decimal("2.0") >>> a / be Decimal('617.06') >>> x = 1234.12 >>> y = 2.0 >>> x/y 617.05999999999995 If you use integers offset and scaled by the number of decimal places you want to keep, then you can do integer arithmetic and round off correctly. This is how Esri handles floating point coordinates in SDE because databases do not have a true floating point type, or it is too slow. The disadvantage is that there are limits to the range (extents) if you only use 32 bit integers and there is still some rounding or 'creep' after processing. With 'high precision' 64 bit it is no longer an issue for users. For more details look up "Choosing the resolution and domain of a low-precision spatial reference" in the help for obsolete 9.1 formats to see how you will need to set up your scale and precision to do your own arithmetic with integers and convert back and forth. The Decimal module looks more suitable. ... View more

You need to install a help authoring tool. These are easy to use, and will compile the help file in CHM and web help formats, plus my chosen one will also output a Word Document or PDF. I used to use HelpBreeze, but it has fallen behind so I switched to HelpSmith which could import all my old help files. HelpBreeze is cheaper but a security upgrade in Windows prevents saving edits and HB have not fixed it. HelpSmith cost $200 but since I purchased it I now write a proper integrated help file for every project I am involved with. I now start early with a specification and mockup of the screens right into the help file and then use it to keep track of specification changes. I find it easy to link the geoprocessing tools to topics in the help files. The new default help editing tool at 10.x are now so difficult to use I now just get on with a HelpSmith project at the start. It has many other benefits, the files can be stand-alone, incorporated in other web resources and are easy to maintain in a single document, whereas the ArcGIS help is separate for each tool. www.helpsmith.com ... View more

OK, well then just use the built-in rounding functions in Python and do the calculation using an UpdateCursor. If the round values are a result of floating point storage rather than proper variability then you could use integers. If you want to do real decimal arithmetic there is a Decimal class in the decimal module. The way to store the data would be to use integers and an implied decimal point. Then you use integer arithmetic or convert the integers to decimal type to do the arithmetic and then convert them back again to store. ... View more

Have you indexed everything? Spatial indexes on shapefiles are not automatic, besides you should have moved decisively to filegeodatabases for analysis like this. File geodatabases will be much faster. Where are you scratch folders? They might be different in ArcMap because environments are different. Do you have a scratch filegeodatabase defined? Not just a temp folder. Is your default temp folders on a file server instead of local drive. Some tools are slower in Python, who knows why, but there are ways of making it faster in both systems. ... View more

All you have to do is change the display to round to one decimal place. Leave the values as a double type. You can right-click on the table column to alter the number of decimal places shown, and the values will be nicely rounded. As others point out floating numbers are always approximately converted from binary to decimal except for some edge cases where they have exact equivalents. Some databases support a decimal type where the digits are stored as characters and decimal arithmetic is performed. Only needed for traditional money calculations. I think its a bit odd because even those are rounded if you allow division, but I suppose accountants don't do division! They don't even allow negative numbers. ... View more

I think that you have broken the unwritten rule: Never name any folder, file, table, script, featureclass, field or variable starting with a digit. Sure to end in tears, especially in raster dataset manipulation. ... View more

Temporary tables are better written to a local temporary file geodatabase. Any analysis output written across a network is subject to timeouts, so avoid it and just copy the final result across. Besides it will be much faster. Even faster if you use" in_memory", a special virtual disk location, for small featureclasses. Set arcpy.env.overwriteOutput = True so that you can rerun the script without tripping over the first aborted run. Even then, if you have a crash inside Pythonwin you may need to cleanup and close open files. This was a bit easier with the gp object, with del gp, but now you might have to reopen the GUI, or even log off to cleanup. ... View more

Because you did not surround your script with CODE tags I cannot see if your indenting was correct. You could also use a slightly more compact for-loop to avoid forgetting the next() at the end of the steps. I have changed 'rows' to 'cur' so that it is not confused with 'row' You may be able to use the strip() function instead of split() if all you are doing is removing surrounding whitespace. This is more general than splitting by a space. import arcpy featureClass = arcpy.GetParameterAsText(0) cur = arcpy.UpdateCursor(featureClass) for row in cur:     subtypeValue = row.SUBTYPE     subtypeValueIndex = subtypeValue.split(" ")     row.SUBTYPE = subtypeValueIndex[0]     materialValue = row.STRMAT     strmatValueIndex = materialValue.split("-")     row.STRMAT = strmatValueIndex[0]     photoynValue = row.PHOTOYN     photoynValueIndex = photoynValue.split("-")     row.PHOTOYN = photoynValueIndex[0]     conditionValue = row.CONDITION     conditionValueIndex = conditionValue.split("-")     row.CONDITION = conditionValueIndex[0]     clidmatValue = row.CLIDMAT     clidmatValueIndex = clidmatValue.split("-")     row.CLIDMAT = clidmatValueIndex[0]     accesstypeValue = row.ACCESSTYPE     accesstypeValueIndex = accesstypeValue.split("-")     row.ACCESSTYPE = accesstypeValueIndex[0]     groundtypeValue = row.GROUNDTYPE     groundtypeValueIndex = groundtypeValue.split("-")     row.GROUNDTYPE = groundtypeValueIndex[0]     stepsValue = row.STEPS     stepsValueIndex = stepsValue.split("-")     row.STEPS = stepsValueIndex[0]     skimmerValue = row.SKIMMER     skimmerValueIndex = skimmerValue.split("-")     row.SKIMMER = skimmerValueIndex[0]     lcyclestatValue = row.LCYCLESTAT     lcyclestatValueIndex = lcyclestatValue.split("-")     row.LCYCLESTAT = lcyclestatValueIndex[0]     collectbyValue = row.COLLECTBY     collectbyValueIndex = collectbyValue.split("-")     row.COLLECTBY = collectbyValueIndex[0]     datasourceValue = row.DATASOURCE     datasourceValueIndex = datasourceValue.split("-")     row.DATASOURCE = datasourceValueIndex[0]     erosionValue = row.EROSION     erosionValueIndex = erosionValue.split("-")     row.EROSION = erosionValueIndex[0]     cur.updateRow(row) del row, rows ... View more

You would think that this would be a property of arcpy.Describe(fc) 8-) Oddly enough it is displayed as a column available in ArcCatalog so it must be possible to retrieve using ArcObjects The only place that dates seem to be accessible to a Python script is by parsing the metadata, but that seems to only be dates for metadata updates or geoprocessing tools, not all edits. Again you would think... that this would be a critical piece of metadata that would be in the automatic update function. I see that if you use ArcCatalog there is a date modified showing on the featureclass properties which corresponds to one of the encrypted files in the filegeodatabase folder. So that is the modifled date, it is just not available to us plebs. A scan of these tables will give a number of dates, so we might be able to get a range which is not much use. At 10.1 there is a new system for logging all edits per feature with a user and datetime, but you would have to turn it on. Meanwhile we need a 'date modified' hack. Perhaps the new API for file geodatabases will give up the value? ... View more

The SQL query using a IN (list) is really successful and is very fast, avoiding having to loop over each item. The trick is to format the SQL as a string and surround it with round brackets " ID IN (1,2,3)" If the item is a string then you have to surround each element with single quotes "CODE IN ('1','2','3')" [Note that if you use "normal" field names without spaces special characters etc etc then you do not need to surround the field name with double quotes, making expression assembly easier] I recommend that you create the expression as a string and print it out to check for debugging. Then the error message has an easy feedback and you will see your logic error easily. # NearFIDlst = [33950, 7700, 4340] sqlQuery = "FID in (" + ','.join([str(x) for x in NearFIDlst]) +")" print sqlQuery codelst = ["abc","xyz","123"] sqlQuery2 = "CODE in (" + ','.join(["'"+x+"'" for x in codelst]) +")" print sqlQuery2 FID in (33950,7700,4340) CODE in ('abc','xyz','123') ... View more

Here is a similar problem I had a long time ago. I needed to draw a right angled line across the road at each drain (cunvert) point that was recorded in an asset management system. ArcIMS could not handle linear referencing. I got the angles by using MakeRouteEventLayer, then it was just a matter of creating a new line. It could be a lot more elegant at 10.x but the maths is the same. I wrote out to a file and used the sample/CreateFeaturesFromTextFile, now deprecated and removed. # culvert.py # create culvert lines from point events # must have loc_angle and length_m fields # default length_m is 10 metres # angle is from event wizard which is parallel to road # which must be wrong at 9.1 # Input # arg1 point event layer # arg2 loc_angle field # arg3 length_m field # arg4 output featureclass # projection from event layer # Kim Ollivier # 26 September 2006 import sys, string, os, win32com.client, math def newculvert(strID,strX,strY,strAngle,strSize) :     """     Create a line array at right angles to angle     and add to new featureclass already defined     """     id = int(strID)     x0 = float(strX)     y0 = float(strY)     radAngle = float(strAngle) / 180.0 * 3.14159     size = float(strSize) / 2.0     if size == 0 :         size = 10     x1 = x0 - size * math.cos(radAngle)     y1 = y0 - size * math.sin(radAngle)     x2 = x0 + size * math.cos(radAngle)     y2 = y0 + size * math.sin(radAngle)             # print "%8d %12.3f %12.3f %8.3f %4.1f" % (id,x0,y0,angle,size)     # print "            %12.3f %12.3f %12.3f %12.3f" % (dx1,dy1,dx2,dy2)     print >> f1,"%d 0" % id     print >> f1,"0 %f %f 0.0 0.0" % (x1,y1)     print >> f1,"1 %f %f 0.0 0.0" % (x2,y2)     return True #----------------- main ----------------- gp = win32com.client.Dispatch("esriGeoprocessing.GpDispatch.1") if len(sys.argv) != 6 :     print "Need 6 parameters, setting defaults for debugging"     strInputFeatureClass = "e:/council/clutha/shape/culvert3.shp"     strIDField = "DRAIN_ID"     strAngleField = "LOC_ANGLE"     strSizeField = "LENGTH_M"     strOutputFeatureClass = "e:/council/clutha/shape/culvert_ln.shp" else :     strInputFeatureClass = sys.argv[1]     strIDField = sys.argv[2]     strAngleField = sys.argv[3]     strSizeField = sys.argv[4]     strOutputFeatureClass = sys.argv[5] if not os.path.exists(strInputFeatureClass) :     print strInputFeatureClass,"not found"     sys.exit() # set the directory and workspace to the source shapefile location os.chdir(os.path.dirname(strInputFeatureClass)) gp.Workspace = os.path.dirname(strInputFeatureClass) f1 = open("culvert.txt","w") print >> f1,"Polyline" # Find the name of the shape field using Describe. strShapeFieldName = gp.Describe(strInputFeatureClass).ShapeFieldName # Check other fields are correct eFlds = gp.ListFields(strInputFeatureClass) fld = eFlds.Next() fldOk = False print gp.GetCount_management(strInputFeatureClass),"gp count" while fld :     if fld.Name.upper() == strIDField.upper() :         fldIDOk = True     if fld.Name.upper() == strAngleField.upper() :         fldAngleOk = True     if fld.Name.upper() == strSizeField.upper() :         fldSizeOk = True     fld = eFlds.Next()     if not fldIDOk:     strErr = "No ID field found "+strIDField     gp.AddError(strErr)     print strErr else :     print strIDField," found" if not fldAngleOk:     strErr = "No Angle field found "+strAngleField     gp.AddError(strErr)     print strErr else :     print strAngleField," found" if not fldSizeOk:     strErr = "No Size field found "+strSizeField     gp.AddError(strErr)     print strErr else :     print strSizeField," found" if not (fldIDOk or fldAngleOk or fldSizeOk) :     print "Fatal error in field names"     sys.exit(2)     nCulvert = 0 # open a Search Cursor objRowList = gp.SearchCursor(strInputFeatureClass) objRow = objRowList.Next() # try :     while objRow:         objGeometry = objRow.GetValue(strShapeFieldName)         objPart = objGeometry.GetPart(0)         strX = objPart.X         strY = objPart.Y         strID = objRow.getValue(strIDField)         strAngle = objRow.GetValue(strAngleField)         strSize = objRow.GetValue(strSizeField)                 # create a new line feature         ok = newculvert(strID,strX,strY,strAngle,strSize)         if ok :             nCulvert += 1         # cleanup or looks like a memory leak         del objRow, objPart, objGeometry         objRow = objRowList.Next() except Exception, e :     print e     gp.AddError(e)     del objRow, objPart, objGeometry     print "Error in Search Cursor"     print gp.GetMessages()     gp.AddError("Error in Search Cursor")     sys.exit(2) print >> f1,"END" f1.close() print nCulvert,"culverts total" gp.AddMessage(str(nCulvert)+" culverts found") # create a shapefile gp.OverwriteOutput = 1 test = "culvert_ln.shp" gp.CreateFeaturesFromTextFile("culvert.txt",".",test,"") print "Created new shapefile",test #================================================== ... View more