% TODO: 
% talk about spectra of different 


\documentclass[11pt,twoside,titlepage]{article} % Document font size and equations flushed left

\usepackage[english]{babel} % Specify a different language here - english by default
\renewcommand{\familydefault}{\sfdefault}
\usepackage{caption}
\usepackage{subcaption}
\usepackage{tikz}
\usepackage{amsmath}
\usepackage[a4paper,left=1.9cm,top=2.0cm,right=1.9cm,bottom=2.0cm,bindingoffset=1.9cm]{geometry}
\renewcommand{\baselinestretch}{1.5}
\usepackage{natbib}
\usepackage{hyperref}
\usepackage{fancyhdr}
\usepackage{mdframed}
\usepackage{pdfpages}
\usepackage{siunitx}

\newmdenv[
  backgroundcolor=gray!20,
  skipabove=\topsep,
  skipbelow=\topsep,
  align=center,
]{skinnybox}


\begin{document}
\sisetup{per-mode = reciprocal}
\pagestyle{fancy}
\fancyhead[HL]{Natalie Ko-Ferrigno}

\renewcommand{\floatpagefraction}{.9}%

\begin{titlepage}
  \begin{center}
    \huge{\textsc{What are the effects of linear mixing and how can they be mitigated for Triple-Exposure PTV}}
    \par
    \vspace{5mm}
    \large{Word count: 8884}\par
    \vspace{20mm}
    \large{Natalie Ko-Ferrigno}

    \vspace{10mm}
    Supervisor: Dr John Lawson
    
    \vspace{\fill}
    This report is submitted in partial fulfillment of the requirements for the MEng Aeronautical and Astronautical Engineering, Faculty of Engineering and Physical Sciences, University of Southampton
  \end{center}
\end{titlepage}

\begin{skinnybox}
    \centering
    \textbf{Declaration}
\end{skinnybox}

I, Natalie Ko-Ferrigno, declare that this thesis and the work presented in it are my own and has been generated by me as the result of my own original research.

\vspace{2em}
I confirm that:
\begin{enumerate}
  \item This work was done wholly or mainly while in candidature for a degree at this University;
  \item Where any part of this thesis has previously been submitted for any other qualification at this University or any other institution, this has been clearly stated;
  \item Where I have consulted the published work of others, this is always clearly attributed;
  \item Where I have quoted from the work of others, the source is always given. With the exception of such quotations, this thesis is entirely my own work;
  \item I have acknowledged all main sources of help;
  \item Where the thesis is based on work done by myself jointly with others, I have made clear exactly what was done by others and what I have contributed myself;
  \item None of this work has been published before submission
\end{enumerate}
\newpage

\begin{skinnybox}
    \centering
    \textbf{Acknowledgements}
\end{skinnybox}

I would like to thank my supervisor, Dr John Lawson for his guidance and insight into this field.
Without his assistance, this project would not have been as successful or enjoyable.

\newpage


\begin{skinnybox}
  \centering
  \textbf{Abbreviations}
\end{skinnybox}

\makebox[1.5cm]{DNS} Direct Numerical Simulation\par
\makebox[1.5cm]{PIV} Particle Image Velocimetry\par
\makebox[1.5cm]{PTV} Particle Tracking Velocimetry\par
\makebox[1.5cm]{LED} Light Emitting Diode\par
\makebox[1.5cm]{RGB} Red, Green, Blue\par
\makebox[1.5cm]{CFD} Computational Fluid Dynamics\par
\makebox[1.5cm]{RANS} Reynolds Averaged Navier Stokes\par
\makebox[1.5cm]{TE-PTV} Triple Exposure Particle Tracking Velocimetry\par
\makebox[1.5cm]{NAS} Network Attached Storage\par
\makebox[1.5cm]{WDM} Wave Division Multiplexing\par

\newpage


\begin{skinnybox}
  \centering
  \textbf{Abstract}
\end{skinnybox}

In this report, a method to allow for multiple colours of light to be used for a Particle Tracking Velocimetry (PTV) system is proposed.
This is through measuring how bright different colours of Light Emitting Diodes (LEDs) appear to each of the red, green and blue sensors on a camera. 
This technique is leveraged on computer-generated images that replicate a PTV setup to test how effectively this correction may operate.
Then, some images from real hardware are used. 

This report summarises past work that has had similar goals, to use LEDs and their ability to generate multiple colours of light to bring down the cost of PTV hardware.
Additionally, it looks at proposed methods to prepare this data for processing in typical PTV software.

The results are that the method used, called linear unmixing in microscopy, is successful in preparing synthetic data for processing.

\newpage


\tableofcontents
\newpage

\section{Introduction}

\subsection{PIV/PTV}
\label{sec:PIV_overview}

Particle Image Velocimetry (PIV) is a technique to determine a flow field in a region of interest of a fluid. 
To do this, the fluid is seeded with small particles, illuminated and photographed. \citet{annual_review}
These particles are selected so that they are close to neutrally buoyant, so that they follow the motion of any given fluid element almost exactly \citet{piv_guide}. 
Additionally, they should be selected so that they are visible to the camera. 
This can be done through reflectivity, or other methods such as fluorescence.
The bulk motion of a number of these particles can be reconstructed computationally for many locations in the flow, and a vector field of velocity can be calculated. 
This is known as Particle Image Velocimetry (PIV).
With three or more captures, acceleration can be found. 

Particle Tracking Velocimetry (PTV) is an extension of this method, where instead of seeing how flow moves in a fixed area, each particle is individually tracked to find the path that it takes between images. 
Because individual particles are followed, rather than seeing how fluid in an area moves, this is also known as Lagrangian Particle Tracking. \citet{STB}

At their most simple, these techniques involve imaging a laser sheet illuminating a fluid. 
This will give the two-dimensional positions of the tracer particles at a given instant. 
With more and more images, more and more flow derivatives can be calculated, assuming that the time interval between images is set. 
A double-frame camera allows for the calculation of velocity and, through Bernoulli's equation, pressure. 
Three frames allows for acceleration and so on. 

Three-dimensional data is also desirable, and there are many ways to achieve this. 
Multiple cameras can image the same flow-field, and with this data, tomographic techniques can be used to give three-dimensional, three component (3d3c) vector fields. 
In this technique, there is obviously the added cost of multiple cameras, but the accuracy of all three directions can be similar. 
Depth can be encoded through the wavelength of incident light, through projecting a spectrum of colours onto a whole volume (3d3c) \citet{RPIV}, or dual-plane PIV can be used to find the out-of-plane component of velocity with multiple laser sheets (2d3c), with an example given in \citet{Schlicke_2002}.

In all of these techniques, strong illumination is necessary so that high-speed particles do not become distorted from a circle on the sensor to a lengthened smear.
This is inherent to photography, lower-powered illumination needs to be on for longer to transfer the same amount of energy, and the longer the illumination is on, the further any fluid element will move.
For this reason, lasers have been the illumination source of choice for many years. 

In this report, a method will be described and tested to allow for greater access to PTV studies.


\subsubsection{Applications}

PTV is extremely useful when the full flow-field around an experiment is needed. 
LaVision, a commercial provider of PIV/PTV hardware and software, gives the example of an insect in flight, where CFD would be impractical due to the complex geometry, in the paper \citet{insect}. 
Helicopter aerodynamics is another area where the complexity of modelling a system makes a CFD study impractical (\citet{turb_book}). 

An accurate PIV/PTV system would allow for the information of a CFD study, but without the inherent errors from reduced-order modelling of any system, from potential flow to Direct Numerical Simulation (DNS).
With a system such as a flying helicopter, the size of the domain, complexity of the geometry and the difficulty modelling a rotor would make RANS difficult and time-consuming, whereas PIV/PTV will be significantly faster.

\subsubsection{Limitations}
\label{sec:PIV_limitations}

Historically, the use of lasers and high-speed cameras has led to costs of PTV being fairly expensive per experiment. 
This also limits the available volume to be studied, as a larger volume at a given accuracy will need more cameras and illumination. 

To give an idea of the costs, in a grant from the Air Force Office for Scientific Research \citet{cost_estimate}, \$192,355 was received for a camera, timing unit, a laser, and some accessories.
This document was found when trying to find the costs of the brand of laser used in \citet{cost_2}, suggesting that this may still be typical. 
A more recent cost estimate was not used, as the cost of these systems is not widely publicised.
This money also covered the installation of this system, meaning that the equipment will cost less than this figure.
This system will clearly be high-quality, however the cost is clearly a factor.

In comparison, the camera system used in this work amounts to less than £90.
In a personal communication, Joe Schaefer who designed the LEDs and driver circuit, said that each driver and LED cost £40.
This puts the total cost of the system at around £130, a dramatic decrease in cost.

While the quality of the data will be orders of magnitude different, not all applications will need to massive quality and quantity of data from Glezer's system.

Prior work has been done to address these costs. 
For example \citet{Aguirre2017} used smartphones and LED panels to bring down the cost of imaging. 
Through doing this, there were issues owing to the fixed optics of a smartphone limiting the distortion-free area, however the velocity field in a vortex ring was successfully reconstructed.




\newpage
\section{Background}

\subsection{Triple exposure PTV}
In the past, double-exposure cameras were used to allow for simple velocity measurements in PIV systems as a single image would contain two time steps. 
The positions were found from reconstruction of the particles, and the timing between these positions was known. 
Because the change in position and the change in time are both known, the velocity could be determined. 

\begin{figure}[htb!]
  \centering
  \includegraphics[width=\textwidth]{images/te-ptv3.png}
  \caption{A synthetic image, representative of triple exposure PTV. Uncompressed file available \href{https://notes.michcode.space/figures}{online}}
  \label{fig:te-ptv_example}
\end{figure}

Triple-Exposure PTV (TE-PTV) allows for a typical colour camera to act like a triple frame camera.
A low-cost consumer camera is used, along with an RGB LED.
The LED is pulsed three different colours to give the three separate exposures.
This allows for the camera to act like a high-speed monochrome camera, if the pulses of light are sufficiently short, while not incurring the costs.
Figure \ref{fig:te-ptv_example} shows an example of one of these triple-exposure images.

When a series of these images are taken, this allows for information to be gained for an arbitrary number of time steps. 
This is known as time-resolved data.

To capture a single image, the timing diagram in Figure \ref{fig:timing} is followed. 
The camera shutter is opened, and some time is left to allow for this to happen. 
Then, the red LED is illuminated for some amount of time, $\Delta t_L$ after which it is extinguished. 
A time period of $\Delta t_o$ is given to allow for the fluid to flow. 
This is repeated for the other colours and then the shutter is closed. 
\begin{figure}[htb!]
  \centering
  \includegraphics[width=0.5\textwidth]{images/timing diagram.png}
  \hfill
  \begin{tikzpicture}
    \draw [black] (0, -1) -- (1.25, -1) -- (1.25, 0) -- (8.5, 0) -- (8.5, -1)-- (10.75, -1);

    \node at (-1.5, -0.75) {Camera shutter};
  \end{tikzpicture}
  \caption{Timing diagram for capturing colour PTV images} %\
  \label{fig:timing}
\end{figure}

$\Delta t_L$ must be sufficient to allow for enough light to be collected that a recognisable image can be formed, but short enough that the particles remain as defined shapes, as was briefly mentioned in Section \ref{sec:PIV_overview}.
Additionally, it is important that the collected light is not bright enough to saturate the image sensor, so that the response of any given pixel is as close to linear as possible.

\begin{figure}[htb!]
  \centering
  \includegraphics[width=\textwidth]{images/motion and timing.png}
  \caption{Graphic to demonstrate particle motion}
  \label{fig:motion}
\end{figure}

Figure \ref{fig:motion} demonstrates how this will function for an experiment, where a single particle is at different positions at different times, and this results in different colours. 
Note that while the particle is shown to be dark in this demonstration to aid visibility, a real particle would be white to allow for more light to be reflected.
Looking at any small cluster of particles in \ref{fig:te-ptv_example}, the same shape can be found in a different colour above or below. 
A TE-PTV algorithm will track each particle in every cluster to come up with a large number of velocity vectors for the whole flow.
To track the particle, there is at first an initialisation phase to establish the initial path. 
Then, when the velocity is known at \(t=n\), the position of the particle is marched forwards in time to \(t=n+1\), assuming that the particle will continue on a similar velocity to at \(t=n\).
Near where the particle is predicted to be at \(t=n+1\), all the nearby particles are found and evaluated as to if they are the same particle.
This process is then repeated. 
Significantly more detail on this procedure is given in \citet{STB}.

When scaling a system to multiple cameras, there will inevitably be delays in capturing an image as any computer network will inherently have delays. 
Through holding all shutters open and illuminating the flow field for a relatively short $\Delta t_L$ relative to the exposure time of the camera, all the images will be more likely to be synchronised.

\subsection {The Bayer mask}
A Bayer mask is a way of arranging colour filters on a camera sensor to accurately acquire colour information \citet{bayer}.
Camera sensors are typically monochromatic, with the wavelength of incident light having little impact on the value read out to the camera processor. 
To resolve colour information from a monochromatic sensor, colour filters are placed over each filter. 
Knowing the arrangement of these filters allows for the raw intensity values from each pixel on the camera sensor to be converted into colour data, allowing for colour images.
Each colour filter will have a particular transmission curve with some overlap across the visible spectrum to allow for a range of hues to be recorded. 
For monochromatic PTV experiments, Bayer masks reduce resolution, as only up to half of the pixels will be sensing to their full capability, so monochromatic cameras can be used. 
For TE-PTV, this colour information is vital to the experiment working.

\begin{figure}[htb!]
  \centering
  \includegraphics[width=0.5\textwidth]{images/bayer.png}
  \caption{Bayer filter representation}
  \label{fig:bayer}
\end{figure}

\subsection{Colour cross correlation}

To help reduce the cost of PTV systems, the use of colour to encode temporal information is investigated here. 
To do this, there will be pulses of different colours of light at different times. 
Because the timing of each pulse is known and the position will be reconstructed, the velocity and acceleration can be calculated. 
Repeating this method will allow for time-resolved data to be acquired at a low cost. 

As discussed, a colour camera will use a Bayer filter to recover colour information. 
However, a the filters used are almost always, to some degree, transparent at all visible wavelengths. 
This means that a green laser would be detected in the red and blue sensors, and that blue and red lasers will also be detected in all the other sensors.
Figure \ref{fig:spectrum} shows how each colour of sensor will respond when illuminated by a light source of the given wavelength. 
From this figure, it is evident that any wavelength will lead to a response in all sensors.
This image is a redrawing of the figure in \citet{sensor_datasheet}.

\begin{figure}[htb!]
  \centering
  \includegraphics[width=0.5\textwidth]{images/sensor_spectrum.png}
  \caption{A representation of the spectral response of the image sensor used} %TODO: change this to the actual spectrum
  \label{fig:spectrum}
\end{figure}

Additionally, LEDs are not monochromatic sources of light, meaning that even if the filters were perfectly opaque outside of their intended operating window, some of the spectrum may be detected by the other-coloured sensors regardless. 
Both of these combined leads to a pulse from a single colour typically being detected on all of the image sensors.

Because a particle imaged with the red LED will lead to a corresponding blue and green image, it may appear that this pixel is stationary. 
This is shown in Figure \ref{fig:spill-graphic}.
This is exactly the same image as Figure \ref{fig:motion}, however now the red and blue particles are detected in the green channel.
Green is picked arbitrarily here, and the same effect is visible in all channels to varying degrees. 



\begin{figure}[htb!]
  \centering
  \includegraphics[width=\textwidth]{images/spill explainer.png}
  \caption{A visual explainer of the origins and effects of coupling}
  \label{fig:spill-graphic}
\end{figure}

Because the particle appears in all channels at all time steps, a likely consequence of this effect is that an algorithm with no correction will detect three particles, where there was only one.
I will be naming these extra particles 'spectral particles' for the rest of this work.
These particles will likely appear near stationary, as the nearest particle in any other colour channel will be the corresponding spectral particle.
This will likely fluctuate, as not all real particles will cause enough mixing to lead to a spectral particle forming.

These are like the ghost particles described in \citet{annual_review}, however with some important distinctions. 
A ghost particle's energy mostly comes from one camera, meaning that they can be fairly easy to detect and remove. 
For a spectral particle, all cameras will show a particle at that position with zero velocity, meaning that additional cameras will not solve the problem.
A ghost particle will appear randomly in the flow field, and will only be present for a single frame, or in this case, a single colour channel.
A spectral particle will always appear very close to a position where a real particle is, and it will only exist for the length of one image, after which the real particle has moved elsewhere.



In the past, a variety of methods have been tried to minimise this effect. %\
The most simplistic method that will be discussed is to simply ignore low intensity particles as presented in \citet{Aguirre2017}. 
There is only a small amount of coupling, and there is a high signal to noise ratio, meaning that any spectral particles will not be anywhere near as bright as a real particle.
For a demonstration of other technology, this is adequate, and their work clearly demonstrates that this method works well for the hardware used. 
However, the circumstances that allow this to work will not necessary extend to most imaging systems.

In \citet{Aruga2020}, an interesting method is used based on the correlation coefficient between two colours at a peak. 
Notably, this only uses blue and green light, and only the blue channel is affected by coupling in this case. 
This method used the spatial correlation between particles, with a manually-selected upper limit of correlation coefficient that will remove a particle. 
A self described limitation with this method is that at high seeding densities, many real particles will also be rejected as two overlapping real particles will look like one particle. 

The method presented in \citet{Charonko} gives a matrix of correlation coefficients to show how illumination of one colour will lead to a detection in the different channels. 
For each pixel, this matrix is used to solve for what the intensity should be ideally. 
To find this matrix, the sensed values of each combination of colour channels were plotted against each other, and then a straight line passing through the origin was fitted to these colours. 
In this experiment, there were two types of particles of different sizes. 
These types of particles had different colouring, however they were both illuminated by white light so coupling was still an issue to be addressed. 
To remove the coupling, this matrix was inverted, and then each pixel of the image was multiplied by this matrix to find the decoupled value. 
This method is significantly more robust, allowing for higher seeing densities than Aruga's method, and applicable to more imaging systems than Aguirre's. 
Other methods were also presented, but will not be discussed here. 


\citet{wang2017} built on Charonko's work by increasing the order of the fitted function to a polynomial. 
Additionally, the experimental method was closer to what will be done here, where white particles are illuminated by different colours of light. 
The extension to a polynomial fit is justified in this paper because there is a curve to the data, however the second-order coefficient is around 5 orders of magnitude smaller than the next term. 


In this work, \citet{Charonko}'s method will be investigated further. 
Wang's additions do not seem necessary, as a linear fit works well for the intensities in this system, as shown in Section \ref{sec:removal}. 
Additionally, hardware changes are inevitable in the future, meaning that even if simply clipping low intensities works for this system, there is no guarantee of future success. 
High seeding densities are desireable, as described in Section \ref{sec:PIV_overview}, meaning that a correlation based method may not work. 
Additionally, red, green and blue light will all be used, and they all affect each other, meaning that this method quickly becomes overcomplicated. 

Results in Charonko do not mention this method, however Wang does present some results. 
For this reason, expected errors can be found from Wang. 

\subsection{Decoupling channels}
\label{sec:removing}

The technique of linear unmixing, shown in \citet{Charonko}, will be used. To do this, the intensities detected from the camera will be represented in vector $I$, with some coefficients in matrix $C$. 
Matrix $S$ gives the values after they have been decoupled.
Each of these coefficients represents how the intensity of one colour channel will be shown in the other channels. 
This gives the equation \(\textbf{S} = (a\underline{\textbf{C}}) \textbf{I}\), which can be inverted and solved.
\(\textbf{S}\) is the decoupled vector of different colours for each pixel, \(\underline{\textbf{C}}\) is a matrix to show the influence that each chanel has on every other channel, \(\textbf{I}\) is a vector to show the detected values of red, green and blue at each sensor on the camera, and \(a\) is a scaling factor that will be used to change how intense the coupling is.
This can then be inverted and solved, giving the image without spill. 
Due to the presence of noise, it is not possible to solve these equations exactly, however the hope is that these equations can be close to solved. 
To calculate the coefficients in $C$, some representative images were taken, keeping all of the camera settings as close to what would be used as possible. 
The peaks were identified on the image, and the coefficients found through a least-squares fitting of a linear function. 
Because the leading diagonal contained values of only 1, and because there was more coupling due to green, the values were deemed as reasonable.

Some sensitivity testing was done in the computational phase of the work to see if an underestimation of the values was likely to affect the results of experiments. 
This was fairly simple, $C$ was just multiplied by some constants between 0 and 1.5 and the solver was re-run.

Charonko's work focussed on microfluidics, rather than channel flows, meaning that some of the issues they faced may not be present, and the flow used here may have issues of its own. 
For example, they report issues with a particle being turned into a checkerboard due to the Bayer filter. 

\subsection{Linear unmixing in other fields}
When viewing organic material from tissues to individual bacteria through a microscope, fluorescent dyes are used to highlight important features or follow material as it travels through the cell. 
Sometimes, multiple dyes are used to track multiple different processes. 
In this case, the images are separated into a number of different images, with one for each layer. 

One of the methods of doing this is through spectral imaging \citet{spectral_imaging}. 
This is a technique where the whole spectrum is captured at every point on an image. 
This has the advantage of allowing many dyes to be used, as long as all the spectra of each dye is linearly independent. 
This technique is not usable in the context of PTV, as the authors note that this is a very slow technique. 
One of the methods described to acquire this spectral data is through using a monochromatic sensor, with bandpass filters to allow certain wavelengths to be measured. 
In principle, this is not too far from how a camera functions, meaning that some of the techniques and insight in this related field will be relevant. 

In \citet{Timo2005}, it is described that using $n$ detector channels, a maximum of $n$ different fluorophores can be separated, else there will be multiple solutions. 
However, if there are more detector channels than fluorophores, the accuracy of the bandpass filter is less important. 

The interpretation of this for TE-PTV is that for a standard camera, only three different LED colours should be used.  
Additionally, because the filters on a camera let through a wide range of wavelengths, using a second, dissimilar sensor may give more accurate decoupled channels. 
Through doing this, more LED colours could potentially be added and the effective frame rate of the cameras could be increased. 
The response of the sensor to a given wavelength cannot be fixed by the end user, meaning that significant improvements may be found by moving to an overdefined architecture.

Linear unmixing is also relevant in remote sensing applications \citet{remote_sensing}. 
In these use cases, a full spectrum of a point in space is known, and the spectrum that it reflects is a linear combination of the reflectances of the materials in the sampled area. 
Linear unmixing here is more used to separate these spectra to determine what the land use may be. 
Because the data here is close to continuous both spectrally and spatially, in comparison to the notably discrete data in TE-PTV, there are unlikely to be useful insights.

Wave-Division Multiplexing (WDM) has a similar problem, where upon demultiplexing a channel of wavelength \(\lambda_i\), the channels \(\lambda_{i\pm 1}\) may affect the signal.
This is called inter-channel crosstalk. 
Some of the concepts found in solving this issue may be applicable here. 
To get around the problem, filtering the signals that go into the multiplexer is suggested \citet{optical_ppt}. 
This means a bandpass filter is put in front of each input of the multiplexer.
To bring this to TE-PTV, this would mean finding a wavelength near each LED's peak that primarily activates one channel and using a bandpass filter for that one channel.
For higher-end systems, this may be feasible.
However, for a low cost system, this would likely be cost prohibitive. 
Because a non-negligible amount of the energy of the LED leads to coupling effects, as evidenced by the brightness of spectral particles, this energy would now be inaccessible to the camera.
This would necessitate higher-powered LEDs or longer exposure times, neither of which are desireable.

\subsection{Shadowgraphy}
In some work, the particle field is illuminated, and then reflected light is detected by the sensor. 
This can potentially lead to a more simple system. 
Because LED illumination is less bright than laser illumination, the field is instead illuminated from behind. 
This is because detecting a shadow is easier than detecting reflected light, because there is a higher signal to noise ratio. 
This will lead to the background also being images, which will need to be corrected for. 
In \citet{Aguirre2017}, this problem was already presented and resolved. 
Images of the system with no particles will be captured, and then these 'light' images can be subtracted from experimental images. 
Additionally, the image will be inverted so that it works with a more traditional image processing pipeline. 
Due to the fact that noise will be present in images, multiple frames of the illuminated background will be captured and averaged, so noise should not be a factor. 


\subsection{Aims and Objectives}

This work aims to reduce the effects of linear mixing in a Triple-Exposure PTV system to a point where data gathered from the system can be processed in typical software.

To do this, linear mixing will have to be simulated when generated synthetic images, to allow it to be studied. 
Then, a series of tests should be made to determine how effective any countermeasures are at reducing its impact.

Additionally, experimental work should be carried out to get a full set of tools to prepare data for use.
This requires a simple example of the full system to be built to capture data.
The data will be used to develop tools. 
\newpage
\section{Method}

\subsection{Measured quantities}
\label{sec:measured}
The quantities presented here are only relevant to the computational section. 
This is because the nature of the work in computational and experimental settings is different. 

Because the source velocities are known, all of these values can be compared to actual values before any processing has taken place, allowing error to be determined. 
This gives another source to compare solutions after decoupling, as well as previous work.

The computational section of this work shows how linear unmixing may work on experimental TE-PTV data, whereas the experimental data only shows how well linear unmixing will apply to any image.
For this reason, the quantities in this subsection are only applicable to the computational testing. 

In this testing, the value of \(a\), as defined in Section \ref{sec:removal}, was set to values of \(0.0\), \(0.5\), \(1.0\) and \(1.5\).
These values were chosen to reflect any errors introduced by the tracking algorithm and the decoupling algorithm, the best, typical, and worst case scenarios respectively.  
The study of turbulent channel flows has produced typical profiles of flow velocities and other statistics, which will be used to evaluate the accuracy of the decoupled data. 

The experimental section was to test how images from the camera may need to be processed to be ready to use in the tracking algorithm.

\subsubsection{Bulk velocity}

As the flow is turbulent and unsteady, it isn't meaningful to measure the velocity at one particular place. 
Instead, the bulk velocity is measured. 
To calculate this, bins of \SI{1}{\mm} are taken in the streamwise direction. 
The flow will look like a typical turbulent flow field for this bin. 
The mean is taken for each velocity in this slice. 
For each \SI{1}{mm} bin in the streamwise direction, this is done. 
This gives a series of velocities, and the mean is taken for these. 
This is repeated for each time step. 

This would give the same number as \(U_b = \dot V \times A\).

The bulk velocity was known for the synthetic data, meaning that this could be used to evaluate if the reconstruction was at all accurate.

\subsubsection{Velocity profile in outer units}
The profile of streamwise velocity against wall-normal distance is also plotted. 
It can be hard to draw anything quantitative from this graph, but it allows for easier interpretation of other parameters. 
The time-averaged profile is well-known, so can be compared to expected values.  
This profile is shown in \ref{fig:textbook_bl}

\begin{figure}[htb!]
  \centering
  \includegraphics[width=0.6\textwidth]{images/textbook_profile.png}
  \caption{A time-averaged plot of velocity along a turbulent channel flow, and a typical instantaneous turbulent channel flow. Taken from \citet{turb_book}}
  \label{fig:textboox_bl}
\end{figure}

The flow under study is turbulent, meaning that it will be unsteady. 
If only one instant of the flow is considered at one position, this will mean that the flow will almost certainly not look like the profile above.
To resolve this, multiple positions in the streamwise direction will be used. 
This will mean that more points are sampled, so the mean should correspond better to the expected profile. 
If this is not adequate, multiple instants in time can also be considered, so that transient structures can dissipate, or move downstream.

\subsubsection{Wall units}
A common method of investigating turbulent boundary layers is through the use of wall units, denoted by a +. 
These are non-dimensionalised and can be used to get meaningful data on the flow, with further details given in \citet{turb_book}. 

\begin{gather}
  u_\tau = \sqrt{\frac{\tau_{w}}{\rho}} \\
   u^+ = \frac{u}{u_\tau} \\ 
  y^+ = \frac{yu_\tau}{\nu}
\end{gather}


These will be another way to determine how effective the system is. 
Because a wall unit is smaller than an outer unit, this may be better for testing small regions of flow. 
Wall units are primarily focused on the interaction between the flow and a wall, and as the wall normal distance increases, there are generally fewer changes in a flow. 
For this reason, variables are typically plotted against \(\log y^+\).
So that at small values of \(y^+\), where there will be the most resolution on the graph, there is enough data to use, cosine spacing is used to sample the data.
This is a method where the width of the bins are not constant.
Their width is based on a cosine curve, so that there are very thin points near \(y^+=0\), and the width of the bin increases 
For bin number \(n\), the highest value of \(y\) in that bin will be \(\cos\left(\frac{n+1}{N}\right) * h\) and the lowest value will be \(\cos \left(\frac{n-1}{N}\right)\).

\subsubsection{Law of the wall}
The velocity profile is well defined, and a typical graph is shown in \ref{fig:wall_profile}.
To the left of the dashed lines is the log region, where viscosity dominates.
This is a very small region of the flow.
To the right of the dashed lines is the log-law region.
Between the lines is the buffer region, where the flow profile is undefined. 

Further information is given in \citet{turb_book}.

\begin{figure}[htb!]
  \centering
  \includegraphics[width=0.6\textwidth]{images/LOTW.png}
  \caption{Plot of the law of the wall, reproduced from \citet{turb_book}}
  \label{fig:wall_profile}
\end{figure}

\subsubsection{Velocity fluctuations}
\label{sec:vel_textbook}
Each $y$ location will have a corresponding time-average velocity, $\bar U(y)$. 
Then, the fluctuation of each particle can be found, $\bar U(y) - U(x, y, x, t)$. 
Finally, the root-mean squared value of \(U(x, y, z, t)\) can be found. 
Using these fluctuations, a variety of plots can be produced in both wall units and external units. 
Some of these are described in \citet{turb_book}. 
One example that will be used to verify the data is shown in Figure \ref{fig:uprmsp_yp}. 
This is taken from DNS data, as measuring a real flow to this level of accuracy may prove difficult, as has been discussed. 
The friction Reynolds number of the TE-PTV images is 1265, so the expected profile can be estimated from Figure \ref{fig:uprmsp_yp}.

\begin{figure}
  \begin{center}
    \includegraphics[width=0.95\textwidth]{images/uprmsp_yp.png}
  \end{center}
  \caption{Plots of $u'^{+}_{rms}$ against $\log y^+$, reproduced from \citet{turb_book} with permission, with source data from \href{https://torroja.dmt.upm.es/channels/data/statistics/}{Universidad Politécnica de Madrid}}
  \label{fig:uprmsp_yp}
\end{figure}


\subsection{Hardware and software used}
\subsubsection{Gigatrack}
Gigatrack is a set of MATLAB scripts to allow for PTV techniques, and it can generate synthetic data, as well as track particles. 
Gigatrack was used due to its low cost (free) in comparison to commercial software, and the ease of adaptability. 
Currently, there is no published work using Gigatrack, but the synthetic image generation and tracking is well developed. 

Gigatrack has the capabilities to define a camera, with variables such as numerical aperture, number of pixels vertically and horizontally and the distance between each pixel (the pitch).
There is a model for how a point source of light will propagate through a camera lens and be detected on the image sensor.
This camera model can then be placed along a simulated flow, with the number of cameras, their separation in the streamwise direction and their field of view height all defined. 
The Bayer mask is also simulated.
Source data from previously-validated flow is used to determine how the particles move, and this can easily be changed.
Figure \ref{fig:Gigatrack_flowchart} demonstrates the files used and what their purpose is.


\begin{figure}[htb!]
  \centering
  \includegraphics[width=\textwidth]{images/Gigatrack.png}
  \caption{Flowchart to demonstrate how Gigatrack functions}
  \label{fig:Gigatrack_flowchart}
\end{figure}

For all of the relevant graphs in this report, the settings used are described in Table \ref{tab:Gigatrack_settings}.
These were kept consistent throughout, as there is no sense in changing them, and this would invalidate the data.

\begin{table}[htb!]
  \begin{tabular}{|l|l|}
    \hline
    \textbf{Variable} & \textbf{Value}\\
    \hline
    \hline
    Sensor width (pixels) & \(4056\)\\
    \hline
    Sensor height (pixels) & \(3040\)\\
    \hline
    Numerical aperture (no unit) & \(4.0\)\\
    \hline
    Focal length (\si{\mm}) & \(20\)\\
    \hline
    Pixel pitch (\si{\mm}) & \(1.55\times 10^{-3}\)\\
    \hline
    Number of cameras & \(4\)\\
    \hline
    Separation in \(x\) direction (\si{\mm}) & \(134\)\\
    \hline
    Field of view height (\si{mm}) & \(304\) \\
    \hline
    Source data & Johns Hopkins University Turbulent Channel Flow data\\
    \hline
  \end{tabular}
  \caption{Settings used for all runs of Gigatrack}
  \label{tab:Gigatrack_settings}
\end{table}

The data source is from DNS data generated at Johns Hopkins University, and properties of the flow are given \citet{JHUTCF}.
It was rescaled so that the viscosity matches water.

This layout is demonstrated in Figure \ref{fig:syn_calib}.
In this figure, the particles are coloured according to their speed.
It is clear that the flow is unsteady, as is to be expected.
Particles are only shown when they appear in the field of view of at least 3 cameras. 
The field of view is displayed as a frustum the same colour as the camera.

\begin{figure}[htb!]
  \centering
  \includegraphics[width=0.6\textwidth]{images/syn_calib.png}
  \caption{Drawing of the experimental setup, with each particle coloured according to its speed}
  \label{fig:syn_calib}
\end{figure}

\subsubsection{Additions to Gigatrack}
A file was made to unmix images. 
It takes in the image with coupling, and then outputs images that have it removed. 
To do this, images are given, as is the \(C\) matrix.
The source image is loaded, and debayered to ensure that the colours are represented by a vector, so that this colour vector can be multiplied by \(C^{-1}\). 
This image cannot be saved as-is, because Gigatrack does not accept colour images, and only accepts monochrome data.
This has the benefit of a significantly smaller file size.
To convert the image back, the debayering process is reversed. 
Selecting the position from the relevant colour channel of the debayered image appears to give identical results. 
There does not appear to be a rebayer method in MATLAB, so this was the only way to do it.

Additionally, some files were written to automatically plot the graphs given later when given a directory.

Because many runs of Gigatrack were done, and each run involved multiple functions, a script was written to handle running all of these.

Finally, \verb|proc_rgblpt.m| outputs a file for \verb|batch_rgblpt.m| to process. 
In all cases, this was placed in the same directory and given the exact same name. 
When running multiple instances of Gigatrack in parallel, this lead to multiple processes trying to read and write data out of turn leading to undesired results, as is all too common in parallel computing. 
This file was placed in a different directory for each run to avoid this becoming an issue.

\subsubsection{Hardware setup}
For this experiment, an optical table was used to securely fasten various components and allow for the images to be taken. 
An LED was attached to a vertical aluminium extrusion. 
Then, a fresnel lens was placed in front of the LED to focus the light onto the camera sensor. 
The camera, a Raspberry Pi HQ camera 2, was placed in line with the light source and fresnel lens.
The camera and driver circuit for the LED were connected to a Raspberry Pi, which was used to trigger the LED to flash, and the camera to capture an image. 
The LED, driver circuit, and software were all developed by Joseph Schaefer, a previous student. 
The region behind the fresnel lens and in front of the camera was the area used to place the tank. %\
This is shown in Figure \ref{fig:bench}.
The tracer particles used were Vestosint 2159 polyamide particles, with a mean diameter of \SI{10}{\mu \m}. 
These were mixed with water, and poured into the tank.

\begin{figure}[htb!]
  \begin{center}
    \includegraphics[width=0.75\textwidth]{images/bench.jpeg}
  \end{center}
  \caption{The lab bench used, with all important objects other than the particles}
  \label{fig:bench}
\end{figure}

\subsection{Experimental images}
\label{sec:experimental_prep}
The main issue for experimental images is quantifying the amount of coupling, and trying to mitigate any coupling that is present. 
To see the coefficients for cross-coupling, each colour channel was pulsed individually. 

The fresnel lens used was visible in the images, as was the edge of the container and the surface of the water.
These would lead to tracking difficulties.
An example of these errors Figure \ref{fig:light}.
\begin{figure}[htb!]
  \begin{center}
    \includegraphics[width=0.65\textwidth]{images/fresnel_rgb.jpg}
  \end{center}
  \caption{An image of the fresnel lens with the LED illuminated}
  \label{fig:light}
\end{figure}

Additionally, the fact that both the camera and the LED were slightly out of alignment, due to the coarse hole spacing of the optical table, manifests in the lower portion being significantly brighter than the rest of the image. 

This issue is similar to in \citet{Aguirre2017}, where inconsistent lighting lead to an inhomogenous background. 
The solution used in that work was to have a representative image of the experimental setup with no particles, and then subtract this background image from subsequent experimental images. 
In this way, only changes from the background will be shown. 
All of these differences should be caused by shadows. 

To do capture an image in this way, multiple images were taken with the particles in motion.
The fresnel lens and the other hardware should remain static, while particles in motion would not.
This would mean that an averaged image from this series should only show the parts of the image that should be removed.
By subtracting this average image from each individual image in the series, only the particles and other objects in motion should be shown.
This has the additional benefit of inverting the image, which is needed.

Due to the low gain, noise was not expected to be a significant factor, meaning that a low number of these background images should be sufficient. 
These will be averaged to minimise the small amounts of noise present.

To access time resolved data, multiple images must be captured. 
The current imaging system presents problems with this, as a significant amount of time is needed to write the image buffer of roughly 25MB to disk. %\
If this amount of time is not given, the image is partially written, with the rest of the image left as black as a new image is captured. 
Figure \ref{fig:toofast} demonstrates this. 

\begin{figure}[htb!]
  \begin{center}
    \includegraphics[width=0.65\textwidth]{images/toofast.jpg}
  \end{center}
  \caption{A result of trying to capture images faster than the SD card could handle}
  \label{fig:toofast}
\end{figure}

This image shows single-channel blue illumination, where the next image was requested to be taken shortly after this one. 
There was not enough time to write the image, so this failed. 
This may be an inherent limitation of using a Raspberry Pi as the temporary storage, and therefore writing these images to an SD card. 
The SD Association has a list of speed classes, ranging from 2MB/s to 90MB/s which will give rough estimations of the time to write an image to the SD card, and it suggests that at best, an image could be captured at around 3Hz, potentially limiting the usefulness of the system for high reynolds number flows. 
Further work could potentially be carried out to cache the images in RAM and allow images to be captured more quickly with a thread reading in camera data and another thread writing the data to the disk. 
Additionally, there is the option of writing the images to a high-speed Network Attached Storage (NAS), which would likely be more convenient as the number of cameras in an array grows. 
While recording video instead of single frames may allow for a higher repetition rate, the compression used to reduce the data rate to manageable levels may mean that the particles, which are small, relatively slow moving objects, my be removed.


\newpage
\section{Results}
\subsection{Computational results}
\label{sec:comp_results}
These will each have a table with the relevant statistics, and my interpretation of how and why this happened.

The values of spill were all clamped to a value of 1 to ensure realistic values. 
The leading diagonal was kept at a value of 1.
The investigated values are 0, representing a system with no coupling, 0.5, representing a best-case scenario, 1 as a representative estimation, and 1.5, to give a worst case scenario.

\subsubsection{Finding the $C$ matrix}
\label{sec:removal}
%\
As was done in \citet{wang2017}, the intensity of each channel from the experimental bench was plotted against each other channel. 
These points then have a linear function fitted to them. 
The points and curves are displayed in Figure \ref{fig:fit}. 

\begin{figure}[htb!]
  \centering
    \begin{subfigure}[b]{0.5\textwidth}
      \includegraphics[width=\textwidth]{images/r_cor.png}
    \end{subfigure}
    \hfill
    \begin{subfigure}[b]{0.5\textwidth}
      \includegraphics[width=\textwidth]{images/g_cor.png}
    \end{subfigure}
    \hfill
    \begin{subfigure}[b]{0.5\textwidth}
      \includegraphics[width=\textwidth]{images/b_cor.png}
    \end{subfigure}
    \caption{plots of the intensities of the different pixels relative to each other, with lines of best fit.}
    \label{fig:fit}
  \end{figure}

Looking at the graph for red and green, the data is quite linear, unlike what is presented in \citet{wang2017}. 
In the blue channel, there is clearly saturation which results in nonlinearities. 
This is fairly easy to deal with, when doing experimental testing, the intensity of the illumination will avoid saturating the sensor where possible. 
This has resulted in the correlation between green and blue to be significantly higher than is representative of the linear region, as seen in the rightmost yellow line. 
For this reason, this value was manually adjusted to be roughly half of what it should be

Additionally, there are some anomalous results on this same channel, where it appears that perhaps two images have significantly higher mixing than the others.

The matrix with each coefficient is listed below. 

\[ C = 
  \begin{bmatrix}
  1 & 0.3201 & 0.2223 \\
  0.1568 & 1 & 0.3700\\
  0.1095 & 0.4000 & 1\\
  \end{bmatrix}
\]
The first index of each direction is red, the next green and the last is blue. 
The columns correspond to the sensed colour vector, and the rows correspond to the decoupled colour vector. 
\(C_{13}\) is how much blue will be caused by red, and \(C_{31}\) is how much red will be caused by blue.

The correlation between red and blue is lower than that between the channels and green. 
Obviously, the difference in wavelengths between green and blue is less than the difference in wavelengths between blue and red, meaning that this seems representative.
Additionally, the leading diagonal is all value of one. 
For this reason, this matrix is reasonable.


\subsubsection{Velocity profiles in outer units}
\label{sec:bl}
In all of the graphs shown in this section, all of the available points from a given time step were plotted.
This was done by iterating over every particle in a reconstruction from Gigatrack, and then scattering its wall normal (\(y\)) position against its streamwise velocity (\(u\)).
\(h\) is defined to be the half of the distance between the two walls. 
In this case, that equates to \SI{300}{\mm}.
\(U_B\) is the bulk velocity of the flow, the volumetric flow rate divided by the cross sectional area of the channel. 
For this flow, it is roughly \SI{168.4}{\mm\per\s}.

Figure \ref{fig:bl_comparison} shows how effective linear unmixing was. 
This used the reconstruction at \(a=1.5\).
From all of the particles at almost zero velocity in Figure \ref{fig:comparison_pre}, it can be concluded that spectral particles are indeed a major problem for cameras with large amounts of coupling.
While Figure \ref{fig:comparison_post} shows that the processing is imperfect, it brings the reconstruction from being completely unusable to being somewhat representative.


\begin{figure}[htb!]
  \centering
    \begin{subfigure}[b]{0.6\textwidth}
      \includegraphics[width=\textwidth]{images/1_dirty.png}
      \caption{}
      \label{fig:comparison_pre}
    \end{subfigure}
    \hfill
    \begin{subfigure}[b]{0.6\textwidth}
      \includegraphics[width=\textwidth]{images/1_clean.png}
      \caption{}
      \label{fig:comparison_post}
    \end{subfigure}

    \caption{Images before (top) and after (below) the correction has been applied.}
    \label{fig:bl_comparison}

\end{figure}

In Figure \ref{fig:bl_all}, it is clear that all of the intensities have a similar profile throughout the flow. 
At higher intensities, the few outliers seem to disappear. 
This is shown in Figure \ref{fig:peak_bl}, corresponding to an intensity of 1.5 times expected coupling, being more packed to the middle, with fewer extremities. 
This may limit the effectiveness of the current system in highly turbulent flows when this temporal data is needed. 
However, this suggests that there may still be uses in time-averaged flows where the bulk behaviour is the useful value. 

\begin{figure}[htb!]
  \centering

  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/boundary_layer_nospill_0.png}
    \caption{$a=0.0$}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/boundary_layer_nospill_0.5.png}
    \caption{$a=0.5$}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/boundary_layer_nospill_1.png}
    \caption{$a=1.0$}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/boundary_layer_nospill_1.5.png}
    \caption{$a=1.5$}
    \label{fig:peak_bl}
  \end{subfigure}

  \caption{Plots of streamwise velocity against wall-normal position for various intensities of coupling}
  \label{fig:bl_all}
\end{figure}

The fact that some data is reconstructed outside of the boundaries of the channel suggests that there are some issues with the reconstruction code. 
As this happens at all intensities of spill, with seemingly no correlation to amount of spill, and that some of these erroneous reconstructions shows that this is not an effect of cross coupling.

In all graphs, the region from \(\pm 1.0\) to \(\pm 0.5\) appears to have a greater interquartile range. 
Because this is consistent in all the plots, and asymmetric, this may be large structures in the streamwise profile as has been documented.
To resolve this, data from more times should be plotted.
As this is not the focus of this report, this was not done.

\subsubsection{Mean velocity}
\label{sec:mean_vel}
Table \ref{tab:mean_vels} shows the bulk velocity of the flow compared to the intensity of the coupling. 
This demonstrates that the more intense the coupling, the less accurate the result will be. 
The cause for this is clear from Section \ref{sec:bl}; there are more and more particles that appear static, bringing down the mean. %\
\begin{table}[htb!]
  \begin{center}
    \begin{tabular}[c]{|c|c|c|c|c|}
      \hline
      \textbf{intensity} & 0.0 & 0.5 & 1 & 1.5\\
      \hline
      \textbf{decoupled value} & 170.1 & 170.1 & 167.9 & 150.6\\
      \hline
      \textbf{coupled value} & 170.1 & 168.2 & 33.4 & 7.4 \\ 
      \hline
    \end{tabular}
  \end{center}
  \caption{Mean velocity for a variety of conditions. Units of \si{\mm\per\s}}
  \label{tab:mean_vels}
\end{table}

While at higher intensities, there are clearly issues, the decoupling method is clearly effective at bringing data from unusable to somewhat meaningful. 
Looking at Figure \ref{fig:bl_all}, the lower bulk velocity is clearly a result of increasing numbers of particles where the correction was not sufficient to remove the spectral particles. 
For this application, these particles could be removed by assuming that at certain wall-normal positions, velocity should be above a certain fraction of the bulk velocity, much like the approach that \citet{Aguirre2017} took to removing the coupling. 
However, this also suffers from the same limitations as that decoupling method. 
It is not resilient to cases where meaningful data is stored in that region. 
This means that if this system were to be used on a flow with fewer well-defined characteristics, an alternate system would have to be used. 
This result suggests that further work is needed to better understand how to reduce the coupling at high intensities. 

\subsubsection{Law of the Wall}
The values of $\kappa$ and $B$ used were those found in \citet{turb_book}, and they were $0.4$ and $5$ respectively. 
These constants may be incorrect for this flow, as discussed in \citet{Bailey} so they should be taken as an indication of the shape of the data, rather than a quantitative measure of error. 

Regardless, these graphs should be compared to Figure \ref{fig:wall_profile}

\begin{figure}
  \centering
  \begin{subfigure}[b]{0.4\textwidth}
    \includegraphics[width=\textwidth]{images/uplusclean_0.png}
    \caption{\(a=0\)}
    \label{fig:uplusclean_0}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.4\textwidth}
    \includegraphics[width=\textwidth]{images/uplusclean_0.5.png}
    \caption{\(a=0.5\)}
    \label{fig:uplusclean_0.5}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.4\textwidth}
    \includegraphics[width=\textwidth]{images/uplusclean_1.png}
    \caption{\(a=1\)}
    \label{fig:uplusclean_1}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.4\textwidth}
    \includegraphics[width=\textwidth]{images/uplusclean_1.5.png}
    \caption{\(a=1.5\)}
    \label{fig:uplusclean_1.5}
  \end{subfigure}
  \hfill

  \caption{$U^+$ plotted against $y^+$ for different intensities}
  \label{fig:uplus}
\end{figure}
Clearly, as the magnitude of the spill increases, so does the amount of noise in the data. 
Once again, the correction significantly improves the data collected, even if this correction is far from perfect. 
Denoising methods may be able to get the data to fit empirical results, however denoising methods are inherently lossy, and would mean data would be less meaningful. 

To generate this data, more points were naturally needed closer to the wall. 
Cosine spacing was used to sample closer to the wall at higher rates than at large distances.
Here, $N$ represents the total number of bins. 
A very low value of bins would mean that there is not enough data to find anything meaningful, and a high value would mean that individual errors would strongly affect the result and give a noisy curve. 
A value of $N=250$ was used for the graphs in Figure \ref{fig:uplus}


As discussed in \ref{sec:mean_vel}, higher intensities tend to lead to more small values. 

These bins each cover a range of wall normal positions. 
This means that for a large bin, there will be some real particles of a velocity near the bulk velocity, and some spectral particles near zero velocity.
Because the value plotted is the mean velocity of all particles in the bin, the end result of this will be that the point shown on the plot is between zero, and where the law of the wall predicts it should be. 
This explains why Figure \ref{fig:uplusclean_1.5} has many particles below the trend line, but so few above. 

As seen in Section \ref{sec:bl}, at \(a=0.5\), there is almost no effect on the result. 
This is likely because Gigatrack detects a particle when any of the RGB values are above a threshold value, and the interpretation of the result here and in Section \ref{sec:bl} is that on aggregate, this low value of coupling is not enough to cause significant issues.
The results in Section \ref{sec:mean_vel} suggests that the reconstructed datasets do differ slightly. 
This is likely due to particles that appear very bright due to their position, or with particles of a more typical brightness that coincide with a lot of noise.

\subsubsection{Velocity fluctuations}
Figure \ref{fig:upprmsp_lyp} shows the fluctuation values for \(a=0.0, 0.5, 1.0, 1.5\). 
It appears that the Figure \ref{fig:upprms_lyp_0} matches the profile shown in Section \ref{sec:vel_textbook}, albeit with some scaling issues. 
The ideal profile is reproduced in Figure \ref{fig:uprmsp_yp}.

As is clear from Figure \ref{fig:uplus}, increasing how coupled each channel of the camera is leads to more spectral particles, which will naturally increase the root-mean squared (rms) value of the velocities in each bin.
This leads to Figures \ref{fig:upprms_lyp_1}, \ref{fig:upprms_lyp_1.5} having clusters of data above \(u^{+'}_{rms} = 3\).
In Figure \ref{fig:upprms_lyp_1.5}, the data out of this upper cluster appears to have lower rms values than is predicted from the typical profile. 
The reason for this was mentioned in Section \ref{sec:bl}, and is visible in Figure \ref{fig:peak_bl}. 
The values appear to get compressed down onto the mean flow profile. 
The reason for this is unclear, however it is demonstrated in Section \ref{sec:overcorrection} to be a persistent issue.

\begin{figure}[htb!]
  \centering
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/turb_intclean_0.png}
    \caption{\(a=0.0\)}
    \label{fig:upprms_lyp_0}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/turb_intclean_0.5.png}
    \caption{\(a=0.5\)}
    \label{fig:upprms_lyp_0.5}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/turb_intclean_1.png}
    \caption{\(a=1.0\)}
    \label{fig:upprms_lyp_1}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/turb_intclean_1.5.png}
    \caption{\(a=1.5\)}
    \label{fig:upprms_lyp_1.5}
  \end{subfigure}
  \hfill
 
  \caption{Plot of $u'^+_{rms}$ against $y^+$ for multiple different levels of intensity}
  \label{fig:upprmsp_lyp}
\end{figure}

\subsubsection{Overcorrection}
\label{sec:overcorrection}
One potential method to better suppress  spectral particles is to use more intense correction. 
This will have an effect, because Gigatrack detects a particle based only on the intensity of a pixel, so further darkening pixels will result in fewer false positives. 

To test this, an intensity of 1.5 is used when generating the images.
As seen above, this value leads to significant errors, so any improvement should hopefully be clear.

To achieve a more aggressive correction, each image was generated with the \(C\) matrix multiplied elementwise by \(a\), as before. 
However, when the images were passed into \verb|janitor.m|, there was an additional constant, so the calculation was \(a*b * C\), where \(b\) is the amount of overcorrection.

In Figure \ref{fig:overcorrection_bl}, more aggressive correction clearly decreases the number of spectral particles, as is to be expected. 
It appears that some real particles come at the cost of these spectral particles, as the region outwards of \(h=\pm 0.5\) appears slightly more sparse in Figure \ref{fig:bl_1.5c} than in Figure \ref{fig:bl_1.1c}.

To gauge what effect this will have at the extreme values, Figure \ref{fig:bl_3c} is included, where \(b\) was set to \(3\). 
The structure of a turbulent channel flow is lost, there are far fewer particles, and this data is essentially meaningless.
This shows that at the extreme, using overcorrection can be harmful when trying to improve the quality of the data. 
For these flow profiles, more work is needed to determine how best to use overcorrection.

\begin{figure}[htb!]
  \centering
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/boundary_layer_nospill_1.1corr.png}
    \caption{\(b=1.1\)}
    \label{fig:bl_1.1c}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/boundary_layer_nospill_1.3corr.png}
    \caption{\(b=1.3\)}
    \label{fig:bl_1.3c}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/boundary_layer_nospill_1.5corr.png}
    \caption{\(b=1.5\)}
    \label{fig:bl_1.5c}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/boundary_layer_nospill_3corr.png}
    \caption{\(b=3\)}
    \label{fig:bl_3c}
  \end{subfigure}
  \hfill

  \caption{Plot of the mean flow profile at 1.5 spill with two levels of correction}
  \label{fig:overcorrection_bl}
\end{figure}



\begin{table}[htb!]
  \begin{center}
    \begin{tabular}[c]{|c|c|c|c|}
      \hline
      \textbf{Overcorrection factor} & 1.0 & 1.1 & 1.2\\
      \hline
      \textbf{Mean flow velocity (\si{\mm\per\s})} & 152.8 & 156.5 & 163.8\\
      \hline
    \end{tabular}
  \end{center}
  \caption{Mean velocity for a variety of conditions. Units of \si{\mm\per\s}}%\
  \label{tab:mean_vel}
\end{table}



\subsection{Experimental results}

\begin{figure}[htb!]
  \centering
  \begin{subfigure}[b]{0.7\textwidth}
      \includegraphics[width=\textwidth]{images/pre_debayer2.png}
      \caption{}
      \label{fig:pre_debayer}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.7\textwidth}
      \includegraphics[width=\textwidth]{images/post_debayer.png}
      \caption{}
      \label{fig:post_debayer}
    \end{subfigure}
    \hfill
    \caption{Images from the sensor, and with light postprocessing, to ease understanding}
    \label{fig:exp_img}
\end{figure}

The aims for testing on images were to have code ready to take data from the image sensor and pass it on to Gigatrack. 
The camera system used was a Raspberry Pi 4B attached to a Raspberry Pi HQ cam 2.
This allowed for the camera to be configured and triggered in software, allowing for maximum flexibility and consistency.

Code for capturing the images and saving the resulting images was already written by Joe Schaefer.
Figure \ref{fig:exp_img} shows the images before and after minor postprocessing. 
Figure \ref{fig:pre_debayer} shows the file as it was recorded from the image sensor, and Figure \ref{fig:post_debayer} has debayered and rotated the image.
The rotation was necessary because the image sensor was flipped upside down, meaning that in the black and white image, gravity is acting up the page.
The debayered image is there for clarity, it is clearly hard to determine colour from an image with the bayer filter still applied.

\begin{figure}[htb!]
  \centering
  \begin{subfigure}[b]{0.7\textwidth}
      \includegraphics[width=\textwidth]{images/process_bayerd2.png}
      \caption{}
      \label{fig:process_bayered}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.7\textwidth}
      \includegraphics[width=\textwidth]{images/process_debayerd.png}
      \caption{}
      \label{fig:process_debayered}
    \end{subfigure}
    \hfill
    \caption{Images from the sensor, and with light postprocessing, to ease understanding}
    \label{fig:exp_processed}
\end{figure}
In this image, the green and blue channels of the LED were activated for \SI{50}{\mu \s}, while the red channel was activated for \SI{150}{\mu\s}.
It may have been useful to have an even longer time on for the red LED, as it is clearly not as bright as the other colours. 

The liquid in the tank was primarily water, with a very small number of polyamide particles. 
Additionally, dish soap was added to reduce the surface tension and allow the polyamide particles to mix with the water.
These were then agitated by lightly shaking the box, in the hopes that the motion would lead to interesting images. 

Many of these images were taken sequentially. 
This was done so that each pixel could be averaged over the image sequence to remove the background.

Figure \ref{fig:exp_processed} shows these images, again bayered and debayered.
The ridges of the fresnel lens are still slightly visible in the background, suggesting that this background removal was not entirely successful.
In Figure \ref{fig:process_debayered}, some small dots are visible, primarily in the green channel.
Around the green is some red light, however this is very close to the particle.
A blue version of each particle is hard to see. 
This is likely because the time gap between each flash was around \SI{500}{\mu \s}, meaning that the particles will have hardly moved. 

\subsubsection{Finding the $C$ matrix for experimental images}
The LED system changed from working on the Gigatrack implementation of linear unmixing to working with real hardware, as the LEDs used before were unnecessarily powerful. 
This necessitated another attempt to find the $C$ matrix.

Images were captured using monochromatic red, green and blue illumination, and then processed in the same way as before. \

This resulted in the graphs visible in Figure \ref{fig:experimental_curves}, and a matrix as given. 
\[ C_{exp} = 
  \begin{bmatrix}
  1 & 0.5957 & 0.4641 \\
  0.3588 & 1 & 0.4479\\
  0.1834 & 0.2621 & 1\\
  \end{bmatrix}
\]

\begin{figure}[bth!]
  \centering
  \begin{subfigure}[b]{0.32\textwidth}
    \includegraphics[width=\textwidth]{images/exp_r_corr.png}
    \caption{}
    \label{fig:exp_r_corr}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.32\textwidth}
    \includegraphics[width=\textwidth]{images/exp_g_corr.png}
    \caption{}
    \label{fig:exp_g_corr}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.32\textwidth}
    \includegraphics[width=\textwidth]{images/exp_b_corr.png}
    \caption{}
    \label{exp_b_corr}
  \end{subfigure}
  \hfill

  \caption{}
  \label{fig:experimental_curves}
\end{figure}

\subsubsection{Testing on images}
\label{sec:exp_test}
With \(C_{exp}\) found and the background of the image removed, the next step of preparing the data for Gigatrack is unmixing. 
Using the same script as before, but with new code to list all of the images in a directory through MATLAB's \verb|dir|, this was made simple. 
Figure \ref{fig:exp_clean} shows how this looks on a real image.
The main effect appears to be that the image is less red. 
Upon close inspection, there is more significant changes. 

\begin{figure}[bht!]
  \centering
  \includegraphics[width=0.6\textwidth]{images/exp_spillless.png}
  \caption{The same image as before, with spill removed}
  \label{fig:exp_clean}
\end{figure}

Along the side of the container, which manifests as an mostly continuous vertical line, the green and blue light is more pure. 
Blue particles even appear to be present in the background of the image. 
A conclusion from the fact that both before and after unmixing, green is the predominant colour, is that the red and blue LEDs need to be turned on for longer. 

Figure \ref{fig:bubble} shows one of the details on the surface of the water, near the edge of the container. 
The unmixing has worked, as the colours are a lot more pure. 
The more significant issue is that the colours are not all in the same location. 
This could be the result of two issues, either chromatic aberration, or the fact that the actual LEDs on the board are not in exactly the same location.
Either way, this will need to be addressed. 
The fact that the LEDs are not the same everywhere is a simple fix, as there will be a constant shift across the whole image. 

Chromatic aberration is due to the fact that a material will have a different refractive index for each wavelength, meaning that the more a beam of light is refracted, the larger an issue this will be. 
The result of this is that the further off-axis a point on the camera lens is, the more it will be shifted. 
Additionally, each channel will have a different shift. 
In this case, it appears that where there was coupling, it was removed anyway with the linear unmixing. 
This will not necessarily be the case with all optical systems.

\begin{figure}[bht!]
  \centering
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/bubble_clean.png}
    \caption{After decoupling}
    \label{fig:bubble_clean}
  \end{subfigure}
  \hfill
  \begin{subfigure}[b]{0.45\textwidth}
    \includegraphics[width=\textwidth]{images/bubble_mix.png}
    \caption{Before decoupling}
    \label{fig:bubble_dirty}
  \end{subfigure}
  \caption{A part of one of the images before and after decoupling to show the effects of decoupling and the effects of chromatic aberration}
  \label{fig:bubble}
\end{figure}


\newpage
\section{Further work}

\subsection{Hardware iterations}
The issues with the camera and LED being misaligned as shown in Section \ref{sec:experimental_prep} will be a reasonably fast fix. 
Manufacturing adjustable, 3- or 6-degree of freedom mounts will significantly improve the ability to centre these optical elements. 
While there was limited capability to move the camera and LED, this was primarily limited to moving these vertically along aluminium extrusion. 
The hole spacing on the optical table was too coarse to allow for fine adjustments in the positions, which would result in better images. 

Additionally, the fresnel lens was hardly constrained, limiting the repeatability of the experiments.

\subsection{Chromatic aberration}
Because the position of a particle will be wrong to different degrees at different time steps from chromatic aberration, there will also be erroneous velocities and accelerations. 
While this may be limited with multiple cameras, this is still an issue that needs to be addressed. 
This will include finding a transfer function from position on image to positional error for each colour channel, and then distorting the resulting image as appropriate to remove this error. 

\subsection{First particle tracks}
The system developed is still far off creating its first point cloud. 
To do this, many steps more steps are needed. 
To start off with, more hardware is needed.

The iterations mentioned above will enable better optical performance, simplifying postprocessing.
The tank currently in use is too small to generate interesting flows.
For 3 dimensional data, more LEDs and cameras are needed, and these will naturally bring significant additional complexity. 
These all naturally need to be calibrated. 

Then, more software will be needed. 
The tools written will need to be better integrated into Gigatrack. 
The unmixing script works in a simple process, but is very inflexible in that it needs a list of filenames passed in. 
A better system would find the files itself. 
There was difficulty with \verb|imread| reading a 16 bit .tiff as an 8 bit .tiff, so the tools for finding the average of an image are currently written in Python.

As the optical system expands, the software to control the electronics will inevitably encounter bugs, which will need to be removed.
Gigatrack's synthetic image generation appears to limit the position of the cameras to one dimension, which would limit the ability to create multi-degree of freedom arrays. 
This would make testing of multiple cameras stacked in the wall-normal direction difficult

Once these are addressed, the system would need to be tested against a known flow.
This has its own difficulties, as seeding densities of the tracer particles will need to be high enough to capture flow structures, but low enough that tracking of the particles is still successful.
This will then need to be validated.

Once these steps have been completed, the system will be able to scale, which will naturally bring about its own issues.

\subsection{POD filtering}
Through discussing this project with others, the idea of using a Proper Orthogonal Decomposition (POD) filter to remove the background was suggested. 
Recent from \citet{PODfilter} suggest that this is a promising method to remove the background from PIV systems. 
With current experimental images, this may be difficult to implement as capturing an image sequence presents challenges (see Section \ref{sec:experimental_prep}).

\cleardoublepage

\bibliographystyle{unsrtnat}
\bibliography{references}

\newpage

\appendix
\section{Appendix A}
\label{app:RA}
\begin{skinnybox}
  \centering
  Risk Assessment
\end{skinnybox}
This shows the risk assessment for beginning lab work as approved.
\includepdf[pages=-]{images/Risk Assessment.pdf}

\section{Appendix B}
\label{app:MS}
\begin{skinnybox}
  \centering
  Method Statement
\end{skinnybox}
This shows the method statement for beginning lab work as approved.
\includepdf[pages=-]{images/Method Statement.pdf}

\section{Appendix C}
\label{app:COSHH}
\begin{skinnybox}
  \centering
  COSHH Assessment
\end{skinnybox}
This shows the COSHH assessment for beginning lab work as approved.
\includepdf[pages=-]{images/COSHH Assessment.pdf}

\end{document}
