Intensity-corrected 4D light-in-flight imaging

Light-in-flight (LIF) imaging is the measurement and reconstruction of light's path as it moves and interacts with objects. It is well known that relativistic effects can result in apparent velocities that differ significantly from the speed of light. However, less well known is that Rayleigh scattering and the effects of imaging optics can lead to observed intensities changing by several orders of magnitude along light's path. We develop a model that enables us to correct for all of these effects, thus we can accurately invert the observed data and reconstruct the true intensity-corrected optical path of a laser pulse as it travels in air. We demonstrate the validity of our model by observing the photon arrival time and intensity distribution obtained from single-photon avalanche detector (SPAD) array data for a laser pulse propagating towards and away from the camera. We can then reconstruct the true intensity-corrected path of the light in four dimensions (three spatial dimensions and time).


I. INTRODUCTION
As light travels and interacts with objects, photons are scattered in all directions. Light-in-flight (LIF) imaging is the process of capturing scattered photons using detectors with high temporal resolution such that light's path can be reconstructed. Three-dimensional LIF imaging was first captured using a holographic plate to record the spherical wavefronts of pulses reflected by mirrors; the technique involved no mechanical processes and achieved a temporal resolution of 800 ps [1]. This demonstrated real-time imaging of light undergoing dynamic processes, previous imaging  1: (a) Light in "real space" is scattered at (x1, y1, z1, t1) in all directions. A proportion of this scattered light travels to the camera, and is recorded in "camera space" as a signal at (x3, y3, t3), where x3 and y3 are pixel positions. (b) The remaining light travels across the field-of-view, and is scattered at position and time (x2, y2, z2, t2). This event is recorded at (x4, y4, t4). (c) Birds-eye-view of the two scattering events where R1 and R2 are the distance between the camera and the first and second scattering events respectively, and α1 and α2 are the scattering angles. The time difference for the two events in "camera space" is t4 − t3 = ∆t + (R2 − R1)/c, whereas the time difference in "real space" is t2 − t1 = ∆t. Rayleigh scattering effects observed by the camera are dependent on α1, α2, R1 and R2, which differ significantly for the two scattering events shown. Focusing effects from the imaging optics also contribute to the intensity signal. An image rendered from the perspective of the camera shows the right side of the beam, which is closer to the camera, is larger compared to the left side of the beam. This corresponds to a lower energy density on the right hand side, and therefore a brighter image on the left. was static and time averaged [2,3]. Further work proposed a mechanism for correcting distortion affects when imaging light [4]. Recent 3D LIF holography techniques use a scattering medium and achieve higher temporal resolutions [5,6].
The field of light-in-flight imaging was recently revolutionised by Velten et al. [7], who imaged femtosecond laser pulses propagating through a scattering media using a streak camera. This new LIF method allowed for scattering light dynamics to be observed at unprecedented temporal and spatial resolutions. However, this method requires a scanning mechanism to build a 2D image resulting in an acquisition time of one hour. Other methods for capturing 3D LIF involve transient imaging using a photonic mixer device (PMD), achieving a nanosecond temporal resolution with a one minute acquisition time [8]. In addition, other 3D LIF imaging methods include time encoded amplified imaging and computer tomography, which achieve nanosecond and picosecond temporal resolutions respectively [9,10].
The type of scattering that is observed is dependent on the medium that the light propagates through. For example, light has been captured propagating through fibre optics [11] and heated rubidium vapor [12]. When light travels through air, Rayleigh scattering is the dominant effect, and this was captured by Gariepy el al. who demonstrated three-dimensional LIF imaging using a single-photon avalanche detector (SPAD) array camera [13]. In this work, the light propagated in one plane that was perpendicular to the axis normal to the detector. Following this, it was recognised that relativistic effects, where the apparent velocity of light would deviate from c, could be observed with LIF [14,15] and these principles have allowed four-dimensional LIF reconstruction to be demonstrated for multiple paths of light over a narrow range of angles [16]. This was generalised, using a megapixel camera and machine learning techniques, to capture 4D LIF imaging of multiple pulses following arbitrary paths in space [17].
These technologies have already been used in fluorescence lifetime imaging [18], light detection and ranging (LIDAR) imaging through scattering media [19] and imaging around corners [20][21][22]. Ultimately, the ability to accurately capture the full scattering dynamics of light could lead to new approaches when imaging deep inside the human body, see Ref. [23] for an overview of light-in-flight research. Recent research tackling the problem of imaging in highly scattering media has shown computational imaging approaches can provide images in two [24] and three dimensions [25].
Our work builds on the recent LIF research by developing a model to compensate for distortions in the the recorded intensity, as well as relativistic effects previously observed, and reconstruct the 4D path of laser pulses. To do this, it is necessary to understand the underlying physics of light scattering in air and the relationship to imaging optics. This is illustrated in Fig. 1, where light scattered at a time t 1 from an object at a location (x 1 , y 1 , z 1 ) propagates a distance R 1 to a camera. The remaining light continues to propagate to position (x 2 , y 2 , z 2 ) where another scattering event occurs at time t 2 , and the scattered light travels a distance R 2 to the camera. The total time taken for the pulse to travel between the two scattering events is t 2 −t 1 = ∆t. Whereas, the two scattering events are recorded by the camera at times t 3 and t 4 respectively, and the difference in arrival time recorded by the camera is t 4 − t 3 = ∆t + (R 2 − R 1 )/c. This means the arrival time data recorded by the camera is different to the true propagation times of light and is ultimately dependent on the propagation angle. For the case of a camera, there is a mapping of an event occurring in three spatial dimensions and time to a camera with two spatial dimensions and time. The third spatial (z) dimension is collapsed and contained within the temporal data of the camera.
Rayleigh scattering effects are also observed by the camera and are dependent on the scattering angles (α 1 and α 2 ) and propagation distances (R 1 and R 2 ). These variables vary along light's path and so the intensity contribution is dependent on the position along the path. Furthermore, focusing effects contribute to the recorded intensity profile and are dependent on the perpendicular distances between the scattering event and camera. This is shown in the camera view image in Fig. 1 (c) which depicts an integrated image of the laser pulse's path across the field-of-view of the camera used to render the image. The depth of field is increased such that the whole path is in focus. The pulse, travelling towards the camera from left to right, is further away from the camera on the left-hand side and is therefore focused to a smaller size than the right-hand side of the pulse. This results in the intensity of the pulse increasing as the distance between the camera and pulse increases.
In this work, we are able to measure and subsequently correct all of the effects mentioned above. That is to say, we can correct both the temporal distortion, arising from relativistic effects, and the intensity distortions, resulting from Rayleigh scattering and the imaging optics. We demonstrate intensity-corrected light-in-flight imaging using a SPAD array, recording data for a laser pulse propagating at large and small angles with respect to the observation axis of the camera. The relativistic effects result in apparent speed-of-light velocities that span several orders of magnitude, and the intensity effects lead to observed intensities changing by at least a factor of two along the pulse's path.

II. THEORY
Consider a pulse of light that travels in three dimensions and is imaged using a camera with high temporal resolution. To develop the theoretical framework, we introduce the concept of the "camera space" to indicate where the data is recorded and the "real space" to indicate the three dimensional space in which the light pulse travels. It is the goal The intensity of a segment of the pulse is dependent on several factors: Rayleigh scattering, focusing effects and the integrated pulse length. These intensity contributions are derived using the variables shown above where θ is the propagation angle of the pulse relative to the optical axis, d is the distance between the centre of the pulse and imaging lens, f is the focal length of the imaging lens, r is the distance between the imaging lens and the nearest edge of the segment, θ1 is given by Equ. (2), θ2 is given by Equ. (3), x is the distance along the SPAD array to a given pixel and ∆ is the active pixel width. Combining these effects, the intensity of scattered photons along the beam in camera space (I(x; θ, A, f, ∆)) is derived and given in Equ. (7). The relativistic effects are explained using the same variables by Equ. (9).
for the work to convert the camera space data to the real space as accurately as possible. The inversion of the camera space data to the real space path enables the true intensity-corrected light path to be reconstructed.
Light-in-flight data is subject to intensity and relativistic effects observed in the camera space. The intensity of scattered photons along the beam in camera space is derived by considering the intensity contribution from one segment of the beam on one pixel. The intensity of a segment of the beam is calculated using the schematic in Fig. 2 where laser pulses traveling across the field-of-view, at propagation angle θ with respect to the observation axis, are imaged using a SPAD array. A proportion of the photons within each pulse are scattered by air molecules and travel through the imaging lens aperture. Different segments of the beam are imaged by different pixels within the SPAD array and the intensity contribution from each segment is dependent on focusing effects from the imaging optics (I f ), Rayleigh scattering (I r ), and integrated path length (I s ).
The intensity of the beam in camera space (I(θ 1 , θ 2 ; θ)) is given by where B is a normalisation constant dependent on integration time and laser power, x is the distance along the SPAD array to a given pixel, A is the sensor width, f is the focal length of the lens, ∆ is the active pixel width, θ is the propagation angle relative to the observation axis and r is the distance between the imaging lens and the nearest edge of the segment, θ 1 satisfies and θ 2 satisfies The first contribution to the intensity of a segment of the beam is from focusing effects in the imaging optics of the system and is given by This contribution is a result of parts of the beam which are further away from the lens focusing to a smaller point on the SPAD array with higher energy density. The second contribution is from photons undergoing Rayleigh scattering with air molecules and is given by where I 0 is the intensity constant, n is the refractive index, d r is the scattering particle diameter and λ is the wavelength of scattered light. Rayleigh scattering is dependent on the scattering angle and distance between the SPAD array and the pulse, which both change along the beam. The final contribution to the intensity of one pixel is from the integrated path length, which is the segment length imaged by each pixel, given by This results in pixels at the edge of the SPAD array seeing a larger length of pulse than pixels in the middle of the SPAD array. By combining these effects and substituting Equ.
(2)-(6) into Equ. (1) the intensity of scattered photons along the beam recorded by the SPAD array (I(x; θ, A, f, ∆))is found to be where C is a normalisation constant which includes the Rayleigh scattering constants, the integration time of the SPAD array and the optical power of the laser. Equation (7) assumes sin(θ 2 − θ 1 ) ≈ θ 2 − θ 1 . Finally, the Rayleigh effect is shown by measuring the central pixel intensity (I c (θ; f, ∆)) for different values of θ. This intensity is independent of focusing effects as d is constant for all θ and is given by which is derived by substituting x = A/2 into Equ. (7). This equation is a modified version of the Rayleigh scattering effect and introduces a normalisation factor that takes into account the length of pulse imaged by the central pixel.
Relativistic effects seen in the camera space result in the pulse appearing to travel at apparent velocities different to the speed of light. The arrival time in camera space is dependent on θ and d as shown in Fig. 2 (a). The arrival time difference between the central pixel and an arbitrary pixel (∆t(x; θ, A, f )) is given by where c is the speed of light in air. From the above equations, the relativistic and intensity effects observed in the camera space can be modelled and compared to experimental data.

III. EXPERIMENTAL SET-UP
The relativistic and intensity effects of LIF imaging are investigated using the experimental set-up shown in Fig. 3. The system includes a SPAD array camera, a 532 nm short pulsed laser (Teem Photonics STG-03E-1x0), and an optical constant fraction discriminator used as a trigger. The impact of intensity and relativistic effects are more pronounced when the light travels at large or small angles with respect to the optical axis of the camera, which corresponds to light travelling towards and away from the camera, see Fig. 3

(b) and (c) respectively.
Laser pulses, which have a pulse width of ≈ 500 ps, are expanded to a beam waist of ≈ 5 mm and collimated via two lenses of focal length 100 mm and 400 mm respectively, resulting in a Rayleigh range of 150 m. This ensures there are no intensity effects due to the beam diverging as it travels across the field-of-view of the sensor. The laser pulses are directed to a constant fraction discriminator acting as a trigger with 200 ps jitter, which sends 4 kHz transistor-transistor logic (TTL) pulses to the SPAD array. The TTL pulse starts the timer for each of the 32 x 32 pixels operated in time-Correlated Single Photon Counting (TCSPC) mode. Histograms of photon counts are  recorded for every pixel over 1024 time bins, each with a width of 55 ps. The pixel area and active area are 50 µm x 50 µm and 6.95 µm x 6.95 µm respectively, giving a fill factor of 1.9 %. An 8 mm focal length c-mount lens is used to image the beam onto the sensor. The aperture of the lens can be stopped down to extend the depth of field, and this is essential to reduce blurring and ensure the entire path of the beam is in focus on the camera. The mirrors used to direct the beam towards the SPAD array are placed outside the field-of-view so only photons scattered by air molecules are collected by the imaging lens, thus avoiding saturation effects and allowing Rayleigh scattering to be observed. Finally, to observe the Rayleigh scattering effects, the laser and trigger were placed on a rotation stage system which allows θ to be easily varied.
When measuring I(x; θ, A, f, ∆), it is important for the whole of the beam to be in focus. This is demonstrated in Fig. 4 (a) and (b) which shows electron multiplying charge-coupled device (EMCCD) intensity images of the beam travelling from right to left away from the SPAD array with an open and closed aperture respectively. When the lens aperture is open, out of focus light contributes to the intensity image resulting in part of the beam being out of focus and less intense than predicted by Equ. (7). When the lens aperture is closed, only in focus light is incident on the detector and the intensity effects predicted are observed. This condition requires longer acquisition times to collect sufficient photon counts to build an intensity image.

IV. RESULTS
In order to achieve intensity-corrected 4D LIF imaging, it is important to remove the noise present in the SPAD array data. This has been achieved by fitting Gaussian functions to each pixel and setting the pixel intensity to zero if the standard deviation of the Gaussian is outside an acceptable range. Furthermore, for noisy pixels within the beam path, interpolation is performed.
Camera space results for laser pulses travelling from left to right towards the SPAD array at θ = 167.0 • ± 0.5 • are shown in Fig. 5. Three frames, each with time duration 55 ps, of laser pulses travelling across row 17 of the SPAD array are shown at 2.1 ns, 2.5 ns and 2.8 ns in Fig. 5 (a) to (c). The time taken for the pulse to travel across the SPAD array is less than the time bin width, resulting in the pulse being present for all pixels of row 17 in a single frame.
The data and fitted model for row 17 as a function of position and time is given in Fig. 5 (d) and (e). The photon intensity decreases as pixel number increases in both the data and the fitted model. This is because the left hand side of the beam is further away from the SPAD array and so focuses to a smaller area on the detector with a higher energy density. The total time taken for the pulse to travel across the SPAD array was measured to be 18 ps, and using equations (9) and (7) to fit a 2D Gaussian function to the data, θ was calculated to be 165.5 • ± 0.1 • . The error in θ was calculated by numerically creating 10 statistically identical data sets and fitting to these. These additional data sets were created by sampling from a Poisson distribution with a mean and variance determined by the initial experimental data. Finally, the apparent velocity of the pulse as a function of propagation time is shown in Fig. 5 (f) and varies from 6.8 c to 7.0 c. This superluminal apparent velocity is entirely due to the pulse travelling toward the camera.
Using the data obtained by the SPAD array, the camera space data is converted to real space data using θ and normalising the data by the fitted intensity values at the appropriate times. Three frames from the real space movie of the pulse traveling towards the SPAD array are shown in Fig. 6 (a) to (c) at 0.0 ns, 0.4 ns and 0.7 ns. The beam diameter used for the real space reconstruction was 5 mm and the pulse width was 15 cm. These values were taken from measurements and known values of the pulse.
Camera space data was also recorded for laser pulses travelling from left to right away from the SPAD array at θ = 13.0 • ± 0.5 • using the set-up shown in Fig. 3 (b). Three frames of laser pulses travelling across row 17 are shown at 1.7 ns, 2.5 ns and 3.3 ns in Fig. 7 (a)-(c). The pulse length present in each frame is shorter for light travelling away from the SPAD array, indicating lower apparent velocities. The data and model used to calculate θ are shown in Fig. 7 (d) and (e) respectively. The total time taken for light to travel in camera space is 1.4 ns, and the curvature of the fitted function indicates the pulses appear to decelerate as they travel away from the SPAD array. Using Equ. (9) and (7) to fit to the data, θ was estimated to be 13.9 Finally, the apparent velocity of the pulse as a function of propagation time is shown in Fig. 7 (f) and varies from 0.19 c to 0.05 c. This results in a ratio of the fastest to slowest apparent velocities equal to 156. This is the largest ratio of super to subluminal apparent velocities in 4D LIF imaging; the previous highest ratio was 17, reported in Ref. [17].
Using the camera space data in Fig. 7, we can then recreate the real space data of the pulse travelling away from the SPAD array. Three frames from a movie showing the real space are shown in Fig. 8 at times of 0.0 ns, 0.4 ns and 0.7 ns.
Our final experiment demonstrates the angle dependence of scattering in air for light-in-flight imaging, i.e. Rayleigh scattering. This is achieved by placing the pulsed laser on a rotation stage, see Fig. 9, allowing θ to be easily altered, and recording the intensity of the central pixel of the camera. The central pixel intensity is only dependent on θ as the distance between the centre of the rotation stage and SPAD array is constant for all values of θ. This removes the effects of focusing and the inverse square dependence, which are both present in the first experiment. Fig. 9 (b) shows the observed experiment data in good agreement with the predictions of Rayleigh scattering, see Equ. (8). It should be noted that the effects of Rayleigh scattering were present in the previous experimental results but were harder to Note that the timescale is the camera space time. This is longer than the real space time, leading to the apparent subluminal velocities. isolate.

V. CONCLUSION
Relativistic effects, focusing, and Rayleigh scattering all play a significant role in the observed signal for light-inflight imaging. By modelling these effects we have been able to invert SPAD array data and reconstruct the true 4D path of laser pulses, showing a strong agreement between experiment and theory. We demonstrate the validity of our model by fitting to data obtained for light travelling towards and away from a SPAD array and comparing the temporal and intensity distributions to the model. The ratio of the apparent velocity of the pulses travelling towards and away is over two orders of magnitude and is the highest ratio observed for light-in-flight imaging.

Disclosures
The authors declare no conflicts of interest.