# Mathematical Estimation of Particulate Air Pollution Levels by Multi-angle Imaging

Or Vernik, Barak Fishbain

► To cite this version:

Or Vernik, Barak Fishbain. Mathematical Estimation of Particulate Air Pollution Levels by Multi-angle Imaging. 13th International Symposium on Environmental Software Systems (ISESS), Feb 2020, Wageningen, Netherlands. pp.249-257, 10.1007/978-3-030-39815-6_24. hal-03361876

HAL Id: hal-03361876

https://inria.hal.science/hal-03361876

Submitted on 1 Oct 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Mathematical Estimation of Particulate Air
Pollution Levels by Multi-angle Imaging

Or Vernik¹ and Barak Fishbain²

¹ Dept. of Applied Mathematics, Faculty of Mathematics, Technion - Israel Institute of Technology,

² Department of Environmental, Water and Agriculture Engineering, Faculty of Civil and Environmental Engineering, Technion - Israel Institute of Technology, Haifa, Israel

**Abstract.** Air pollution control and mitigation are important factors in wellbeing and sustainability. To this end, air pollution monitoring has a significant role. Today, air pollution monitoring is mainly done by standardized stations. The spread of those stations is sparse and their cost hinders the option of adding more. Thus, arises the need for cheaper and available means to assess air pollution. In this article, a method for assessing air pollution levels by means of multi angle imaging is presented. Specifically, the focus is on estimating images' blur as an indication for PM (Particulate Matter) ambient levels. The suggested method applies back-projection Radon transform. By back projection methodology, particles' concentration at each voxel in a 3D space is reconstructed from photos taken from a few different angles.

**Keywords:** Air pollution monitoring · Multi-angle imaging · Filtered Back Projection.

# 1 Introduction

Exposure to PM is known to be one of the predominant factors in morbidity and mortality, causing the premature death of millions of people a year. This is especially true in developing countries with growing industry [1]. Thus, measuring and monitoring PM pollution is a very important task the world needs to deal with in order to ensure public health.

Currently, PM monitoring is standardly done by Air Quality Monitoring (AQM) stations that are considered accurate [2]. Due to their size and cost, the stations are sparsely dispersed and so we get a low-resolution concentration map. This limitation is typically addressed by interpolation schemes. However, the interpolation is a complicated task as PM concentrations are characterized by high spatial variability and may present different behaviors for different fractions and geographical areas [3]. These drawbacks have lead researchers to seek other approaches such as Micro Sensing Units (MSUs) for air pollution monitoring [4–7], either high-resolution methods or less costly and portable monitors. Those methods can be used for taking measurements, validating new technologies or giving hints in interpolating existing in-situ measurements.
Here we focus on *visual* means for assessing PM levels in the atmosphere. To this end, an image processing technique for evaluating the extinction coefficient using visibility cameras was suggested by Graves and Newsam [8]. The extinction coefficient is a measure that quantifies local radiance attenuation and is used as a standard for measuring atmospheric visibility. Using this measure, an image-processing method to extract a quantitative measurement of atmosphere transmission (the ability of radiation to pass through the atmosphere) from standard photos was presented. The transmittance has a closed form relationship with the extinction coefficient (inversely related), that then can be computed. The research here aims at improving this notion by looking at a 3D volume with changing concentrations at each voxel, as opposed to the method above, which assumes a homogeneous area.

Methods for retrieving PM by multi-angle imaging are in use by satellites. Instruments such as Multi-angle Imaging SpectroRadiometer (MISR) and Airborne Multi-angle SpectroPolarimetric Imager (AirMSPI) are currently applied for remote aerosols sensing [9, 10]. However, these methods have low temporal resolution, as it takes a satellite few hours (or few days) to complete their orbit (latitude dependent), but more importantly, these devices measure the full vertical atmospheric column, which makes it less relevant to the amount of pollution at ground level that can indicate health hazards [11]. On top of these, this method is expensive and cannot be easily deployed.

The underlying assumption here is that due to the particles' optical properties of scattering, different concentrations of particles will cause the light field to scatter differently [11]. The higher the concentration of scatterers is, the more refractions the light will go through. Refracted light impairs human visibility and expression for this phenomenon is evident in standard photos [8]. A photo acquires a snapshot of the light field coming from all angles, refracting by the objects and particles in the scenery and eventually hits the camera lens. Particles suspended in the air cause a level of blurriness in the acquired photo. Having this relationship between blurriness and amount of scatterers in the air leads us to the approach of estimating the particles' concentration by measuring image blurriness.

The blur of an image is positively correlated with the integral of the attenuation (extinction) coefficient over the Line Of Sight (LOS) to the object. Thus, by measuring the blurriness of one image we cannot infer the concentration at each voxel. The intensity of the light's direct transmission is attenuated by the scattering (absorption by particulate matter is negligible) and is computed by:

$$I = I_0 e^{-\int_{\text{LOS}} \beta^{(l)} dl}. \qquad (1)$$

Where, $I_0$ is the initial intensity of the radiation source. Our objective then is to find $\beta$, which is strongly related with the PM concentration denoted n, by:

$$\beta(\bar{x}) = \sigma \cdot n(\bar{x}). \qquad (2)$$

Where, $n(\bar{x})$ is the PM concentration at position $\bar{x}$ in space and $\sigma$ is the extinction cross section. $\sigma$ is depended upon the light's wavelength and the particle's
shape and size, thus for our purposes it will be assumed constant. Assuming we know the values of the extinction coefficient integral from a large number of angles (having a correlation with the blurriness measure), the inverse problem is mathematically solvable by the Radon transform [12].

Radon transform is widely used in Computational Tomography (CT) applications, allowing us to reconstruct the object being scanned from its measured projections. The classically used form of the Radon transform and that we are going to describe, is the 2D Radon. A projection is the integral on the extinction coefficient which equals to:

$$p(y) = \int \underbrace{\beta(x, y)}_{\text{our objective}} dy = -\underbrace{\ln\left(\frac{I(y)}{I_0}\right)}_{\text{assumed known}}. \qquad (3)$$

The most effective way for reconstructing $\beta$ from the projections is by using Filtered Back Projections (FBP) that is equivalent to the Inverse Radon Transform (IRT) but less computationally costly.

The Radon Transform:

$$p_{\theta}(r) = \int \beta(r \cdot \cos \theta - s \cdot \sin \theta, r \cdot \sin \theta + s \cdot \cos \theta) ds. \qquad (4)$$

FBP:

$$I(x, y) = \int_{0}^{\pi} p(r, \theta) * q(r) d\theta. \qquad (5)$$

$q(r)$ is the LPF, for example Ram-Lak filter [13].

Relying on these principals, this research will define a method for three-dimensional PM concentrations reconstruction by multi-angle imaging. The proposed method will be visibility cameras based, aimed to be a simple, cheap and portable technique for air pollution assessment.

## 2 Monte-Carlo Simulation

To illustrate the potential of the suggested method, a simulation of the image acquisition process of the light field in a 3D volume, taken from few different angles is presented. The simulation allows for the reconstruction of PM suspended at each voxel in the volume.

We relate to the physical behavior of light propagating through a volume containing PM in different concentrations as a random process [14]. The light, taken in its particulate sense of photons, passes a medium, which has spatially variable optical depth (a measure of the light ability to propagate through the medium) denoted by $\tau$. At each stage, $\tau$ is sampled from the optical depth Cumulative Density Function (CDF):

$$F(\tau) = \int_{0}^{\tau} e^{-\tau'} d\tau' = 1 - e^{-\tau}. \qquad (6)$$
Using Monte-Carlo method, we get a random optical depth sampled from its CDF. For this, we use the uniform distribution U[0,1] CDF for sampling a random number u, and get $\tau$ by:

$$ \tau = F^{-1}(u), \quad u = rand(). \qquad (7) $$

By using the random $\tau$ we sampled, we can determine $l$, the distance the ray propagates until the next diffraction, by the relationship:

$$ \tau = \int_{0}^{l} \beta(x, y, z)dl = \int_{0}^{l} \sigma \cdot n(x, y, z)dl. \qquad (8) $$

Where $\beta$ is the extinction coefficient, $\sigma$ is the extinction cross section and $n$ is the PM concentration at each voxel with its coordinates denoted by $\{X, Y, Z\}$. $l$ is found by numeric integration.

Once $l$ is found, the scattering angle after the collision has to be determined. The angle is computed in a similar fashion, assuming randomness in the process and relying on Mie scattering theory under the assumption of spherical particles [15]. The Mie Theorem provides a physical solution for the scattering of an electromagnetic wave by spherical uniform particles about the size of the light wavelength. Using Monte-Carlo, we will sample a random number from the CDF of the Mie scattering angle on the intersection plane [14] and a fully random [0,2$\pi$) zenith angle. The found angles, will give us the new direction of the ray's propagation.

We repeat this process of finding new distance to the next collision, new direction and the intensity at each stage, until the ray exits our volume at some point. This will be done for each ray in the packet of rays entering from each radiated voxel. For now, the image acquisition process is done by summing the values of the exiting rays' intensities for each pixel at the volume's boundaries (visualization of the simulation shown at figure 1).

The received photos (each received from different angle) are used to estimate the level of pollution by applying the Blur Metric (BM) suggested by Frederique Crete [16]. The BM estimation of the image blurriness is based upon the notion that once an image is blurred, blurring it again will result in smaller differences than blurring a sharp image. We thus receive a value between zero and one indicating the effective blurriness of the image from which we can deduce the general level of pollution (amount of scatters).

## 3 Preliminary Results

As for now, we have created a volumetric grid containing random, uniform distributed PM concentration at each voxel (Figure 2a). In order to create a more realistic distribution we smoothed the outcome with a 3D box filter (Figure 2b).

We then simulated with Monte-Carlo method the rays enter the volume at specific points and angle and diffracting by the particles. The ray's entrance angle determines which faces of the cube are going to be affected by the radiation. The
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Sun rays hit the medium and scattered until exiting the grid boundaries. Image is received on the bottom plane as a result of rays intersecting the plane.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Initial uniform random values in the range [0,10] (a) and the smoother version with smaller spatial gradients in concentration (b).</figcaption></figure>
code is designed to plot random number of rays track as demonstrated in figure 3. For each voxel in the face that the radiation intersects with, we simulate a beam that hits the voxel and propagating through the volume until exiting from the grid boundaries.

<figure><img src="image_4.png" /><figcaption>Fig. 3. on the left, simulation of random rays entering the experiment volume with same angle, and scattering until exiting the volume. On the right, same simulation zoomed in; we can see the red arrow indicating entrance point and angle and the blue arrows indicating new directions after scattering.</figcaption></figure>

All rays start with some initial intensity, which will be reduced by the collisions with the aerosols in the volume by the following connection:

$$I_{after} = I_{before} \cdot \bar{\omega}. \qquad (9)$$

Where, $\bar{\omega}$ is the single scattering albedo (the attenuation in intensity) (SSA) of PM and $I$ is the intensity.

<figure><img src="image_5.png" /><figcaption>Fig. 4. received images at the X-Y plane (the floor of the volume). Sun light is coming from the top left corner direction to all three affected faces. The image on the left was acquired with 100 rays per pixel and presented with low resolution. The image on the right has higher spatial resolution (x5 than the left image) and 500 rays per pixel.</figcaption></figure>

For better understanding of the results we repeated our experiment, this time instead of all rays having the same initial intensity, we used a black and white
image of 'Lena' [17] that will filter our rays when entering the upper face of the volume. The rays are now coming in a straight angle from above. The results are a highly noisy image of the original Lena at the bottom plane (figure 5). Repeating the same initial conditions only with growing levels of PM concentrations and measuring the blurriness of the image using our blur metric, gives us a strong correlation between the effective image blur and the level of pollution as we can see in figure 1. Those results will be better examined, evaluated and validated as we continue our work.

<figure><img src="image_6.png" /><figcaption>Fig. 5. from left to right. 1. Lena's original B&W image. 2. Bottom plane received image. 3. Using LPF on the received image.</figcaption></figure>

<table><caption>**Table 1.** the measure of blurriness increases with the rise in PM pollution level as we can see in different ranges we examined and different scattering simulations.</caption><thead><tr><th>scattering angle<br>Range of<br>random<br>concentrations</th><th>random<br>scattering</th><th>Mie<br>scattering</th></tr></thead><tbody><tr><td>[0,5]</td><td>0.4961</td><td></td></tr><tr><td>[4,9]</td><td>0.5185</td><td></td></tr><tr><td>[8,13]</td><td>0.6396</td><td></td></tr><tr><td>[12,17]</td><td>0.7141</td><td></td></tr><tr><td>[0,0.5]</td><td></td><td>0.4763</td></tr><tr><td>[0,0.25]</td><td></td><td>0.4689</td></tr><tr><td>[0,0.1]</td><td></td><td>0.4533</td></tr></tbody></table>

## 4 Conclusions

The preliminary results indicate a certain connection between image blur and pollution level. As for now, this study's concept is theory based. We have reason
to believe that the results will be highly correlated with real-life measurements. The system discussed here is integral based. Hence, it tends to be stable and presents small perturbation at the system's output as a result of small perturbation at its inputs. Therefore, we expect this theoretical exercise to show similar results in real-life applications. Albeit, this is still a work in progress, the study needs to continue and establish the found correlation. The model is not yet fully built and there are no concrete reconstruction results to test at this stage.

At next stage, a full back projection (FBP) based reconstruction scheme will be used for finding the original concentrations in space, having only the images as input. Currently, the most simplified model is used. We aim to gradually add more layers of accuracy, better modelling of the physical phenomena, using enhanced air pollution dispersion modelling tools (e.g. GRAL) [18], adding effects of light polarisation and effects of ground and background reflections. Once the model is complete, sensitivity analysis will be done as well as validating concentrations against ground truth values and cross validation. Future work will also aim at testing the scenario simulation in real life. Our proposed method has the potential to be cheaper, more feasible and with higher spatial and temporal resolution than existing standard methods (remote sensing, AQM stations, etc.). Those advantages add up to the simplicity of the method which may be deployed by anyone with the possession of commodity cameras.

## Acknowledgments

The authors would like to thank Holodovsky Vadim and Schechner Y. Yoav, for sharing their code for Mie scattering, that has helped to push the paper forward in terms of accuracy.

## References

1. J. S. Apte, J. D. Marshall, A. J. Cohen, and M. Brauer, "Addressing global mortality from ambient pm2.5," *Environmental Science & Technology*, vol. 49, pp. 8057-8066, jul 2015.

2. S. Moltchanov, I. Levy, Y. Etzion, U Lerner, D.M. Broday and B. Fishbain, "On the feasibility of measuring urban air pollution by wireless distributed sensor networks," *Sci Total Environ*, vol. 502, pp. 537-547, jan 2015.

3. J. P. Pinto, A. S. Lefohn, and D. S. Shadwick, "Spatial variability of pm2.5 in urban areas in the united states," *Journal of the Air & Waste Management Association*, vol. 54, pp. 440-449, feb 2012.

4. N. Castell, F. R. F. Dauge, P. Schneider, M. Vogt, U. Lerner, B. Fishbain, D. Broday, and A. Bartonova, "Can commercial low-cost sensor platforms contribute to air quality monitoring and exposure estimates?," *Environment International*, vol. 99, pp. 293-302, feb 2017.

5. B. Fishbain, U. Lerner, N. Castell, T. Cole-Hunter, O. Popoola, D. M. Broday, T. M. Íñiguez, M. Nieuwenhuijsen, M. Jovasevic-Stojanovic, D. Topalovic, R. L. Jones, K. S. Galea, Y. Etzion, F. Kizel, Y. N. Golumbic, AyeletBaram-Tsabari, T. Yacobi, D. Drahler, J. A.Robinson, D. Kocman, M. Horvat, V. Svecova, A. Arpaci, A. Bartonova, N. Castell-balaguer, O. Popoola, D. M. Broday,
T. Martinez Iñiguez, M. Nieuwenhuijsen, M. Jovasevic-Stojanovic, D. Topalovic, R. L. Jones, K. S. Galea, Y. Etzion, F. Kizel, Y. N. Golumbic, A. Baram-Tsabari, J. Robinson, D. Kocman, M. Horvat, V. Svecova, A. Arpaci, and B. Alena, "An evaluation tool kit of air quality micro-sensing units," *Science of the Total Environment*, vol. 575, pp. 639-648, jan 2017.

6. U. Lerner, O. Hirshfeld, and B. Fishbain, "Optimal Deployment of a Heterogeneous Environmental Sensor Network," *Jorunal of Environmental Informatics*, 2018.

7. S. Moltchanov, I. Levy, Y. Etzion, U. Lerner, D. M. Broday, and B. Fishbain, "On the feasibility of measuring urban air pollution by wireless distributed sensor networks," *Science of the Total Environment*, vol. 502, pp. 537-547, 2015.

8. N. Graves and S. Newsam, "Using visibility cameras to estimate atmospheric light extinction," in *2011 IEEE Workshop on Applications of Computer Vision (WACV)*, pp. 577-584, Jan 2011.

9. D. Diner, J. Beckert, T. Reilly, C. Bruegge, J. Conel, and R. Kahn, "Multi-angle imaging spectroradiometer (misr) instrument description and experiment overview," *IEEE Transactions on Geoscience and Remote Sensing*, vol. 36, pp. 1072-1087, jul 1998.

10. D. Diner, F. Xu, M. J. Garay, J. V. Martonchik, B. E. Rheingans, and S. Geier, "The airborne multiangle spectropolarimetric imager (airmspi): a new tool for aerosol and cloud remote sensing," *Atmospheric Measurement Techniques*, vol. 36, aug 2013.

11. Y. Etzion, D. Broday, and B. Fishbain, "Analysis of image color and effective bandwidth as a tool for assessing air pollution at urban spatiotemporal scale," *Proceedings of SPIE - The International Society for Optical Engineering*, vol. 8657, 02 2013.

12. L. Yaroslavsky, *Digital holography and digital image processing: principles, methods, algorithms*. 2013.

13. J. Bernal and J. Sanchez, "Use of filtered back-projection methods to improve ct image reconstruction," 09 2009.

14. A. Levis, Y. Schechner, A. Aides, and A. Davis, "Airborne three-dimensional cloud tomography," pp. 3379-3387, 12 2015.

15. J. E. Hansen and L. D. Travis, "Light scattering in planetary atmospheres," *Space Science Reviews*, vol. 16, pp. 527-610, Oct 1974.

16. F. Crété-Roffet, T. Dolmière, P. Ladret, and M. Nicolas, "The blur effect: Perception and estimation with a new no-reference perceptual blur metric," *Human Vision and Electronic Imaging*, vol. 12, 03 2007.

17. D. C. Munson, "A note on lena," *IEEE Transactions on Image Processing*, vol. 5, pp. 3-3, Jan 1996.

18. D. Ottl, P.-J. Sturm, and R. Almbauer, "Application of the lagrangian dispersion model gral to assess the air quality in the city of graz," in *Transport and air pollution, 12 th Symposium*, pp. 75-79, ,. 2003.