-
Notifications
You must be signed in to change notification settings - Fork 7
/
canny.fut
146 lines (126 loc) · 4.77 KB
/
canny.fut
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
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
-- Port of the canny example from Accelerate:
-- https://github.com/AccelerateHS/accelerate-examples/blob/master/examples/canny/src-acc
--
-- Only implements the first seven[sic] data-parallel stages, not the final
-- recursive algorithm (which in the Accelerate implementation is also
-- done in Repa).
--
-- ==
-- nobench input {
-- 50f32
-- 100f32
-- [[2092830683i32, 1394728708i32, -33326269i32],
-- [1986857019i32, 1313352304i32, -424163480i32],
-- [1813961220i32, -967167319i32, -503802777i32],
-- [-2133263610i32, -1327310708i32, -1863700019i32],
-- [-873352195i32, -1591419765i32, 378033750i32]]
-- }
-- output {
-- [1, 4, 6, 7]
-- }
--
-- compiled input @ data/lena256.in
-- output @ data/lena256.out
--
-- compiled input @ data/lena512.in
-- output @ data/lena512.out
def luminanceOfRGBA32(p: i32): f32 =
let r = i8.i32(p >> 24)
let g = i8.i32(p >> 16)
let b = i8.i32(p >> 8)
let r' = 0.3 * f32.i8(r)
let g' = 0.59 * f32.i8(g)
let b' = 0.11 * f32.i8(b)
in (r' + g' + b') / 255.0
def clamp(lower: i64, x: i64, upper: i64): i64 =
i64.max lower (i64.min upper x)
def orientUndef: i32 = 0
def orientPosD: i32 = 64
def orientVert: i32 = 128
def orientNegD: i32 = 192
def orientHoriz: i32 = 255
def edgeNone: f32 = 0.0
def edgeWeak: f32 = 0.5
def edgeStrong: f32 = 1.0
def toGreyscale [h][w] (img: [h][w]i32): [h][w]f32 =
map (\row -> map (255.0*) (map luminanceOfRGBA32 row)) img
def gaussianX [h][w] (img: [h][w]f32): [h][w]f32 =
tabulate_2d h w (\x y ->
#[unsafe]
let a = img[clamp(0,x-2,h-1),y] * (1.0 / 16.0)
let b = img[clamp(0,x-1,h-1),y] * (4.0 / 16.0)
let c = img[clamp(0,x+0,h-1),y] * (6.0 / 16.0)
let d = img[clamp(0,x+1,h-1),y] * (4.0 / 16.0)
let e = img[clamp(0,x+2,h-1),y] * (1.0 / 16.0)
in a + b + c + d + e)
def gaussianY [h][w] (img: [h][w]f32): [h][w]f32 =
tabulate_2d h w (\x y ->
#[unsafe]
let a = img[x,clamp(0,y-2,w-1)] * (1.0 / 16.0)
let b = img[x,clamp(0,y-1,w-1)] * (4.0 / 16.0)
let c = img[x,clamp(0,y+0,w-1)] * (6.0 / 16.0)
let d = img[x,clamp(0,y+1,w-1)] * (4.0 / 16.0)
let e = img[x,clamp(0,y+2,w-1)] * (1.0 / 16.0)
in a + b + c + d + e)
def gradiantMagDir [h][w] (low: f32) (img: [h][w]f32): [h][w](f32,i32) =
tabulate_2d h w (\x y ->
#[unsafe]
let v0 = img[clamp(0, x-1, h-1), clamp(0, y-1, w-1)]
let v1 = img[clamp(0, x+0, h-1), clamp(0, y-1, w-1)]
let v2 = img[clamp(0, x+1, h-1), clamp(0, y-1, w-1)]
let v3 = img[clamp(0, x-1, h-1), clamp(0, y+0, w-1)]
let v4 = img[clamp(0, x+1, h-1), clamp(0, y+0, w-1)]
let v5 = img[clamp(0, x-1, h-1), clamp(0, y+1, w-1)]
let v6 = img[clamp(0, x+0, h-1), clamp(0, y+1, w-1)]
let v7 = img[clamp(0, x+1, h-1), clamp(0, y+1, w-1)]
let dx = v2 + (2.0*v4) + v7 - v0 - (2.0*v3) - v5
let dy = v0 + (2.0*v1) + v2 - v5 - (2.0*v6) - v7
let mag = f32.sqrt(dx * dx + dy * dy)
let theta = f32.atan2 dy dx
let alpha = (theta - (f32.pi/8.0)) * (4.0/f32.pi)
-- Normalise the angle to between [0..8)
let norm = alpha + if alpha <= 0.0 then 8.0 else 0.0
let dir = if f32.abs(dx) <= low && f32.abs(dy) <= low
then 0
else i32.min (64 * (1 + i32.f32(norm) % 4)) 255
in (mag, dir))
def nonMaximumSuppression [h][w] (low: f32) (high: f32) (magdir: [h][w](f32,i32)): [h][w]f32 =
tabulate_2d h w (\x y ->
let (mag, dir) = magdir[x,y]
let offsetx = if dir > orientVert then -1
else if dir < orientVert then 1
else 0
let offsety = if dir < orientHoriz then -1 else 0
let (fwd, _) = #[unsafe]
magdir[clamp(0, x+offsetx, h-1),
clamp(0, y+offsety, w-1)]
let (rev, _) = #[unsafe]
magdir[clamp(0, x-offsetx, h-1),
clamp(0, y-offsety, w-1)]
in if dir == orientUndef || mag < low || mag < fwd || mag < rev
then 0.0
else if mag >= high then 1.0 else 0.5)
def selectStrong [h][w] (img: [h][w]f32): []i32 =
let strong =
map (\x -> if x == edgeStrong then 1 else 0)
(flatten img)
-- The original Accelerate implementation used an exclusive scan
-- here, so we have to play with the indices.
let targetIdxAndLen = scan (+) 0 strong
let targetIdx = init targetIdxAndLen
let len = i64.i32 (last targetIdxAndLen)
let (indices', values) =
unzip(map3 (\i target_i strong_x ->
if strong_x == 0
then (-1, 0)
else (i64.i32 target_i, i32.i64 i+1))
(indices targetIdx) targetIdx (strong[1:]))
in spread len 0 indices' values
def main [h][w] (low: f32) (high: f32) (img: [h][w]i32): []i32 =
img
|> toGreyscale
|> gaussianX
|> gaussianY
|> gradiantMagDir low
|> nonMaximumSuppression low high
|> selectStrong