Skip to content

Commit 8f3e639

Browse files
committed
Adding Convolution Chapter
1 parent 782b13b commit 8f3e639

31 files changed

+943
-3
lines changed

.gitignore

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -515,4 +515,5 @@ vscode/
515515
# Coconut compilation files
516516
**/coconut/*.py
517517

518+
# aspell
518519
*.bak

SUMMARY.md

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,10 @@
1313
* [Affine Transformations](contents/affine_transformations/affine_transformations.md)
1414
* [Bit Logic](contents/bitlogic/bitlogic.md)
1515
* [Taylor Series](contents/taylor_series_expansion/taylor_series_expansion.md)
16+
* [Convolutions](contents/convolutions/convolutions.md)
17+
* [Convolutions in 1D](contents/convolutions/1d/1d.md)
18+
* [Convolutions of images (2D)](contents/convolutions/2d/2d.md)
19+
* [Convolutional Theorem](contents/convolutions/convolutional_theorem/convolutional_theorem.md)
1620
* [Tree Traversal](contents/tree_traversal/tree_traversal.md)
1721
* [Euclidean Algorithm](contents/euclidean_algorithm/euclidean_algorithm.md)
1822
* [Monte Carlo](contents/monte_carlo_integration/monte_carlo_integration.md)

contents/IFS/IFS.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -225,7 +225,7 @@ The code examples are licensed under the MIT license (found in [LICENSE.md](http
225225

226226
##### Text
227227

228-
The text of this chapter was written by [James Schloss](https://github.com/leio) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
228+
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).
229229

230230
[<p><img class="center" src="../cc/CC-BY-SA_icon.svg" /></p>](https://creativecommons.org/licenses/by-sa/4.0/)
231231

contents/computus/computus.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -318,7 +318,7 @@ The code examples are licensed under the MIT license (found in [LICENSE.md](http
318318

319319
##### Text
320320

321-
The text of this chapter was written by [James Schloss](https://github.com/leio) and is licensed under the [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by-sa/4.0/legalcode).
321+
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).
322322

323323
[<p><img class="center" src="../cc/CC-BY-SA_icon.svg" /></p>](https://creativecommons.org/licenses/by-sa/4.0/)
324324

contents/convolutions/1d/1d.md

Lines changed: 400 additions & 0 deletions
Large diffs are not rendered by default.

contents/convolutions/2d/2d.md

Lines changed: 185 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,185 @@
1+
# Convolutions on Images
2+
3+
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.
4+
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).
6+
7+
The extension of one-dimensional convolutions to two dimensions requires a little thought about indexing and the like, but is ultimately the same operation.
8+
Here is an animation of a convolution for a two-dimensional image:
9+
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>
16+
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 discussion in the previous section.
19+
In code, using simple boundaries, a two-dimensional convolution might look like this:
20+
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 axis of the output domain *and* the filter.
27+
Because this indexing is slightly non-trivial, we felt it was worth writing a separate chapter for two-dimensional convolutions.
28+
29+
It is worth highlighting common filters used for convolutions of images.
30+
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 naïve edge detection or image derivatives.
31+
32+
## The Gaussian kernel
33+
34+
The Gaussian kernel serves as an effective *blurring* operation for images.
35+
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.
36+
37+
As a reminder, the formula for any Gaussian distribution is
38+
39+
$$
40+
g(x,y) = \frac{1}{2\pi\sigma^2}e^{-\frac{x^2+y^2}{2\sigma^2}},
41+
$$
42+
43+
where $$\sigma$$ is the standard deviation and is a measure of the width of the Gaussian.
44+
A larger $$\sigma$$ means a larger Gaussian; however, remember that the Gaussian must fit onto the filter, otherwise it will be cut off!
45+
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.
46+
As a general rule of thumb, the larger the filter and standard deviation, the more "smeared" the final convolution will be.
47+
48+
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.
49+
50+
{% method %}
51+
{% sample lang="jl" %}
52+
[import:30-47, lang:"julia"](../code/julia/2d_convolution.jl)
53+
{% endmethod %}
54+
55+
Though it is entirely possible to create a Gaussian kernel whose standard deviation is independent on the kernel size, we have decided to enforce a relation between the two in this chapter.
56+
As always, we encourage you to play with the code and create your own Gaussian kernels any way you want!
57+
As a note, all the kernels will be scaled at the end by the sum of all internal elements.
58+
This ensures that the output of the convolution does not have an obnoxious scale factor associated with it.
59+
60+
Below are a few images generated by applying a kernel generated with the code above to a black and white image of a circle.
61+
62+
<p>
63+
<img class="center" src="../res/circle_blur.png" style="width:100%" />
64+
</p>
65+
66+
67+
In (a), we show the original image, which is just a white circle at the center of a $$50\times 50$$ grid.
68+
In (b), we show the image after convolution with a $$3\times 3$$ kernel.
69+
In (c), we show the image after convolution with a $$20\times 20$$ kernel.
70+
Here, we see that (c) is significantly fuzzier than (b), which is a direct consequence of the kernel size.
71+
72+
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.
73+
74+
## The Sobel operator
75+
76+
The Sobel operator effectively performs a gradient operation on an image by highlighting areas where a large change has been made.
77+
In essence, this means that this operation can be thought of as a naïve edge detector.
78+
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 a final output array.
79+
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.
80+
Each gradient will be created by convolving our image with their corresponding Sobel operator:
81+
82+
$$
83+
\begin{align}
84+
S_x &= \left(\begin{bmatrix}
85+
1 \\
86+
2 \\
87+
1 \\
88+
\end{bmatrix} \otimes [1~0~-1]
89+
\right) = \begin{bmatrix}
90+
1 & 0 & -1 \\
91+
2 & 0 & -2 \\
92+
1 & 0 & -1 \\
93+
\end{bmatrix}\\
94+
95+
S_y &= \left(
96+
\begin{bmatrix}
97+
1 \\
98+
0 \\
99+
-1 \\
100+
\end{bmatrix} \otimes [1~2~1]
101+
\right) = \begin{bmatrix}
102+
1 & 2 & 1 \\
103+
0 & 0 & 0 \\
104+
-1 & -2 & -1 \\
105+
\end{bmatrix}.
106+
\end{align}
107+
$$
108+
109+
The gradients can then be found with a convolution, such that:
110+
111+
$$
112+
\begin{align}
113+
G_x &= S_x*A \\
114+
G_y &= S_y*A.
115+
\end{align}
116+
$$
117+
118+
Here, $$A$$ is out input array or image.
119+
Finally, these gradients can be summed in quadrature to find the total Sobel operator or image gradient:
120+
121+
$$
122+
G_{\text{total}} = \sqrt{G_x^2 + G_y^2}
123+
$$
124+
125+
So let us now show what it does in practice:
126+
127+
<p>
128+
<img class="center" src="../res/sobel_filters.png" style="width:100%" />
129+
</p>
130+
131+
In this diagram, we start with the circle image on the right, and then convolve them with the $$S_x$$ and $$S_y$$ operators to find the gradients along $$x$$ and $$y$$ before summing them in quadrature to get the final image gradient.
132+
Here, we see that the edges of our input image have been highlighted, showing outline of our circle.
133+
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 }}.
134+
135+
In code, the Sobel operator involves the same operations: first finding the operator in $$x$$ and $$y$$ and then applying them with a traditional convolution:
136+
137+
{% method %}
138+
{% sample lang="jl" %}
139+
[import:49-63, lang:"julia"](../code/julia/2d_convolution.jl)
140+
{% endmethod %}
141+
142+
With that, I believe we are at a good place to stop discussions on two-dimensional convolutions.
143+
We will definitely return to this topic in the future as new algorithms require more information.
144+
145+
## Example Code
146+
147+
For the full code in this section, we have modified the visualizations from the [one-dimensional convolution chapter](../1d/1d.md) to add a two-dimensional variant for an image of random white noise.
148+
We have also added code to create the Gaussian kernel and Sobel operator and apply it to the circle, as shown in the text.
149+
150+
{% method %}
151+
{% sample lang="jl" %}
152+
[import, lang:"julia"](../code/julia/2d_convolution.jl)
153+
{% endmethod %}
154+
155+
<script>
156+
MathJax.Hub.Queue(["Typeset",MathJax.Hub]);
157+
</script>
158+
159+
### Bibliography
160+
161+
{% references %} {% endreferences %}
162+
163+
## License
164+
165+
##### Code Examples
166+
167+
The code examples are licensed under the MIT license (found in [LICENSE.md](https://github.com/algorithm-archivists/algorithm-archive/blob/master/LICENSE.md)).
168+
169+
##### Images/Graphics
170+
- 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).
171+
- 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).
172+
- 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).
173+
- 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).
174+
175+
##### Text
176+
177+
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).
178+
179+
[<p><img class="center" src="../../cc/CC-BY-SA_icon.svg" /></p>](https://creativecommons.org/licenses/by-sa/4.0/)
180+
181+
##### Pull Requests
182+
183+
After initial licensing ([#560](https://github.com/algorithm-archivists/algorithm-archive/pull/560)), the following pull requests have modified the text or graphics of this chapter:
184+
- none
185+
Lines changed: 74 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,74 @@
1+
using DelimitedFiles
2+
using LinearAlgebra
3+
4+
function convolve_cyclic(signal::Array{T, 1},
5+
filter::Array{T, 1}) where {T <: Number}
6+
7+
# output size will be the size of sign
8+
output_size = max(length(signal), length(filter))
9+
10+
# convolutional output
11+
out = Array{Float64,1}(undef,output_size)
12+
sum = 0
13+
14+
for i = 1:output_size
15+
for j = 1:output_size
16+
if (mod1(i-j, output_size) <= length(filter))
17+
sum += signal[mod1(j, output_size)] * filter[mod1(i-j, output_size)]
18+
end
19+
end
20+
21+
out[i] = sum
22+
sum = 0
23+
24+
end
25+
26+
return out
27+
end
28+
29+
function convolve_linear(signal::Array{T, 1}, filter::Array{T, 1},
30+
output_size) where {T <: Number}
31+
32+
# convolutional output
33+
out = Array{Float64,1}(undef, output_size)
34+
sum = 0
35+
36+
for i = 1:output_size
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]
40+
end
41+
end
42+
43+
out[i] = sum
44+
sum = 0
45+
end
46+
47+
return out
48+
end
49+
50+
function main()
51+
52+
# sawtooth functions for x and y
53+
x = [float(i)/200 for i = 1:200]
54+
y = [float(i)/200 for i = 1:200]
55+
56+
# Normalization is not strictly necessary, but good practice
57+
normalize!(x)
58+
normalize!(y)
59+
60+
# full convolution, output will be the size of x + y
61+
full_linear_output = convolve_linear(x, y, length(x) + length(y))
62+
63+
# simple boundaries
64+
simple_linear_output = convolve_linear(x, y, length(x))
65+
66+
# cyclic convolution
67+
cyclic_output = convolve_cyclic(x, y)
68+
69+
# outputting convolutions to different files for plotting in external code
70+
writedlm("full_linear.dat", full_linear_output)
71+
writedlm("simple_linear.dat", simple_linear_output)
72+
writedlm("cyclic.dat", cyclic_output)
73+
74+
end

0 commit comments

Comments
 (0)