Skip to content

Commit 2a2574a

Browse files
committed
adding draft of convolution chapters
1 parent 098c5b6 commit 2a2574a

18 files changed

+380
-741
lines changed

contents/convolutions/1d/1d.md

Lines changed: 13 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -1,9 +1,3 @@
1-
# TODO:
2-
3-
* Make sure cyclic convolution works with larger filter
4-
* LICENSING
5-
6-
71
# Convolutions in 1D
82
As mentioned in the [introductory section for convolutions](../convolutions.md), convolutions are methods that allow mathematicians to "blend" two seemingly unrelated functions; however, this definition is not very rigorouos, so it might be better to think of a convolution as a method to apply a filter to a signal or image.
93
This, of course, brings up more questions: what is a filter? What is a signal? How is this all related to images?
@@ -349,6 +343,19 @@ MathJax.Hub.Queue(["Typeset",MathJax.Hub]);
349343

350344
The code examples are licensed under the MIT license (found in [LICENSE.md](https://github.com/algorithm-archivists/algorithm-archive/blob/master/LICENSE.md)).
351345

346+
##### Images/Graphics
347+
- The image "[Square Wave](../res/square_wave.png)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
348+
- The image "[Triangle Wave](../res/triangle_wave.png)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
349+
- The video "[Triangle Square Convolution](../res/triangle_square_conv.mp4)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
350+
- The video "[Gaussian Square Convolution](../res/1d_gaussian.mp4)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
351+
- The video "[Gaussian Random Convolution](../res/1d_rand_gaussian_cyclic.mp4)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
352+
- The video "[Double Convolution](../res/double_gaussian.mp4)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
353+
- The image "[Sawtooth Wave](../res/sawtooth.png)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
354+
- The video "[Sawtooth Square Convolution](../res/1d_sawtooth.mp4)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
355+
- The video "[Full Random Convolution](../res/1d_rand_gaussian_full.mp4)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
356+
- The video "[Simple Random Convolution](../res/1d_rand_gaussian_simple.mp4)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
357+
358+
352359
##### Text
353360

354361
The text of this chapter was written by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).

contents/convolutions/2d/2d.md

Lines changed: 72 additions & 25 deletions
Original file line numberDiff line numberDiff line change
@@ -2,16 +2,29 @@
22

33
For this section, we will no longer be focusing on signals, but instead images, where an image is an array filled with elements of red, green, and blue values.
44
That said, for the code examples, greyscale images may be used such that each array element is simply composed of some floating-point value.
5+
In addition, we will not be discussing boundary conditions too much in this chapter and will instead be using the simple boundaries introduced in the section on [one-dimensional convolutions](../1d/1d.md).
56

67
The extension of one-dimensional convolutions to two dimensions requires a little thought about indexing and the like, but is ultimately the same operation.
78
Here is an animation of a convolution for a two-dimensional image:
89

9-
ADD ANIMATION
10+
<div style="text-align:center">
11+
<video style="width:90%" controls loop>
12+
<source src="../res/2d.mp4" type="video/mp4">
13+
Your browser does not support the video tag.
14+
</video>
15+
</div>
1016

11-
In this case, we convolved the image with a 3x3 square filter.
12-
In code, using simple boundaries, it would look like this:
17+
In this case, we convolved the image with a 3x3 square filter, all filled with values of $$\frac{1}{9}$$.
18+
This created a simple blurring effect, which is somewhat expected from the previous section.
19+
In code, using simple boundaries, a two-dimensional convolution might look like this:
1320

14-
ADD CODE
21+
{% method %}
22+
{% sample lang="jl" %}
23+
[import:4-28, lang:"julia"](../code/julia/2d_convolution.jl)
24+
{% endmethod %}
25+
26+
This code is very similar to what we have shown in previous sections; however, it essentially requires 4 iterable dimensions because we need to iterate through each element of the output domain *and* the filter.
27+
Because the indexing is slightly non-trivial, we felt it was worth writing a separate chapter for two-dimensional convolutions.
1528

1629
It is worth highlighting common filters used for convolutions of images.
1730
In particular, we will further discuss the Gaussian filter introduced in the section on [one-dimensional convolutions](../1d/1d.md), and then introduce another set of kernels known as Sobel operators, which are used for naive edge detection or image derivatives.
@@ -21,40 +34,47 @@ In particular, we will further discuss the Gaussian filter introduced in the sec
2134
The Gaussian kernel serves as an effective *blurring* operation for images.
2235
Like we showed with the triangle filter in the section on [one-dimensional convolutions](../1d/1d.md), the Gaussian kernel effectively adds a small amount of the neighboring elements to each pixel in an image.
2336

24-
As a reminder, the foruma for any Gaussian distribution is
37+
As a reminder, the formula for any Gaussian distribution is
2538

2639
$$
2740
g(x,y) = \frac{1}{2\pi\sigma^2}e^{-\frac{x^2+y^2}{2\sigma^2}},
2841
$$
2942

3043
where $$\sigma$$ is the standard deviation and is a measure of the width of the Gaussian.
3144
A larger $$\sigma$$ means a larger Gaussian; however, Remember that the Gaussian must fit onto the filter!
45+
As a note, some definitions of $$\sigma$$ allow users to have a separate deviation in $$x$$ and $$y$$, but for the purposes of this chapter, we will assume they are the same.
3246
As a general rule of thumb, the larger the filter and standard deviation, the more "smeared" the final convolution will be.
3347

3448
At this stage, it is important to write some code, so we will generate a simple function to return a Gaussian kernel with a specified standard deviation and filter size.
35-
As a note, all the kernels will be scaled at the end by the sum of all internal elements.
36-
This simply ensures that the output of the convolution does not have an obnoxious scale factor associated with it.
3749

38-
ADD CODE
50+
{% method %}
51+
{% sample lang="jl" %}
52+
[import:30-47, lang:"julia"](../code/julia/2d_convolution.jl)
53+
{% endmethod %}
3954

55+
Though it is entirely possible to create a Gaussian kernel whose standard deviation is independent on the kernel size, we have decided to keep them the same for this chapter.
56+
As a note, all the kernels will be scaled at the end by the sum of all internal elements.
57+
This simply ensures that the output of the convolution does not have an obnoxious scale factor associated with it.
4058

41-
Below are a few kernels generated with the above code along with their application to a ___ image.
59+
Below are a few kernels generated with the above code along with their application to a black and white heart image.
4260

43-
Add 1 + 2x4 image with original image going into each lane and then showing the convolution result for small, large, and somewhat in-between.
61+
<p>
62+
<img class="center" src="../res/circle_blur.png" style="width:100%" />
63+
</p>
4464

45-
In (a), we show a 3x3 kernel with a standard deviation of 1.
46-
In (b), we show a 10x10 kernel with a standard deviation of 5.
47-
In (c), we show a 3x3 kernel with a standard deviation of 5.
48-
Finally, in (d), we show a 10x10 kernel with a standard deviation of 1.
4965

50-
We see that (a) and (d) are somewhat similar, (b) has the most blurred boundary, and (c) is partially blurred.
66+
In (a), we show the original image, which is just a black circle at the center of a $$50\times 50$$ grid.
67+
In (b), we show the image after convolution with a $$3\times 3$$ kernel.
68+
In (c), we show the image after convolution with a $$20\times 20$$ kernel.
69+
Here, we see that (c) is significantly fuzzier than (b), which is a direct consequence of the kernel size.
5170

52-
There is a bit more to this, but I think this is a good place to wrap up the discussion on the Gaussian blurring operation before moving onto a slightly more complex convolutional method: the Sobel operator.
71+
There is a lot more that we could talk about, but I think this is a good place to wrap up the discussion on the Gaussian blurring operation before moving onto a slightly more complex convolutional method: the Sobel operator.
5372

5473
## The Sobel operator
5574

56-
The Sobel operator effectively performs a gradient operation on an image by highlighting areas where a large change has been made and can be considered a naive edge detector.
57-
That is to say that the $$n$$-dimensional Sobel operator is composed of $$n$$ separate gradient convolutions that are then combined together into one, output array.
75+
The Sobel operator effectively performs a gradient operation on an image by highlighting areas where a large change has been made.
76+
In essence, this means that this operation can be considered to be a naïve edge detector.
77+
That is to say that the $$n$$-dimensional Sobel operator is composed of $$n$$ separate gradient convolutions (one for each dimension) that are then combined together into one, output array.
5878
Again, for the purposes of this chapter, we will stick to two dimensions, which will be composed of two separate gradients along the $$x$$ and $$y$$ directions.
5979
Each gradient will be created by convolving our image with their corresponding Sobel operator:
6080

@@ -94,36 +114,63 @@ G_y &= S_y*A.
94114
\end{align}
95115
$$
96116

117+
Here, $$A$$ is out input array or image.
97118
Finally, these gradients can be summed in quadrature to find the total Sobel operator or image gradient:
98119

99120
$$
100121
G_{\text{total}} = \sqrt{G_x^2 + G_y^2}
101122
$$
102123

103-
At this point, it does not likely make sense, so let us now show what it does in practice:
124+
So let us now show what it does in practice:
104125

105-
ADD IMAGE
126+
<p>
127+
<img class="center" src="../res/sobel_filters.png" style="width:100%" />
128+
</p>
106129

107-
Here, we see that the edges in the image are highlighted.
108-
This is why the Sobel operators are also known as naive edge detectors and are integral components to many more sophisticated edge detection methods like one proposed by Canny [CITE].
130+
In this diagram, we start with the circle image on the right, and then use the $$S_x$$ and $$S_y$$ operators on them to find the gradients along $$x$$ and $$y$$ before summing them in quadrature to get the final image gradient.
131+
Here, we see that the edges in the image are highlighted, showing the outline of our circle image.
132+
This is why the Sobel operator is also known as naïve edge detection and is an integral component to many more sophisticated edge detection methods like one proposed by Canny {{ "canny1986computational" | cite }}.
109133

110-
For now, it is a good idea to show what this operation might look like in code:
134+
In code, applying the Sobel operator involves first finding the operator in $$x$$ and $$y$$ and then applying them with a traditional convolution:
111135

112-
ADD CODE
136+
{% method %}
137+
{% sample lang="jl" %}
138+
[import:49-63, lang:"julia"](../code/julia/2d_convolution.jl)
139+
{% endmethod %}
113140

114-
With that, I believe we are ate a good place to stop discussions on two-dimensional convolutions.
141+
With that, I believe we are at a good place to stop discussions on two-dimensional convolutions.
115142
We will definitely return to this topic in the future as new algorithms require more information.
116143

144+
## Example Code
145+
146+
For the full code in this section, we have modified the code from the [one-dimensional convolution chapter](../1d/1d.md) to add a two-dimensional variant for an image of random white noise.
147+
We have also added code to create the Gaussian kernel and Sobel operator.
148+
149+
{% method %}
150+
{% sample lang="jl" %}
151+
[import, lang:"julia"](../code/julia/2d_convolution.jl)
152+
{% endmethod %}
153+
117154
<script>
118155
MathJax.Hub.Queue(["Typeset",MathJax.Hub]);
119156
</script>
120157

158+
### Bibliography
159+
160+
{% references %} {% endreferences %}
161+
121162
## License
122163

123164
##### Code Examples
124165

125166
The code examples are licensed under the MIT license (found in [LICENSE.md](https://github.com/algorithm-archivists/algorithm-archive/blob/master/LICENSE.md)).
126167

168+
##### Images/Graphics
169+
- The image "[8bit Heart](../res/heart_8bit.png)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
170+
- The image "[Circle Blur](../res/circle_blur.png)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
171+
- The image "[Sobel Filters](../res/sobel_filters.png)" was created by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
172+
- The video "[2D Convolution](../res/2d.mp4)" was created by [James Schloss](https://github.com/leios) and [Grant Sanderson](https://github.com/3b1b) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
173+
127174
##### Text
128175

129176
The text of this chapter was written by [James Schloss](https://github.com/leios) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).

contents/convolutions/code/julia/1d_convolution.jl

Lines changed: 9 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -1,20 +1,20 @@
11
using DelimitedFiles
22
using LinearAlgebra
33

4-
function convolve_cyclic(signal1::Array{T, 1},
5-
signal2::Array{T, 1}) where {T <: Number}
4+
function convolve_cyclic(signal::Array{T, 1},
5+
filter::Array{T, 1}) where {T <: Number}
66

77
# output size will be the size of sign
8-
output_size = max(length(signal1), length(signal2))
8+
output_size = max(length(signal), length(filter))
99

1010
# convolutional output
1111
out = Array{Float64,1}(undef,output_size)
1212
sum = 0
1313

1414
for i = 1:output_size
1515
for j = 1:output_size
16-
if (mod1(i-j, output_size) <= length(signal2))
17-
sum += signal1[mod1(j, output_size)] * signal2[mod1(i-j, output_size)]
16+
if (mod1(i-j, output_size) <= length(filter))
17+
sum += signal[mod1(j, output_size)] * filter[mod1(i-j, output_size)]
1818
end
1919
end
2020

@@ -26,17 +26,17 @@ function convolve_cyclic(signal1::Array{T, 1},
2626
return out
2727
end
2828

29-
function convolve_linear(signal1::Array{T, 1}, signal2::Array{T, 1},
29+
function convolve_linear(signal::Array{T, 1}, filter::Array{T, 1},
3030
output_size) where {T <: Number}
3131

3232
# convolutional output
3333
out = Array{Float64,1}(undef, output_size)
3434
sum = 0
3535

3636
for i = 1:output_size
37-
for j = 1:i
38-
if j < length(signal1) && i-j+1 < length(signal2)
39-
sum += signal1[j] * signal2[i-j+1]
37+
for j = max(1, i-length(filter)):i
38+
if j <= length(signal) && i-j+1 <= length(filter)
39+
sum += signal[j] * filter[i-j+1]
4040
end
4141
end
4242

Lines changed: 125 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,125 @@
1+
using DelimitedFiles
2+
using LinearAlgebra
3+
4+
function convolve_linear(signal::Array{T, 2}, filter::Array{T, 2},
5+
output_size) where {T <: Number}
6+
7+
# convolutional output
8+
out = Array{Float64,2}(undef, output_size)
9+
sum = 0
10+
11+
for i = 1:output_size[1]
12+
for j = 1:output_size[2]
13+
for k = max(1, i-size(filter)[1]):i
14+
for l = max(1, j-size(filter)[2]):j
15+
if k <= size(signal)[1] && i-k+1 <= size(filter)[1] &&
16+
l <= size(signal)[2] && j-l+1 <= size(filter)[2]
17+
sum += signal[k,l] * filter[i-k+1, j-l+1]
18+
end
19+
end
20+
end
21+
22+
out[i,j] = sum
23+
sum = 0
24+
end
25+
end
26+
27+
return out
28+
end
29+
30+
function create_gaussian_kernel(kernel_size)
31+
32+
kernel = zeros(kernel_size, kernel_size)
33+
34+
# The center must be offset by 0.5 to find the correct index
35+
center = kernel_size * 0.5 + 0.5
36+
37+
sigma = sqrt(0.1*kernel_size)
38+
39+
for i = 1:kernel_size
40+
for j = 1:kernel_size
41+
kernel[i,j] = exp(-((i-center)^2 + (j-center)^2) / (2*sigma^2))
42+
end
43+
end
44+
45+
return normalize(kernel)
46+
47+
end
48+
49+
function create_sobel_operators()
50+
Sx = [1.0, 2.0, 1.0]*[-1.0 0.0 1.0] / 9
51+
Sy = [-1.0, 0.0, 1.0]*[1.0 2.0 1.0] / 9
52+
53+
return Sx, Sy
54+
end
55+
56+
function compute_sobel(signal)
57+
Sx, Sy = create_sobel_operators()
58+
59+
Gx = convolve_linear(signal, Sx, size(signal) .+ size(Sx))
60+
Gy = convolve_linear(signal, Sy, size(signal) .+ size(Sy))
61+
62+
return sqrt.(Gx.^2 .+ Gy.^2)
63+
end
64+
65+
# Simple function to create a square grid with a circle embedded inside of it
66+
function create_circle(image_resolution, grid_extents, radius)
67+
out = zeros(image_resolution, image_resolution)
68+
69+
for i = 1:image_resolution
70+
x_position = ((i-1)*grid_extents/image_resolution)-0.5*grid_extents
71+
for j = 1:image_resolution
72+
y_position = ((j-1)*grid_extents/image_resolution)-0.5*grid_extents
73+
if x_position^2 + y_position^2 <= radius^2
74+
out[i,j] = 1.0
75+
end
76+
end
77+
end
78+
79+
return out
80+
end
81+
82+
function main()
83+
84+
# Random distribution in x
85+
x = rand(100, 100)
86+
87+
# Gaussian signals
88+
y = [exp(-(((i-50)/100)^2 + ((j-50)/100)^2)/.01) for i = 1:100, j=1:100]
89+
90+
# Normalization is not strictly necessary, but good practice
91+
normalize!(x)
92+
normalize!(y)
93+
94+
# full convolution, output will be the size of x + y
95+
full_linear_output = convolve_linear(x, y, size(x) .+ size(y))
96+
97+
# simple boundaries
98+
simple_linear_output = convolve_linear(x, y, size(x))
99+
100+
# outputting convolutions to different files for plotting in external code
101+
writedlm("full_linear.txt", full_linear_output)
102+
writedlm("simple_linear.txt", simple_linear_output)
103+
104+
# creating simple circle and 2 different Gaussian kernels
105+
circle = create_circle(50,2,0.5)
106+
107+
normalize!(circle)
108+
109+
small_kernel = create_gaussian_kernel(3)
110+
large_kernel = create_gaussian_kernel(25)
111+
112+
small_kernel_output = convolve_linear(circle, small_kernel,
113+
size(circle).+size(small_kernel))
114+
large_kernel_output = convolve_linear(circle, large_kernel,
115+
size(circle).+size(large_kernel))
116+
117+
writedlm("small_kernel.txt", small_kernel_output)
118+
writedlm("large_kernel.txt", large_kernel_output)
119+
120+
# Using the circle for Sobel operations as well
121+
sobel_output = compute_sobel(circle)
122+
123+
writedlm("sobel_output.txt", sobel_output)
124+
125+
end

0 commit comments

Comments
 (0)