Latest Contributions by NeilAyres

Welcome back Curtis, It's good to see that the crisis has been averted and that we can again get your insightful contributions here. Cheers, Neil ... View more

If I am understanding you correctly... When you add a data layer to a new blank map document, the data frame will adopt the coordinate system of the first data layer. For the online base maps, this is web Mercator. Then, if you want to be in a different coordinate system, simply right click on the data frame, select the coordinate system tab and select a new coord sys. All data added to the frame will be reprojected on-the-fly into whatever system the data frame is using. The rider to this is when the input data has "no projection" or has an incorrectly specified system. But that's another debate. Cheers, Neil ... View more

Aside from the tiff renaming problems pointed out by others.. Note you are going from PSAD_1956_UTM_Zone_19S (based on D_Provisional_S_American_1956) & WKID 24879 to Peru96_UTM_Zone_19S (based on D_SIRGAS-Chile) & WKID 5389. So the coordinate system is not changing but the datum is. You need to specify a datum transformation for this operation, otherwise your output and input will be no different. As I do not have any data available which is in either of these systems, I cannot comment here as to which transformation will be appropriate for your data. But I suggest that you do some by hand tests first and get the transformation issue properly sorted out before you write a script to reproject all of your imagery. Cheers, Neil ... View more

Why not use the Calculate Geometry option on the column, then select Z. Cheers, Neil ... View more

jp, Looking at the python code : try:     # enable the overwrite option...     arcpy.env.overwriteOutput = 1             # check out spatial analyst license...     arcpy.CheckOutExtension("spatial")                         # script parameters...     landcover = sys.argv[1]         # input raster dataset     edgeWidth = float(sys.argv[2])         # number     wksp = sys.argv[3]              # output raster location     FragMap_lyr = sys.argv[4]       # display layer ##    landcover = r"C:\Projects\ForestFragModel\test_landcover\lc_forest_nonforest.img" ##    edgeWidth = 305        # string ##    wksp = r"C:\Temp"              # output raster location ##    FragMap_lyr = "ff_test"     s = time.clock()     #-----------------------------------------------------------     # CREATE OUTPUT RASTER NAME FROM WORKSPACE AND LAYER NAME...                 # if workspace is a geodatabase, no extension needed...     if ".gdb" in wksp or ".mdb" in wksp:         FragMap = "%s\\%s" % (wksp,FragMap_lyr)     # otherwise use ".img" extension (no length / character restrictions)...     else:            FragMap = "%s\\%s.img" % (wksp,FragMap_lyr)     # provide user feedback...     arcpy.AddMessage("\nOUTPUT DATASET NAME IS:")     print "\nOUTPUT DATASET NAME IS:"     arcpy.AddMessage("\n\t%s\n" % FragMap)     print "\n\t%s\n" % FragMap This is the first part, taking in the inputs and setting the outputs. The pathing is not handled very elegantly imho, but perhaps you are selecting a dataset rather than a workspace (sys.argv(3)). I see the your successful run you used an *.img image file instead of outputting to a geodatabase. Perhaps try that. Otherwise, a general point, I would not have any spaces in the path or the filenames. Never used this tool, so that's about it from me. Good luck, Neil ... View more

Running script Landscape_Fragmentation_Tool... OUTPUT DATASET NAME IS:     C:\Users\Beth\Documents\ArcGIS\lftmod\landfragmod.gdb\C:\Users\Beth\Documents\ArcGIS\lftmod\landfragmod.gdb\lfttestsat Well you can see that this path is invalid. It looks like the prog is not resolving the difference between the dataset name "lfttestsat" and the path the the gdb, but is concatenating them together. Never used the Landscape Fragmentation tool. But if it worked before this, then what are you now doing differently? In this example you are storing your outputs in a fgdb. Previously were you using shapefiles perhaps? Good luck, Neil ... View more

And as an additional tip when it comes to extracting stream lines from dem's... I always create a curvature surface and use that as a weighting into the flow accumulation step. Because of the math of the curvature function, upwardly concave areas (like a stream) are negative, so invert the curvature output first, (X -1)  so that +ve values are associated with upwardly concave areas and -ve values indicate the top of the hill (for example). This idea comes from a paper by Tarboten & Ames from about 1992. Otherwise follow the steps in the hydrology section of the help files as suggested by Jim. Cheers, Neil ... View more

I would stay with rasters and spatial analyst. But extracting geological features automatically from a dem might be quite tricky. But you must explore the various filters and other surface operators to enhance abrupt changes of slope. sa.Curvature might be one example. There are high and low pass filters in there as well. But again, sorting out which feature is a "lineament" from say a stream line, will be a challenge. If you want a good book on analyzing this sort of stuff my bible is : "The Image Processing Handbook" by John C. Russ. Highly recommended. Good luck, Neil ... View more

Perhaps something like this : sysFieldList = [] attFields = arcpy.ListFields(fc_name) for attField in attFields:     if attField.required:         sysFieldList.append(attField) print "System Fields", sysFieldList Cheers, Neil ... View more

rai = arcpy.mapping.ListLayers(mxd, "Railroads", df)[0] EWmove = float(arcpy.GetParameter(0)) NSmove = float(arcpy.GetParameter(1)) def shift_features(in_features, EWmove=None, NSmove=None): in_features = rai point.X += EWmove point.Y += NSmove with arcpy.da.UpdateCursor(rai, ['SHAPE@XY']) as cursor:   for row in cursor:    cursor.updateRow([[row[0][0] + (EWmove or 0), row[0][1] + (NSmove or 0)]]) return arcpy.RefreshActiveView() shouldn't this be... rai = arcpy.mapping.ListLayers(mxd, "Railroads", df)[0] EWmove = float(arcpy.GetParameter(0)) NSmove = float(arcpy.GetParameter(1)) def shift_features(in_features, EWmove=None, NSmove=None):     rai = in_features         with arcpy.da.UpdateCursor(rai, ['SHAPE@XY']) as cursor:             for row in cursor:                 cursor.updateRow([[row[0][0] + (EWmove), row[0][1] + (NSmove)]]) return arcpy.RefreshActiveView() Cheers, Neil ... View more

You don't actually have to update each vertex according to the Café Python page.. http://arcpy.wordpress.com/page/2/ just update the SHAPE@XY token instead, then the whole feature moves. like : with arcpy.da.UpdateCursor(in_features, ['SHAPE@XY']) as cursor:     for row in cursor:         cursor.updateRow([[row[0][0] + (x_shift or 0), row[0][1] + (y_shift or 0)]]) Cheers, Neil ... View more

The append tool allows you to input a whole bunch of feature at once into the append routine. If your variable shpList is a real python list like ["Shape1.shp", "Shape2.shp", ...] Then just pass the entire list into Append_management instead of doing it one at a time. Cheers, Neil ... View more

The geometry methods are very useful. like : dist = geom1.distanceTo(geom2) See the Café Python page for some useful tips. http://arcpy.wordpress.com/ Cheers, Neil ... View more