1. Introduction

Standard imaging techniques rely on photosensitive elements with spatial resolution, such as photographic film or electronic matrix sensors. Single-pixel imaging (SPI) is an alternative imaging method in which the photosensitive element has no spatial resolution, that is, it is represented by only one pixel. The image in this method is not registered by a photosensitive element, but is calculated, or reconstructed [1].

Using just one pixel has a number of advantages over the matrix-based imaging. Firstly, the spectral sensitivity of a single-pixel detector can be selected from a very wide range, in contrast to a matrix sensor. If we consider only the visible and near-infrared range of light, then the existing technology of mass production of silicon matrices makes them outstanding in terms of price, resolution, and performance. However, if we are interested in the spectral range beyond the spectral sensitivity of silicon, then the characteristics of existing matrices deteriorate significantly compared to silicon ones. For example, relatively inexpensive silicon matrices with a resolution of tens of megapixels are currently available, while in the short-wave infrared range (around a wavelength of 1.5 microns) or long-wave infrared range (around a wavelength of 10 microns), the resolution of the best matrices is about 1 megapixel, the price increases by orders of magnitude, and their practical availability is severely limited.

Secondly, the mode of operation of a single-pixel detector is also much more flexible. For example, it can be a detector operating in the single-photon mode, which is important for registering very weak signals at the single-photon level in quantum information science. Detecting such signals in the telecommunications shortwave infrared range is a technically difficult task. One of the modern solutions is a matrix based on an array of superconducting detectors with a resolution of 1 kilopixel (32 by 32 pixels) [2]. With the single-pixel imaging, it is much easier to obtain the same resolution, using only one superconducting detector [3]. In addition, it is much easier to implement new types of detection with a single detector, for example, to discriminate weak signals with minimal error [4-7].

A distinctive feature of single-pixel imaging is that the image in this approach is calculated based on the use of a sequence of different light patterns. In this case, we measure not the spatial distribution of light intensity (which is what the matrix sensor does), but its integral amount corresponding to each pattern. Therefore, an important element of a single-pixel imaging system is computational resources, without which the method cannot work. Modern computational algorithms can restore an image even with an incomplete set of measurements using compressed sampling methods [9,7] or machine learning [10-13]. A disadvantages of the single-pixel imaging method is that the resolution of the calculated image is relatively low. The capabilities of the SPI technology are usually demonstrated at low resolutions, for example, 32 by 32 pixels, and rarely exceed 128 by 128 pixels. This is due to the fact that increasing the resolution requires an increase in data acquisition time (more light patterns are required) and more computational resources for image reconstruction algorithms to work.

In this paper, we investigate the possibilities of reducing the number of patterns and increasing the resolution of images in the SPI method via the use of light patterns based on modifications of Hadamard matrices. In particular, we show what limitations arise related to the computing resources of the SPI, and how they can be overcome using modified Hadamard matrices. We also propose new approaches to increasing the resolution of the SPI, taking into account the realistic limitations of computational resources.

2. Scaling the image resolution in the single-pixel imaging method

Let’s start with the basic definitions of the SPI method. The desired image of an object is traditionally represented as a matrix , each element of which represents a pixel. In the case of color images, a pixel consists of three brightness values of the red, green and blue channels; in the case of gray images, a pixel is represented by one value of total brightness, normalized from 0 to 1, where 0 corresponds to black and 1 to white. We will consider the latter case (gray images), and for simplicity, we will consider the image to be square with a size of pixels. Thus, in the case we are considering, the total number of elements in the matrix is real numbers from 0 to 1.

The object is illuminated by spatially inhomogeneous light (the object is sampled by a light pattern), which has a certain brightness distribution over the object, which can also be represented by a matrix , the elements of which are normalized from 0 to 1, where 0 is the minimum (zero) illumination, and 1 — maximum illumination. We refer to a pixel with the maximum brightness (one) as a light pixel, and with the minimum (zero) — as a dark pixel.

Next, we measure the total amount of light reflected from the object (or a value proportional to it), which is determined by the component-by-component product of the matrices . The signal measured in this way is a single real number, and the greater it is, the better the coincidence of the dark and light pixels of the object and the light pattern.

By sequentially changing patterns and measuring the corresponding signals, we can create a system of linear equations of the form

where — this is a matrix, where each row of length consists of a one-dimensional representation of a light pattern with sequential numbering of all pixels in the pattern, — a column of elements (pixels of the desired image, numbered by one index line by line), and — a column of measured signals. Then the calculation, or reconstruction, of the image of the object is the solution of a system of linear equations (1). The standard way to solve a system of linear equations (1) in matrix form is represented as

This requires that the matrix is a square matrix, that is, the number of patterns must be equal to the number of pixels in the desired image, and all patterns must be linearly independent (the determinant of is nonzero).

An unobvious point is related to the calculation of the inverse matrix In general, the matrix represents illumination of each pixel in each pattern. Spatially inhomogeneous light patterns can be physically implemented in different ways. At the quantum level, it is possible to control the spatio-temporal shape of the wave function of individual photons [14-16]. At the level of classical light, one can use a random speckle pattern [17], or create an arbitrary light distribution using a spatial light modulator [18]. There are two main technologies for digital spatial light modulators: liquid crystal and micro mirrors. In both cases, initially homogeneous light falls on a matrix of elements, the reflection or transmission of which can be controlled independently of each other. In the case of liquid crystal modulators, the light transmission coefficient can continuously change from 0 to 1, and in the case of micro-mirror modulators, the reflection coefficient is discrete: either 0 or 1. Accordingly, in the first case, one can create a matrix with real coefficients, and in the second case — from integers.

 

Fig. 1: The calculation time of the inverse matrix (in seconds), depending on the linear size of the image. Statistical averages and standard deviation are shown for a set of 30 random matrices of each size.

 

As we have already noted, the size of the matrix is determined by the desired resolution of the final image. To obtain an image with a size of pixels, it is required to reverse the matrix with a size of , which leads to a rapid increase in the required computational resources. Let’s consider scaling the time to invert a random matrix consisting of 0 and 1 (we obtained similar results for random matrices consisting of real numbers from 0 to 1) when using a conventional medium-sized laptop (quad-core Core i5 processor, 1.4 GHz). Calculations were performed in Wolfram Mathematica using the standard function Inverse[P]. For different , 30 random matrices of size were created and the time to invert each matrix was measured. In Fig. 1 we show the calculation results for in the range from 10 to 20, from which it can be seen that adding two pixels to the linear resolution of the image actually doubles the calculation time of the inverse matrix. For 16 by 16 images (that is, the matrix of a size 256 by 256), the calculation time is about 1 second, for 32 by 32 images — a few minutes, and for 50 by 50 images — about a day. Despite the fact that the speed of calculations can in principle be increased by choosing other software or hardware, the scaling of the calculation time means that it is impossible to significantly increase the image resolution in this approach in practice. Even considering that the inverse matrix needs to be calculated only once, and then only matrix multiplication of the inverse matrix by the measured sequence of signals is required, the calculation time of the inverse matrix is still a limiting factor for obtaining a sufficiently large resolution, for example 100 by 100 pixels or more.

3. Hadamard matrices and their modifications for single-pixel imaging

One of the possible solutions to the problem of obtaining high-resolution images by the SPI method is to use matrices for which an analytical expression for the inverse matrix is known. The obvious case is the choice of , where is a unit matrix of a size . The undoubted advantage of this choice is an extremely simple way to calculate the image . Each pattern consists of only one light pixel, and in different patterns this light pixel changes its index, that is, with this pattern, the object is sequentially “scanned” by one “running” pixel, and the amount of reflected light directly corresponds to the image of the object. The disadvantage of this choice of patterns is the small amount of reflected light, since only one pixel out of is illuminated, which in practice leads to a small signal-to-noise ratio [3,18,19].

To increase the signal-to-noise ratio, it is advisable to choose the matrix such that the number of light pixels is greater. In practice, Hadamard matrices are a frequent choice of matrices of the pattern system for SPI. According to the standard definition, the Hadamard matrix is a square matrix consisting of and , the rows of which are orthogonal to each other. One of the properties of Hadamard matrices is the ability to quickly calculate the inverse matrix , based on the orthogonality property: . According to Hadamard’s hypothesis, Hadamard matrices exist for all multiples of 4, although for the purposes of SPI, the question of constructing Hadamard matrices of arbitrary size is not very important. As we have already noted, we can consider square images with a linear size equal to a power of two (16 by 16, 32 by 32, etc.): . In this case, there is an explicit construction of Hadamard matrices in a recursive form (Sylvester construction):

Here, the index for the Hadamard matrix means a matrix of a size , and in the right part, parentheses show a matrix composed of four matrices, that is, the linear size of is twice the linear size of . At the first step of such a recursive construction, there will be a matrix of one element . Note that a similar recursive construction is found in other areas, for example, when constructing a complete set of generalized Bell states [20].

An important property of Hadamard matrices is that these matrices have the greatest determinant among all matrices with elements of . This means that the system (1) is well-conditioned and its solution will be stable in numerical calculations under the presence of noise in the measured signal.

To implement inversion of the matrix through its transposition, it is possible to normalize its elements so that the determinant is equal to 1. Then the elements of the matrix are equal to instead of . Inversion of such a modified Hadamard matrix is a transpose operation , which is computationally simple even for high image resolution.

The problem with implementing the matrix of the pattern system in the form of Hadamard matrices is that it is a matrix of positive and negative numbers ( for or for ), and the light patterns should consist of of non-negative numbers, since they are related to illumination (or the transmission of light through a spatial light modulator), which cannot be negative.

For this reason, a pair of patterns are used in SPI instead of one pattern formed out of a Hadamard matrix. The first pattern from this pair consists of a row of the Hadamard matrix, where negative elements are replaced by 0, and the second pattern from this pair is obtained from the first pattern by replacing 0 by 1 and 1 by 0, that is, the second pattern is a “negative” of the first one. Next, two signals and corresponding to the first and the second patterns are measured, and the difference between these signals is taken. This differential signal effectively corresponds to a pattern consisting of and , that is, the original Hadamard matrix. Indeed, if we consider the pixel-by-pixel difference of such patterns, then in the case when a given pixel in the first pattern is 1, in the second pattern it is 0, and the difference is . If this pixel in the first pattern is 0, then in the second pattern it is 1, and the difference is . Thus to implement the Hadamard pattern system in the SPI, a doubled set of patterns has to be shown, which doubles the exposure time of the object.

 

Fig. 2: A set of 4 by 4 pixels square patterns (left), based on a 16 by 16 modified Hadamard matrix (right). For clarity, the eighth pattern is graphically highlighted, and its correspondence to the eighth row in the Hadamard matrix is shown: line-by-line numbered pixels in the pattern make up a row in the Hadamard matrix.

 

In Fig. 2 the patterns formed from the rows of the Hadamard matrix with the replacement of by are shown on the left. For clarity, the image size is 4 by 4 pixels, and the Hadamard matrix (shown on the right) has a size of 16 by 16. As an example, the eighth pattern is graphically highlighted, and its correspondence to the eighth row in the Hadamard matrix is shown.

Let’s pay attention to a number of peculiarities of patterns based on Hadamard matrices. The first pattern (for an arbitrary image resolution) is always completely light, since the first row in the Hadamard matrix consists of only . Accordingly, the first pixel in each pattern is also always light (the first column in the Hadamard matrix consists of only ). Further, for each pattern in the set, there is a pattern transposed to it. In Fig. 2 the patterns are arranged in a 4 by 4 table, and the pairs “pattern —the pattern transposed to it” are noticeable as symmetrical pairs relative to the diagonal of the table. Some patterns match the transposed ones (shown diagonally in the table). An important advantage of patterns based on Hadamard matrices is that the number of light and dark pixels in all patterns is the same (except for the completely light first pattern), which greatly increases the signal-to-noise ratio in the practical implementation of an SPI system.

As we noted above, the standard use of Hadamard matrices in SPI requires a doubled set of patterns. Let’s show a way to reduce the required set and leave the number of patterns equal to the number of pixels in the image. Note that if we take not the difference between the signals and , but the sum of , then the result corresponds to a pattern consisting only of light pixels ( and ), which is the first pattern from the Hadamard matrix. Then we do not measure the signal , but consider it as , and the desired differential signal is . In this case, there is no need to show a doubled set of patterns and the exposure time is halved compared to the standard method.

Note that the inverse matrix in this approach is no longer equal to the transposed one. To modify the Hadamard matrix with the normalization of the determinant by 1, we have the matrix

then the inverse matrix is equal to

with the replacement of the first element by .

The resulting relationship between the Hadamard matrices defined by recursive construction (3) and the matrix makes it possible to significantly simplify the calculation of the image compared to the method of direct matrix multiplication (2). The fact is that matrix multiplication requires storing the entire matrix in the computer memory, which requires a lot of computing resources for a large number of pixels in an image. For example, for a 100 by 100 pixels image, the matrix has a size of 10000 by 10000, which is elements that need to be calculated and stored in the memory. However, a fairly simple linear relationship between the matrices (3), (4) and (5) allows us to calculate the coefficients of the matrix element by element independently of each other. From the construction (3), one can get a direct expression for the elements of the normalized Hadamard matrix in the form

where the sequences and are the binary representation of the numbers and (row and column number), respectively. In this way we avoid calculating and storing the entire matrix at once, which is extremely convenient when implementing image calculation in a practical SPI device at the level of microcontrollers and field programmable gate arrays.

4 Block sampling

Obtaining high-resolution images by selecting very large Hadamard matrices is not the most optimal solution from a practical point of view. The problem is related to the fact that when the image resolution is, for example, 256 by 256 pixels, 65536 patterns consisting of 65536 bits are needed, which requires about half a gigabyte of memory, with appropriate scaling when choosing even higher resolutions. For the electronics that control digital spatial light modulators, this can cause difficulties with implementation.

Let’s consider another approach to scaling the image resolution by splitting the original image into blocks. The idea is that we can sample the original image not entirely at once, but in separate parts, for example, in rectangular blocks. Note that when obtaining a one-dimensional representation of the light pattern (1), we used square matrices in connection with the use of Hadamard matrices with a linear size equal to the power of two. With the reverse transformation, we can split the column into a matrix that is not necessarily square. For example, an 8 by 8 pixel image requires 64 patterns, each of which is represented by an 8 by 8 binary matrix. However, we can use a block-wise partition of the image into 4 parts, each measuring 4 by 4 or 8 by 2, and sample each part with patterns based on a 16 by 16 Hadamard matrix, similar to how a 4 by 4 square image is sampled. The limiting case of such a partition is a set of 16 one-dimensional images of a size 16 by 1.

 

Fig. 3: A set of rectangular 8 by 2 pixels patterns (left), based on a 16 by 16 modified Hadamard matrix (right). For clarity, the eighth pattern is graphically highlighted, and its correspondence to the eighth row in the Hadamard matrix is shown: line-by-line numbered pixels in the pattern make up a row in the Hadamard matrix.

 

For illustrative purposes, in Fig. 3 rectangular patterns of a size 8 by 2 are shown on the left, formed from rows of a Hadamard matrix of a size 16 by 16 with the replacement of by . As an example, the eighth pattern is graphically highlighted, and its correspondence to the eighth row in the Hadamard matrix is shown, similarly to Fig. 2.

The time of complete sampling of the image is determined by the sum of the sampling times of all blocks of patterns, and the time of image calculation is the sum of the calculation times of all blocks of the image. In the case of Hadamard matrices, we use a complete set of patterns, so the total sampling time does not change compared to sampling with a larger Hadamard matrix, but the image calculation time decreases, since smaller matrix operations are required. For example, calculating a 1024 by 1024 image requires matrix multiplication with a matrix size of , which is unrealistic. In our approach, to obtain such a resolution, block sampling can be performed with 64 blocks of pixels, which is quite easy to implement on a medium-sized computer. Note that the general calculation of all individual blocks can be implemented more efficiently using specialized hardware and software, for example, graphics cards and parallel programming.

This approach can be implemented in an experimental setup, the schematic diagram of which is shown in Fig. 4. To obtain a high resolution of the object (schematically depicted tree), the spatial light modulator creates light patterns with unequal vertical and horizontal resolution (in this example, the vertical resolution is higher than the horizontal one). Next, the part of the object that overlaps with such patterns is sampled and reconstructed. After that, the sequence of patterns is repeated again, but shifted along the coordinate with a lower resolution. In this way the image of the entire object can be reconstructed, the final resolution of which may be significantly higher than the resolution of the system of patterns used.

Note that in addition to a simpler implementation, taking into account the realistic limitations of the control electronics, the proposed method opens up a number of new possibilities. For example, in this way it is possible to obtain a non-uniform resolution over the entire image field. In particular, when restoring an image, we can analyze which parts of the image have small-scale structures (neighbouring pixels differ greatly in intensity), and further increase the resolution of patterns only for these parts of the image. As a result, sampling time can be significantly reduced, while obtaining highly detailed image. The proposed method can also be used in combination with other image reconstruction methods, for example, compressed sampling [21].

 

Fig. 4: A schematic representation of a possible experimental implementation of single-pixel imaging with block sampling. To obtain an image of an object, the spatial light modulator (SLM) creates a high-resolution light pattern along one of the coordinates (in this case, the vertical), while the resolution decreases along the other coordinate. Such a pattern covers only a part of the object, the photodetector (PD) measures the amount of light from only a given illuminated part of the object, and the image is reconstructed only for this part. Next, the pattern is shifted along the horizontal coordinate, and the sampling process is repeated.

Conclusion

We have considered the method of single-pixel imaging using light patterns obtained from the modified Hadamard matrices. We have shown that it is possible to reduce the sampling time of an object by reducing the number of patterns compared to the standard approach used in most modern works. To do this, we modified both the matrix of the patterns and the corresponding inverse matrix, and also changed the procedure for measuring the signal from a single-pixel detector. If the sampling time is fixed, then the proposed method improves the image quality by increasing the signal-to-noise ratio. For the proposed modification of the pattern matrix, we have shown the procedure for calculating the image without directly calculating the entire inverse matrix at once. Further, we have shown that using single-pixel imaging with Hadamard matrices, provided realistically limited computing resources, it is possible to obtain high-resolution images via the use of block sampling, with appropriate modification of the applied light patterns.

The software developed and used in this article has no license restrictions.

Acknowledgements

The study was supported by a grant from the Russian Science Foundation No. 23-22-00381, https://rscf.ru/project/23-22-00381/.

References

1. M. F. Duarte, M. A. Davenport, D. Takhar, J. N. Laska, T. Sun, K. F. Kelly, and R. G. Baraniuk, “Single-pixel imaging via compressive sampling,” IEEE Signal Processing Magazine, vol. 25, no. 2, pp. 83–91, 2008.

2. E. E. Wollman, V. B. Verma, A. E. Lita, W. H. Farr, M. D. Shaw, R. P. Mirin, and S. W. Nam, “Kilopixel array of superconducting nanowire single-photon detectors,” Opt. Express, vol. 27, pp. 35279–35289, 2019.

3. M. Shcherbatenko, M. Elezov, N. Manova, K. Sedykh, A. Korneev, Y. Korneeva, M. Dryazgov, N. Simonov, A. Feimov, G. Goltsman, and D. Sych, “Single-pixel camera with a large-area microstrip superconducting single photon detector on a multimode fiber,” Applied Physics Letters, vol. 118, p. 181103, 2021.

4. D. Sych and G. Leuchs, “Practical receiver for optimal discrimination of binary coherent signals,” Phys. Rev. Lett., vol. 117, p. 200501, 2016.

5. M. Elezov, M. Scherbatenko, D. Sych, and G. Goltsman, “Active and passive phase stabilization for the all-fiber michelson interferometer,” Journal of Physics: Conference Series, vol. 1124, p. 051014, 2018.

6. M. Shcherbatenko, M. Elezov, D. Sych, and G. Goltsman, “Towards the fiber-optic Kennedy quantum receiver,” EPJ Web Conf., vol. 220, p. 03011, 2019.

7. M. L. Shcherbatenko, M. S. Elezov, G. N. Goltsman, and D. V. Sych, “Sub-shot-noise-limited fiber-optic quantum receiver,” Phys. Rev. A, vol. 101, p. 032306, 2020.

8. E. J. Candes, J. K. Romberg, and T. Tao, “Stable signal recovery from incomplete and inaccurate measurements,” Communications on Pure and Applied Mathematics, vol. 59, no. 8, pp. 1207–1223, 2006.

9. D. Donoho, “Compressed sensing,” IEEE Transactions on Information Theory, vol. 52, no. 4, pp. 1289–1306, 2006.

10. T. Shimobaba, Y. Endo, T. Nishitsuji, T. Takahashi, Y. Nagahama, S. Hasegawa, M. Sano, R. Hirayama, T. Kakue, A. Shiraki, and T. Ito, “Computational ghost imaging using deep learning,” Optics Communications, vol. 413, pp. 147–151, 2018.

11. C. F. Higham, R. Murray-Smith, M. J. Padgett, and M. P. Edgar, “Deep learning for real-time single-pixel video,” Scientific Reports, vol. 8, 2369, 2018.

12. A. L. Mur, P. Leclerc, F. Peyrin, and N. Ducros, “Single-pixel image reconstruction from experimental data using neural networks,” Opt. Express, vol. 29, pp. 17097–17110, 2021.

13. S. Xie, L. Peng, and L. Bian, “Large-scale single-pixel imaging via deep learning,” in Optoelectronic Imaging and Multimedia Technology IX (Q. Dai, T. Shimura, and Z. Zheng, eds.), vol. 12317, p. 1231703, International Society for Optics and Photonics, SPIE, 2023.

14. D. Sych, V. Averchenko, and G. Leuchs, “Generic method for lossless generation of arbitrarily shaped photons,” Phys. Rev. A, vol. 96, p. 053847, 2017.

15. V. Averchenko, D. Sych, G. Schunk, U. Vogl, C. Marquardt, and G. Leuchs, “Temporal shaping of single photons enabled by entanglement,” Phys. Rev. A, vol. 96, p. 043822, 2017.

16. V. Averchenko, D. Sych, C. Marquardt, and G. Leuchs, “Efficient generation of temporally shaped photons using nonlocal spectral filtering,” Phys. Rev. A, vol. 101, p. 013808, 2020.

17. D. V. Strekalov, B. I. Erkmen, and N. Yu, “Ghost imaging of space objects,” Journal of Physics: Conference Series, vol. 414, p. 012037, 2013.

18. M. D. Aksenov and D. V. Sych, “Optimal data acquisition methods for single-pixel imaging,” Journal of Russian Laser Research, vol. 39, no. 5, p. 492–498, 2018.

19. D. Sych and M. Aksenov, “Computational imaging with a single-pixel detector and a consumer video projector,” AIP Conference Proceedings, vol. 1936, p. 020016, 2018.

20. D. Sych and G. Leuchs, “A complete basis of generalized Bell states,” New Journal of Physics, vol. 11, p. 013006, 2009.

21. D. V. Sych, “Optimization of compressed sampling in single-pixel imaging,” Bull. Lebedev Phys. Inst., vol. 51, pp. 202–205, 2024.