Estimation method, image processing method, and program

Abstract

To reduce processing time in a method for estimating a Fried length.SOLUTION: Provided is a method for estimating a Fried length from a motion image including a plurality of frames, comprising the steps of: acquiring first information by averaging the plurality of frames S101; acquiring second information by Fourier transforming the first information S102; acquiring third information by squaring the second information S103; acquiring fourth information by differentiating the third information twice with respect to frequency S104; and estimating the Fried length based on the third information, the fourth information, and an optical transfer function S106.SELECTED DRAWING: Figure 1

Description

[Technical Field]

[0001] The disclosure of this specification relates to a method for estimating the Fried length from a moving image affected by fluctuations. [Background technology]

[0002] To reduce the degradation of image quality of moving images due to atmospheric turbulence, a method for estimating the Fried length from moving images is known. The Fried length contains information on image quality degradation such as blurring and positional shift of the subject image.

[0003] Non-Patent Document 1 discloses a method for estimating the Fried length based on the Fourier transform of a moving image and the frame average of the moving image. [Prior art documents] [Non-patent literature]

[0004] [Non-Patent Document 1] von der Luhe, Estimating Fried's parameter from a time series of an arbitrary resolved object imaged through atmospheric turbulence, J. Opt. Soc. Am. A 1984, 1, 510-519. Summary of the Invention [Problem to be solved by the invention]

[0005] In the Fried length estimation method of Non-Patent Document 1, a Fourier transform is performed for each frame of a set of frames (frame group) that make up a moving image, and then Fourier root mean square data is calculated using absolute squaring. The process of calculating Fourier root mean square data in Non-Patent Document 1 requires a large amount of calculation and takes a long time to process because a Fourier transform is performed for each frame. [Means for solving the problem]

[0006] An estimation method according to one embodiment of the present invention is a method for estimating the Fried length from a moving image including a plurality of frames, and is characterized by comprising the steps of: acquiring first information by averaging the plurality of frames; acquiring second information by Fourier transforming the first information; acquiring third information by squaring the second information; acquiring fourth information by second-order differentiating the third information with respect to frequency; and estimating the Fried length based on the third information, the fourth information, and an optical transfer function. [Effects of the Invention]

[0007] According to the above-described estimation method, the processing time for estimating the Fried length can be reduced. [Brief explanation of the drawings]

[0008] [Figure 1] 3 is a flowchart of an estimation method according to the first embodiment. [Figure 2] 1 is a block diagram of an image processing system 1 according to a first embodiment. [Figure 3] 1 is an external view of an image processing system 1 according to a first embodiment. DETAILED DESCRIPTION OF THE INVENTION

[0009] DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS Hereinafter, embodiments of the present invention will be described in detail with reference to the accompanying drawings. In the drawings, the same reference numerals are used to designate the same components, and redundant explanations will be omitted.

[0010] It is assumed that light emitted from each point on the subject passes through a space with a fluctuating refractive index, enters the imaging system, and is incoherently imaged on the imaging surface. This results in a video in which image degradation fluctuates over time. In the following, it is assumed that the subject is stationary. The value I at spatial frequency ν of the Fourier transform of frame number t of a video image is given by Equation 1, where OTF (Optical Transfer Function) is the optical transfer function due to the imaging system and fluctuations, O is the Fourier transform of the subject, and N is the Fourier transform of noise.

[0011]

number

[0012] For simplicity of explanation, we will omit the symbols ν and t, which represent spatial frequency and frame, respectively. Considering the frame average of both sides of Equation 1, we obtain Equation 2 below.

[0013]

number

[0014] Here, the frame average of A is In general, frame averaging and Fourier transform are interchangeable, so Equation 2 can be obtained by averaging frames of real-space image data and then performing a Fourier transform. The former is referred to as frame average data, and the latter as Fourier average data. In the following, unless otherwise specified, frame average data and Fourier average data may both be referred to as frame average.

[0015] In addition, in this embodiment, it is assumed that the Fourier transform and inverse Fourier transform are performed by fast Fourier transform (FFT), but a general discrete Fourier transform may also be used. The frame average is generally considered to be equal to the number of frames captured by continuously exposing the image plane for a long period of time. Therefore, may be a single frame of A taken under long exposure, where the mean noise <n>is usually considered to be zero, so Equation 2 becomes Equation 3.

[0016]

number

[0017] Furthermore, taking the absolute squares of both sides of Equation 3 gives Equation 4. The absolute squares on the left side are expressed as Fourier mean absolute square data.

[0018]

number

[0019] Next, when both sides of Equation 4 are second-order differentiated with respect to frequency ν, Equation 5 is obtained.

[0020]

number

[0021] Here, vh and vw represent the components of the Cartesian coordinates of the ν plane. Note that vh and vw may be the same. The symbol for the sum of vx and vy in the third term represents the sum of overlapping combinations of vh and vw.

[0022] Regarding Equation 5, in the low frequency region near the origin of frequency ν, the magnitude of OTF generally has a maximum value at the origin, and therefore its derivative is 0 as shown in Equation 6.

[0023]

number

[0024] Therefore, the third term in Equation 5 becomes 0.

[0025] In this embodiment, we consider a situation in which the low-frequency structure of the object can be seen to some extent in the fluctuation. In other words, since the object spectrum is not 0 at the frequency origin, we assume a situation in which the change in the object's frequency components is more gradual than the fluctuation components. In other words, we assume Equation 7.

[0026]

number

[0027] Considering the above equations 6 and 7, equation 5 can be approximately estimated as equation 8 near the frequency origin.

[0028]

number

[0029] Dividing both sides of equation 4 by the corresponding sides of equation 8 gives us equation 9, which uses the second derivative data with respect to frequency.

[0030]

number

[0031] From Equation 9, the left side is calculated from the group of video frames, and the optical transfer function on the right side is <otf>is expressed by the statistical theory of light as a function of the Fried length r0, or is determined experimentally. In this embodiment, the statistical theory of light refers to the theory that uses Kolmogorov's turbulence theory for the fluctuation of the refractive index in a medium, or a theory based on this theory.

[0032] That is, by determining the wave structure function that represents the deterioration of light collection due to fluctuations from the power spectral density of the refractive index fluctuations, <otf>It is possible to find <otf>is also called the long-exposure OTF. Alternatively, the spectral density may be determined using non-Kolmogorov theory, and the long-exposure OTF may be determined in the same manner as above. However, this is not limiting, and the long-exposure OTF may be obtained by Fourier transforming a single frame captured by actually exposing a point image for a long time, or by transforming the average of multiple frames captured with short exposures.

[0033] In Kolmogorov's theory, there are multiple spectral density shapes, but the specific long-exposure OTF derived from the most standard one is given by the following Equation 10.

[0034]

number

[0035] Here, r0 is the Fried length, v is the frequency, f is the focal length, λ is the wavelength, and pow(A,B) is the Bth power of A. Furthermore, alpha and beta are constants. λ is also called the wavelength of interest in this embodiment. In this embodiment, alpha = 5 / 3 and beta = 3.44 are used. However, this is not limiting. If the estimated result of the Fried length r0 described below can be considered to be equal, values ​​changed from these numerical values ​​may be used.

[0036] The constant alpha in this embodiment may be 1.00 or more and 2.25 or less. It is more preferable that it is 1.25 or more and 2.00 or less. It is even more preferable that it is 1.50 or more and 1.80 or less. The constant beta in this embodiment may be 2.50 or more and 4.50 or less. It is even more preferable that it is 2.75 or more and 4.25 or less. It is even more preferable that it is 3.00 or more and 4.00 or less.

[0037] Next, we substitute equation 10 into the right-hand side of equation 9, and set vh and vw equal to the magnitude of the frequency ν. Furthermore, if we set ν to a small value near the origin, we can derive equation 11, the fitting that forms the basis for estimating the Fried length r0. Here, we use the fact that alpha is considered to be less than 2 in Kolmogorov's theory.

[0038]

number

[0039] Since the terms other than r0 in Equation 11 are known, r0 is determined. Specifically, a sampling position several pixels away from the frequency origin is set as v in Equation 11, and the Fried length r0 is calculated. Alternatively, the Fried lengths r0 estimated at multiple positions near the frequency origin may be averaged. Furthermore, Equation 11 may be transformed into a general optimization problem to solve r0. For example, the least squares method of Equation 12 as follows may be used.

[0040]

number

[0041] Furthermore, optimization may be performed by setting constraints within the range in which r0 can be taken. In this case, optimization may also be performed by adding a regularization term. In the above, the Fried length r0 was estimated using the long-exposure OTF as Equation 10. However, this is not limited to this. For example, the Fried length r0 can be estimated in the same manner as above by defining the long-exposure OTF as Equation 10 based on a general statistical theory of light, such as an equation other than Equation 10 based on Kolmogorov's theory or a non-Kolmogorov theory.

[0042] As mentioned above, the Fourier mean absolute square data | By using the square of | and its second derivative, the Fried length r0 can be estimated.

[0043] Here, with reference to Non-Patent Document 1, a conventional method for estimating the Fried length r0 will be explained. Conventionally, not only Fourier mean absolute square data but also Fourier mean absolute square data is used. Fourier mean absolute square data is obtained by calculating the absolute square from the Fourier transform I of a frame and then averaging the frames. In this case, each frame is generally captured with an exposure time that is short compared to the fluctuation of the fluctuation. The conventional method calculates the ratio between this Fourier mean absolute square data and the Fourier mean absolute square data. Focusing on the low-frequency region, this ratio can be expressed as Equation 13 using the corresponding OTF term.

[0044]

number

[0045] According to Non-Patent Document 1, the denominator of the right-hand side of Equation 13, like the numerator, is a function of r0. Therefore, the entire right-hand side of Equation 13 is an equation of r0, and as above, the Fried length r0 can be estimated by finding the left-hand side of Equation 13 from video frames.

[0046] In the conventional method of Equation 13, the Fourier transform is performed on the left side. When calculating the Fourier mean absolute square data, a Fourier transform is performed for each frame, requiring multiple Fourier transforms for the number of frames. In contrast, in this embodiment, only the Fourier mean absolute square data is used from Equation 9. With this configuration, when estimating the Fried length, the Fourier transform only needs to be performed on one frame obtained by taking the frame average, thereby reducing the processing time in the Fried length estimation method. Furthermore, in this embodiment, in addition to the method of performing second-order differentiation of the Fourier mean absolute square data in the denominator of Equation 9 at all positions in frequency space, a method of performing difference calculation twice in frequency space can also achieve the same effect.

[0047] In floating-point arithmetic used in digital image processing, the above difference calculation is a subtraction process, so the computational load is less than that of a Fourier transform, which uses multiplication and division. Alternatively, as a variant, the second-order differentiation can be performed by multiplying coordinate values ​​in real space. For simplicity of symbols, differentiation with respect to vh and vw is represented by a subscript. The symbol "*" at the top right of the character represents a complex conjugate. In this case, the second-order differentiation of the denominator in Equation 9 is expressed by Equation 14.

[0048]

number

[0049] The method for determining each term in Equation 14 will be described in detail later. The frame average at the pixel position (x, y) of the frame before the Fourier transform is Ir. At this time, can be expressed as Equation 15.

[0050]

number

[0051] Differentiating Equation 15 with respect to vh results in the Fourier transform of the multiplication of the frame average Ir in the brackets of the integrand function of Equation 16 and the coordinate value x in real space.

[0052]

number

[0053] Next, differentiating equation 16 with respect to vw similarly results in multiplication by y.

[0054]

number

[0055] The second derivative of Equation 15 can be obtained by further differentiating Equations 16 and 17.

[0056]

number

[0057] When the second derivative is repeated with respect to the same coordinates, it can be expressed by Equations 19 and 20.

[0058]

number

[0059] Also, The complex conjugate of is given by: Equation 21 can be obtained based on Equation 15 if Ir is a real number.

[0060]

number

[0061] In this case, equation 21 is differentiated to obtain equations 22 and 23.

[0062]

number

[0063] Equation 21 is differentiated twice to obtain equation 24.

[0064]

number

[0065] When the second derivative is repeated with respect to the same coordinates, it can be expressed by Equations 25 and 26.

[0066]

number

[0067] A feature of Equations 15 to 26 regarding the calculation of the right-hand side of Equation 14 is that the differential term in Equation 14 is calculated by multiplying the frame average Ir by the coordinate value. Also, the complex conjugate term in Equation 14 can be calculated by a forward complex conjugate operation or an inverse Fourier transform.

[0068] With this configuration, Equation 14 can be found by integrating coordinate values, which requires less calculation than the Fourier transform. Also, assuming that vh and vw are in the same direction is preferable because it further reduces the amount of calculation. Furthermore, performing a complex conjugate calculation in the forward direction is preferable because it further reduces the amount of calculation.

[0069] In this manner, in this embodiment, by actually measuring the Fourier mean absolute square data and differentiating it twice according to Equation 9, the Fried length can be estimated with a smaller amount of calculation than in the conventional method.

[0070] Furthermore, an image processing method can be provided that corrects image degradation and estimates the subject using the Fried length r0 estimated as described above. A specific example of this method will be described in detail below. As shown in Equation 27, the frame average Substituting the estimated r0 into Equation 9, we obtain <otf>Then, Wiener filtering is performed to obtain the estimated object image Oest.

[0071]

number

[0072] Here, IFT is the inverse Fourier transform. Γ is the regularization parameter of the Wiener filter. Also, the inverse Fourier transform is <otf>, the correction kernel is (F-1[ <otf>) As an example, the application of the Richardson-Lucy method is expressed as Equation 28.

[0073]

number

[0074] Here, the arrow in the formula indicates updating of the object Oest. The estimation result of the Fried length according to the embodiment can be used in any image processing method that corrects image degradation by characterizing it using the Fried length r0, without being limited to the above example.

[0075] The Fried length is the statistical, time-averaged Fried length of light in the pupil of an imaging system under fluctuation, and as long as it can be defined, the medium is not limited to a specific one. The frame group in this embodiment may be composed of one or more frames. For example, it may be one frame of a video captured with a long exposure, or a captured image. Alternatively, one or more frames may be extracted from multiple frames and used as the frame group for estimating the Fried length. [Example]

[0076] Next, a flow of estimating the Fried length in the first embodiment will be described.

[0077] In the flow of estimating the Fried length in the first embodiment, effective values ​​of moving images (frames), the number of frames to be processed, information on the imaging system (focal length, aperture diameter, wavelength of interest, pixel interval of the sensor), etc. are used as necessary. In addition, the method of estimating r0 in this embodiment may use digital data expressed by coefficients such as alpha and beta as Equation 11.

[0078] The value used in the estimation flow of each embodiment may be a value that is the distribution center or peak of the sensitivity characteristic of the sensor mounted in the imaging system, or may be any value within the wavelength range handled by the imaging system. The value used in the estimation flow of each embodiment may be stored in a storage device of an image processing device that handles digital data, such as a personal computer. In this case, any value (executed value) may be selected during execution of the Fried length estimation flow and used in the Fried length estimation flow, etc.

[0079] Here, each step of the flow of estimating the Fried length will be described with reference to FIG.

[0080] In step S101, frame average data (first information) is obtained for a selected number of frames of the selected moving image.

[0081] In step S102, the frame average data obtained in S101 is Fourier transformed to obtain Fourier average data (second information). In this embodiment, the Fourier transform uses FFT, although it is not limited to this.

[0082] In step S103, the frame Fourier mean data obtained in S102 is subjected to absolute squaring to obtain Fourier mean absolute square data (third information). This data corresponds to the numerator of Equation 9.

[0083] In step S104, the Fourier mean absolute square data obtained in S103 is second-order differentiated with respect to frequency to obtain second-order differential data (fourth information). At this time, the second-order differentiation is performed by repeating the difference calculation twice at each frequency position. This data corresponds to the denominator in Equation 9. Alternatively, step S104 may be performed by multiplying the frames of the video by position coordinates as described above.

[0084] In step S105, the second-order differential data obtained in S104 is divided by the Fourier mean absolute square data to obtain a ratio, thereby obtaining differential data (fifth information). This data corresponds to the value on the left side of Equation 9.

[0085] In step S106, the differential data obtained in S105 is fitted using Equation 11 to determine the Fried length r0.

[0086] In step S107, image processing of the moving image is performed based on the Fried length r0. At this time, a restoration filter is generated based on the Fried length r0 that represents blurring due to the influence of fluctuations in the moving image, and processing is performed on the moving image. Note that this is not limited to this, and a sharpening filter may be generated instead of the restoration filter, and processing of the moving image may be performed.

[0087] In this embodiment, the fitting in step S106 is performed using the value corresponding to the right side of Equation 11 for the value corresponding to the right side of Equation 9, but this is not limited to this. In this case, for example, the value may be replaced with a value generally given in the statistical theory of light, such as the representative Kolmogorov theory or non-Kolmogorov theory, or a value derived separately. Note that the fitting of r0 in this embodiment may be performed using a general optimization method, such as the least squares method expressed by Equation 12, as described above. In this case, various coefficients or constants corresponding to the coefficients beta and alpha may be used as needed. Note that the above-mentioned effective values ​​may be configured to be selectable by the user as needed.

[0088] Next, an image processing system 1 according to a first embodiment will be described with reference to FIGS.

[0089] The image processing system 1 includes an image processing device 100, an input device 111, and a display device 112, and estimates the Fried length.

[0090] The image processing device 100 includes a storage unit 101 , a processing unit 102 , and a connection unit 103 .

[0091] The storage unit 101 is a component for storing electronic data, such as a semiconductor memory, a hard disk, or a server on a network, and stores a program capable of executing the processing of each step of estimating the Freed length. It may also store an execution value used for estimating the Freed length. It may also store a program capable of executing the processing of correcting image degradation, which will be described later, as needed. Furthermore, the storage unit 101 may be configured to store the results obtained by each process.

[0092] The processing unit 102 is configured with one or more processors such as CPUs, and can execute the processing of each step of estimating the Fried length.

[0093] The connection unit 103 connects the image processing device 100 to the input device 111 and the display device 112, and transmits and receives various commands and data.

[0094] In the image processing system 1, an image being processed can be checked via the display device 112, and instructions for executing processing can be given via the input device 111. Here, the input device 111 is, for example, a keyboard or a mouse. The display device 112 is, for example, a liquid crystal display or a projector.

[0095] It should be noted that the image processing system 1 according to this embodiment is not limited to this. For example, a program capable of executing the processing of each step of estimating the Freed length may be stored in the storage medium 113 in Fig. 3. In this case, the processing unit 102 is configured to be able to read and execute the program from the storage medium 113. It should be noted that the storage medium 113 is, for example, a semiconductor memory, a hard disk, a server on a network, etc.

[0096] Although the embodiments and examples of the present invention have been described above, the present invention is not limited to these embodiments and examples, and various combinations, modifications, and changes are possible within the scope of the gist of the present invention.

[0097] The embodiments of the present invention include the following methods and configurations.

[0098] [Method 1] 1. A method for estimating a Fried length from a video including a plurality of frames, comprising: averaging the plurality of frames to obtain first information; obtaining second information by Fourier transforming the first information; obtaining third information by squaring the second information; obtaining fourth information by second-order differentiation of the third information with respect to frequency; an estimation method comprising a step of estimating the Fried length based on the third information, the fourth information, and an optical transfer function.

[0099] [Method 2] 2. The estimation method according to Method 1, comprising the step of obtaining the optical transfer function based on Kolmogorov's theory.

[0100] [Method 3] When the Fried length is r0, the frequency is ν, the focal length is f, the wavelength is λ, and the constants related to the optical transfer function are alpha and beta, The estimation method according to Method 1, wherein the optical transfer function is expressed as exp(-beta*pow((fλν / r0),alpha)).

[0101] [Method 4] The estimation method according to Method 1, characterized in that the constant alpha for the optical transfer function is greater than or equal to 1.00 and less than or equal to 2.25.

[0102] [Method 5] The estimation method according to Method 1, characterized in that the constant beta related to the optical transfer function is greater than or equal to 2.50 and less than or equal to 4.50.

[0103] [Method 6] 6. An image processing method for correcting degradation in the moving image using the Fried length estimated by the estimation method according to any one of Methods 1 to 5.

[0104] [Program 1] A program causing a computer to execute the image processing method according to any one of Methods 1 to 6.

[0105] [Configuration 1] A storage medium storing the program described in Program 1. [Explanation of symbols]

[0106] S101: Acquiring first information S102: Acquiring second information S103: Acquiring third information S104: Acquiring fourth information S106: Step of estimating Fried length< / otf> < / otf> < / otf> < / otf> < / otf> < / otf> < / n>

Claims

1. 1. A method for estimating a Fried length from a video including a plurality of frames, comprising: averaging the plurality of frames to obtain first information; obtaining second information by Fourier transforming the first information; obtaining third information by squaring the second information; obtaining fourth information by second-order differentiation of the third information with respect to frequency; an estimation method comprising a step of estimating the Fried length based on the third information, the fourth information, and an optical transfer function.

2. 2. The method of claim 1, further comprising the step of obtaining the optical transfer function based on Kolmogorov's theory.

3. When the Fried length is r0, the frequency is v, the focal length is f, the wavelength is λ, and the constants related to the optical transfer function are alpha and beta, 2. The estimation method according to claim 1, wherein the optical transfer function is expressed as exp(-beta*pow((fλν / r0), alpha)).

4. 2. The estimation method according to claim 1, wherein the constant alpha relating to the optical transfer function is equal to or greater than 1.00 and equal to or less than 2.

25.

5. 2. The estimation method according to claim 1, wherein the constant beta relating to the optical transfer function is equal to or greater than 2.50 and equal to or less than 4.

50.

6. 6. An image processing method for correcting degradation in the moving image by using the Fried length estimated by the estimation method according to claim 1.

7. A program causing a computer to execute the image processing method according to claim 6.

8. A storage medium storing the program according to claim 7.

Images