METHOD AND DEVICE FOR IMAGE CORRECTION

Description

[0001] This application relates to methods and devices for image correction and corresponding computer programs. In particular, this application relates to such methods and devices which can be used to correct the effects of air streaks and / or for image stabilization, especially for videos, i.e., image sequences. Such image corrections, in which the effects of disturbances such as air streaks are compensated for, are also referred to as image restoration or video restoration.

[0002] When optical systems, such as telescopes or cameras, are used in the atmosphere, atmospheric turbulence can disrupt the optical image. This turbulence can be caused, for example, by solar radiation, which in turn causes atmospheric turbulence. Such optical systems are often coupled with cameras or image sensors to enable image acquisition.

[0003] Such air streaks can lead to flickering, i.e., movement, in the image and / or blurring in the image that can be viewed through the optical system or recorded via the optical system.

[0004] Unwanted movement in the image can also occur if the optical system is unintentionally moved during image acquisition, for example, due to vibrations of a tripod. With scanning (probing) image acquisition devices, such as a line sensor that captures an object line by line, or a single-pixel sensor that scans an object, additional disturbances can arise from scanning movements. An example of a line sensor is a slit-scan laser scanning microscope (LSM), while an example of a single-pixel sensor combined with object scanning is an electron microscope or a conventional laser scanning microscope. Such movements can have translational and / or rotational components; that is, both displacement and rotation can occur.

[0005] Image disturbances that cause movement, such as the flickering mentioned above or unintentional movements of the optical system, are particularly disruptive in video recordings, i.e., in the case where a large number of successive images are taken.

[0006] With digital image capture using the image sensor commonly used today, it's possible to improve captured images through subsequent computer-aided processing. For example, there are various approaches to correcting unwanted motion or flickering through image processing. Some approaches reduce flickering by registering distinctive patterns, such as prominent edges, in the image. In other words, prominent edges are identified in a sequence of images, such as a video, and the images are corrected so that the patterns align in successive frames. This approach works well when the images contain such distinctive patterns, like edges, for example, when there are many houses with windows. It's also suitable for nature shots with bushes, grass, etc.This becomes difficult because distinctive edges or other suitable patterns are less common or nonexistent. If registering distinctive patterns fails in such a case, it can also negatively impact subsequent image sharpening if a sharpening algorithm is applied to multiple consecutive images.

[0007] Another approach is so-called "speckle imaging," which is described, for example, in C.J. Carrano, "Speckle Imaging over Horizontal Paths," Lawrence Livermore National Laboratory, July 8, 2002. Further information can be found, for example, in G. Weigelt and B. Wirnitzer, OPTICS LETTERS Vol. 8, No. 7, 1983, pp. 389ff., or Taylor W. Lawrence et al., OPTICAL ENGINEERING Vol. 31, No. 3, pp. 627ff., 1992. A major disadvantage of speckle imaging is that it requires significant computing time or correspondingly powerful computers. Speckle imaging was originally used primarily in astronomy, where, at least in many applications, sufficient time for image post-processing and / or corresponding computing power is available.For terrestrial imaging, which often needs to be done in real time, implementing such an approach based on speckle imaging is at least very complex, as it requires specially designed hardware (e.g., ASICs, application-specific integrated circuits), or even impossible. Furthermore, the speckle imaging algorithm processes moving objects, such as walking people or moving vehicles, as well as air currents. This can lead to artifacts in which these moving objects are typically blurred or even partially disappear.

[0008] German patent application DE 10 2015 120 967.5 discloses a method in which images are divided into tiles and the tiles are transformed into the spatial frequency domain. Shifts are identified by evaluating the arguments of the spectral density function of the tiles in the spatial frequency domain, and at least some of the identified shifts are compensated, in particular by averaging. The tiles are then back-transformed into the spatial domain, and the back-transformed tiles are reassembled. The averaging is performed separately for the argument and magnitude of the spectral density function in the spatial frequency domain.

[0009] This approach works well with relatively low levels of interference and noise. However, limitations may arise with stronger turbulence and / or higher noise levels. In particular, larger, more complex fluctuations in the images at higher frequencies, and correspondingly fine spatial structures, may sometimes not be reliably detected and restored.

[0010] US patent 2012 / 288212 A1 discloses a method for image processing that uses an optical transfer function across multiple images. Based on this optical transfer function, an image from a sequence of images is then restored. For example, air streaks are removed or reduced.

[0011] US patent 2004 / 005098 A1 also discloses such a method, where an atmospheric telescope point spread function is used, where "telescope" refers to an image acquisition optic.

[0012] US 2014 / 270565 A1 discloses another such method, which uses a power spectral density of atmospheric turbulence.

[0013] It is therefore an object of the present invention to provide further improved methods, devices and corresponding computer programs for image correction, which deliver improved results compared to the above-mentioned prior art, particularly in the case of larger, more complex fluctuations and / or higher noise.

[0014] For this purpose, a method according to claim 1, a device according to claim 9, and a computer program according to claim 12 are provided. The dependent claims define advantageous and preferred embodiments.

[0015] Typical fluctuation signatures can be determined through experiments and assigned to specific image corrections (restoration information), which are then applied to the images when carrying out the above procedure.

[0016] In this way, unwanted motion can be removed from a sequence of images with comparatively little computational effort. Because the process operates in the spatial frequency domain, for example, compared to methods working in spatial space, no registration of prominent image parts such as edges is necessary. This allows the method to operate robustly even in situations where such prominent image parts are absent or barely present. Furthermore, if the Nyquist-Shannon sampling theorem is observed, motion compensation can be performed with sub-pixel accuracy, for example, in the range of 0.01 pixels. This is often not possible with methods operating in spatial space.

[0017] The sequence of images can be captured with a short exposure time, which can be so short that virtually no motion blur occurs, and a high frame rate. The exposure time can be less than or equal to 1 / 100 s, less than or equal to 1 / 500 s, or less than or equal to 1 / 1000 s, and / or the frame rate can be greater than or equal to 5 frames per second, greater than or equal to 10 frames per second, greater than or equal to 20 frames per second, and / or greater than or equal to 50 frames per second. With such images, lower frequencies in the spatial frequency domain, which have a high amplitude compared to higher frequencies, can be evaluated.

[0018] The use of fluctuation signatures makes it possible to determine larger, stronger, or more complex fluctuations.

[0019] According to the invention, the fluctuation signature is identified based on a comparison of the fluctuation signature in the spatial frequency space based on the sequence of images with several stored fluctuation signatures.

[0020] The restoration is then carried out based on restoration information assigned to the stored fluctuation signatures.

[0021] Calibration allows for the prior determination of which fluctuation signals will trigger the calibration.

[0022] The fluctuation signature consists of a plurality of complex-valued points, where each point is assigned to an image within the sequence of images. Thus, the fluctuation signature can describe the behavior of fluctuations over the time of the sequence of images.

[0023] The process can further include sharpening over corresponding tiles of the sequence of images before or after assembly.

[0024] The procedure may include further: Subdividing the images into one or more tiles, transforming the tiles into the spatial frequency space, identifying linear shifts in the tiles, compensating at least some of the identified linear shifts, and back-transforming the tiles into the spatial space, wherein the subdivision, transformation, identification, compensation, and back-transformation are repeated N · N ≥ 1, storing the identified linear shifts, with the fluctuation signature being based on the shifts.

[0025] The tile size when dividing images into tiles can be determined depending on the size of air streaks during image capture and / or depending on the degree of camera shake during image capture.

[0026] By gradually reducing the size, smaller air turbulence cells can also be taken into account and smaller shifts can be gradually removed, which further improves the image quality.

[0027] The factor by which the tiles are reduced in size can differ between at least two of the repetitions.

[0028] This can further reduce the computational effort by dynamically adapting the tiles used to the requirements.

[0029] The images can be color images. Subdivision, transformation, and inverse transformation can be performed separately for each color channel.

[0030] Identifying the linear shift can involve determining a proportionality factor in the spatial frequency space by fitting a plane in the spatial frequency space.

[0031] The fitting can be done based on differences from image to image in the sequence of images.

[0032] The method can involve averaging across multiple images in the spatial frequency domain. This allows noise to be suppressed.

[0033] The procedure can further include the application of a Wiener filter if the fluctuation signature is not identifiable, e.g., due to noise. This allows cases of high noise to be taken into account.

[0034] One or more steps of the process described above can be dynamically adjusted depending on the tiles and / or the content of the images. This allows for a flexible approach depending on the requirements.

[0035] One or more steps of the process described above can be performed in parallel to increase speed. For example, the transformation and / or other described calculations, which can be relatively time-consuming depending on the image size, can be performed in parallel. For example, the calculation for each tile or for groups of tiles can be performed in separate computing units, e.g., separate kernels. In other embodiments, the transformation can also be performed in parallel for different rows and / or columns of the image.

[0036] According to the invention, a device for image processing according to claim 9 is provided.

[0037] The device can be configured to perform one of the procedures described above.

[0038] The device can further comprise an optical system and a camera system coupled to the optical system for capturing the sequence of images. The camera system can, for example, have a 2D image sensor in which n x m pixels, where n, m > 1, are captured simultaneously. However, it can also have a 1D image sensor (line sensor) or a 0D image sensor (a single pixel) and operate in a scanning mode to, for example, capture a desired object.

[0039] In addition, a computer program with program code is provided which, when executed on a computing device, performs the procedure as described above.

[0040] The computer program can, for example, be stored on a data carrier or in memory, or be provided via a network.

[0041] The invention is explained in more detail below with reference to embodiments and the accompanying drawings. Figure 1 shows a schematic representation of a device according to an exemplary embodiment, Figure 2 shows a flowchart to illustrate a process according to an exemplary embodiment, and Figures 3 and 4 The diagrams illustrate some steps of the exemplary embodiment of the Figure 2 . Figure 5 shows a diagram illustrating a function in the spatial frequency domain. Figure 6 illustrates a shift in image content for two images in one dimension. Figure 7 shows a diagram illustrating the determination of a proportionality factor. Figure 8 shows data storage in a ring buffer using the method of Figure 2 . Figure 9 shows a diagram illustrating a fluctuation signature.

[0042] Several exemplary embodiments are explained in detail below. It should be noted that these embodiments serve only for illustration and are not to be interpreted as limiting. For example, embodiments with a large number of elements or components are presented. This should not be interpreted to mean that all of these elements or components are essential for the implementation of the embodiments. Rather, other embodiments may have fewer elements or components than those shown or described, and / or alternative elements or components. In still other embodiments, additional elements or components may be provided, such as elements used in conventional optical systems or computing devices. Components or elements from different embodiments may be combined with one another unless otherwise specified.Variations and modifications described for one embodiment can also be applied to other embodiments. The scope of protection of the invention and its embodiments is defined by the appended claims.

[0043] The present invention relates to image correction. For the purposes of this application, an image is understood to be a single image or an image of a video composed of a plurality of successive images. In the case of videos, such an image is also referred to as a frame. For the purposes of this application, a sequence of images refers to a series of images recorded one after the other. Such a series of images can, for example, be recorded or played back as a video, or as a series of single images recorded in quick succession.

[0044] Figure 1shows a device according to an exemplary embodiment. The device of Figure 1The system comprises a camera device 10 coupled to an optical system 16. The optical system 16 can be, for example, a telescope, a binocular telescope, or a microscope to which the camera device 10 is coupled. In another embodiment, the optical system 16 can simply be a camera lens of the camera device 10. One or more objects to be recorded are imaged onto an image sensor of the camera device 10 via the optical system 16, for example, a CCD sensor (charge-coupled device) or a CMOS sensor (complementary metal-oxide semiconductor). As already explained, the camera device can include a 2D image sensor. However, images can also be acquired line by line with a 1D image sensor or point by point (0D image sensor) by scanning.The terms camera device and image sensor are to be understood broadly here and refer to any type of device capable of capturing images of one or more objects, including, for example, electron microscopes and laser scanning microscopes (LSMs). The camera device 10 can thus capture a sequence of digital images via the optical system 16. This sequence of images can be captured with a short illumination time, which can be so short that virtually no motion blur occurs, and a high frame rate. For example, the illumination time can be less than or equal to 1 / 100 s, less than or equal to 1 / 500 s, or less than or equal to 1 / 1000 s, and / or the frame rate can be greater than or equal to 5 frames per second, greater than or equal to 10 frames per second, greater than or equal to 20 frames per second, and / or greater than or equal to 50 frames per second.In such images, lower frequencies in the spatial frequency domain, which exhibit a high amplitude compared to higher frequencies, can be evaluated. Images acquired with short illumination times may contain a comparatively high amount of noise, but this does not significantly impede the application of the techniques described here.

[0045] The device of Figure 1The system further comprises a computing unit 11, which may, for example, have one or more processors 13 and associated memory 14, such as non-volatile memory (e.g., flash memory), random access memory (RAM), or read-only memory (ROM). In particular, program code may be stored in the memory 14, which enables the execution of the image corrections described in more detail below when the program code is executed on the processor 13 in the computing unit 11. The program code can be supplied to the computing unit 11, for example, via a data carrier (memory card, CD or DVD, etc.) or via a network.

[0046] In order to process the images received from the camera unit 10, the computer unit 11 can receive the images via an interface 12 over a connection 15. This connection can be established directly via a wired or wireless connection, for example, if the camera unit 10 and the computer unit 11 are located in close proximity, such as within a microscope system. However, it is also possible for the images to be written to a memory card or other storage media in the camera unit 10 and then transferred from this storage medium to the computer unit 11. Transmission via computer networks, such as the internet, possibly with one or more intermediate stations, is also possible.The processing of the images in the computing unit can take place immediately after the images are captured, particularly in real time, but also at a later time. In other embodiments, the computing unit 11 can be integrated into the camera unit 10, for example, within a video camera. Thus, the representation of the... Figure 1 no special spatial relationship between the camera device 10 with the optical system 16 on the one hand and the computing device 11 on the other.

[0047] The computing device can perform a series of processing steps on a recorded sequence of images, which are described below with reference to the Figure 2These steps can be implemented, as already explained, by means of appropriate programming of the processor 13. The processor 13 can be a general-purpose processor such as a CPU, but it can also be, for example, a special graphics processor, another type of processor, or a combination thereof. Special graphics processors are frequently used within cameras, for example, while general-purpose processors are often used for image processing on a computer. The implementation of the present invention is not limited to any particular arrangement in this respect. In principle, it is also possible to implement the image processing described below using special hardware components, such as application-specific integrated circuits (ASICs) or field-programmable gate arrays (FPGAs).

[0048] Figure 2 shows a flowchart of a process according to an exemplary embodiment, which, for example, is explained above in the Figure 1 The discussed device, in particular the computing unit 11, may be implemented.

[0049] The process involves using a sequence of images. To process a specific image within this sequence, a certain number of neighboring images are used, for example, 15 images before and 15 images after the specific image. The specific image can also be referred to as the image to be restored. In this example, restoring a specific image would therefore require a block of 31 images, preferably taken consecutively. The number of 15 images before and after the image to be restored is merely an example, and any other number can be used, depending on the available images and computing power. The number of images used before and after the image to be restored does not have to be the same. The images are preferably taken in quick succession with a shorter aperture time.

[0050] In a sequence of images, various events can lead to changes in the image content over time, i.e., from image to image. Some of these events are undesirable and can be corrected by the method according to the invention.

[0051] The camera can perform translational and / or rotational movements, which are undesirable in the case of vibrations and the like. This movement leads to a corresponding movement of the image content. An entire recorded object can also move. Furthermore, as already explained, flickering can be caused by atmospheric turbulence. Finally, only parts of the recorded image may move, for example, when a moving object such as a car or truck is being filmed.

[0052] The procedure of Figure 2As will be explained first, this method is particularly suitable for black and white or grayscale images. Its application to color images will be explained in more detail later.

[0053] In step 20 of the Figure 2 The images in the sequence of images are divided into tiles, i.e., into sub-areas.

[0054] The subdivision into tiles from step 20 is described with reference to the Figures 3 and 4 explained. Preferably, the tiles have one or more of the following properties: The size of the tiles should be as small as possible compared to the largest areas of disturbance in the image; for example, in the case of compensating for air turbulence, they should be smaller than the largest cells of the air turbulence. Image shifts caused by camera shake and similar factors, such as unstable telescope tripod mounts, should be significantly smaller than the size of the tiles—for example, at most half, at most one-third, or at most one-quarter of the tile size. The tiles should be easily reassembled at the end of the process.

[0055] A tile size for step 20 could, for example, have an edge length of one-quarter or one-eighth of the image size. The tiles can be square. For a typical image with 2048 x 1024 pixels, a tile size of 256 x 256 pixels could be chosen, although other sizes are also possible. Preferably, the tiles overlap so that each pixel not located at the edge is present in four tiles. For illustration, [the following is shown]. Fig. 3 such a division. For example, the Fig. 3 Image 30 is divided into nine tiles (3 x 3) 31A-31J. Each tile comprises four of the small squares, with the reference symbol in the center of each. To illustrate this, tile 31A is shown in Fig. 3 Shown in bold.

[0056] As can be seen, in this division only the pixels in the corners belong to one tile, other edge pixels belong to two tiles each, and pixels in the small middle squares belong to four tiles each. Figure 3 This is for illustrative purposes only, and depending on the image size, the number of tiles may be different, especially larger, for example as in the pixel example above.

[0057] Preferably, the tiles are multiplied by a weighting function. This suppresses artifacts during the transformation in step 21. Several suitable window functions for this weighting are known individually. The weighting in two dimensions for the case of Fig. 3 is schematically in Fig. 4The three tiles are represented as follows: each is multiplied in every direction by a weighting function 40, which decreases towards the edge of the tile. In the overlapping areas, the weighting functions preferably behave such that their sum is constant, e.g., 1 (indicated by a line 41), so that the tiles can later be recombined into a complete image by simply adding the values. A suitable window function is the Hanning window function, although other window functions can also be used.

[0058] The tile decomposition is applied to each image in the sequence of images (for example, the 31 images in the numerical example above).

[0059] It should be noted that if the movement to be compensated is the same or occurs throughout the entire image, then only a single tile can be used, i.e., this single tile then encompasses the entire image.

[0060] In step 21, a transformation from the spatial domain to the spatial frequency domain is performed on each of the tiles, whereby this transformation is a two-dimensional transformation corresponding to the two-dimensional image. In a preferred embodiment, a Fourier transform, in particular a fast Fourier transform (FFT), is used. However, other types of transformation can also be used, for example, other spectral transformations, on which the exact calculation is then performed in the subsequent step 22. For example, a wavelet transform can be used. Each image in spatial space then corresponds to a corresponding spectral density in spatial frequency space.

[0061] This transformation and / or the calculations described below, which can be relatively time-consuming depending on the image size, can be performed in parallel in some embodiments. For example, the calculation for each tile or for groups of tiles can be performed in separate computing units, e.g., separate kernels. In other embodiments, the transformation can also be performed in parallel for different rows and / or columns of the image. For example, a fast 2D Fourier transform (FFT) is performed in two sections, e.g., first in a first direction (x-direction, row-wise) and then in a second direction (y-direction, column-wise). In each of these two sections, such an FFT can be performed for the different rows or columns in different kernels or otherwise with different computing units, e.g.,with 1024 lines in 1024 kernels.

[0062] In step 22, a proportionality factor for a frequency-proportional shift of the tiles in the spatial frequency domain is determined. This procedure is also referred to as "plane fitting" within the scope of this application and is explained in detail below. To simplify the following description, a one-dimensional grayscale image is assumed, i.e., the image in the spatial domain can be described by a function i(x). Two-dimensional images are treated accordingly, except that the calculations are then performed in two dimensions. Transforming this into the spatial frequency domain results in a frequency-dependent function i(f).

[0063] In the following discussion, h(x) denotes a Dirac impulse in the spatial domain, which, when applied to the image (i(x)), determines a displacement. The Dirac impulse h(x) corresponds to a function H(f) in the spatial frequency domain.

[0064] An image shifted, for example by disturbances, is denoted by g(x), and the correspondingly shifted image in the spatial frequency domain by G(f). Thus, g(x) = h(x)*i(x) in the spatial domain, where * denotes the convolution operator, or H(f) · I(f) = G(f) in the spatial frequency domain. For a linear shift, the following shift theorem applies:

[0065] Here, F denotes the spectral transformation (e.g., Fourier transform) that was used for the transformation to the spatial frequency domain in step 21. x₀ describes a shift, and e is Euler's number of the exponential function. f denotes the frequency. From equation (1) and the definitions above, we obtain:

[0066] As can be seen, the displacement in spatial space by x₀ corresponds to a change in the spectral density argument in spatial frequency space. Comparing the two equations above (1) and (2) also yields |H(f)| = 1, i.e., the magnitude of the function H(f) = 1, and the argument arg {H(f)} is proportional to -f. This is called frequency proportionality. As can also be seen, the proportionality factor is a measure of the displacement x₀. Thus, by determining the argument, the displacement x₀ for the respective tile can be determined.

[0067] Fig. 5 illustrates the function H(f) with magnitude 1 (circle with radius 1 in the in Fig. 5 (represented complex plane) and an argument that determines the position of a function value on the circle.

[0068] Examples of shifts are in Fig. 6 graphically represented. In the Fig. 6Curves 61 and 62 show examples of the shift in multiples of the pixel spacing over the coordinate x for two different images i(x) of a series of images, which are used for the procedure of Fig. 2 was recorded. Image content can be shifted in the spatial frequency domain by a multiple of a pixel spacing of the image (or by π at the highest frequency). The image ends are denoted by 63, and 60 is an example of a region in which a plane fit is performed on a tile. Starting from the total image size between the image ends 63, the images are processed in step 20 of the procedure. Fig. 2 The area is divided into tiles. In the first decomposition, a tile in an example is as large as area 60. A layer fit is then performed for this area (as well as for other tiles corresponding to other areas).

[0069] The Fig. 7Figure 61 shows the argument of the function H(f) over the frequency. The slope of curve 61 corresponds to the proportionality factor x₀, i.e., a shift. For more complex shifts, an average slope is determined for each tile size in one dimension, for example, by fitting a straight line, which in the two-dimensional case (i.e., for two-dimensional images) corresponds to fitting a plane. Therefore, this determination of the shift can also be called a plane fit.

[0070] For each image in the sequence of images, the spatial frequency-transformed tiles G(f) are defined, i.e., G1(f), G2(f)..., where the index denotes the image number in the sequence. Now, Hn(f) must be determined for each image (where n again denotes an image number). Initially, it is assumed that I(f) remains constant from image to image, meaning the image is static apart from disturbances. The consideration of moving components in the images will be explained later.

[0071] Differences arg{G n (f)}-arg{G n-1 (f)} are now calculated, representing the fluctuations from image to image. As explained above, H(f) · I(f) = G(f), and multiplying the complex values ​​H(f) and I(f) corresponds to multiplying their magnitudes and summing their arguments. Since the magnitude of H(f) is equal to 1, as explained, this difference calculation yields the function H(f), which characterizes the shift from one image to the next.

[0072] In other words, a linear shift is determined for each tile in step 22 of the procedure. This linear shift is then subtracted from the tile's value.

[0073] In step 23, a reverse transformation into spatial space is then performed, for example by an inverse Fourier transformation.

[0074] Step 24 checks whether a desired number of calculations has been completed. If not, smaller tiles are used in step 25 (for example, the tile size is halved), and steps 21 to 23 are repeated with the smaller tiles.

[0075] The plane fits are therefore performed with progressively smaller tiles. The determined shift is passed from tile level to tile level. For each iteration, the calculated shift, i.e., the determined proportionality factor in the plane fit, is stored. By subtracting the calculated linear shift in each plane fit, the residual fluctuation from image to image, even at the highest frequencies, remains below ±π in the spatial frequency domain for the progressively smaller tiles. The processing and storage of the images can be done in a ring buffer, as schematically represented as a ring buffer in Fig. 8is shown. Fig. 8Figure 80 shows a ring buffer 80. A pointer 83 points to an image position (i.e., to an image within a sequence of images) that is currently being updated. A pointer 82 points to a position of the result of the image restoration, for example, opposite the position that is currently being updated. Memory locations 21 contain everything necessary for the execution of the procedure for this image, in particular shifts, which are included in each iteration in step 22. In some embodiments, an old portion can also be subtracted from the mean during the update, and a new portion added, to avoid a complete recalculation of the mean. Such a mean can generally be used to reduce noise. For example, |G(f)| can be averaged over all or several images of the respective sequence of images, for example, over 31 images in the numerical example above.In an implementation using the ring buffer 80 of the . Fig. 8 In one implementation, an average can be performed for each target (pointer 82) across all images in the ring buffer 80. In this case, the entire average is performed for each target. Alternatively, an old value |G k (f)| / 31 (for 31 images) is subtracted from the average, and a new value |G i (f)| / 31 is added, where i is the index of the image that falls out of the average when the pointer "moves forward," and k is the index of the image that is added. This reduces the computation time.

[0076] Regarding the layer fit, the following should also be noted: The magnitudes of I(f) typically decrease significantly from the origin by several orders of magnitude as the magnitude of the frequency f increases, up to a highest detected, upper frequency. Therefore, in exemplary implementations for the plane fit, only the values ​​of G(f) at low spatial frequencies are used, since noise can dominate at high values. At medium and high frequency magnitudes, it can happen that arg{Gn(f)}-arg{Gn-1 (f)} is no longer within the range of ±π. Multiplication by the weighting function ( Fig. 4 ) in the spatial domain, in the frequency domain, there is a convolution that changes the argument of G(f) and, with regard to the plane fit, can especially distort H(f) at small frequency values.

[0077] The above aspects can complicate the plane fit. To solve or mitigate these problems, various measures can be taken in exemplary embodiments of the invention: For each tile size, the function G(f) is only used at lower frequencies below a frequency threshold that can be determined experimentally during calibration. Very low spatial frequencies below a further frequency threshold, which can also be determined during calibration, are not used. The weighting function of Fig. 4 is chosen appropriately.

[0078] In some embodiments, if components at lower spatial frequencies of the fluctuation change only slowly over time, differences from one image to the next in a sequence of images can be used for the plane fit, which results in a number of stages in tiling, i.e., a number of iterations of the loop. Fig. 2 , can reduce.

[0079] After the last iteration of steps 21-23 of the procedure Fig. 2The remaining fluctuation, plus rotational components not captured by the linear shifts compensated by the plane fits, leaves a range of values ​​for the argument of less than ±π. To this range, the fluctuations determined for the larger tiles with the plane fits, along with the shifts of a few ±π, are added. Then, in step 26, a fluctuation signature is detected and compensated based on this signature. An example of such a fluctuation signature for a frequency f is shown in Fig. 9Given a spatial frequency f, the values ​​Gn(f) can be calculated from I(f) by multiplying it by the fluctuation signature (where Gn(f) is defined more generally here than in equations (1) and (2). Here, Gn(f) is the result of the fluctuation signature, not just the linear shift as in equations (1) and (2), and thus denotes the image in the spatial frequency space with all image disturbances that contribute to the fluctuation signature). The abscissa represents the magnitude and the ordinate the argument of the complex fluctuation signature, with the arguments for higher frequencies having a range from - several π to + several π, where one π at the highest frequency corresponds to a shift of one pixel in the spatial domain. The individual points of the fluctuation signature are therefore assigned to different Gn(f). The smallest amounts are typically in the range of 0.3 to 0.1; the stronger the air turbulence, the lower the amounts.

[0080] To identify the fluctuation signature, a variety of typical fluctuation signatures can be stored through calibration measurements, which are briefly explained below. By comparing an existing fluctuation signature with the stored fluctuation signatures, a fluctuation signature can be identified. After identifying the fluctuation signature, the spectral density value corresponding to value 1 of the fluctuation signature is entered into the restored spectral density. This value corresponds to the mean value obtained by averaging, as mentioned in the aforementioned patent application DE 10 2015 120 967.5. As explained above, the values ​​Gn(f) can be calculated from I(f) by multiplying it by the fluctuation signature.Conversely, if the values ​​G n (f) are known for a given frequency (by spectral transformation of the recorded images) and a suitable fluctuation signature has been identified, the value I(f) at the point with magnitude 1 and argument 0 of the fluctuation signature is obtained (cf. . Fig. 9 This is the value we are looking for and is stored as the restoration result, e.g., in ring buffer 80 of the Fig. 8In general, the fluctuation signature extends like a cloud over a certain number of images, over many or even all spatial frequencies of a current smallest tile used in the process, or even over all smallest tiles. For an even better implementation, the contents of the entire ring buffer can be used (i.e., not just a certain number of images) to further improve the restoration at a frequency f. For the sake of simplicity, the restoration for a frequency f has been explained here in the one-dimensional case. For two-dimensional images, the frequency f is to be replaced by a two-dimensional frequency fx, fy for each smallest tile i,j (i, j being the coordinates of the tile) of a given image.Each fluctuation signature is therefore assigned a spatial frequency fx, fy of a tile and an image, at which the restored spectral density value is entered into the restoration result.

[0081] The fluctuation signature can also take into account the distance in the complex plane between images. In the example of the Fig. 9 Each "point" can be assigned a number corresponding to a specific image in the sequence of images. Based on this, in one embodiment, each or almost each point can be assigned a distance in the complex plane. Fig. 9 This can be assigned to a point in a previous image, either as a complex difference or as a difference of arguments. Based on this, the degree of fluctuation from image to image can be taken into account, thereby further improving the restoration.

[0082] For shorter time intervals, such as minutes, the fluctuation characteristics typically remain unchanged in the presence of image disturbances, and a limited number of images are sufficient to capture the fluctuation signature, for example, images taken during one second. A weighting function can be used to give less weight to the values ​​from images taken at the beginning and end of the time interval. After a longer period, such as half an hour, the fluctuation signature may change, for example, due to a change in heat radiation and the associated turbulence.

[0083] The fluctuation signature is initially multidimensional. For practical handling and identification, there are various ways to reduce the complexity of this multidimensional fluctuation signature and thus facilitate processing. These methods may also depend on the intensity and characteristics of the air turbulence. In principle, approaches generally known from pattern recognition can be used.

[0084] An example of this is a decomposition into basis functions and scaling values, where different scaling values ​​reflect air turbulence of varying intensity. In many cases, it may be sufficient to consider only a few dimensions (e.g., dependencies), as in Fig. 9 depicted. In the Fig. 9Specifically, only one dependency between images is considered, which generally yields good results. However, other dependencies, dimensions, and / or information can also be taken into account, such as those resulting from all the data in ring buffer 80.

[0085] As briefly mentioned earlier, separate experiments and calibrations can clarify which complex value should ultimately be used as the restoration result for each detected fluctuation signature, instead of the conventional averaging. This approach includes compensating for contrast losses due to air turbulence, as well as those caused by optical imaging, using a Viennese inverse filtering approach. The fluctuation signature typically depends on the intensity of air turbulence, the distance between the object and the lens, its position in the image, and the scene being captured.

[0086] Using fluctuation signatures instead of conventional averaging over arguments (and, if applicable, magnitudes) is more robust and less susceptible to noise. Furthermore, fluctuation signatures facilitate the separation of effects caused by air turbulence from those caused by moving objects. In addition to air turbulence, camera and lens movements caused, for example, by wind or small, unintentional camera pans, can also be compensated for.

[0087] The fluctuation signatures can be obtained in various ways during the calibration process mentioned above. For example, a typical scene can be captured under varying intensities of air turbulence, ranging from no turbulence to the strongest turbulence present. The image captured without air turbulence then essentially becomes the target for the restoration. As the turbulence intensifies, the various fluctuation signatures appear and are related to the restoration goal. This means that the necessary restoration procedures are determined to transform the image captured with turbulence into the target image. Both the fluctuation signature and its relationship to the restoration goal are stored, and the fluctuation signatures can be decomposed into basis functions and scaling values, as mentioned above.

[0088] This process also determines the steps involved in later decomposing the images into progressively smaller tiles during restoration. The intensity of air turbulence is taken into account. As a result, in some implementation examples, the tile size is gradually reduced so that fluctuation signatures extending over many πs can still be reliably captured during restoration.

[0089] It should be noted that different fluctuation signatures can occur in different image areas.

[0090] In some embodiments, a plane fit with respect to the change in the argument according to the displacement theorem from image to image in the spatial frequency domain can be performed to compensate for the interfering motion and also to determine a motion to be separated (i.e., a motion that actually exists, e.g., due to the movement of a recorded object). This can be used particularly with 2D image sensors. A motion to be separated might occur, for example, when video recordings are made from an aircraft. The motion to be separated is then caused by the movement of the aircraft. The same applies to other vehicles. Such motion can also occur in other applications, e.g., with a surgical microscope, for example, the movement of a surgeon with their instruments might be related to the movement of a microscope head above the stationary patient (which might need to be compensated for).

[0091] If movement occurs in only some parts of the image, multiple tiles are preferably used. This allows image areas with moving objects to be treated differently than image areas without such objects. This preserves the movement of objects as much as possible while compensating for unwanted movements. The number and size of the tiles can be dynamically adjusted, for example, by using pattern recognition to identify where and in which tiles moving objects occur.

[0092] To account for movement of a portion of the captured image (for example, a car or truck), some implementations allow tiles where such an event is detected—for instance, based on a relatively large movement in the respective image area—to be left unchanged (unrestored). To account for such effects even more precisely, tiles can be used that are small in relation to the size of air turbulence. In this case, spatially and temporally adjacent tiles are used to detect movements of only parts of the image and subsequently perform compensation if necessary.

[0093] If it is detected that a tile only exhibits flickering due to air streaks, this can be compensated as described above, i.e., by plane fits and fluctuation signatures.

[0094] If a flicker is detected in a tile and a moving object is also fully contained within the tile, the two can be separated by examining the dependency over time, since flicker causes fluctuations, while, for example, the displacement in the case of movement increases linearly. This separation then allows only the fluctuation over time caused by flicker to be corrected.

[0095] If it is detected that a moving object in a tile is only partially present in the area of ​​that tile, in addition to flickering, the time-shifting effect of flickering in neighboring tiles can be determined by analyzing them. This shift can then be interpolated for the current tile and subsequently used for compensation. For example, the shift from neighboring tiles can be averaged. Compensation is again based on the shift theorem with a frequency-proportional factor.

[0096] In addition to the operations described above, other conventional methods can also be applied to the image, such as image sharpening. Furthermore, the described method can also be corrected using the methods described in DE 10 2015 120 967.5.

[0097] A threshold can be set for detecting the fluctuation signature. This threshold determines the point at which a signature is considered recognized, i.e., how similar a present signature must be to a stored signature. If no signature can be detected within this threshold, for example due to high noise levels, the amount to be restored is set to 0 to remove informationless noise from the image. Essentially, a Wiener filter is applied to the image.

[0098] The procedure of Figure 2 As mentioned, this can initially be applied to black and white or grayscale videos. However, it can be extended to color videos or other colored sequences of images. In one embodiment, this involves dividing the video into tiles and performing the transformation or inverse transformation (steps 20, 21, 23, and optionally 25 of the process). Figure 2) are performed separately for different color channels. For example, in a sequence of RGB (red, green, blue) images, these steps are performed separately for the red channel, the green channel, and the blue channel. Steps 22 and 26 can then be performed on grayscale values ​​calculated based on the red, green, and blue components. For example, in step 22, a grayscale spectral density is calculated by averaging the spectral densities for the red, green, and blue channels, and based on this grayscale spectral density, the necessary changes to the spectral density argument are determined as described above. These changes can then be applied to the spectral densities for the red, green, and blue channels.

[0099] The averaging can also be performed as a weighted averaging. For example, it can be taken into account that the human eye is more sensitive to green light than to red or blue light, so the green channel can be weighted more heavily than the red or blue channel. A similar approach can be used for averaging and / or sharpening step 26.

[0100] Furthermore, in some embodiments, as mentioned, the method of Figure 2 This can be expanded to allow for the differentiation of "real" movements from unwanted movements (for example, those caused by the aforementioned air streaks). Such "real" movements could be objects moving through the recorded images, such as a moving vehicle or a person.

[0101] In some implementations, a moving object can be detected by observing that the spectral density argument varies more in one tile than in neighboring tiles during the same time interval, for example, by more than a predefined threshold. For instance, the spectral density argument in one tile might vary by + / - 90 degrees in a given time interval, while in neighboring tiles, the argument might only vary by + / - 30 degrees. This suggests that the tile in question is moving, whereas the change in the argument in neighboring tiles could be caused by atmospheric turbulence or similar phenomena.

[0102] If it is determined that a moving object is present in a particular tile, several approaches are possible. In some embodiments, such a tile can simply be displayed without restoration (i.e., without performing, for example, step 22 and / or step 26 for that tile). In other embodiments, the components of the changes, for example, in the spectral density argument for the moving object and for the flicker, can be separated to compensate for only the flicker component without distorting the moving object component. In some embodiments, this can be done based on the size of the spatial change. For example, the area of ​​large-scale flicker caused by large-scale air streaks will often be larger than that caused by small, moving objects.

[0103] If necessary, this can be done based on the tile sizes in the different passes.

[0104] In other examples, a distinction is made based on the direction of movement. A moving object often moves uniformly or at least always in the same direction over a period of time, whereas, for example, shimmering caused by air currents has a fluctuating character, so that the displacements can change direction. Other properties that distinguish disturbances caused by air currents or other unwanted movements from normally moving objects can also be used as a basis for this distinction.

[0105] As already explained, the above-described examples are applied to a sequence of images, for example, frames from a video. This allows for the compensation of disturbances such as air turbulence or vibrations. However, it is also possible to output only a single image, as is common with conventional single-frame cameras. In this case, a camera device such as the one described in Figure 1The camera setup depicted captures a sequence of images, including the final output image as well as a series of images before and after it, enabling the process described above. These images can be captured at a high frame rate, for example, 300 frames per second. This allows for the compensation of distortion or shifts caused by flicker in the final output image, as described above, and optionally, the correction of moving objects. Furthermore, sharpening can be applied, adjusted to compensate for detected disturbances, such as atmospheric turbulence.

[0106] As already mentioned, the procedure of Figure 2It can also be used for image stabilization, i.e., to compensate for vibrations or camera shake. In this case, it may be sufficient to use only one tile size (i.e., one pass of steps 21 to 23 of the Figure 2 (to use). In this application, typically only the translational displacement, i.e., a displacement in a specific direction, is compensated. For this purpose, an analysis of the displacement (i.e., the change in the spectral density argument after the transformation at 21) over time can be performed as described, distinguishing displacements due to vibrations from object movements (as described above). Thus, in some implementations, the displacement can be compensated very accurately, sometimes with sub-pixel precision.

Claims

1. Computer-implemented method for image processing, comprising: identifying a fluctuation signature in the spatial frequency domain on the basis of a sequence of images, wherein the fluctuation signature is identified on the basis of a comparison of the fluctuation signature in the spatial frequency domain on the basis of the sequence of images with multiple stored fluctuation signatures, and wherein the fluctuation signature is a characteristic distribution of the spectral density function over multiple images in the spatial frequency domain, and restoring at least one image of the sequence of images on the basis of the detected fluctuation signature, wherein the restoring takes place on the basis of restoration information assigned to the stored fluctuation signatures.

2. Method according to Claim 1, wherein the fluctuation signature has a plurality of complex-valued points, wherein each point is assigned to an image of the sequence of images.

3. Method according to either of Claims 1-2, further comprising: subdividing the images into one or more tiles, transforming the tiles into the spatial frequency domain, identifying linear shifts in the tiles, compensating for at least a portion of the identified linear shifts, and inverse transforming the tiles into the space domain, wherein the subdividing, the transforming, the identifying, the compensating and the inverse transforming are repeated N · N ≥ 1, wherein the identified linear shifts are stored, wherein the fluctuation signature is based on the shifts.

4. Method according to Claim 3, wherein identifying the linear shift comprises ascertaining a proportionality factor in the spatial frequency domain by fitting a plane in the spatial frequency domain.

5. Method according to Claim 4, wherein the fitting is carried out on the basis of differences from image to image of the sequence of images.

6. Method according to any of Claims 1-5, further comprising applying a Wiener filter if the fluctuation signature is not identifiable.

7. Method according to any of Claims 1-6, wherein one or more method steps of the method are adapted dynamically depending on the tiles and / or a content of the images.

8. Method according to any of Claims 1-7, comprising averaging in the spatial frequency domain over multiple images.

9. Apparatus for image processing, comprising: a computing device having at least one processor (13) and a memory (14), wherein a sequence of images is storable in the memory (14), wherein the processor is configured to carry out the following steps: identifying a fluctuation signature in the spatial frequency domain on the basis of the sequence of images, wherein the fluctuation signature is identified on the basis of a comparison of the fluctuation signature in the spatial frequency domain on the basis of the sequence of images with multiple fluctuation signatures stored in the memory (14), and wherein the fluctuation signature is a characteristic distribution of the spectral density function over multiple images in the spatial frequency domain, and restoring at least one image of the sequence of images on the basis of the detected fluctuation signature, wherein the restoring takes place on the basis of restoration information assigned to the stored fluctuation signatures.

10. Apparatus according to Claim 9, wherein the apparatus is configured to carry out the method according to any of Claims 2-8.

11. Apparatus according to Claim 9 or 10, further comprising an optical system (16) and a camera device (10), coupled to the optical system (16), for recording the sequence of images.

12. Computer program having a program code which, when executed on a computing device (11), causes said computing device to carry out the method according to any of Claims 1-8.