BookNotesSlidesP-sets Syllabusdraft · v0.2.49.1
◀▲⌂▶
💬Comments welcome. To leave a note, select any text and click the note / highlight button that pops up — or open the panel with the tab at the top-right (‹). Notes are visible only inside our private review group.

45.5 Problem Set 2 — Convolution and the bilateral filter⧉

Problem Set 2: Convolution and the Bilateral Filter. Due Wednesday September 30, 9:00pm.

45.5.1 Summary⧉

Important: This problem set is longer than the previous two. Reading the slides carefully can easily save you some headaches!

Image processing functions from the previous pset have been moved to basicImageManipulation.cpp and Image.cpp.

45.5.2 Smart Accessor⧉

When resampling images or performing neighborhood operations, we might end up accessing a pixel at x,y coordinates that are outside the image. To handle such cases gracefully, it is good practice to write accessors that take x,y as inputs and check them against the bounds of the image before retrieving a value.

Multiple options are possible when a pixel is requested outside the bounds; but the most common are to return a black pixel or the value of the closest valid pixel (i.e.\ of the closest edge). For the latter, clamp x to [0..width-1] and y to [0..height-1], and perform the lookup there.

Black padding (a.k.a.\ zero-padding) is more appropriate for applications such as scale and rotation (which we'll see in the next assignment) whereas edge-pixel padding generally looks better for convolution.

45.5.3 Blurring⧉

In the following problems, you will implement several functions in which you will convolve an image with a kernel. This will require that you index out of the bounds of the image. Handle these boundary effects with the smart accessor from Section Smart Accessor. Also, process each of the three color channels independently.

We have provided you with a function Image impulseImg(int k) that generates a $k \times k \times 1$ grayscale image that is black everywhere except for one pixel in the center that is completely white. If you convolve a kernel with this image, you should get a copy of the kernel in the center of the image. An example of this can be seen in in a2_main.cpp. It can also be useful to have simple images, such as an all white Image e.g. Image im(w,h,c); im.set_color(1.0, 1.0, 1.0);.

Box blur⧉

Problem set 2, Cambridge2
Problem set 2, Cambridge2
Problem set 2, boxBlurred1
Problem set 2, boxBlurred1

Figure 1. Result of blurring an image with a box width of k=9 and clamp=true

General kernels⧉

Now, let us implement a more general convolution function that uses an arbitary kernel. From now on, we'll perform this using the Filter class. This class contains a buffer of floats that represent the filtering kernel, together with its spatial footprint.

To create a filter use the constructor Filter(const vector<float> &fData, int fWidth, int fHeight). This takes in a row-major vector containing the values of the kernel (just like in the Image class), and the width and height of the kernel respectively (to avoid alignment and coordinate round-off issues, we'll use odd width and height). See a2_main.cpp for an example.

Note that the kernels are not flipped ahead of time. You may find the original format of convolution helpful, i.e. $(I\otimes g)(x)=\int_{x'}I(x')g(x-x')\mathrm{d}x'$. Tip: To access the (x,y) location in a kernel from within the class type operator()(x, y).

Your test images.

Put the images you make in Input/debug/ and name them <part>_<n>_<what>.png, lowercase, at most $256\times256$ pixels: conv_1_… png to conv_3_… png here, and bila_1_… png to bila_3_… png later in the set. make prepare copies Input/debug/ into your submission, so they travel with your code. Write the filtered versions into Output/ under the same name with _out added, and keep the code that generates the images in a2_main.cpp.

Pay attention to indexing, $(0,0)$ denotes the upper left corner of the Image, but for our kernels we want the center to be in the middle. This means you might need to shift indices by half the kernel size. Test your function with impulseImg, a constant image and real images of your choice.

Gradients⧉

Image gradients measure the variation in intensity between neighboring pixels. We estimate gradients at each of pixel in the image using finite differences. Since this operation is identical regardless of the spatial location, we can write it as a convolution. A kernel often used to estimate the gradient is the Sobel kernel. Here are their weights(Note that the Sobel kernel is not scaled so that the output keeps the same range as the input, i.e. $[0,1]$. For visualizing image gradients, e.g. filtered results using a horizontal Sobel kernel as in Fig. 2(b), you might want to rescale the convolution output to [0,1] considering the possible output range convolved by this filter ($[-4,4]$ for the Sobel filter).) for the horizontal and vertical components of the gradient are respectively:

$$ \left[\begin{array}{ccc} -1 & 0 & 1 \\ -2 & 0 & 2 \\ -1 & 0 & 1 \\ \end{array}\right]\text{ and } \left[\begin{array}{rrr} -1 & -2 & -1 \\ 0 & 0 & 0 \\ 1 & 2 & 1 \\ \end{array}\right] $$

The gradient magnitude is defined as the square root of the sum of the squares of the two components(No need to rescale it to $[0,1]$, just keep it as the original range.) :

$$ ||\nabla I|| = \sqrt{\left|\frac{\partial I}{\partial x}\right|^2 + \left|\frac{\partial I}{\partial y}\right|^2} \;. $$
Problem set 2, lounge view
Problem set 2, lounge view
Problem set 2, sobelFiltered
Problem set 2, sobelFiltered
Problem set 2, imageMagnitude
Problem set 2, imageMagnitude

Figure 2. Result of filtering with a horizontal Sobel kernel and computing the gradient magnitude.

Gaussian Filtering⧉

The Gaussian filter is ubiquitous in image processing. It is a key building block of many algorithms. We'll start by implementing a one-dimensional version of the Gaussian filter, then the 2D extension. The mathematical expression for a continuous normalized 1D Gaussian is:

$$ f(x) = \frac{1}{\sqrt{2\pi \sigma^2}}\exp\left(-\frac{x^2}{2\sigma^2}\right) $$

Notice that $\frac{1}{\sqrt{2\pi \sigma^2}}$ is the normalization factor.

In practice, we will not care about this exact normalization factor since we are working with truncated and discrete Gaussians. Instead, we'll make sure our discrete filter is normalized, i.e. its weights sum to one.

1D Horizontal Gaussian filtering

2D Gaussian Filtering

Let us now extend our skinny 1D Gaussian to a the full-blown 2D version(Similar to the 1D case, $\propto$ means “proportional to” considering a normalization factor, which should be determined after truncation and discretization.), e.g.

$$ f(x,y) \propto \exp\left(-\frac{x^2 + y^2}{2\sigma^2} \right) \;. $$

Separable 2D Gaussian Filtering

Filtering our images with a 2D Gaussian kernel gives us a nice blur effect, but it is computationally costly. Fortunately we can exploit the separability of the Gaussian kernel, and make the 2D convolution almost as fast as its 1D counterpart.

Verify that you get the same result with the full 2D filtering as with the separable Gaussian filtering.

You can compare the running time of your two implementations following the example in the starter code in a2_main.cpp.

Problem set 2, gaussBlurHorizontal
Problem set 2, gaussBlurHorizontal
Problem set 2, gaussBlur2D seperable
Problem set 2, gaussBlur2D seperable

Figure 3. Result of blurring an image with a horizontal and 2D Gaussian kernel for sigma = 3.0, truncate=3.0, clamp=true

Sharpening⧉

$$ I_\mathtt{sharpen}=I+\mathtt{strength}\times(I-I\otimes g)\,. $$

45.5.4 Denoising using Bilateral Filtering⧉

The Bilateral filter is defined as:

$$ \begin{aligned}\notag I_{\text{out}}(x,y) = & \frac{1}{Z} \sum_{x',y'} G([x - x', y - y'], \sigma_D) \times G(\|I(x,y) - I(x',y')\|,\sigma_R) \times I(x',y') \\ \notag \text{with }Z = & \sum_{x',y'}G([x - x', y - y'], \sigma_D) \times G(\|I(x,y) - I(x',y')\|,\sigma_R)\end{aligned} $$

where $I$ is the input image, $Z$ is a normalization factor and $G$ is a Gaussian kernel. The notation $G([x-x', y-y'], \sigma_D)$ is shorthand for $G(x-x', \sigma_D) \times G(y-y', \sigma_D)$. The bilateral filter is very similar to a convolution, but the “kernel” varies spatially and depends on the color difference between a pixel and its neighbors. The different standard deviations trade-off the spatial (domain), $\sigma_D$, and color differences (range), $\sigma_R$. It is important to note that unlike other kernels, the normalization here depends on the pixel values and cannot be computed ahead of time. The range Gaussian should be computed using the 3D Euclidean distance in RGB, that is $\|I(x,y) - I(x',y')\|$ above is a shorthand for $\sqrt{\sum_c (I(x,y,c) - I(x',y',c))^2}$ (loop over three color channels).

Warning: Because this is not a straightforward convolution, you will not be able to use the ${\tt Filter}$ class implementation for this part. The normalization process should proceed only once, that is ignoring it for each individual Gaussian kernels but applying it at the last step for all three color channels.

Try your filter on the provided noisy image lens as well as on simple test cases (e.g. an image with a rectangle).

Problem set 2, lens
Problem set 2, lens
Problem set 2, bilaterial RGB
Problem set 2, bilaterial RGB
Problem set 2, bilaterial YUV
Problem set 2, bilaterial YUV

Figure 4. Results of denoising using a bilateral filter. RGB Bilateral Filtering: bilateral(im, 0.1,1.0,3.0,true); YUV Bilateral Filtering: bilaYUV(im, 0.1,1.0,4.0,3.0,true);

6.8370 only (or 5% Extra Credit): YUV version⧉

We want to avoid chromatic artifacts by filtering chrominance more than luminance. This is because the human visual system is more sensitive to low frequencies in the chrominance components.

In all cases, make sure you compute the range Gaussian with respect to the full YUV coordinates, and not just for the channel you are filtering. We recommend a spatial sigma four times bigger for U and V as for Y.

45.5.5 Extra credit⧉

Here are ideas for extensions you could attempt, for 5% each. At most, on the entire assignment, you can get 10% of extra credit:

45.5.6 Submission⧉

Turn in your files to the online submission system and make sure all your files are in the asst directory under the root of the zip file. If your code compiles on the submission system, it is organized correctly. The submission system will run code in your main function, but we will not use this code for grading. The submission system should also show you the image your code writes to the ./Output directory

In the submission system, there will be a form in which you should answer the following questions: