forked from AFM-SPM/TopoStats
-
Notifications
You must be signed in to change notification settings - Fork 0
/
heightthresholding.py
123 lines (103 loc) · 5.13 KB
/
heightthresholding.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
#!/usr/bin/env python2
import gwy
# Set the settings for each function from the saved settings file (~/.gwyddion/settings)
s = gwy.gwy_app_settings_get()
# Turn colour bar off
s["/module/pixmap/ztype"] = 0
# Define the settings for image processing functions e.g. align rows here
# s['/module/linematch/method'] = 2
# 'align_rows' function
s["/module/linematch/direction"] = 0
s["/module/linematch/do_extract"] = False
s["/module/linematch/do_plot"] = False
s["/module/linematch/masking"] = 1
s["/module/linematch/max_degree"] = 0
s["/module/linematch/method"] = 1
s["/module/linematch/trim_fraction"] = float(0.05)
def otsuthresholdgrainfinding(data, k):
# 'align_rows' function
# s["/module/linematch/direction"] = 0
# s["/module/linematch/do_extract"] = False
# s["/module/linematch/do_plot"] = False
# s["/module/linematch/masking"] = 2
# s["/module/linematch/max_degree"] = 0
# s["/module/linematch/method"] = 0 # uses median
# s["/module/linematch/trim_fraction"] = float(0.05)
# Select channel 'k' of the file
gwy.gwy_app_data_browser_select_data_field(data, k)
datafield = gwy.gwy_app_data_browser_get_current(gwy.APP_DATA_FIELD)
mask = gwy.DataField.new_alike(datafield, False)
# Apply a 1.5 pixel gaussian filter
data_field = gwy.gwy_app_data_browser_get_current(gwy.APP_DATA_FIELD)
data_field.filter_gaussian(1.5)
# # Mask data that are above thresh*sigma from average height.
# # Sigma denotes root-mean square deviation of heights.
# # This criterium corresponds to the usual Gaussian distribution outliers detection if thresh is 3.
# datafield.mask_outliers(mask, 1)
# # Shift contrast - equivalent to 'fix zero' - essential for next step to work
datafield.add(-datafield.get_min())
# Calculate min, max and range of data to allow calculation of relative value for grain thresholding
min_datarange = datafield.get_min()
max_datarange = datafield.get_max()
datarange = max_datarange + min_datarange
# Calculate Otsu threshold for data
o_threshold = datafield.otsu_threshold()
o_threshold = o_threshold + min_datarange
# Calculate relative threshold for grain determination
rel_threshold = 100 * (o_threshold / datarange)
# Mask grains using either relative threshold of the data, i.e. a percentage
# this can be set manually i.e. 35 or
# as the otsu threshold expressed as a percentage which is rel_threshold
# this will fail unless the min_datarange is negative
datafield.grains_mark_height(mask, rel_threshold, False)
# # Invert mask to maks things below the membrane
mask.grains_invert()
gwy.gwy_process_func_run("align_rows", data, gwy.RUN_IMMEDIATE)
# # gwy.gwy_process_func_run("level", data, gwy.RUN_IMMEDIATE)
# gwy.gwy_process_func_run('flatten_base', data, gwy.RUN_IMMEDIATE)
# Select channel 'k' of the file
gwy.gwy_app_data_browser_select_data_field(data, k)
datafield = gwy.gwy_app_data_browser_get_current(gwy.APP_DATA_FIELD)
mask = gwy.DataField.new_alike(datafield, False)
# # Mask data that are above thresh*sigma from average height.
# # Sigma denotes root-mean square deviation of heights.
# # This criterium corresponds to the usual Gaussian distribution outliers detection if thresh is 3.
# datafield.mask_outliers(mask, 1)
# # Shift contrast - equivalent to 'fix zero' - essential for next step to work
datafield.add(-datafield.get_min())
# Calculate min, max and range of data to allow calculation of relative value for grain thresholding
min_datarange = datafield.get_min()
max_datarange = datafield.get_max()
datarange = max_datarange + min_datarange
# Calculate Otsu threshold for data
o_threshold = datafield.otsu_threshold()
o_threshold = o_threshold + min_datarange
# Calculate relative threshold for grain determination
rel_threshold = 100 * (o_threshold / datarange)
# Mask grains using either relative threshold of the data, i.e. a percentage
# this can be set manually i.e. 35 or
# as the otsu threshold expressed as a percentage which is rel_threshold
# this will fail unless the min_datarange is negative
datafield.grains_mark_height(mask, rel_threshold, False)
gwy.gwy_process_func_run('zero_mean', data, gwy.RUN_IMMEDIATE)
# # Invert mask to make things below the membrane
mask.grains_invert()
# Calculate pixel width in nm
dx = datafield.get_dx()
# Calculate minimum feature size in pixels (integer)
minsize = int(2000e-9 / dx)
# Remove grains smaller than (size) in integer pixels
mask.grains_remove_by_size(minsize)
mask.grains_invert()
gwy.gwy_process_func_run('zero_mean', data, gwy.RUN_IMMEDIATE)
mask.grains_invert()
# Numbering grains for grain analysis
grains = mask.number_grains()
print max(grains)
# Update data to show mask, comment out to remove mask
s['/module/pixmap/draw_mask'] = True
data['/0/mask/red'] = 0.1234
#excluding mask, zero mean
stats = datafield.area_get_stats_mask(mask, gwy.MASK_EXCLUDE, 0, 0, datafield.get_xres(), datafield.get_yres())
datafield.add(-stats[0])
return data, mask, datafield, grains