I though this was an interesting problem, and, waffle-head that I am, went at it with an entirely raster approach.
- This conveniently avoids issues of topology!
- There's a lot of variability in these ice extent lines so I took the approach of averaging distances for equal length reaches of coastline, instead of thousands of tangent lines. You could simplify the coastline just for the allocation processing to straighten it enough to allow more bands if you wanted to split the coastline more.
I think this is a neat example of a problem that can be expressed and solved very simply using map algebra. That said, I am in awe of Xander's code up-thread. Xander probably provided a more accurate solution than this raster approach, pretty as it is.
from arcpy.sa import *
env.extent = "coast_simple"
env.cellSize = 10000
# create a raster of coastline cells
arcpy.PolylineToRaster_conversion("coast_simple", "OBJECTID", "ccell")
# Create coastline cell zones
# find the southernmost cell (this approach would vary with coastline config)
arcpy.RasterToPoint_conversion("ccell", "ccellpt")
arcpy.AddXY_management("ccellpt")
where = "ABS(POINT_Y - {}) < .001".format(arcpy.Describe("ccellpt").extent.YMin)
arcpy.MakeFeatureLayer_management("ccellpt", "cmin", where)
# Create 20 equal-length reaches along coast line
cc = ExtractByMask(CostDistance("cmin", Con(IsNull("ccell"), 1e6, 1)), "ccell")
cc1 = Slice(cc, 20, "EQUAL_AREA")
# separate ocean into areas corresponding to each coastline reach
calloc = EucAllocation(cc1)
# calculate distance from (any) coastline
cdist = EucDistance(cc1)
# get statistics for ice lines along coast for each year
for year in [2008, 2009, 2010, 2011]:
arcpy.SelectLayerByAttribute_management("ice_edges", "", "YEAR = {}".format(year))
ZonalStatisticsAsTable(ExtractByMask(calloc, "ice_edges"), "VALUE", cdist, "stats_{}".format(year), "DATA")