- Implement a bilateral filter in Python.
- Given a noisy image, students will be able to adjust the parameters of a bilateral filter to achieve maximum noise reduction,
- Students will combine flash and no-flash photos using the cross-bilateral filter to generate high quality images in low-light conditions.

- Basic familiarity with the concept of filtering, and how it can be implemented using a convolution.
- Basic knowledge of programming in Python.

As our noble digital photographer ventures out into the wilderness
armed with their trusted cameras, they are repeatedly hounded in their
quest for truth by the vile and iniquitous *Noise Monster*. In this lab
we will cover one of the tools the lionhearted photographer
can use to combat this fiendish interloper: the bilateral filter.

If you have past experience with filtering, you may be aware of how it can be used to remove unwanted frequencies (a fancy term for noise) from an input signal. What, then, makes the bilateral filter so special? To appreciate the advantages offered by a bilateral filter, it would be useful to first understand why an old-fashioned Gaussian filter proves insufficient when it comes to processing images.

Let's load up a noisy image in Python and convolve it with a Gaussian filter:

While the amount of noise has certainly been reduced, the astute viewer will notice that so has most of the detail around object edges. Everything appears blurry now. This is why a Gaussian—or any other convolutional filter, for that matter—is not very useful as a noise supression tool for images: it filters indiscriminately, leading to loss of resolution. What we want, instead, is a filter that works on the `internal' pixels of an object, but leaves the edges alone.

**Q1: The following two figures show a crude decomposition of our
image into edge (left) and non-edge regions (right). Do you notice anything different about the color composition of
the regions in each case? How many different colors can you count in each region?**

In class, we saw how the *median filter* can handle edge-aware filtering. Let us now see how a bilateral filter handles the same image.

The result, in this case, is much closer to what we want.

** Q2: Based on the color composition of the edge and non-edge regions
shown above, what do you think the bilateral filter is doing?**

The bilateral filter is a Gaussian that acts strongly on regions of uniform color, and lightly on regions with high color variance. Since we expect edges to have high color variance, the bilateral filter acts as an edge-preserving or edge-aware filter.

Let us dive into the details of how the bilateral filter works. Once again, the Gaussian distribution provides a good launching pad. First, let's remind ourselves of its 2D form—here unnormalized. Let's assume the Gaussian of variance \(\sigma\) is centered at a point \(p\) and we're sampling the PDF for a point at \(q\): \(q\) is \(d\) away from \(p\).

The Gaussian filter applied at a pixel index \(p\) in image \(I\) can be written as:

where \(I_q\) is the value at pixel index \(q\), and \(S\) denotes a small neighborhood
of pixels around \(p\). The \(\mathcal{N}\) denotes a Gaussian (which is also called a \(\mathcal{N}\)ormal
distribution), and \(w\) is the normalization factor—this preserves image brightness during filtering.
A simple Python implementation of this equation is provided in Listing 2. While this isn't the most efficient
implementation of Gaussian filtering—we would typically use a `numpy`

implementation—it is helpful for understanding
the Gaussian filter and the changes we need to make to a Gaussian filter to obtain a bilateral filter.

*Note:* If you have trouble mentally mapping this code to your conceptual understanding of filtering, feel free to review the course slides on filtering.

Let's try to speed the code up a bit. We will do this using a python library
called *Cython*. The main source of latency in our Python code of Listing 2 is the nested for loops.
Cython allows us to run these as compiled C code, which is very fast.
You can install Cython using Pip:

The changes that we need to make to the python code from Listing 2 are shown in Listing 4:

- We import C versions of specific functions
`exp`

and`sqrt`

, and - We add the keyword
`cdef`

to our variable definitions, and we strongly type our variables with`int`

and`double`

.

Every time we make changes to our code, and before we run our code, we must ask Cython to compile the changes we made.
We've set this up for you: download `setup.py`

and execute the following command:

Try running `$python our_gaussian.py`

now. You should see a *significant* improvement in speed.
You can do this for any Python code. Sometimes it's impossible to vectorize code into `numpy`

functions. Cython gives us another option to make Python code faster.

In Equation 1, the contribution of each *I _{q}* to the
filtered value of

It turns out that one *controllable* way to modulate filtering in edge regions is to
weight the contribution of each *I _{q}* to

- Modify the code in Listing 4 to include an additional weight based on the color difference between pixels, as shown in Equation 2.
- Modify the normalization factor \(w\) in your code such that the filter weights still sum to 1. How can we achieve this across the three color channels?

- The weight in Equation 2 decreases linearly with respect to the difference in color. Modify your code so that weight is exponential with respect to the difference in color. Does this make the results better or worse?
- Make sure to correct your normalization factor \(w\).

Now, let's go a step farther, and weight the color difference based on a *Gaussian*.

- Modify the code so that the weight is Gaussian with respect to the color difference.
- Introduce a new sigma to control this second Gaussian on the color information. You may need to increase sigma depending on your kernel size!

Congratulations! You have successfully implemented a bilateral filter:

Since the standard definition uses a Gaussian as the weight decay
function, bilateral filters are commonly defined by the variance
values of the two Gaussians that determine the
weights: \(\textrm{BF}(\sigma_1, \sigma_2)\). This is what
makes them powerful *controllable* edge-aware filters through varying
these sigma values. We will see this in the next section.

**Freedom!** Now that you understand what role each term in Equation 3 plays, you need
not be constrained by convention. You may choose to use any
combination of filters—box, tent, or even something novel
you've created yourself—depending on the task at hand, and how
the filter affects the result. The important thing to keep in mind is
that the weight is now a function of both color, and space. In fact,
you could even come up with a *trilateral* filter that used a
third domain, say, the image gradients, in addition to the color and
spatial distance. The last task of this lab encourages you to think
about designing such techniques.

In the last section we saw how a bilateral filter is controlled by two parameters: the variance of the Gaussian filter in the spatial domain, and the variance of the Gaussian in the color domain. We will refer to these parameters as the spatial and range sigmas respectively. Moreover, as the two Gaussians that make up the bilateral filter are being multiplied together, a value close to zero for either one of the sigmas makes the whole filter ineffective. On the other hand, If the two values are too large, the bilateral filter no longer preserves edges well.

**Q3: Based on the image content, we may choose to assign a relatively
higher value to the spatial sigma, or to the range sigma. Given our
understanding of the bilateral filter, what kind of images do you
think would require a higher spatial sigma, and what kind a higher
range sigma?**

To put our intuition about the effect of each parameter into practice, we will denoise the images in Dataset I. In each case, our goal will be to remove as much noise from the image as possible, while maintaining the original image edges. Our measure of success will be the L2 distance from a ground truth image (also included in the dataset).

- Write a script that varies the spatial and range sigma values for an image in Dataset 1, and outputs a 5x5 2D grid of filtered images showing the variation. An example is here: parameter_grid.png (GB = simple Gaussian blur). The
`pyplot.subplot`

function in the`matplotlib`

library will help lay out the figure. You might need to experiment with your parameter values. What do you notice about the result? - Apply the bilateral filter you wrote for Task I to each image in Dataset I. Determine the approximate range of the spatial and range sigma values that minimize the L2 distance from the ground truth.
- List a feature of each image that determines the unique spatial/range sigmas for each image.

*Note: This is an optional component of the lab that demonstrates a useful application of the bilateral filter beyond denoising.*

Task II demonstrated how image content determines the effectiveness of
bilateral filtering. In some instances, however, a variant called
the *cross bilateral filter* can be used to double-cross the
aforementioned rule. The idea behind cross bilateral filtering is that
the range and spatial Gaussian filters can act on two *different*
images! This proves especially useful when we have multiple
photographs of a single scene captured with different camera
settings.

The cross bilateral filter is not restricted to the task of denoising; it can be used to improve the overall quality of a photo. In fact, your next task explores this very feature.

Images captured in low-light situations tend to be noisy, and lack sharp edge resolution. A common solution to this problem is to use a flash. However, while the flash image does have sharp edge definition, it suffers from an unpleasing direct-lighting effect.

- Capture a pair of images of the same scene, in low-light conditions, with and without a flash. You may want to use a tripod to make sure your camera does not move between captures.
- Modify your bilateral filter code so that the range and spatial Gaussians are applied to two different input images.
- Use your implementation of the cross bilateral filteral to combine the flash/no-flash image pairs into a high quality photo. An example is shown in Fig 7.
- What parameters did you use, and why does this approach work?

We started off with the problem of filtering an image to remove noise while maintaining sharp detail in the edges. We saw how the bilateral filter solves this problem by combining two different sources of information: the weights of the filter were determined by color proximity, in addition to spatial proximity. The cross bilateral filter took the idea of drawing on multiple sources of information a step further. The RGB values for the filter came from two altogether different photos of the same scene.

In fact, this idea recurs throughout the field of computational photography. An imaging problem is solved by drawing on information from different representations, and even modalities. Together, the combined sources make up for their individual shortcomings, and provide a more comprehensive (or perhaps simply a more beautiful) picture of the world.

Please upload your Python code, input/result images, and any notes of interest as a PDF to Gradescope. Please use writeup.tex for your submission.

Bilateral filters are useful for *so many* image processing applications. Sylvain Paris has a gentle introduction to the topic and its many applications here.

This lab was developed by Numair Khan and the 1290 course staff. Numair created the logos and drew Sir Cameralot, with inspiration from HergĂ© and Tintin. Thanks to Sylvain Paris for the parameter grid example.