JP2008117416A - Multi-parameter highly-accurate simultaneous estimation method in image sub-pixel matching - Google Patents

Abstract

<P>PROBLEM TO BE SOLVED: To provide a multi-parameter highly-accurate simultaneous estimation method in image sub-pixel matching, which can estimate correspondence parameters between images stably, highly precisely, concurrently, at high speed with a small amount of computation. <P>SOLUTION: The multi-parameter highly-accurate simultaneous estimation method comprises: a step of determining a sub-sampling position where an N-dimensional similarity value between images obtained at discrete positions is maximum or minimum on a line in parallel with a certain parameter axis, and determining an N-dimensional hyperplane that most approximates the determined sub-sampling position; a step of determining N of the N-dimensional hyperplanes with respect to each parameter axis; a step of determining an intersection point of N of the N-dimensional hyperplanes; and a step of setting the determined intersection point as a sub-sampling grid estimation position for the correspondence parameter between images that gives a maximum value or a minimum value of N-dimensional similarity in an N-dimensional similarity space. <P>COPYRIGHT: (C)2008,JPO&INPIT

Description

  The present invention can be used in many image processing fields such as position measurement by image, remote sensing, cartography using aerial photography, stereo vision, image stitching (mosaicing), moving image registration, and three-dimensional modeling. More particularly, the present invention relates to a multi-parameter high-precision simultaneous estimation method in image sub-pixel matching.

  In many image processing fields such as stereo image processing, tracking, image measurement, remote sensing, and image registration, region-based matching is adopted as a basic process to obtain inter-image displacement.

  Region-based matching is characterized in that the size and shape of the region of interest can be arbitrarily set, the region of interest can be offset with respect to the pixel of interest, and the calculation is intuitive and can be implemented directly. When the displacement between images is obtained by sub-pixel by area-based matching, a method of estimating sub-pixel displacement by using similarity evaluation values obtained discretely in units of pixels is generally used.

  In order to obtain correspondence between images, an evaluation value called similarity is used. The similarity is a function that receives the relative positional relationship between the images and outputs the similarity between the images at the positional relationship. Since the degree of similarity is determined at discrete positions in units of pixels, when the corresponding positional relationship between images is obtained based only on the degree of similarity, the resolution is limited to the unit of pixels. Therefore, the correspondence between images with subpixel resolution is obtained by interpolation based on the similarity value.

  Conventionally, in order to increase the position resolution, subpixel estimation is often performed by fitting interpolation of similarity evaluation values. However, there has been a problem that sufficient accuracy could not be obtained because the horizontal direction and the vertical direction of the image were independently inserted when the amount of parallel movement of the image was determined with subpixel resolution.

  In short, conventionally, in region-based matching using a digital image, sub-pixel estimation using a similarity evaluation value is performed independently in the horizontal / vertical direction of the image in order to improve the resolution of displacement estimation.

  Conventionally, when the displacement between images is obtained using the density gradient method, the subpixel displacement can be obtained from the beginning. In the concentration gradient method, the horizontal and vertical directions are combined. However, when the inter-image displacement is 1 pixel or more, it is necessary to reduce the image. Normally, it is implemented at the same time as changing the image scale. Since the implementation of the concentration gradient method is an iterative calculation, there are drawbacks that it is difficult to estimate the calculation time and to implement it in hardware.

Furthermore, in the conventional method in which an image is DFT and then a complex conjugate product is obtained and inverse DFT is performed, the size of the region of interest must be 2 ⁿ , and various techniques are required. It is often used in the field of fluid measurement, and is characterized in that the amount of calculation can be reduced. However, it is rarely used in the vision field. In addition, since the assumption that the measurement target is a flat surface is used, there is a drawback that it is difficult to apply to a measurement target with unevenness.

  Conventionally, as described in detail below, the horizontal direction and the vertical direction of an image are independently estimated by interpolation. However, according to such a conventional method, there is a problem that an estimation error occurs as shown in FIG.

  Here, a conventional method in which an image is considered to be independent in the horizontal direction and the vertical direction and subpixel estimation of displacement is also performed independently is referred to as a “one-dimensional estimation method”.

In order to explain the problems of the conventional one-dimensional estimation method, first, the problem of obtaining the horizontal displacement between images is considered. The true inter-image displacement is (d _s , d _t ), and the rotation angle of similarity with anisotropy is θ _g . In the one-dimensional estimation method, R (-1, 0), R (0, 0), R (1, 0) (□ in FIG. 1) are used as the discretization similarity and
Estimated. FIG. 1 is a diagram for explaining a conventional one-dimensional subpixel estimation method. In the figure
Is the position estimated by the conventional one-dimensional subpixel estimation method using three similarity values at s = -1, 0, 1, and when d _t ≠ 0 and θ _g ≠ 0, The estimate has an error with respect to the true peak position.

Even if the maximum similarity position on the straight line t = 0 can be correctly estimated in the above equation 1, as shown in FIG. 1, the true horizontal image displacement d _s (the horizontal component of ▲ in FIG. 1) is obtained. On the other hand, a large estimation error is included. That is, when all of the following conditions are true, a horizontal direction estimation error occurs,
It becomes.
Vertical displacement d _t ≠ 0.
-Two-dimensional similarity has anisotropy.
The rotation angle θ _g ≠ 0 of the two-dimensional similarity having anisotropy.
Most images meet the above conditions.

The same applies to the vertical direction. For example, when the SSD self-similarity of the corner area in the image shown in FIG. 2A is obtained, it is as shown in FIG. This self-similarity has anisotropy and θ _g ≠ 0, which indicates that the conventional one-dimensional estimation method may cause a subpixel estimation error. In addition, there is a possibility that an estimation error occurs in the sub-pixel estimation using the texture image shown in FIG.

  In addition, depending on the field of computer vision, the search range can be limited to one dimension by using constraint conditions such as epipolar constraint, so it may be sufficient to perform high-precision matching of one-dimensional signals. However, there are many applications where it is necessary to obtain a two-dimensional displacement with high accuracy using a two-dimensional image. For example, motion estimation, area segmentation by motion estimation, target tracking, mosaicing, remote sensing, super-resolution, and the like.

  In the area-based matching of images, when the displacement between images can be expressed only by parallel movement, many methods have been proposed and used in practice. For example, as a typical method, a subpixel estimation method combining region-based matching and a similarity interpolation method estimates subpixel displacement between images as follows.

When the images for which the corresponding positions are to be obtained are f (i, j) and g (i, j), the SAD between the images with respect to the displacement (t _x , t _y ) in integer units is calculated as follows.
However, SAD represents the sum of absolute values of luminance differences (Sum of Absolute Differences), and W represents a region of interest for which correspondence is sought. SAD becomes 0 when the images completely match, and takes a larger value as the mismatch increases. Since SAD represents a dissimilar evaluation value, it is dissimilarity.

Instead of SAD, the following SSD is often used.
However, SSD represents the sum of squared differences of luminance differences. The SSD also becomes 0 when the images completely match, and takes a larger value as the mismatch increases. Since SSD also represents an evaluation value that is not similar, it is dissimilarity.

The following correlation coefficient (ZNCC: Zero-mean Normalized Cross-Correlation) may be used instead of SAD or SSD. ZNCC considers the pixel density in the attention area as statistical data and calculates a statistical correlation value between the attention areas.
However,
Is the concentration average within each region. R _ZNCC takes a value from −1 to 1, and becomes 1 when the images completely match, and takes a smaller value as the mismatch increases. A feature of ZNCC is that almost the same degree of similarity can be obtained even if there is a change in density or contrast between images.

  ZNCC is called similarity because it represents an evaluation value of similarity, but SAD and SSD represent dissimilarity. Here, in order to simplify the description, hereinafter, SAD, SSD, and ZNCC are collectively described as “similarity”.

Displacement in integer units where similarity is minimum (when SAD or SSD) or maximum (when ZNCC)
The integer unit displacement between images can be obtained by searching for. After obtaining the integer unit displacement, the subpixel displacement is estimated as follows.
Using Equation 5 above, the inter-image subpixel displacement that is finally obtained is
It becomes.

  In such a conventional method, there is a problem that a correct response cannot be obtained when there is a rotation or a scale change between images. Even if the rotation angle between the images is small, it is determined that the two images are different from each other, and there is a problem that the corresponding position cannot be found or the correspondence is detected at an incorrect position.

  On the other hand, when the correspondence between images cannot be expressed only by translation, that is, when there is a rotation or scale change between images, it is necessary to estimate the rotation and scale change at the same time. Conventionally, at this time, one of the images (template image) is translated, rotated, and scaled using interpolation processing, then matching is performed, and each parameter is changed while calculating the similarity. Or it was necessary to search for a parameter that would be the maximum value. Next, this parameter search method will be described.

When f (i, j) is a template image and g (i, j) is an input image, the parallelism (t _x , t _y ), the rotation angle θ, and the similarity R _SSD with the scale change s as parameters are Calculate as follows.
However,
Represents an interpolated image of the image g (i, j) at the position (i, j) when the translation (t _x , t _y ), the rotation angle θ and the scale change s are given. When creating an interpolation image, bilinear interpolation (Bi-Linear interpolation), bicubic interpolation (Bi-Cubic interpolation), or the like is used. Also, SAD or ZNCC may be used instead of SSD.

An update parameter that calculates the initial parameter (t _x , t _y , θ, s) ^<0> by some method, calculates R _SSD (t _x , t _y , θ, s) ^<0> , and gives a smaller R _SSD (t _x , t _y , θ, s) ^<1> is searched. Continue updating until the parameter is updated and the results no longer change. At this time, numerical solutions such as Newton-Raphson method, Steepest (Gradent) Descent method (steepest descent method), and Levenberg-Marquadt method are used.

The initial parameters (t _x , t _y , θ, s) ^<0> were estimated by normal region-based matching, assuming that the correspondence between images can be expressed only by translation, for example.
Used
Etc.
Masao Shimizu and Masatoshi Okutomi, "Precise Subpixel Estimation Method for Image Matching", IEICE Transactions D-II, July 2001, J84-D-II, No. 7, p.1409 -1418 Masao Shimizu and Masatoshi Okutomi, “The Meaning and Properties of Precise Subpixel Estimation in Image Matching”, IEICE Transactions D-II, December 2002, J85-D-II, Vol. 12, p.1791-1800 Masao Shimizu and Masatoshi Okutomi, “Precise Sub-Pixel Estimation on Area-Based Matching”, Proc. 8th IEEE International Conference on Computer Vision (ICCV2001) (Proc. 8th IEEE International Conference on Computer Vision (ICCV2001)), (Vancouver, Canada), July 2001, p.90-97 Masao Shimizu and Masatoshi Okutomi, “An Analysis of Sub-Pixel Estimation Error on Area-Based Image Matching”, Prok. 14th International Conference on Digital Signal Processing (DSP2002) (Proc. 14th International Conference on Digital Signal Processing (DSP2002)), (Greece, Santorini), July 2002, Volume II, p.1239-1242 (W3B.4) C. Quentin Davis, Zohar Z.・ C. Quentin Davis, Zoher Z. Karu and Dennis M. Freeman, “Equivalence of subpixel motion estimators based on Equivalence of subpixel motion estimators based on optical flow and block matching), IEEE International Symposium for Computer Vision 1995, Florida, USA, November 1995, p.7-12. Samson J. Timoner and Dennis M. Freeman, "Multi-image gradient-based algorithms for motion estimation", Optical Engineering, (USA), September 2001, Vol.40, No.9, p.2003-2016 Jonathan W. Brandt, “Improved Accuracy in Gradient-Based Optical Flow Estimation”, International Journal of Computer Vision , (USA), 1997, Vol. 25, No. 1, p. 5-22 Q. Ten, M.C.・ N.・ Q. Tian and MN Huhns, "Algorithms for Subpixel Registration", Computer Vision, Graphics and Image Processing, (USA), 1986 35, p.220-233 Sean Borman, Mark A. Robertson, Robert L. Stevenson, “Block Matching Subpixel Motion Estimation From Noisy, Under-Sampled Frame Block-Matching Sub-Pixel Motion Estimation from Noisy, Under-Sampled Frames --- An Empirical Performance Evaluation ", SPIE Visual Communications and Image Processing 1999 (SPIE Visual Communications and Image Processing 1999) Processing 1999), (California, USA)

  Such multi-parameter optimization has a problem that a correct result cannot be obtained unless the initial value is appropriate. There is also a problem that a solution cannot be obtained due to the influence of image pattern, noise, optical distortion of the taking lens, and the like. In addition, since iterative calculation is used, there is a problem that the calculation time cannot be estimated. For this reason, in implementing the algorithm, the total response time of the completed system cannot be accurately estimated, and it has been extremely difficult to implement the hardware. Further, in terms of calculation time, it is necessary to interpolate an image at each stage of repeated calculation. However, since the image interpolation calculation requires a large amount of calculation, a lot of calculation time is required. Furthermore, in terms of estimation accuracy, the properties of the interpolated image differ depending on the image interpolation method, so that the final parameter estimation accuracy and the nature of the estimation error differ depending on the interpolation method employed.

  The present invention has been made under the circumstances as described above, and an object of the present invention is to solve the problems of the conventional method by considering the N-dimensional similarity and to reduce the inter-image correspondence parameter with a small amount of calculation. An object of the present invention is to provide a multi-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image, which enables simultaneous estimation at high speed with high stability.

  The present invention relates to a multi-parameter high-accuracy simultaneous estimation method in image sub-pixel matching, and the object of the present invention is to provide N-dimensional similarity with N (N is an integer of 4 or more) inter-image correspondence parameters as axes. Multi-parameter high accuracy in sub-pixel matching of images that simultaneously estimate the N inter-image correspondence parameters with high accuracy using N-dimensional similarity values between images obtained at discrete positions in a degree space In the simultaneous estimation method, a sub-sampling position where an N-dimensional similarity value between the images is maximum or minimum on a straight line parallel to a certain parameter axis is obtained, and the obtained sub-sampling position is determined. Obtaining the closest N-dimensional hyperplane; obtaining N N-dimensional hyperplanes for each parameter axis; and intersecting the N N-dimensional hyperplanes. And the step of setting the intersection as the sub-sampling grid estimated position of the inter-image correspondence parameter that gives the maximum value or the minimum value of the N-dimensional similarity in the N-dimensional similarity space. Is achieved.

  The present invention also relates to a three-parameter high-accuracy simultaneous estimation method in image sub-pixel matching. The object of the present invention is to provide discrete positions in a three-dimensional similarity space with three inter-image correspondence parameters as axes. A three-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image for simultaneously estimating the three inter-image correspondence parameters with high accuracy using a three-dimensional similarity value between images obtained by Obtaining a sub-sampling position at which a three-dimensional similarity value between images is maximized or minimized on a straight line parallel to a certain parameter axis, and obtaining a plane that most closely approximates the obtained sub-sampling position; Obtaining three planes for each parameter axis; obtaining intersections of the three planes; and calculating the intersection points in the three-dimensional class. Effectively be achieved as a step of a sub-sampling grid estimated position of the image between the corresponding parameter giving the maximum or minimum value of the three-dimensional similarities at degrees space.

  Furthermore, the present invention relates to a two-parameter high-accuracy simultaneous estimation method in image sub-pixel matching, and the above object of the present invention is to use a two-dimensional similarity value between discretely obtained images in a continuous region. A two-parameter high-accuracy simultaneous estimation method in sub-pixel matching of images for accurately estimating the position of the maximum value or the minimum value of two-dimensional similarity, wherein the two-dimensional similarity value between the images is relative to a horizontal axis A step of obtaining a subsampling position that is maximum or minimum on a parallel straight line, obtaining a straight line (horizontal extreme value line HEL) that most closely approximates the obtained subsampling position, and a two-dimensional similarity value between the images. The subsampling position that is maximum or minimum on a straight line parallel to the vertical axis is obtained, and a straight line (vertical) that most closely approximates the obtained subsampling position An extreme value line VEL), a step of obtaining an intersection point of the horizontal extreme value line HEL and the vertical extreme value line VEL, and a sub-pixel estimated position that gives the intersection a maximum value or minimum value of the two-dimensional similarity. This step is effectively achieved.

  According to the present invention, translation, rotation, and scale change between images can be simultaneously estimated with high accuracy. That is, according to the multi-parameter high-accuracy simultaneous estimation method for sub-pixel matching of images according to the present invention, it is repeatedly calculated not only when there is a parallel movement but also when there is a rotation or a scale change in the correspondence between images. The inter-image correspondence parameters can be estimated simultaneously with high accuracy without using. In the present invention, since iterative calculation is not used, there is no need to interpolate an image for each calculation as in the prior art.

  Further, according to the present invention, a small number of interpolated images are created using template images, and parallel movement, rotation, and scale change are performed using similarity values between images using these interpolated images. At the same time, it can be estimated with high accuracy. Unlike the parameter optimization methods that are generally used, the present invention does not use iterative calculations, so the calculation is easy and the calculation time can be predicted, and the convergence that becomes a problem when performing iterative calculations It is not necessary to guarantee the instability, and instability due to the initial value does not occur. That is, according to the present invention, it is possible to eliminate the convergence problem and the initial value problem that occur in the conventional method.

  Hereinafter, preferred embodiments of the present invention will be described with reference to the drawings and mathematical expressions.

Example 1:
In the first embodiment of the multi-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image according to the present invention, the maximum value or the minimum value of the two-dimensional similarity is obtained using the two-dimensional similarity determined in the image sampling period. The two-dimensional sub-pixel displacement (that is, the displacement between images in the horizontal and vertical directions) that gives the same is simultaneously estimated with high accuracy.

  Here, the conventional method in which the image is considered to be independent in the horizontal direction and the vertical direction and the subpixel estimation of the displacement is also independently called “one-dimensional estimation method”, whereas the image according to the first embodiment of the present invention is The multi-parameter high-accuracy simultaneous estimation method in sub-pixel matching performs two-dimensional sub-pixel estimation using two-dimensional similarity at discrete positions calculated from a two-dimensional image. It is referred to as “two-dimensional estimation method” or “two-dimensional subpixel estimation method”, and is further referred to as “two-dimensional estimation method” for short.

  The two-dimensional estimation method according to Example 1 of the present invention to be described later is a high-precision simultaneous estimation method using region-based matching, and the calculation amount is slightly increased as compared with the conventional “one-dimensional estimation method”. An overwhelmingly high accuracy two-dimensional sub-pixel estimation is possible.

(1) Two-dimensional image model and two-dimensional similarity model Here, several types of image models are taken and the two-dimensional similarity is considered. It shows that the two-dimensional similarity obtained from different image models can be approximated by the same two-dimensional similarity model.

  Region-based matching is to evaluate the similarity between two images and determine the maximum or minimum position as a displacement between the two images. This similarity is usually obtained in accordance with the sampling period of the image, and a value is obtained for the discrete displacement position. Matching in pixel units corresponds to searching for the maximum value or the minimum value from these similarities. Subpixel estimation is to interpolate these similarities to estimate the subpixel displacement position that gives the maximum or minimum value at the subpixel. In the following description, “maximum or minimum” is simply expressed as “maximum”. When using SAD or SSD similarity, it means “minimum”, and when using correlation similarity, it means “maximum”.

(1-1) Corner model First, as a two-dimensional image model capable of region-based matching and sub-pixel estimation, a one-dimensional image model
Is simply expanded to consider a corner image represented by the following equation (8). Here, σ _image is the sharpness of the density edge. This corner image is shown in FIG.
However, the image coordinate system is uv, and the similarity coordinate system (displacement coordinate system) is st.

SSD similarity when two-dimensional matching is performed using the image represented by Equation 8 as a reference image and the exact same image as an input image (in the present invention, SSD similarity is sometimes referred to as “self-similarity”). ) Can be obtained by the following equation (9). The SSD similarity indicates that the smaller the value, the more similar, and in this case, the darker the location, the higher the similarity. Since the self-similarity is a similarity using exactly the same image, the similarity at the position of displacement (s, t) = (0, 0) is the smallest. This similarity is shown in FIG.
However, although the total range is a range of interest of an arbitrary shape, in FIG. 4B, the calculation was performed using an 11 × 11 rectangular area.

  In FIG. 4A, since no shear component is included in the corners, that is, the change in density in the horizontal direction does not depend on the vertical position, the horizontal direction and the vertical direction of the two-dimensional image can be handled independently. At this time, the degree of similarity of Equation 9 can also be handled independently in the horizontal direction and the vertical direction. Therefore, when the subpixel displacement position is estimated, the horizontal direction and the vertical direction can be handled independently.

  However, an actual image is not necessarily composed of 90 degree corners. For example, when three-dimensional information is reconstructed using a stereo image obtained by photographing a building from the ground, the corner area of the building is requested to correspond to the other image, but the corner is photographed at 90 degrees. It is rare (see FIG. 2 (a)).

Therefore, two two-dimensional images having density edges at v = 0
Are images rotated by θ _{a on} the left and θ _{b on} the right, respectively.

These θ _a and θ _b are defined as follows using the rotation angle θ _g and the corner angle θ _c of the entire image.

At this time, _a two-dimensional image model I _c (u, v) expressed using I _a (u, v) and I _b (u, v) is defined by the following equation (13).

As shown in FIG. 5, this two-dimensional image model has as parameters the rotation angle θ _g and corner angle θ _{c of the} entire _image and the sharpness σ _image of the two-dimensional density edge. The range of 0 ≦ θ _g ≦ π / 4 and 0 <θ _c ≦ π / 2 is considered. The range other than this can be expressed by inverting the left and right sides of the image and inverting the shade. When θ _g = π / 4 and θ _c = π / 2, the image model is the same as the simple image model of the above equation (8). 6 (a), (b), (c), (g), (h), and (i) show examples of the image model represented by the above equation (13). 6 (d), (e), (f), (j), (k), (l), FIG. 6 (a), (b), (c), (g), (h), ( The similarity corresponding to the two-dimensional image of i) is shown.

(1-2) Repetitive Texture Model When similarity is not directional-dependent, that is, if the similarity is isotropic, no error occurs even if subpixel estimation is performed independently in the horizontal and vertical directions. However, when the degree of similarity is direction-dependent (referred to as “anisotropic” in the present invention), that is, when the differential value differs depending on the direction, an error will occur if subpixel estimation is performed independently in the horizontal and vertical directions. May occur.
In the two-dimensional image model I _c (u, v) shown in the above equation 13, anisotropy appeared in the self-similarity when θ _c ≠ π / 2. However, as can be understood intuitively, there are other textures that cause anisotropy in self-similarity. As an example of this, the SSD self-similarity of the image shown in FIG. 3A is shown in FIG. Although the corner included in the texture is θ _c = π / 2, anisotropy appears in the similarity.

As a two-dimensional image model that causes anisotropy in self-similarity even though θ _c = π / 2, the following equation 14 is considered.

As shown in FIG. 7, the two-dimensional image model I _t (u, v) of Equation 14 has u-direction spatial frequency f _u , v-direction spatial frequency f _v , and image rotation angle θ _g as parameters. Consider a range of 0 ≦ θ _g ≦ π / 2. When f _u ≠ f _v , anisotropy occurs in the self-similarity. The anisotropy of the self-similarity corresponds to the anisotropy of the spatial frequency in the texture. The two-dimensional image model I _c (u, v) introduced in the previous section can be considered to have spatial frequency anisotropy when θ _c ≠ π / 2.

  FIG. 8 shows the two-dimensional image model and the SSD self-similarity.

(1-3) Gaussian function model As a two-dimensional image model that causes anisotropy in the self-similarity, the following two-dimensional Gaussian function can be considered.
However, as shown in FIG. 9, σ is a standard deviation of a Gaussian function, k is an anisotropy coefficient (k> 0), and θ _g is a rotation angle. Consider a range of 0 ≦ θ _g ≦ π / 2. The two-dimensional image model and the SSD self-similarity are shown in FIG.

(1-4) Two-dimensional similarity model Unlike the one-dimensional subpixel estimation, the two-dimensional subpixel estimation has a large number of parameters for the two-dimensional image model, so the SSD similarity obtained from the two-dimensional image model should be examined directly. Is not a good idea. The image model picked up here includes many of the properties of the image, but images handled in actual use include a combination of a plurality of the image models picked up here and higher-order images.

Therefore, in the subsequent study, the two-dimensional similarity model that approximates the two-dimensional similarity obtained from these image models is examined instead of directly using the several types of image models studied so far. As a two-dimensional similarity model, a two-dimensional Gaussian function expressed by the following equation 16 is adopted.
However, (d _s , d _t ) is a true displacement between images, σ is a standard deviation of a Gaussian function, k is an anisotropy coefficient (k> 0), and θ _g is a rotation angle. Consider a range of 0 ≦ θ _g ≦ π / 2.

For comparison, FIG. 11 shows an example of the SSD self-similarity obtained from the two-dimensional image model taken up here and an example of the similarity of the above equation (16). However, when (d _s , d _t ) = (0,0), θ _g = 0, k = 1, so as to best approximate the discrete self-similarity calculated from the two-dimensional image model, The above 16 parameters were obtained by numerical calculation.

  Although the above description has been performed in continuous regions, the similarity obtained from actual image data is sampled in units of pixels. This is shown in FIG. The two-dimensional estimation method according to the first embodiment of the present invention estimates the two-dimensional displacement that gives the maximum similarity value in a continuous region with high accuracy using the two-dimensional similarity obtained in discrete units as described above. To do.

(2) Two-dimensional estimation method according to Embodiment 1 of the present invention (2-1) Principle of the two-dimensional estimation method according to Embodiment 1 of the present invention <Description in a continuous region>
Here, the principle of the two-dimensional estimation method according to the first embodiment of the present invention will be described using the above-described two-dimensional similarity model. The present invention is characterized in that a similarity maximum value position in a continuous region is estimated with high accuracy using a similarity value obtained in a discretization unit. Here, the principle of the first embodiment of the present invention will be described in a continuous region.

FIG. 13 shows a two-dimensional similarity model displayed with contour lines. The parameters of this two-dimensional similarity model are (d _s , d _t ) = (0.2,0.4), σ = 1.0, k = 0.7, and θ _g = π / 6. The contour lines indicate levels from 10% to 90% with respect to the maximum value of the two-dimensional similarity model. The contour line of the two-dimensional similarity model is an ellipse. The center of this ellipse is the maximum value position of the two-dimensional similarity model. The FIG. 13 (a), the is shown a major axis of the two-dimensional similarity model, this angle corresponds to the rotation angle theta _g.

When the two-dimensional similarity model of Equation 16 is partially differentiated by s and set to 0, the relationship representing the straight line s = at + b in FIG. 13B can be obtained.

Since the condition to consider is k> 0, the denominator of each coefficient of the above equation 17 is not equal to 0. Further, when the two-dimensional similarity model is isotropic, that is, when k = 1, the above equation 17 becomes s = _ds , and no estimation error occurs. In the first embodiment of the present invention, this straight line is referred to as a horizontal extreme value line HEL (Horizontal Extremal Line). HEL passes through a contact point between an ellipse representing a contour line in FIG. 13B and a horizontal straight line.

On the other hand, when the two-dimensional similarity model of Equation 16 is partially differentiated with respect to t and set to 0, the relationship representing the straight line t = As + B in FIG. 13B can be obtained.

Similar to Equation 17, the denominator of each coefficient is not zero. When k = 1, Equation 18 becomes t = _dt , and no estimation error occurs. In Embodiment 1 of the present invention, this straight line is referred to as a vertical extreme value line VEL (Vertical Extremal Line). VEL passes through a contact point between an ellipse representing a contour line in FIG. 13B and a vertical straight line.

The intersection of HEL and VEL is a point where when the two-dimensional similarity model of Equation 16 is partially differentiated in different directions, it becomes 0 in both directions. This is the maximum similarity value of the two-dimensional similarity model. It is nothing but the position to make. Actually, when calculating this intersection,
As a matter of course, the displacement (d _s , d _t ) of the two-dimensional similarity model is expressed.

In this case, the denominator ≠ 0 is a condition for the intersection of HEL and VEL. The denominator is calculated as follows.
Since the range of k is k> 0, this denominator 1−aA ≠ 0.

  In summary, the s-direction differential and the t-direction differential of the two-dimensional similarity can be expressed by two lines of HEL and VEL, and the intersection of the two lines is the maximum value position of the two-dimensional similarity. Therefore, in Embodiment 1 of the present invention, two-dimensional subpixel estimation can be performed by calculating HEL and VEL using similarity values obtained in discretization units and calculating their intersections.

(2-2) Implementation Method of Two-Dimensional Estimation Method According to Embodiment 1 of the Present Invention <Application to Discrete Domain>
Hereinafter, the two-dimensional estimation method according to the first embodiment of the present invention that estimates the similarity maximum value position in the continuous region using the discretely calculated similarity values between images will be specifically described.

  First, HEL is obtained using discrete similarity values. Since HEL is a straight line that passes through the maximum similarity value in the horizontal direction, HEL can be obtained by determining the subpixel similarity maximum value position on two or more horizontal lines. This is shown in FIG. In the figure, concentric ellipses represent contour lines of a two-dimensional similarity model. The HEL indicated by a straight line is obtained using the similarity values obtained discretely at the square positions in the figure.

In FIG. 14A, the position giving the maximum similarity on the straight line t = 0 is
If the similarity at the displacement position (s, t) is R (s, t),

This estimation result includes an estimation error as described in the previous chapter, which will be described in the next section. The position giving the maximum similarity on the straight line t = -1 and the straight line t = 1 respectively
Then,
However, i _{t = −1} and i _{t = 1} are integer position offsets that give the maximum similarity on the straight line t = −1 and the straight line t = 1 (described later).

A straight line passing through the maximum similarity positions of these three points is HEL. Actually, these three points may not completely lie on a straight line due to image patterns, noise included in the image, differences between images, and the like, so an approximate straight line is obtained with the least squares. A straight line approximating these three points with the least squares can be obtained by the following equation (24).

  Next, VEL is obtained using the similarity value obtained in the discretization unit. Since VEL is a straight line that passes through the maximum similarity in the vertical direction, VEL can be obtained by determining the subpixel similarity maximum value position on two or more vertical lines. This situation is shown in FIG. In the figure, concentric ellipses represent contour lines of a two-dimensional similarity model. Using the similarity value obtained in the discretization unit at the □ position in the figure, VEL indicated by a straight line is obtained.

In FIG. 14B, the position giving the maximum similarity on the straight line s = 0 is
Then,

The position giving the maximum similarity on the straight line s = -1 and on the straight line s = 1 respectively
Then,
However, i _{s = −1} and i _{s = 1} are integer position offsets that give the maximum similarity on the straight line s = −1 and on the straight line s = 1.

A straight line passing through the maximum similarity positions of these three points is VEL. A straight line approximating these three points with the least squares can be obtained by the following equation (28).

The intersection of HEL and VEL, that is, the intersection of Equations 24 and 28 is the subpixel estimated position of the maximum two-dimensional similarity in the continuous region.
It has become.

By the way, the integer position offsets i _{t = −1} and i _{t = 1} that give the maximum similarity on the straight line t = −1 and the straight line t = 1 in the equations 26 and 27 are not necessarily guaranteed to be zero. For example, FIG. 15B shows a two-dimensional similarity of (d _s , d _t ) = (0.2,0.4), σ = 1, k = 0.3, and θ _g = π / 8. Although it is a degree model, the maximum similarity on the straight line t = −1 for obtaining HEL for this model is in the vicinity of s = −2. Therefore,
When calculating, it is necessary to re-search the integer position that gives the maximum similarity on the straight line t = −1 and the straight line t = 1. The same applies to the vertical direction.
When calculating, it is necessary to re-search the integer position that gives the maximum similarity on the straight line s = −1 and the straight line s = 1.

At this time, the similarity in the range of ± 3 shown in FIG. 15B is used, and −2 ≦ i ≦ + 2 is considered. This search range is a range determined based on a large number of realistic images. If the parameter range of the two-dimensional similarity model is not limited, the re-search range at this time becomes infinite. The reason for this is that when (d _s , d _t ) = (0,0), considering D (k, θ _g ) of the following _equation 30 for the HEL of _equation 17,

The reason why D (k, θ _g ) approaches the maximum value is that when θ _g → 0, k → 0, because D (k, θ _g ) → cos θ _g / sin θ _g → ∞. By the way, when k → 0, the anisotropy of the two-dimensional similarity becomes infinite, and HEL and VEL almost coincide. Even if there is an image that creates such a two-dimensional similarity, it is only necessary to search the maximum similarity on HEL and VEL that are almost identical. For example, when D (k, θ _g ) → ∞ in Equation 30, 3 points on VEL
The similarity value at s is obtained and parabolic fitting is performed.

  In the conventional one-dimensional estimation, the subpixel position is estimated using the five similarity values shown in FIG. At this time, subpixel estimation was performed using the similarity of three points on the straight line t = 0 (□ mark) to obtain a horizontal subpixel displacement (● mark). In addition, subpixel estimation was performed using the similarity of three points on the straight line s = 0 (marked with □) to obtain a vertical subpixel displacement (marked with ●). The intersection of the two straight lines shown in FIG. 15A is a two-dimensional subpixel estimated value obtained by the conventional method, and it can be seen that a large estimation error is included.

  On the other hand, the two-dimensional estimation method according to the first embodiment of the present invention estimates the two-dimensional subpixel position using the 25 similarity values shown in FIG. Since the intersection of HEL and VEL accurately passes through the maximum subpixel position of two-dimensional similarity, the estimation is extremely accurate. Furthermore, the increase in the calculation amount of the two-dimensional estimation method according to the first embodiment of the present invention is small compared to the conventional “one-dimensional estimation method”. In region-based matching, the similarity value is calculated for most of the calculation time. However, since the two-dimensional estimation method according to the first embodiment of the present invention uses the similarity value already obtained, There is little increase.

(2-3) Combination of the two-dimensional estimation method and the subpixel estimation error reduction method according to the first embodiment of the present invention Parabolic fitting using the similarity of three positions in Equations 21 to 23 and Equations 25 to 27 When the subpixel position is estimated by the estimation error, an estimation error occurs. This estimation error can be canceled by using an interpolation image.

  Since the subpixel positions obtained in Equations 21 to 23 and Equations 25 to 27 are merely one-dimensional estimation as in the past, the subpixel estimation method (described in Non-Patent Document 1) that can reduce the estimation error ( Hereinafter, in order to distinguish from the two-dimensional estimation method according to the first embodiment of the present invention, the subpixel estimation method described in Non-Patent Document 1 is referred to as “subpixel estimation error reduction method” as it is. it can. Here, the application method will be briefly described.

When the two-dimensional image function used for matching is I ₁ (u, v) and I ₂ (u, v), the horizontally interpolated image I _1iu (u, v) is
It can be expressed as. However, (u, v) at this time is an integer. It can be considered that the horizontally interpolated image I _1iu (u, v) is displaced by ( _0.5,0 ) with respect to I ₁ (u, v). When subpixel estimation is performed with the SSD similarity shown below using this interpolated image, it can be considered that the estimation result is also displaced by (0.5, 0). Can be calculated by the following equations 32 and 33.

Equation 33 also includes an estimation error as in Equation 21, but since the phase is reversed, the estimation error can be canceled by Equation 34 below.

Similarly, for horizontal lines t = −1 and t = 1,

However, i _{t = −1} and i _{t = 1} are integer position offsets that give the maximum similarity of R _iu (s, t) on the straight line t = −1 and on the straight line t = 1. These values may be different from Equations 22 and 23. Using the equations 35 and 36, the estimation error can be canceled by the following equations 37 and 38.

The same applies to the vertical direction. The vertically interpolated image I _1iv (u, v) is
It can be expressed as. However, (u, v) at this time is an integer. It can be considered that the vertically interpolated image I _1iv (u, v) is displaced by (0,0.5) with respect to I ₁ (u, v). When subpixel estimation is performed with the SSD similarity shown below using this interpolated image, it can be considered that the estimation result is also displaced by (0, 0.5). Can be obtained by the following equations 40 and 41.

Equation 41 also includes an estimation error as in Equation 25, but since the phase is reversed, the estimation error can be canceled by Equation 42 below.

Similarly, for vertical lines s = −1 and s = 1,

However, i _{s = −1} and i _{s = 1} are integer position offsets that give the maximum similarity of R _iv (s, t) on the straight line s = −1 and on the straight line s = 1. These values may be different from Equations 26 and 27. Using the equations 43 and 44, the estimation error can be canceled by the following equations 45 and 46.

The subpixel estimation values on the respective lines calculated in Equations 34, 37, 38, 42, 45, and 46 are estimated by comparing with the subpixel estimation values estimated by the conventional one-dimensional estimation method. The error has been reduced. Therefore, by using these estimation results estimated by the subpixel estimation error reduction method and further performing two-dimensional subpixel estimation by the two-dimensional estimation method according to the first embodiment of the present invention, more accurate two-dimensional Sub-pixel estimation is possible.
Since the interpolation image only needs to be created once in the horizontal direction and once in the vertical direction, the increase in the amount of calculation is similar to that in the one-dimensional estimation.

(3) Comparison between the estimation error by the two-dimensional estimation method and the estimation error by the conventional estimation method according to the first embodiment of the present invention. The estimation error due to the estimation method is compared with the estimation error due to the two-dimensional estimation method according to Embodiment 1 of the present invention. By comparing the estimation errors using the above-described two-dimensional image, the same image is compared.

  As the two-dimensional subpixel estimation error, it is necessary to consider an estimation error in the horizontal direction and the vertical direction. However, the estimation error is affected by the total of 5 parameters including the 3 parameters of the 2D image and the 2D input displacement between the images, so it is difficult to cover all of them. Therefore, the estimation error by each method is compared using a representative two-dimensional image.

(3-1) Two-dimensional estimation method according to Embodiment 1 of the present invention The above-described corner image is used as a representative two-dimensional image, and θ _g = 0, π / 8, π / 4, θ _c as parameters. = Π / 4, σ _image = 1.0. Two-dimensional images other than 0 ≦ θ _g ≦ π / 4 can be expressed by reversing the left and right sides of the image and by reversing the grayscale, but the actual examination is more two-dimensional images, rotation The angle and corner angle are used. In other two-dimensional images, when the two-dimensional similarity has the same shape, the subpixel estimation error shows the same tendency.

The horizontal direction estimation error error _s and the vertical direction estimation error error _t in the two-dimensional subpixel estimation method according to the first embodiment of the present invention to which the subpixel estimation error reduction method is applied are expressed by the following Expression 48 using Expression 47: be able to.

The horizontal direction estimation error error _s and the vertical direction estimation error error _t in the two-dimensional subpixel estimation method according to the first embodiment of the present invention when the subpixel estimation error reduction method is not applied, Can be expressed as

The two-dimensional subpixel estimation errors of Equations 48 and 49 for three images of θ _g = 0, π / 8, π / 4 are respectively shown in FIGS. 16, 17, and 18 (c), (d), and (e). Shown in (f). (c) and (d) respectively show the horizontal direction estimation error error _s and the vertical direction estimation error error _t in the two-dimensional subpixel estimation method according to the first embodiment of the present invention to which the subpixel estimation error reduction method is applied. (e) and (f) are respectively the horizontal direction estimation error error _s and the vertical direction estimation error error _t in the two-dimensional subpixel estimation method according to the first embodiment of the present invention when the subpixel estimation error reduction method is not applied. Show.

FIG. 16 shows the result using an image of θ _g = 0, but no anisotropy occurs in the two-dimensional similarity at this time. FIG. 17 and FIG. 18 show the results using images of θ _g = π / 8 and π / 4, respectively, and in this case, the subpixel estimation error despite the occurrence of anisotropy in the two-dimensional similarity. Is extremely small.

(3-2) Conventional One-Dimensional Subpixel Estimation Method The conventional one-dimensional estimation method can be obtained by Equation 21 and Equation 25. In addition, an estimated value to which the subpixel estimation error reduction method described in Non-Patent Document 1 is applied can be obtained by Equations 34 and 42. Here, the estimation error when using Equations 21 and 25 is compared.

The above-mentioned subpixel estimation errors of the number 50 for the three images θ _g = 0, π / 8, and π / 4 are shown in FIGS. 16, 17, and 18 (g) and (h), respectively. (g) and (h) show a horizontal direction estimation error error _s and a vertical direction estimation error error _t , respectively, in the conventional one-dimensional estimation method.

FIG. 16 shows the result using an image of θ _g = 0, but no anisotropy occurs in the two-dimensional similarity at this time. At this time, an estimation error does not occur even with the conventional one-dimensional estimation method. FIGS. 17 and 18 show results using images of θ _g = π / 8 and π / 4, respectively. At this time, anisotropy occurs in the two-dimensional similarity, but at this time, a large estimation error occurs. The magnitude of the estimation error depends on the rotation angle of the two-dimensional similarity having anisotropy and the true inter-image displacement.

(3-3) Subpixel estimation by the conventional density gradient method The density gradient method includes two unknowns, a horizontal movement amount and a vertical movement amount, assuming that there is no density change at the corresponding position between images. A constraint condition is obtained for each pixel. Therefore, since these two unknowns cannot be obtained as they are, the attention area is set, and the unknown quantity is obtained by iterative calculation so as to minimize the square sum of the constraint conditions in the attention area. The amount of displacement obtained in this way is not limited in pixel units, and can be obtained in sub-pixel units.

  On the other hand, in the image interpolation method, the region of interest and the surrounding search region are interpolated more densely than the discretization unit. The maximum similarity value is searched in a closely interpolated unit. At this time, unlike the density gradient method, there is no restriction on the amount of movement between images.

  These two types of subpixel displacement estimation methods seem to be completely different implementations and treated as different methods, but are shown to be essentially the same (see Non-Patent Document 5). Here, it is confirmed that the subpixel displacement estimation error changes depending on the order of the interpolation function used when creating the interpolation image. Therefore, the subpixel estimation result by the interpolation image using the primary interpolation is regarded as the result by the density gradient method. Furthermore, the results of subpixel estimation using an interpolated image using a higher-order interpolation function are also examined.

(3-4) Subpixel Estimation by Conventional Bilinear Image Interpolation Assume that two images for obtaining subpixel displacement are I ₁ (u, v) and I ₂ (u, v). Assume that these images are discretized and sampled. That is, u and v are integers. Further, when the displacement between the images is (d _s , d _t ), 0 ≦ d _s ≦ 1 and 0 ≦ d _t ≦ 1. When the true displacement between images is not within this range, the entire image is moved by pixel-by-pixel matching.

Sub-pixel estimation value using bilinear interpolation (bidirectional primary interpolation)
Can be expressed by the following equation (51).
However, the total range is an arbitrary-shaped attention area, and argmin is the smallest value.
Represents an operation for obtaining a set of

The subpixel estimation errors of Formula 51 for the three images of θ _g = 0, π / 8, and π / 4 are shown in FIGS. 16, 17, and 18 (i) and (j), respectively. (i) and (j) show the horizontal direction estimation error error _s and the vertical direction estimation error error _t , respectively.

FIG. 16 shows the result using an image of θ _g = 0, but no anisotropy occurs in the two-dimensional similarity at this time. At this time, no estimation error occurs. FIGS. 17 and 18 show results using images of θ _g = π / 8 and π / 4, respectively. At this time, anisotropy occurs in the two-dimensional similarity, but an estimation error occurs at this time. The magnitude of the estimation error depends on the rotation angle of the two-dimensional similarity having anisotropy and the true inter-image displacement.

  This estimation error is larger than the results of (c) and (d) of FIGS. That is, the two-dimensional estimation method according to the first embodiment of the present invention can estimate the displacement with higher accuracy than the conventional subpixel estimation method based on bilinear interpolation. Since the subpixel estimation method by bilinear interpolation is equivalent to the density gradient method, it can be said that the two-dimensional estimation method according to the first embodiment of the present invention can estimate the subpixel displacement with higher accuracy than the conventional density gradient method. .

  Further, the estimation error corresponding to the density gradient method is much smaller than the results of the one-dimensional estimation methods shown in FIGS. 16, 17 and 18 (g) and (h). In other words, the density gradient method performs two-dimensional sub-pixel displacement with much higher accuracy than the conventional method of performing region-based matching on a pixel-by-pixel basis and sub-pixel displacement using one-dimensional estimation using similarity. Can be estimated. For this reason, it has been recognized that the density gradient method is more accurate than the region-based matching. However, by using the two-dimensional estimation method according to the first embodiment of the present invention, the sub-gradation method can be performed with higher accuracy than the density gradient method. Pixel displacement can be estimated.

(3-5) Subpixel estimation by conventional bicubic image interpolation When an interpolated image is created, a higher-order interpolated image can be created by considering even higher-order terms. Bicubic interpolation (bidirectional cubic interpolation) is an interpolation method that uses 4 × 4 pixel values around the subpixel pixel position to be interpolated. Sub-pixel estimation using bicubic interpolation
Can be expressed by the following formula 52.

However, the total range is an arbitrary-shaped attention area, and argmin is the smallest value.
Represents an operation for obtaining a set of At this time, various weighting functions h (x) have been proposed based on the sinc function, and the following equation 53 is often used.

However, a is an adjustment parameter, and a = −1 and a = −0.5 are often used.

Also, the following number 54 is often used.
B and C are adjustment parameters. When B = 1 and C = 0, a cubic B-spline is approximated. When B = 0 and C = 0.5, a Catmull-Rom cubic characteristic is obtained. As described above, the interpolation characteristics differ depending on how the weighting factor is selected. Here, the parameter a = −0.55 is adopted in the weighting function of Formula 53. However, a = -1 is most frequently introduced in textbooks for image processing. The difference between the two characteristics will be described later, but a = −0.55 was determined numerically so as to reduce the subpixel estimation error.

The subpixel estimation errors of Formula 52 for three images of θ _g = 0, π / 8, and π / 4 are shown in FIGS. 16, 17, and 18 (k) and (l), respectively. (k) and (l) show the horizontal direction estimation error error _s and the vertical direction estimation error error _t , respectively.

FIG. 16 shows the result using an image of θ _g = 0, but no anisotropy occurs in the two-dimensional similarity at this time. At this time, no estimation error occurs. FIGS. 17 and 18 show results using images of θ _g = π / 8 and π / 4, respectively. At this time, anisotropy occurs in the two-dimensional similarity, but an estimation error occurs at this time. The magnitude of the estimation error depends on the rotation angle of the two-dimensional similarity having anisotropy and the true inter-image displacement.

  This estimation error is comparable to the results of the two-dimensional estimations (c) and (d) in FIGS. That is, it can be said that the two-dimensional sub-pixel estimation method according to the first embodiment of the present invention is as accurate as sub-pixel estimation using a higher-order interpolated image when an optimal parameter is selected.

By the way, the parameter a in Formula 53 will be described with a specific calculation example. Equation 53 can be considered as a weighting coefficient obtained by approximating the sinc function within a finite range. The parameter a is for adjusting the approximate characteristic. FIGS. 19A and 19B show the interpolation results of h (x) when a = −0.55 and a one-dimensional function f (i) having a value at an integer position. f (i) is expressed by Equation 7 with σ _image = 0.7. FIGS. 19C and 19D show interpolation results of h (x) when a = −1.0 and a one-dimensional function f (i) having a value at an integer position.

  a = -1.0 is a parameter introduced as bicubic interpolation in many image processing textbooks. Although the sinc function is well approximated, it is difficult to say that at least the function f (i) examined here is well interpolated. Although a = −0.55 is not so good as an approximation of the sinc function, the function f (i) is well interpolated.

  In FIG. 16 to FIG. 18, the parameter a = −0.55 is used, but this value is a result obtained numerically so as to reduce the subpixel estimation error. When a = −0.5 is used, the error is about the same as the result of bilinear interpolation, and when a = −1.0 is used, the error is larger than the subpixel estimation result using bilinear interpolation.

(3-6) Comprehensive comparison of estimation errors by each method In the following, the parameters of the two-dimensional image are changed over a wider range, and the subpixel estimation errors of each method are compared and examined. The two-dimensional image uses the aforementioned corner image. As parameters to be considered, consider a rotation angle θ _g = 0, π / 16, 2π / 16, 3π / 16, 4π / 16, and a corner angle θ _c = π / 4, π / 3. Consider estimation errors for a total of 10 types of two-dimensional images. Also consider σ _image = 1.0.

The inter-image displacement is given for each 0.1 [pixel] in the range of −0.5 ≦ d _s ≦ + 0.5 and −0.5 ≦ d _t ≦ + 0.5. At this time, the magnitude of the two-dimensional estimation error,
RMS error was evaluated.

The evaluated subpixel estimation method evaluated the estimation result which applied the subpixel estimation error reduction method with respect to the conventional one-dimensional estimation method other than each method mentioned above. Therefore,
Two-dimensional estimation (2D-EEC) according to Embodiment 1 of the present invention using the subpixel estimation error reduction method
Two-dimensional estimation (2D) according to Embodiment 1 of the present invention that does not use the subpixel estimation error reduction method
-One-dimensional estimation using subpixel estimation error reduction method (1D-EEC)
-One-dimensional estimation without subpixel estimation error reduction method (1D)
-Two-dimensional estimation using bilinear interpolation (Bi-Linear)
Two-dimensional estimation using bicubic interpolation (a = -0.55) (Bi-Cubic)
The estimation errors of the six estimation methods were evaluated. The evaluation results are shown in Table 1.

Table 1 shows the following.
First, the one-dimensional estimation (1D-EEC) using the subpixel estimation error reduction method can reduce the estimation error with respect to the conventional one-dimensional estimation (1D) when the rotation angle θ _g = 0. However, no difference appears when the rotation angle θ _g ≠ 0. The rotation angle θ _g corresponds to a rotation angle of similarity having anisotropy.

  Next, in the two-dimensional estimation using bilinear interpolation (Bi-Linear) or bicubic interpolation (Bi-Cubic), the estimation error is smaller by about two digits than the one-dimensional estimation (1D). This difference is obvious. Two-dimensional estimation using bilinear interpolation (Bi-Linear) corresponds to the density gradient method. This is why the density gradient method has been considered to be able to estimate subpixels with high accuracy.

  Then, the two-dimensional estimation using bicubic interpolation (Bi-Cubic) is more accurate than the two-dimensional estimation using bilinear interpolation (Bi-Linear). However, the results vary depending on the parameters of the weight function. Adjusting the parameters of the weight function corresponds to adjusting the characteristics of the interpolation filter. That is, by selecting an appropriate interpolation filter for an image, the accuracy of the subpixel estimation method using image interpolation can be improved.

  Finally, the two-dimensional estimation (2D-EEC) according to the first embodiment of the present invention using the subpixel estimation error reduction method has a smaller estimation error than the two-dimensional estimation using bicubic interpolation (Bi-Cubic). Moreover, there is no need for parameter adjustment. The two-dimensional estimation (2D-EEC) according to the first embodiment of the present invention using the subpixel estimation error reduction method is the two-dimensional estimation (2D) according to the first embodiment of the present invention that does not use the subpixel estimation error reduction method. The calculation cost is greater than. In view of this, when the calculation result is emphasized with a focus on the calculation cost, the two-dimensional estimation (2D) according to the first embodiment of the present invention which does not use the subpixel estimation error reduction method is used when the estimation result is similar to the density gradient method. When it is desired to obtain an accurate estimation result, it is possible to selectively use the two-dimensional estimation (2D-EEC) according to the first embodiment of the present invention using the subpixel estimation error reduction method.

(4) An estimation experiment by a computer program that implements the two-dimensional estimation method according to the first embodiment of the present invention. A two-dimensional subpixel estimation experiment was performed using a composite image and a real image by executing the program on a terminal (for example, a computer such as a personal computer or a workstation). The computer program referred to in the present invention may be abbreviated as a program for convenience.

(4-1) Estimation experiment using composite image The above-mentioned corner image was created and a two-dimensional subpixel estimation experiment was conducted. The created composite image is an 8-bit monochrome image with an image size of 64 × 64 [pixels] and an attention range of 11 × 11 [pixels]. Noise is not included in the image. The image parameters are σ _image = 1.0, corner angle θ _c = π / 4, and rotation angle θ _g = 0, π / 8, π / 4. The inter-image displacement was given for each 0.1 [pixel] in the range of −0.5 ≦ d _s ≦ + 0.5 and −0.5 ≦ d _t ≦ + 0.5.

  FIG. 20 shows the results of the conventional one-dimensional estimation and the two-dimensional estimation according to the first embodiment of the present invention to which the subpixel estimation error reduction method is applied. In FIG. 20, the calculation result is represented by a mesh, and the experiment result is represented by ●. The cause of the difference between the experimental result and the calculation result can be considered that the image has 8-bit quantization gradation (256 gradations). In addition, it can be considered that the part having a large step is that the result of searching for the maximum similarity in the integer has moved to the adjacent pixel position. Such movement to the adjacent position is a fundamental drawback of the conventional one-dimensional estimation method.

(4-2) Estimation Experiment Using Real Image When the rotation angle θ _g ≈0 of the two-dimensional similarity obtained from the image used for matching, an estimation error does not occur even with the conventional one-dimensional estimation method, but θ _g When ≠ 0, a large estimation error may occur. Therefore, a target tracking experiment using two types of patterns as matching patterns was performed.

  The matching pattern used is a circular pattern and an inclined edge pattern. These patterns are shown in FIGS. 21 (a) and 21 (e). The size of the circular pattern is about 41 [pixel] in diameter. In addition, the SSD self-similarity of these patterns is shown in FIGS.

  Since the self-similarity of the circular pattern is isotropic, it can be expected that the estimation error is small in both the conventional one-dimensional estimation method and the two-dimensional estimation method according to the first embodiment of the present invention. On the other hand, since the self-similarity corresponding to the inclined edge pattern is anisotropic and inclined, the estimation error is large in the conventional one-dimensional estimation method, and the estimation error is in the two-dimensional estimation method according to the first embodiment of the present invention. Can be expected to be small.

  The matching pattern was printed on cardboard to make a tracking target. The target moves linearly on the screw feed linear stage. The movement of the target was photographed with about 250 time-series images, the first image was used as a reference pattern, and the target was tracked from the subsequent images. The attention area size is 60 × 60 [pixels] (black rectangular areas in FIGS. 21A and 21E), and the search area is ± 8 (17 × 17) with respect to the movement predicted position. The target moves slowly in the lower right direction on the image. If matching is correct and there is no error in subpixel estimation, the measurement trajectory is a straight line. Both the conventional one-dimensional estimation method and the two-dimensional estimation method according to the first embodiment of the present invention use SSD self-similarity and perform subpixel estimation by parabolic fitting to which the subpixel estimation error reduction method is applied.

  FIGS. 21C and 21D show measurement results when a circular pattern is used as the tracking pattern. FIG. 21C shows the result of using the conventional one-dimensional estimation method, and FIG. 21D shows the result of using the two-dimensional estimation method according to Example 1 of the present invention. In both cases, the subpixel estimation error is small and the locus is a good straight line.

  FIGS. 21G and 21H show the measurement results when using an inclined edge pattern as the tracking pattern. FIG. 21G shows the result of using the conventional one-dimensional estimation method, and FIG. 21H shows the result of using the two-dimensional estimation method according to Example 1 of the present invention. In the conventional one-dimensional estimation method, it is clear that a very large subpixel estimation error occurs.

  In the two-dimensional estimation method according to the first embodiment of the present invention, the estimation error is small regardless of the tracking pattern, but the measurement trajectory deviates from the straight line as compared with the circular pattern. This is because the two-dimensional similarity is different from that of the two-dimensional Gaussian model, and therefore it is considered that an error is included when HEL and VEL are approximated by a straight line.

  In reality, it is rare to use such a tracking pattern. However, this experiment shows that an unexpectedly large subpixel estimation error may occur depending on the matching pattern when the conventional one-dimensional estimation method is used. In this experiment, it was found that the object was moving on a straight line, so the subpixel estimation error could be confirmed. However, for example, when reconstructing three-dimensional information of a building using a plurality of images taken from different directions, an inclined edge pattern region may be used as a matching region. The conventional one-dimensional estimation method has a problem that a large estimation error occurs in such a case.

  In Example 1 of the present invention, when obtaining the horizontal extreme value line HEL and the vertical extreme value line VEL, instead of obtaining a straight line approximating the three points with the least squares, it is obtained as a quadratic curve passing through the three points. Further, the two-dimensional estimation method according to the first embodiment of the present invention can be extended so as to reduce errors due to linear approximation.

Example 2 and Example 3:
Embodiments 2 and 3 of the multi-parameter high-accuracy simultaneous estimation method in image sub-pixel matching according to the present invention will be described as follows. The focus points of Example 2 and Example 3 of the present invention are as follows.

(A) Parallel movement and rotation, or parallel movement, rotation, and scale change are parameters corresponding to the image (in the present invention, they are also called parameters for short), and a multidimensional space (this book) with these parameters as axes. In the invention, similarity obtained at discrete positions in the parameter space or multidimensional similarity space) is used. Since the similarity may be calculated at discrete positions, there is no need for repeated calculation.

(B) The similarity obtained at discrete positions is considered to be a sampling of the similarity that is a continuous function in a multidimensional space (ie, parameter space).

(C) When the similarity is modeled under realistic assumptions, the degree of similarity obtained by partially differentiating the similarity with each axis in the multidimensional space (that is, each parameter axis in the parameter space) is set to zero. The plane passes through the maximum or minimum similarity. For example, when SAD or SSD is adopted as the similarity evaluation function, the hyperplane passes through the minimum value of similarity. When ZNCC is adopted as the similarity evaluation function, the hyperplane passes through the maximum value of similarity.

(D) Such hyperplanes can be obtained by the number of axes in the multidimensional space (that is, the number of parameters). Since the intersection coordinates of these hyperplanes coincide with the parameters that give the maximum or minimum similarity in multidimensional space, the maximum value of similarity in multidimensional space can be obtained by obtaining the intersection coordinates of these hyperplanes. Alternatively, the parameter that gives the minimum value can be obtained simultaneously.

  The multi-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image according to the present invention uses the evaluation values obtained discretely by the sampling grid to determine the evaluation value maximum position or minimum position of the high-precision sub-sampling grid accuracy. It can be considered as an estimation method.

  When the correspondence between images is limited to translation only (that is, in the case of Embodiment 1 of the present invention), this sampling grid corresponds to the pixels constituting the image. For this reason, the expression “sub-pixel” was appropriate as an expression of the amount of parallel movement with higher accuracy than the pixel unit. Therefore, as described above, the multi-parameter high-accuracy simultaneous estimation method in image sub-pixel matching according to the first embodiment of the present invention is also called a “two-dimensional sub-pixel estimation method”.

  However, in the present invention in which parallel movement and rotation, or parallel movement, rotation, and scale change are used as inter-image correspondence parameters, that is, the multi-parameter high-accuracy simultaneous estimation method according to the second and third embodiments of the present invention. Since expressions such as pixels and sub-pixels lack generality, the expressions “sampling grid” and “sub-sampling grid” will be used hereinafter.

<1> Three Parameter Simultaneous Estimation Method According to Second Embodiment <1-1> Principle of Three Parameter Simultaneous Estimation Method Here, the principle of the three parameter simultaneous estimation method according to the second embodiment of the present invention will be described. In the present invention, when searching for correspondence between images, a minimum or maximum similarity evaluation function is searched using some similarity evaluation function. The similarity obtained at this time varies depending on the similarity evaluation function to be used, but further varies depending on the image treated as an input.

Here, a specific input image and similarity evaluation function are not touched, and the similarity between images is approximated by the following three-dimensional Gaussian function, and this is used as a three-dimensional similarity model.
Where (s ₁ , s ₂ , s ₃ ) is a search parameter (that is, an inter-image correspondence parameter), σ is a standard deviation of a Gaussian function, and (k ₂ , k ₃ ) is an anisotropy coefficient (k ₂ , k ₃ > 0). Further, (t ₁ , t ₂ , t ₃ ) is a result of rotating and translating (s ₁ , s ₂ , s ₃ ) as follows.

However, (d ₁ , d ₂ , d ₃ ) is the correct value of the search parameter (s ₁ , s ₂ , s ₃ ), and (θ ₁ , θ ₂ , θ ₃ ) is (s ₁ , s ₂ , s), respectively. ₃ ) Consider the range of 0 ≦ θ ₁ , θ ₂ , θ ₃ ≦ π / 2 at the rotation angle around each axis.

If this three-dimensional similarity model is partially differentiated by s ₁ , s ₂ , s ₃ and set to 0, the following three planes can be obtained.

Π ₁ , Π ₂ , and 構成₃ constitute three planes in a three-dimensional space (that is, a parameter space with s ₁ , s ₂ , and s ₃ as axes). The intersection of these three planes is clearly (d ₁ , d ₂ , d ₃ ) and coincides with the correct parameter. The condition that the intersection of the three planes cannot be obtained is when the three-dimensional similarity model is degenerated to less than three dimensions. Under such conditions, k ₂ , k ₃ = 0 or ∞, which does not match the conditions of the three-dimensional similarity model.

<1-2> Implementation of 3-parameter Simultaneous Estimation Method A specific method for obtaining a sub-sampling grid estimated value in the 3-parameter simultaneous estimation method will be described below. For example, these three parameters can be considered as a parallel movement amount (d _s , d _t ) and a rotation angle θ.

  As shown in the previous section “Principle of the 3-parameter simultaneous estimation method”, the similarity value in the parameter space is differentiated with respect to each parameter axis and set to 0, so that the three planes that pass through the maximum value position or the minimum value position of each parameter. Can be made. Therefore, if the three planes can be estimated using the similarity values obtained in the sampling grid, the maximum position or the minimum position of the similarity values of each parameter can be obtained in the sub-sampling grid.

First, how to obtain the plane plane ₁ will be described. The point on the plane _{１ 1} is the point where the result of differentiating the similarity along the s ₁ axis is 0, so the subsampling position is estimated by using the similarity value on the straight line parallel to the s ₁ axis. it can be obtained several points on a plane [pi _1. FIG. 22 is a schematic diagram for explaining the plane _{１ 1} that passes through the maximum value with respect to the parameter axis s ₁ . FIG. 22 shows a sampling position similarity value is obtained in a square, sub-sampling the estimated position calculated on a straight line parallel to s ₁ axis by a black circle. The grid unit maximum value position of the similarity value is (r ₁ , r ₂ , r ₃ ).

_{(s 1, s 2, s} 3) If the similarity value of the _{R (s 1, s 2,} s 3), R (r 1 + i 1, r 2 + i 2, r 3 + i 3) (i 1, Nine points (p _{1 (r2 + i2, r3 + i3)} , r ₂ + i ₂ , r ₃ + i on the plane _{１ 1} using the similarity value of 27 points of i ₂ , i ₃ = -1,0,1) ₃ ) can be estimated. The following parabolic fitting can be used for the estimation.

  Further, as will be described later, it is possible to apply a technique for reducing an estimation error due to parabolic fitting.

When performing parabolic fitting of Formula 60, it is assumed that the similarity value of R (r ₁ , r ₂ + i ₂ , r ₃ + i ₃ ) is maximum, but this assumption may not be satisfied depending on the image pattern. is there. For example, when the shear deformation of the texture included in the image is extreme, R (r ₁ , r ₂ + i ₂ , r ₃ + i ₃ ) may not reach the maximum value at a position of i ₂ = i ₃ = ± 1. Also in this case, it is necessary to obtain the sub-sampling grid estimated position on a straight line that passes through (r ₂ + i ₂ , r ₃ + i ₃ ) and is parallel to the parameter axis s ₁ . In this case, the maximum value position of similarity values obtained on the sampling grid on this straight line is searched, and parabolic fitting of Equation 60 may be performed with respect to the corrected maximum value position r ₁ ′.

In the following description, in order to simplify the mathematical expression notation, along each parameter axis as s ₁ −r ₁ → s ₁ , s ₂ −r ₂ → s ₂ , s ₃ −r ₃ → s ₃ ( (r ₁ , r ₂ , r ₃ ) is translated and (r ₁ , r ₂ , r ₃ ) is considered as the origin position of the parameter axis. Such parallel movement is equivalent to moving the maximum similarity position in grid units to the origin, and does not lose the generality of the present invention. The sub-sampling grid position obtained by the multi-parameter high-accuracy simultaneous estimation method according to the present invention is finally added to (r ₁ , r ₂ , r ₃ ), so that the sub-sampling grid position before the parallel movement can be easily obtained. Can be calculated.

Subsampling grid position of the nine points that are expected to be on the plane [pi ₁ is a real image by the influence of variations between noise and image completely is likely not on the plane. For this reason, these nine points are applied to the plane ₁ by the least square method. Specifically, the sum of squares of distances along the s ₁ axis with nine points and the plane [pi ₁ may be determined the coefficient to minimize. The plane Π ₁ can be obtained as follows.

It becomes. Similarly, the planes Π ₂ and Π ₃ can be obtained as follows.

Intersection of planes _{１ 1} , Π ₂ , Π ₃
However, at the sub-sampling estimated position, it can be obtained as follows.

<2> N Parameter Simultaneous Estimation Method According to Third Embodiment <2-1> Principle of N Parameter Simultaneous Estimation Method Here, the principle of the N parameter simultaneous estimation method according to the third embodiment of the present invention will be described. Even when the number of parameters is greater than 3, the present invention can perform sub-sampling grid high-precision simultaneous estimation based on the same concept. When the similarity value is determined by N parameters from s ₁ to s _N , the similarity value is expressed as n multiplications of a Gaussian function as follows.

However,
Is an N-dimensional vector representing coordinates in the N-dimensional similarity space,
Is a matrix representing rotation and translation in N-dimensional space,
Is an anisotropy coefficient vector (N-1 dimension) of a Gaussian function.

At this time, the N-dimensional similarity model expressed by Equation 66 is partially differentiated by s _i (i = 1, 2,..., N) and set to 0, so that the following N N-dimensional relational expressions are obtained. Can be obtained.

  By solving these relational expressions simultaneously, a sub-sampling grid estimated value of N-dimensional similarity can be obtained. That is, the solution of these relational expressions is an N-dimensional similarity sub-sampling grid estimated value. If the N-dimensional similarity model deviates from approximation of the Gaussian function multiplication model of Equation 66, Equation 67 does not become a complete hyperplane, but at this time, it may be approximated to the hyperplane by the least square method or the like.

<2-2> Implementation of N Parameter Simultaneous Estimation Method Hereinafter, a specific method for actually performing N parameter sub-sampling grid simultaneous estimation will be described. First, when performing sub-sampling grid simultaneous estimation, the maximum position of similarity in the sampling grid is known as a precondition, and its (displacement) position is
And

On the N-dimensional hyperplane Π _i where the result of partial differentiation of the similarity value R with respect to the parameter axis s _i (i = 1, 2,..., N) is 0, the similarity on the straight line parallel to the s _i parameter axis Degree value
3 ^N−1 positions where the subsampling grid is estimated by parabolic fitting
Exists. However,
It is.

here,
Is an N-dimensional vector in which only the i-th element has the value k (k = -1, 0, +1) and all other elements have the value 0. For example,
It is. Also,
Is the j-th (j = 1, 2,..., 3 ^N-1 ) N-dimension in which the i-th element is 0 and each of the other elements takes a value of -1, 0, 1 Is a vector. For example,
It is.

These 3 ^N-1 positions,
Exists on the N-dimensional hyperplane _{ｉ i} or fits to the N-dimensional hyperplane Π _i in the least-square sense. Therefore, the N-dimensional hyperplane Π _i
3 ^N-1 positions,
Against
The relationship is obtained. This is rewritten as follows.

However,
Is a 3 ^N-1 row N matrix,
Is a vector with N elements. Coefficient representing the N-dimensional hyperplane Π _i
Can be obtained by the least square method using a generalized inverse matrix as follows.

As an intersection of N N-dimensional hyperplanes _i (i = 1, 2,..., N) obtained in this way, N parameter sub-sampling grid estimated values can be obtained. This intersection can be determined as follows.

Although it is possible to perform the simultaneous estimation of N parameter in the above way, actually the formula N-dimensional hyperplane [pi _i may be expressed in the following format including the same N unknowns.

  The reason is that the intersection of the hyperplanes can exist stably even if the N-dimensional similarity function is N-dimensionally symmetric.

Moreover, although the coefficient of the N-dimensional hyperplane was obtained by obtaining 3 ^N-1 positions on the ^N- dimensional hyperplane _{ｉ i} , this number can be reduced. For example, 3 ^N−1 = 9 when N = 3, but specifically when j = 3, for example.
(-1, -1,0) (0, -1,0) (+ 1, -1,0)
(-1, 0, 0) (0, 0, 0) (+1, 0, 0)
(-1, + 1,0) (0, + 1,0) (+ 1, + 1,0)

By adopting the following set from these, the number can be made 2 (N−1) +1.
(0, -1,0)
(-1, 0, 0) (0, 0, 0) (+1, 0, 0)
(0, + 1,0)

<3> Comparison between the experimental results of the present invention and the experimental results of the conventional method A time series image shown in FIG. 23 was artificially created. This time-series image is composed of a two-dimensional Gaussian function expressed by shading. The two-dimensional Gaussian function used in this experiment has directionality (anisotropy), and the directionality is set to be a direction other than 0 or 90 degrees. For this reason, the similarity calculated from this image is also anisotropic, and the rotation angle is other than 0 or 90 degrees.

  The image shown in FIG. 23 is the first frame of the time-series image, and the position of the subsequent frames does not change with respect to the first frame, and only the rotation angle changes by 0.1 degree for each frame. There are 31 time-series images in total, and the last frame is an image rotated three times with respect to the first frame. On the other hand, the position of the image does not move with respect to the rotation center.

  FIG. 24 shows the result of measuring the rotation angle of the subsequent frames with respect to the first frame. In the conventional method indicated by +, the image is divided into 6 × 6 = 36 regions, the parallel movement amount is measured for each small region by matching, and the rotation angle is estimated using the obtained parallel movement amount. It is a result. The method according to the present invention indicated by ◯ is the result of applying the present invention with the rotational angle measured by the conventional method as an initial estimated value and measuring the rotational angle with higher accuracy. The situation in which the rotation angle changes by 0.1 degree for each frame can be accurately measured.

  FIG. 25 shows the result of measuring the positions of subsequent frames with respect to the first frame. In the conventional method indicated by +, the image is divided into 6 × 6 = 36 regions, the parallel movement amount is measured for each small region by matching, and the position is estimated as an average of the obtained parallel movement amounts. It is. The method according to the present invention indicated by ◯ is the result of applying the present invention with the position measured by the conventional method as an initial estimated value and measuring the position with higher accuracy. The situation where the position does not change can be measured accurately.

<4> Other Embodiments of the Present Invention Another embodiment of the multi-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image according to the present invention will be described as follows.

<4-1> Combination with Subpixel Estimation Error Reduction Method Since parabolic fitting is performed using the similarity of three positions when the subsampling grid estimated position is obtained by Equation 60, the estimated value is unique to this estimated value. Includes estimation error. This estimation error can be canceled by using an interpolation image. Here, an example will be described in which the sub-pixel position relating to the horizontal displacement between images is estimated with two parameters.

When the two-dimensional image function used for matching is I ₁ (u, v), I ₂ (u, v), the horizontal subpixel position estimated by SSD and parabolic fitting can be expressed as follows.

On the other hand, the horizontal direction interpolation image _{I 1iu} (u, v) of the image _{I 1} (u, v) is
It can be expressed as. However, (u, v) at this time is an integer.

It can be considered that the horizontally interpolated image I _1iu (u, v) is displaced by ( _0.5,0 ) with respect to I ₁ (u, v). When subpixel estimation is performed with the SSD similarity shown below using this interpolated image, it can be considered that the estimation result is also displaced by (0.5, 0). Can be obtained by the following equation.

Equation 78 also includes an estimation error as in Equation 60, but since the phase is reversed, the estimation error can be canceled by the following equation.

<4-2> Affine parameter estimation When the correspondence between images can be expressed by affine transformation, it is necessary to determine six parameters in the following deformation model of 82.

However, (u1, v1) and (u2, v2) are coordinates representing the images I1 and I2, respectively, and the image I2 is deformed to match the image I1.

  Using the similarity values obtained at discrete positions, these six parameters can be simultaneously estimated with high accuracy.

<4-3> Projection Parameter Estimation When the correspondence between images can be expressed by projective transformation, it is necessary to determine eight parameters in the following deformation model of 83.

However, (u1, v1) and (u2, v2) are coordinates representing the images I1 and I2, respectively, and the image I2 is deformed to match the image I1.

  Using the similarity values obtained at discrete positions, these eight parameters can be estimated simultaneously with high accuracy.

<4-4> Combination with High-Speed Search Method The sub-sampling grid estimation method according to the present invention can be considered as a high-precision processing after the correspondence in the sampling grid unit is correctly obtained. Therefore, the various speed-up methods already proposed can be used for the search in the sampling grid unit.

  For example, as a high-speed search method using image density information, there are a hierarchical structure search method using an image pyramid, a method of setting an upper limit value for SSD or SAD, and stopping subsequent integration. Further, as a high-speed search method using image features, there is a method using a density histogram of a region of interest.

  By causing a computer program that implements the multi-parameter high-accuracy simultaneous estimation method in image sub-pixel matching according to the present invention to be executed by an information processing terminal equipped with a CPU (for example, a computer such as a personal computer or a workstation). The multi-parameters can be simultaneously estimated with high accuracy using the composite image or the real image.

  As described below, the present invention can be applied to various image processing fields. For example, in remote sensing and map creation using aerial photographs, multiple images are often connected, that is, mosaiced. At this time, region-based matching is performed to obtain corresponding positions between images. At this time, it is important not only that the correspondence is correct, but also that the correspondence is obtained with a resolution finer than the pixel unit. Furthermore, since the correspondence is obtained using a large number of images, the high speed of calculation is also important. When using region-based matching to obtain a response that exceeds the pixel resolution, sub-pixel estimation is performed using the fitting function using the similarity evaluation value obtained at discrete positions in units of pixels. At this time, the present invention is applied. By doing so, the increase in calculation time is slight, and it is possible to obtain a corresponding position with much higher accuracy.

  In particular, when the image pattern includes an oblique component, that is, when a pattern such as a road intersection is photographed while being tilted with respect to the image coordinate system in making a map or the like, the effect of the present invention is further exhibited.

  Further, in stereo vision, obtaining a corresponding position between a plurality of images is a basic process, and the corresponding position accuracy at this time greatly affects the final result accuracy. When obtaining the corresponding position between images, a constraint condition such as epipolar constraint is used. However, it is necessary to examine the correspondence between images when obtaining the epipolar constraint condition itself, and at this time, the epipolar constraint cannot be used. Therefore, by applying the present invention, it is possible to obtain the epipolar constraint condition with high accuracy, and as a result, high-accuracy stereo processing can be performed.

  Furthermore, when performing MPEG compression, it is necessary to accurately obtain the corresponding position between adjacent frames in order to achieve a high compression rate with high quality. Corresponding positions between adjacent frames are usually obtained in units of ½ pixel. At this time, a large error is included because a one-dimensional simple method using SAD or SSD similarity evaluation values is used. There is a possibility. However, the present invention can obtain a more reliable and highly accurate correspondence between images with a slight increase in the amount of calculation compared to the conventional one-dimensional estimation method, and therefore improves the compression rate of MPEG image generation. Is possible.

FIG. 1 is a diagram for explaining a conventional one-dimensional subpixel estimation method. FIG. 2 is a diagram showing an image example (a) generally used for a three-dimensional reconstruction application and its self-similarity (b). The self-similarity map indicates that the sub-pixel estimation error is generated by using the conventional one-dimensional estimation method because the self-similarity map has the properties of k ≠ 1 and θg ≠ 0. FIG. 3 is a diagram illustrating an example of a texture image having θ _c = π / 2 and its self-similarity. FIG. 4 is a diagram showing a two-dimensional image model and its two-dimensional similarity. (A) shows an image model of σ _image = 0.7, and the size of the image is 11 × 11 [pixel]. (B) shows the self-similarity of the image model, and the range of the similarity is 11 × 11. The dark position of (s, t) corresponds to the displacement of two images having high similarity. The darkest position (s, t) is (0, 0) due to the nature of self-similarity. FIG. 5 is a diagram illustrating a two-dimensional corner image model having a rotation angle θ _g and a corner angle θ _c of the entire image. FIG. 6 is a diagram showing a two-dimensional corner image model and its two-dimensional self-similarity. (a), (b), and (c) are image models of θ _c = π / 2 and σ _image = 1, and θ _g = 0, π / 8, and π / 4, respectively. (d), (e), and (f) are similarities corresponding to (a), (b), and (c), respectively. (g), (h), and (i) are image models of θ _c = π / 4 and σ _image = 1, and θ _g = 0, π / 8, and π / 4, respectively. (j), (k), and (l) are similarities corresponding to (g), (h), and (i), respectively. FIG. 7 is a diagram illustrating a two-dimensional repetitive texture image model having a rotation angle θ _{g of the} entire image, a horizontal spatial frequency f _u, and a vertical spatial frequency f _v . FIG. 8 is a diagram showing a two-dimensional repetitive texture image model and its two-dimensional self-similarity. (a), (b), and (c) are image models of f _u = 0.1 and f _v = 0.1, and θ _g = 0, π / 8, and π / 4, respectively. (d), (e), and (f) are similarities corresponding to (a), (b), and (c), respectively. (g), (h), and (i) are image models of f _u = 0.05 and f _v = 0.1, and θ _g = 0, π / 8, and π / 4, respectively. (j), (k), and (l) are similarities corresponding to (g), (h), and (i), respectively. The self-similarity is obtained from two similar images having a certain displacement. The size of the image is 11 × 11 [pixel], and the range of similarity is from −5 to 5. Figure 9 is a diagram showing a two-dimensional Gaussian image model having the rotation angle theta _g and the horizontal standard deviation σ and the vertical standard deviation kσ of the entire image. FIG. 10 is a diagram showing a two-dimensional Gaussian image model and its two-dimensional self-similarity. (a), (b), and (c) are image models with σ = 2 and k = 1, and θ _g = 0, π / 8, and π / 4, respectively. (d), (e), and (f) are similarities corresponding to (a), (b), and (c), respectively. (g), (h), and (i) are image models of σ = 2 and k = 0.5, and θ _g = 0, π / 8, and π / 4, respectively. (j), (k), and (l) are similarities corresponding to (g), (h), and (i), respectively. The self-similarity is obtained from two similar images having a certain displacement. The size of the image is 11 × 11 [pixel], and the range of similarity is from −5 to 5. FIG. 11 is a diagram showing three types of two-dimensional image models, their two-dimensional self-similarity, and a partial diagram of two-dimensional self-similarity compared to the two-dimensional similarity model. (a) shows a two-dimensional corner image model having θ _c = π / 4, σ _image = 1, and θ _g = 0. (b) shows a two-dimensional repetitive texture image model having f _u = 0.051, f _v = 0.1, and θ _g = 0. (c) shows a two-dimensional Gaussian image model having σ = 2, k = 0.5, and θ _g = 0. (d), (e), and (f) are similarities corresponding to (a), (b), and (c), respectively. (g), (h), (i) are (d), (e), (f) horizontal partial diagrams sampled at integer positions (represented by ◯), and t = 0 horizontal sub-views, and two-dimensional similarity models (Represented by a dotted line). FIG. 12 is a diagram illustrating an example of two-dimensional continuous similarity calculation and discrete similarity calculation. FIG. 13 is a diagram showing a two-dimensional similarity model having (d _s , d _t ) = (0.2,0.4), σ = 1, k = 0.7, and θ _g = π / 6. . The contour lines correspond to 10% to 90% of the function value. In (a), the long axis of the two-dimensional similarity model indicates θ _g = π / 6. In (b), both HEL (horizontal extreme value line) and VEL (vertical extreme value line) pass through the peak value position of the two-dimensional similarity model (d _s , d _t ). Different from axis or short axis. FIG. 14 is a diagram for explaining an estimation process of HEL and VEL. In (a), ● indicates the estimated sub-pixel position obtained from the similarity values (□) on the three horizontal lines t = −1, t = 0, t = 1, and HEL is three in the least squares. Can be estimated from the position of. In (b), VEL can be estimated in a similar process. FIG. 15 is a diagram illustrating positions of similarity values necessary for the conventional one-dimensional estimation method and the two-dimensional estimation method according to Embodiment 1 of the present invention. The two-dimensional similarity model has (d _s , d _t ) = (0.2,0.4), σ = 1, k = 0.3, and θ _g = π / 8. FIG. 16 is a diagram for comparing estimation errors. (a) is an image model used for comparison, that is, a two-dimensional corner model of θ _g = 0, θ _c = π / 4, and σ _image = 1.0. (b) shows the self-similarity of the image model of (a). FIG. 17 is a diagram for comparing estimation errors. (a) is an image model used for comparison, that is, a two-dimensional corner model of θ _g = π / 8, θ _c = π / 4, and σ _image = 1.0. (b) shows the self-similarity of the image model of (a). FIG. 18 is a diagram for comparing estimation errors. (a) is an image model used for comparison, that is, a two-dimensional corner model of θ _g = π / 4, θ _c = π / 4, and σ _image = 1.0. (b) shows the self-similarity of the image model of (a). FIG. 19 is a diagram for explaining an interpolation characteristic depending on a weight function and a parameter. FIG. 20 is a diagram illustrating an experimental result using a composite image. FIG. 21 is a diagram illustrating a result of target tracking experiment using an actual image. FIG. 22 is a schematic diagram for explaining the plane _{１ 1} that passes through the maximum value for the parameter s ₁ in the three-parameter simultaneous estimation method according to the second embodiment of the present invention. FIG. 23 is a composite image used in the experiment. FIG. 24 is a diagram showing the result of image rotation measurement. FIG. 25 is a diagram illustrating a result of image position measurement.

Claims (3)

In an N-dimensional similarity space with N (N is an integer of 4 or more) image correspondence parameters as axes, using the N-dimensional similarity values between images obtained at discrete positions, A multi-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image for simultaneously estimating N inter-image correspondence parameters with high accuracy,
An N-dimensional hyperplane that obtains a sub-sampling position where an N-dimensional similarity value between the images is maximum or minimum on a straight line parallel to a certain parameter axis, and that most closely approximates the obtained sub-sampling position. A step of seeking
Obtaining N N-dimensional hyperplanes for each parameter axis;
Obtaining an intersection of the N N-dimensional hyperplanes;
Setting the intersection as a sub-sampling grid estimated position of the inter-image correspondence parameter that gives the maximum or minimum value of the N-dimensional similarity in the N-dimensional similarity space;
A multi-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image characterized by comprising: In the three-dimensional similarity space around the three image correspondence parameters, the three image correspondence parameters are obtained with high accuracy by using the three-dimensional similarity values between images obtained at discrete positions. A three-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image to be simultaneously estimated,
Obtaining a sub-sampling position at which a three-dimensional similarity value between the images becomes maximum or minimum on a straight line parallel to a certain parameter axis, and obtaining a plane that most closely approximates the obtained sub-sampling position When,
Obtaining three planes for each parameter axis;
Obtaining an intersection of the three planes;
Setting the intersection as a sub-sampling grid estimated position of the inter-image correspondence parameter that gives a maximum value or a minimum value of the three-dimensional similarity in the three-dimensional similarity space;
A three-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image characterized by comprising: Two-parameter high-precision simultaneous in sub-pixel matching of images that accurately estimate the position of the maximum or minimum two-dimensional similarity in a continuous region using two-dimensional similarity values between discrete images An estimation method,
A sub-sampling position where the two-dimensional similarity value between the images is maximum or minimum on a straight line parallel to a horizontal axis is obtained, and a straight line (horizontal extreme value line HEL that most closely approximates the obtained sub-sampling position is obtained. )
A sub-sampling position where the two-dimensional similarity value between the images is maximum or minimum on a straight line parallel to the vertical axis is obtained, and a straight line (vertical extreme value line VEL that most closely approximates the obtained sub-sampling position is obtained. )
Obtaining an intersection of the horizontal extreme value line HEL and the vertical extreme value line VEL;
Making the intersection point a sub-pixel estimation position that gives the maximum or minimum value of the two-dimensional similarity;
A two-parameter high-accuracy simultaneous estimation method in sub-pixel matching of an image characterized by comprising: