1
+ use crate :: core:: PixelBound ;
1
2
use crate :: core:: { ColourModel , Image } ;
2
- use crate :: core:: { PixelBound } ;
3
3
use crate :: processing:: * ;
4
4
use ndarray:: prelude:: * ;
5
- use num_traits:: cast:: { FromPrimitive } ;
6
- use num_traits:: cast:: { ToPrimitive } ;
7
- use ndarray_stats:: QuantileExt ;
5
+ use ndarray_stats:: histogram:: { Bins , Edges , Grid } ;
8
6
use ndarray_stats:: HistogramExt ;
9
- use ndarray_stats:: histogram:: { Grid , Bins , Edges } ;
7
+ use ndarray_stats:: QuantileExt ;
8
+ use num_traits:: cast:: FromPrimitive ;
9
+ use num_traits:: cast:: ToPrimitive ;
10
10
use num_traits:: { Num , NumAssignOps } ;
11
+ use std:: iter:: FromIterator ;
11
12
use std:: marker:: PhantomData ;
12
13
13
14
// Development
@@ -16,14 +17,13 @@ use assert_approx_eq::assert_approx_eq;
16
17
#[ cfg( test) ]
17
18
use noisy_float:: types:: n64;
18
19
19
-
20
20
/// Runs the Otsu Thresholding algorithm on a type T
21
21
pub trait ThresholdOtsuExt < T > {
22
22
/// Output type, this is different as Otsu outputs a binary image
23
23
type Output ;
24
24
25
- /// Run the Otsu threshold detection algorithm with the
26
- /// given parameters. Due to Otsu being specified as working
25
+ /// Run the Otsu threshold detection algorithm with the
26
+ /// given parameters. Due to Otsu being specified as working
27
27
/// on greyscale images all current implementations
28
28
/// assume a single channel image returning an error otherwise.
29
29
fn threshold_otsu ( & self ) -> Result < Self :: Output , Error > ;
@@ -34,8 +34,8 @@ pub trait ThresholdMeanExt<T> {
34
34
/// Output type, this is different as Otsu outputs a binary image
35
35
type Output ;
36
36
37
- /// Run the Otsu threshold detection algorithm with the
38
- /// given parameters. Due to Otsu being specified as working
37
+ /// Run the Otsu threshold detection algorithm with the
38
+ /// given parameters. Due to Otsu being specified as working
39
39
/// on greyscale images all current implementations
40
40
/// assume a single channel image returning an error otherwise.
41
41
fn threshold_mean ( & self ) -> Result < Self :: Output , Error > ;
60
60
61
61
impl < T > ThresholdOtsuExt < T > for Array3 < T >
62
62
where
63
- T : Copy + Clone + Ord + Num + NumAssignOps + ToPrimitive + FromPrimitive
63
+ T : Copy + Clone + Ord + Num + NumAssignOps + ToPrimitive + FromPrimitive ,
64
64
{
65
65
type Output = Array3 < bool > ;
66
66
@@ -76,60 +76,58 @@ where
76
76
}
77
77
78
78
///
79
- /// Calculates Otsu's threshold
80
- /// Works per channel, but currently
79
+ /// Calculates Otsu's threshold
80
+ /// Works per channel, but currently
81
81
/// assumes grayscale (see the error above if number of channels is > 1
82
- /// i.e. single channel; otherwise we need to output all 3 threshold values).
82
+ /// i.e. single channel; otherwise we need to output all 3 threshold values).
83
83
/// Todo: Add optional nbins
84
84
///
85
85
fn calculate_threshold_otsu < T > ( mat : & Array3 < T > ) -> Result < f64 , Error >
86
86
where
87
- T : Copy + Clone + Ord + Num + NumAssignOps + ToPrimitive + FromPrimitive
87
+ T : Copy + Clone + Ord + Num + NumAssignOps + ToPrimitive + FromPrimitive ,
88
88
{
89
89
let mut threshold = 0.0 ;
90
90
let n_bins = 255 ;
91
91
for c in mat. axis_iter ( Axis ( 2 ) ) {
92
- let scale_factor = ( n_bins) as f64
93
- /( c. max ( ) . unwrap ( ) . to_f64 ( ) . unwrap ( ) ) ;
92
+ let scale_factor = ( n_bins) as f64 / ( c. max ( ) . unwrap ( ) . to_f64 ( ) . unwrap ( ) ) ;
94
93
let edges_vec: Vec < u8 > = ( 0 ..n_bins) . collect ( ) ;
95
94
let grid = Grid :: from ( vec ! [ Bins :: new( Edges :: from( edges_vec) ) ] ) ;
96
-
95
+
97
96
// get the histogram
98
97
let flat = Array :: from_iter ( c. iter ( ) ) . insert_axis ( Axis ( 1 ) ) ;
99
- let flat2 = flat. mapv (
100
- |x| ( ( * x) . to_f64 ( ) . unwrap ( ) * scale_factor) . to_u8 ( ) . unwrap ( )
101
- ) ;
98
+ let flat2 = flat. mapv ( |x| ( ( * x) . to_f64 ( ) . unwrap ( ) * scale_factor) . to_u8 ( ) . unwrap ( ) ) ;
102
99
let hist = flat2. histogram ( grid) ;
103
- // Straight out of wikipedia:
104
- let counts = hist. counts ( ) ;
100
+ // Straight out of wikipedia:
101
+ let counts = hist. counts ( ) ;
105
102
let total = counts. sum ( ) . to_f64 ( ) . unwrap ( ) ;
106
103
let counts = Array :: from_iter ( counts. iter ( ) ) ;
107
- // NOTE: Could use the cdf generation for skimage-esque implementation
108
- // which entails a cumulative sum of the standard histogram
104
+ // NOTE: Could use the cdf generation for skimage-esque implementation
105
+ // which entails a cumulative sum of the standard histogram
109
106
let mut sum_b = 0.0 ;
110
107
let mut weight_b = 0.0 ;
111
108
let mut maximum = 0.0 ;
112
109
let mut level = 0.0 ;
113
110
let mut sum_intensity = 0.0 ;
114
- for ( index, count) in counts. indexed_iter ( ) {
111
+ for ( index, count) in counts. indexed_iter ( ) {
115
112
sum_intensity += ( index as f64 ) * ( * count) . to_f64 ( ) . unwrap ( ) ;
116
113
}
117
- for ( index, count) in counts. indexed_iter ( ) {
114
+ for ( index, count) in counts. indexed_iter ( ) {
118
115
weight_b = weight_b + count. to_f64 ( ) . unwrap ( ) ;
119
- sum_b = sum_b + ( index as f64 ) * count. to_f64 ( ) . unwrap ( ) ;
116
+ sum_b = sum_b + ( index as f64 ) * count. to_f64 ( ) . unwrap ( ) ;
120
117
let weight_f = total - weight_b;
121
- if ( weight_b > 0.0 ) && ( weight_f > 0.0 ) {
118
+ if ( weight_b > 0.0 ) && ( weight_f > 0.0 ) {
122
119
let mean_f = ( sum_intensity - sum_b) / weight_f;
123
- let val = weight_b * weight_f
124
- * ( ( sum_b / weight_b) - mean_f)
125
- * ( ( sum_b / weight_b) - mean_f) ;
126
- if val > maximum{
120
+ let val = weight_b
121
+ * weight_f
122
+ * ( ( sum_b / weight_b) - mean_f)
123
+ * ( ( sum_b / weight_b) - mean_f) ;
124
+ if val > maximum {
127
125
level = 1.0 + ( index as f64 ) ;
128
126
maximum = val;
129
127
}
130
128
}
131
129
}
132
- threshold = level as f64 / scale_factor;
130
+ threshold = level as f64 / scale_factor;
133
131
}
134
132
Ok ( threshold)
135
133
}
@@ -153,7 +151,7 @@ where
153
151
154
152
impl < T > ThresholdMeanExt < T > for Array3 < T >
155
153
where
156
- T : Copy + Clone + Ord + Num + NumAssignOps + ToPrimitive + FromPrimitive
154
+ T : Copy + Clone + Ord + Num + NumAssignOps + ToPrimitive + FromPrimitive ,
157
155
{
158
156
type Output = Array3 < bool > ;
159
157
@@ -170,12 +168,11 @@ where
170
168
171
169
fn calculate_threshold_mean < T > ( array : & Array3 < T > ) -> Result < f64 , Error >
172
170
where
173
- T : Copy + Clone + Num + NumAssignOps + ToPrimitive + FromPrimitive
171
+ T : Copy + Clone + Num + NumAssignOps + ToPrimitive + FromPrimitive ,
174
172
{
175
173
Ok ( array. sum ( ) . to_f64 ( ) . unwrap ( ) / array. len ( ) as f64 )
176
174
}
177
175
178
-
179
176
fn apply_threshold < T > ( data : & Array3 < T > , threshold : f64 ) -> Array3 < bool >
180
177
where
181
178
T : Copy + Clone + Num + NumAssignOps + ToPrimitive + FromPrimitive ,
@@ -184,7 +181,6 @@ where
184
181
result
185
182
}
186
183
187
-
188
184
#[ cfg( test) ]
189
185
mod tests {
190
186
use super :: * ;
@@ -208,24 +204,20 @@ mod tests {
208
204
209
205
assert_eq ! ( result, expected) ;
210
206
}
211
-
207
+
212
208
#[ test]
213
209
fn threshold_calculate_threshold_otsu_ints ( ) {
214
- let data = arr3 ( & [
215
- [ [ 2 ] , [ 4 ] , [ 0 ] ] ,
216
- [ [ 7 ] , [ 5 ] , [ 8 ] ] ,
217
- [ [ 1 ] , [ 6 ] , [ 0 ] ] ,
218
- ] ) ;
210
+ let data = arr3 ( & [ [ [ 2 ] , [ 4 ] , [ 0 ] ] , [ [ 7 ] , [ 5 ] , [ 8 ] ] , [ [ 1 ] , [ 6 ] , [ 0 ] ] ] ) ;
219
211
let result = calculate_threshold_otsu ( & data) . unwrap ( ) ;
220
212
println ! ( "Done {}" , result) ;
221
213
222
- // Calculated using Python's skimage.filters.threshold_otsu
223
- // on int input array. Float array returns 2.0156...
214
+ // Calculated using Python's skimage.filters.threshold_otsu
215
+ // on int input array. Float array returns 2.0156...
224
216
let expected = 2.0 ;
225
217
226
218
assert_approx_eq ! ( result, expected, 5e-1 ) ;
227
219
}
228
-
220
+
229
221
#[ test]
230
222
fn threshold_calculate_threshold_otsu_floats ( ) {
231
223
let data = arr3 ( & [
@@ -236,27 +228,23 @@ mod tests {
236
228
237
229
let result = calculate_threshold_otsu ( & data) . unwrap ( ) ;
238
230
239
- // Calculated using Python's skimage.filters.threshold_otsu
240
- // on int input array. Float array returns 2.0156...
231
+ // Calculated using Python's skimage.filters.threshold_otsu
232
+ // on int input array. Float array returns 2.0156...
241
233
let expected = 2.0156 ;
242
234
243
235
assert_approx_eq ! ( result, expected, 5e-1 ) ;
244
236
}
245
-
237
+
246
238
#[ test]
247
239
fn threshold_calculate_threshold_mean_ints ( ) {
248
- let data = arr3 ( & [
249
- [ [ 4 ] , [ 4 ] , [ 4 ] ] ,
250
- [ [ 5 ] , [ 5 ] , [ 5 ] ] ,
251
- [ [ 6 ] , [ 6 ] , [ 6 ] ] ,
252
- ] ) ;
240
+ let data = arr3 ( & [ [ [ 4 ] , [ 4 ] , [ 4 ] ] , [ [ 5 ] , [ 5 ] , [ 5 ] ] , [ [ 6 ] , [ 6 ] , [ 6 ] ] ] ) ;
253
241
254
242
let result = calculate_threshold_mean ( & data) . unwrap ( ) ;
255
243
let expected = 5.0 ;
256
-
244
+
257
245
assert_approx_eq ! ( result, expected, 1e-16 ) ;
258
246
}
259
-
247
+
260
248
#[ test]
261
249
fn threshold_calculate_threshold_mean_floats ( ) {
262
250
let data = arr3 ( & [
@@ -267,7 +255,7 @@ mod tests {
267
255
268
256
let result = calculate_threshold_mean ( & data) . unwrap ( ) ;
269
257
let expected = 5.0 ;
270
-
258
+
271
259
assert_approx_eq ! ( result, expected, 1e-16 ) ;
272
260
}
273
261
}
0 commit comments