Valley_Fun.py 文件源码

def neg_bond(curvature_file, valley_bond_file , ridge_bond_file , covergent_sig_file , divergent_sig_file):

    curve_Array = arcpy.RasterToNumPyArray(curvature_file  , nodata_to_value=0)
    corner = arcpy.Point(arcpy.Describe(curvature_file).Extent.XMin,arcpy.Describe(curvature_file).Extent.YMin)
    dx = arcpy.Describe(curvature_file).meanCellWidth

    valley_thresh , ridge_thresh = segment_size_curvature(curve_Array)


    valley = np.where(curve_Array > valley_thresh , 1 , 0)
    ridge = np.where(curve_Array < ridge_thresh , 1 , 0)

    Lab_Val_Array, num_label_val = ndimage.label(valley , structure = np.ones((3 , 3)))
    Lab_Rid_Array, num_label_rid = ndimage.label(ridge , structure = np.ones((3 , 3)))

    perc90_valley = ndimage.labeled_comprehension(curve_Array, Lab_Val_Array, np.arange(1, num_label_val + 1) , percentile, float, 0)
    perc90_ridge = ndimage.labeled_comprehension(-1 * curve_Array, Lab_Rid_Array, np.arange(1, num_label_rid + 1) , percentile, float, 0)

    area_valley = ndimage.labeled_comprehension(valley, Lab_Val_Array, np.arange(1, num_label_val + 1) , np.sum, float, 0)
    area_ridge = ndimage.labeled_comprehension(ridge, Lab_Rid_Array, np.arange(1, num_label_rid + 1) , np.sum, float, 0)

    main_patch_val = perc90_valley[area_valley == np.max(area_valley)][0]
    main_patch_ridge = perc90_ridge[area_ridge == np.max(area_ridge)][0]


    sig_valley_thresh = otsu(perc90_valley[perc90_valley < main_patch_val])
    sig_ridge_thresh = otsu(perc90_ridge[perc90_ridge < main_patch_ridge])

    perc90_valley = perc90_valley[Lab_Val_Array - 1]
    perc90_valley = np.where(valley == 1 , perc90_valley , 0)

    perc90_ridge = perc90_ridge[Lab_Rid_Array - 1]
    perc90_ridge = np.where(ridge == 1 , perc90_ridge , 0)

    sig_valley_patch = np.where(perc90_valley > sig_valley_thresh , 1 , 0)
    sig_ridge_patch = np.where(perc90_ridge > sig_ridge_thresh , 1 , 0)

    arcpy.NumPyArrayToRaster(perc90_valley , corner,dx,dx , value_to_nodata=0).save(covergent_sig_file)
    arcpy.NumPyArrayToRaster(perc90_ridge , corner,dx,dx , value_to_nodata=0).save(divergent_sig_file)

    arcpy.NumPyArrayToRaster(sig_valley_patch , corner,dx,dx , value_to_nodata=0).save(valley_bond_file)
    arcpy.NumPyArrayToRaster(sig_ridge_patch , corner,dx,dx , value_to_nodata=0).save(ridge_bond_file)

    return valley_thresh