Latest Contributions by XanderBakker

I am attempting to create a series of bivariate maps of an index of habitat suitability (point data) with elevation, slope, and aspect. Because there is no attribute table for any of the data, they cannot be joined. However, I am wondering how (if it is possible) to join and map  the zonal statistics of the two datasets. I have zonal statistics of each one in the same format. Can they be joined in one file to show the unique values of both sets of data? the Join Field tool doesn't seem to work. Zonal Statistics is normally carried out on raster data. Slope, elevation and aspect are normally managed in raster format. You mention that the habitat suitability point data has elevation, slope and aspect, but no attributes. This confuses me a bit. If point data contains slope, aspect they will probably be stored as an attribute. Elevation could be stored as the z-coordinate of the point.  Is it possible to include a screenshot of the situation (table of content and data)? That might explain what you are dealing with. The output would be a bivariate map of habitat suitability and elevation, another of habitat suitability and slope, and one of habitat suitability and aspect. So you want to investigate the relationship between habitat suitability and slope, elevation and aspect. If you split up the island in the areas, would these areas be representative for a habitat? You can use zonal statistics to obtain the stats per area of the rasters (slope, aspect and elevation). Another method could be to classify slope, aspect and elevation and use the combine tool (Spatial Analyst) to create a raster with the unique combination of the three classified rasters. If you analyze which combination is suitable you have access to the areas where this combination occurs on the island. ... View more

Just because it is fun... ...using this code: import arcpy import math def main():     from collections import OrderedDict     # create ordered list of distance vs id + pop     fc = r"C:\Forum\NearEnhanced\gdb\pop_counties.gdb\counties"     fld_pop = "PST025181D"     fld_fips = "FIPS"     fips_start = "29189" # "05021"     max_pop = 5000000 # 1000000     # get start pnt and pop     pnt_start, pop_start = getPointStart(fc, fld_fips, fld_pop, fips_start)     # create ordered dict with distance, fips and pop     dct_data = createDictDistancesAndData(fc, fld_fips, fld_pop, pnt_start, fips_start)     # let's sort on distance to true centroid     dct_sorted = OrderedDict(sorted(dct_data.items(), key=lambda x: x[1][0]))     lst_fips = [fips_start]     pop_tot = pop_start     for oid, lst in dct_sorted.items():         pop = lst[1]         pop_tot += pop         if pop_tot > max_pop:             print pop_tot             break         else:             lst_fips.append(lst[2])     print "{0} in ('{1}')".format(arcpy.AddFieldDelimiters(fc, fld_fips), "','".join(lst_fips)) def getPointStart(fc, fld_fips, fld_pop, fips_start):     where = "{0} = '{1}'".format(arcpy.AddFieldDelimiters(fc, fld_fips), fips_start)     row = arcpy.da.SearchCursor(fc, ("SHAPE@TRUECENTROID",fld_pop), where_clause=where).next()     return row[0], row[1] def createDictDistancesAndData(fc, fld_fips, fld_pop, pnt_start, fips_start):     sr = arcpy.Describe(fc).spatialReference     dct_data = {}     where = "{0} <> '{1}'".format(arcpy.AddFieldDelimiters(fc, fld_fips), fips_start)     with arcpy.da.SearchCursor(fc, ("SHAPE@TRUECENTROID",fld_pop, fld_fips, "OID@"), where_clause=where) as curs:         for row in curs:             pnt = row[0]             pop = row[1]             fips = row[2]             oid = row[3]             dist = getDistanceArcpy(pnt_start, pnt, sr)             dct_data[oid] = [dist, pop, fips]     return dct_data def getDistance(pnt_start, pnt):     return math.hypot(pnt[0] - pnt_start[0], pnt[1] - pnt_start[1]) def getDistanceArcpy(pnt1, pnt2, sr):     pnt_geo1 = arcpy.PointGeometry(arcpy.Point(pnt1[0], pnt1[1]), sr)     pnt_geo2 = arcpy.PointGeometry(arcpy.Point(pnt2[0], pnt2[1]), sr)     return pnt_geo1.distanceTo(pnt_geo2) if __name__ == '__main__':     main() The example is using county with fips 29189 and a max population of 5 million. Change lines 7 to 9 to indicate the featureclass with the counties and the names of the fips and populations fields. Change the settings on line 11 and 12 to indicate the fips to start and the maximum populations. The result is a print of the where clause you can apply to show the counties included in the selection. The code can be adapted to create a feature class with the counties. The code stops with the value before it passes the maximum population. This happens on line 27 to 29. ... View more