WO2024033998A1 - Data processing method, measurement system, and program - Google Patents

Abstract

This data processing method is for analyzing scattered waves from waves emitted to an object and involves: emitting the waves to an object from a plurality of transmission points p₁(x'₁, y'₁, z'₁) arrayed on a y-axis; receiving scattered waves, having reflected at a reflection point (x, y, z) on the object at a reflection factor f(x, y, z), at a plurality of reception points p₂(x'₂, y'₂, z'₂) arrayed on the y-axis, as measurement values s(x'₁, x'₂, y'₁, y'₂, z'₁, z'₂, k); subjecting the measurement values s(x'₁, x'₂, y'₁, y'₂, z'₁, z'₂, k) to triple Fourier transformation to obtain S(k'_x1, k'_x2, k'_y1, k'_y2, z'₁, z'₂, k); performing a variable substitution process in the x-direction in a range of -k_x1≤k_x≤k_x1 (where k_x,nyq≤k_x1) when the Nyquist wavenumber in the x-direction determined by a measurement interval in the x-direction is k_x,nyq; performing a variable substitution process in the y-direction in ranges of -k_y1≤k'_y1≤k_y1 and -k_y2≤k'_y2≤k_y2 (where k'_y1,nyq≤k_y1 and k'_y2,nyq≤k_y2) when the Nyquist wavenumber in a y'₁-direction and the Nyquist wavenumber in a y'₂-direction, which are determined by a measurement interval in the y-direction, are k'_y1,nyq and k'_y2,nyq, respectively; performing triple inverse Fourier transformation; and obtaining a reflection factor f(x, y, z).

Description

Data processing method, measurement system, and program

The present invention relates to a data processing method, a measurement system, and a program for using a computer to process wave measurement data whose value is determined by the frequency of waves such as electromagnetic waves generated in space and the spatial coordinates of space.

Radar devices that non-destructively inspect the inside of non-metallic structures such as concrete and wood have been known. A conventional radar device has an array antenna in which a plurality of antennas are arranged on a plane. The array antenna has, for example, a configuration in which antennas such as planar antennas are lined up in one direction, and a transmitting array antenna and a receiving array antenna are arranged close to each other. Furthermore, in order to accurately measure the inside of a structure, the radar device measures the object to be measured using a wide band of frequencies while changing the frequency of electromagnetic waves at set frequency intervals.

Regarding a radar device having an array antenna, for example, a radar device in which a transmitting array antenna and a receiving array antenna composed of a plurality of planar antennas are formed on a common dielectric substrate is known (Japanese Patent Laid-Open No. 2015-095840). (hereinafter referred to as "Patent Document 1"). In a conventional radar device, the direction in which the planar antennas of the transmitting array antenna are arranged is parallel to the direction in which the planar antennas of the receiving array antenna are arranged. The position of the receiving array antenna in the arrangement direction of the planar antennas is between the two positions of adjacent planar antennas of the transmitting array antenna.

In addition, a scattering tomography method is known that can perform general and high-speed analysis of inverse problems and easily visualize information inside an object (Japanese Patent No. 6557747, hereinafter referred to as "Patent Document 2"). ).

Synthetic aperture processing is used to visualize the inside of a structure from measured data. There are two main types of synthetic aperture processing: addition methods such as the diffraction stacking method, and methods that utilize Fourier transform such as the FK migration method.
In order to realize a practical calculation time, synthetic aperture processing using Fourier transform is practical.

In order to accurately inspect the inside of a structure, radar equipment is desired to have high spatial resolution in measurement. In general, the spatial resolution of data obtained by radiation of waves having a frequency such as electromagnetic waves is determined by the fact that the distance between the structure to be measured and the measurement surfaces of the transmitting array antenna and the receiving array antenna is relatively close, In addition, when the measurement interval of the measurement data is small, it is determined by the wavelength of the center frequency of the wave. Here, the distance between the structure to be measured and the measurement surfaces of the transmitting array antenna and the receiving array antenna is, for example, one-fourth or less of the arrangement length of the array antenna. Moreover, the spatial resolution is the resolution within the plane in which each antenna of the array antenna is arranged.

For example, when the measurement interval of measurement data along the plane where each antenna of an array antenna is arranged is sufficiently small, the theoretical spatial resolution is determined by considering the round trip path of electromagnetic waves, and the frequency of the emitted waves has a frequency band. If the wave is swept at a frequency of 1/4, it will be one quarter of the wavelength of the wave at the center frequency of the frequency band. However, when the measurement interval of measurement data becomes large and exceeds one-fourth of the minimum wavelength of electromagnetic waves, the spatial resolution in actual measurement becomes larger than the ideal spatial resolution. In some cases, the spatial resolution in actual measurements is the measurement interval.

Since conventional radar equipment measures electromagnetic waves in a wide frequency band from low frequencies to high frequencies, the maximum wavelength of electromagnetic waves becomes long. For this reason, each antenna constituting the array antenna becomes large, and the length of the array antenna in the direction in which the antennas are arranged becomes long. As a result, the spacing between the antennas in the transmitting and receiving array antenna becomes longer, the measurement interval of the measured data tends to exceed one-fourth of the minimum wavelength of the electromagnetic waves emitted, and the spatial resolution is lower than the theoretical resolution. Aliasing components are likely to occur in measurement data. In order to reduce the spatial resolution to the ideal theoretical spatial resolution, that is, one quarter of the wavelength of the electromagnetic wave at the center frequency, we must increase the number of antennas in the transmitting/receiving array antenna and shorten the spacing between them. No. However, as described above, it is difficult to reduce the size of a broadband antenna, and therefore it is also difficult to shorten the arrangement interval.

An object of the present invention is to provide a data processing method, a measurement system, and a program that can improve the spatial resolution in measurement while maintaining a constant number of antennas.

The first aspect of the present invention is
A data processing method for analyzing scattered waves of waves emitted to an object,
radiating the wave to the object from a plurality of transmission points p ₁ (x' ₁ , y' ₁ , z' ₁ ) arranged on the y-axis;
The scattered waves reflected at the reflection point (x, y, z) on the object with a reflectance f(x, y, z) are received at a plurality of receiving points p ₂ (x' ₂ , receive the measured value s(x _{' 1 ,} x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) at y' ₂ , z' 2 );
The measured value s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) is subjected to triple Fourier transform using equation (1) to obtain S(k' _x1 , k ' _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k),

If the Nyquist wave number in the x direction determined by the measurement interval in the x direction is k _x,nyq , then variable substitution processing in the x direction is performed in the range -k _x1 ≦k _x ≦k _x1 (however, k _x,nyq ≦k _x1 ). conduct,
If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 ), perform variable substitution processing in the y direction,
The reflectance f(x, y, z) is determined by triple inverse Fourier transform using equation (2).

This is a data processing method.
however,
x' = x' ₁ = x' ₂
z' ₁ = z' ₂ = 0
k is the wave number of the propagating wave,
k' _x1 , k' _y1 , k' _z1 are the waves propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z). Components of the wave vector of a spherical wave,
k' _x2 , k' _y2 , k' _z2 are the waves propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ). Components of the wave vector of a spherical wave,
k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
It is.

The second aspect of the present invention is
A measurement system that analyzes scattered waves of waves emitted to an object,
A transmitting/receiving section,
a transmitter that radiates the wave to the object from a plurality of transmission points p ₁ (x' ₁ , y' ₁ , z' ₁ ) arranged on the y-axis;
The scattered waves reflected at the reflection point (x, y, z) on the object with a reflectance f(x, y, z) are received at a plurality of receiving points p ₂ (x' ₂ , a receiving unit that receives the measured value s(x _{' 1 ,} x' ₂ , y' ₁ , y' ₂ , z _{' 1 ,} z' ₂ , k) at y' 2 , z' 2 );
a transmitting/receiving unit having;
A processing device,
The measured value s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) is subjected to triple Fourier transform using equation (1) to obtain S(k' _x1 , k ' _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k);

If the Nyquist wave number in the x direction determined by the measurement interval in the x direction is k _x,nyq , then the variable substitution process in the x direction is performed in the range -k _x1 ≦k _x ≦k _x1 (k _x,nyq ≦k _x1 ). The steps to take and
If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 );
A step of calculating the reflectance f(x, y, z) by performing triple inverse Fourier transform from equation (2);

a processing device that executes;
It is a measurement system with
however,
x' = x' ₁ = x' ₂
z' ₁ = z' ₂ = 0
k is the wave number of the propagating wave,
k' _x1 , k' _y1 , k' _z1 are the waves propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z). Components of the wave vector of a spherical wave,
k' _x2 , k' _y2 , k' _z2 are the waves propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ). Components of the wave vector of a spherical wave,
k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
It is.

The third aspect of the present invention is
A program that analyzes scattered waves of waves emitted to an object,
_{The measured value} s _{( x ' 1} , _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k);

If the Nyquist wave number in the x direction determined by the measurement interval in the x direction is k _x,nyq , then the variable substitution process in the x direction is performed in the range -k _x1 ≦k _x ≦k _x1 (k _x,nyq ≦k _x1 ). The steps to take and
If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 );
Steps to calculate the reflectance f(x, y, z) by performing triple inverse Fourier transform from equation (2);

This is a program that causes a computer to execute.
however,
x' = x' ₁ = x' ₂
z' ₁ = z' ₂ = 0
k is the wave number of the propagating wave,
k' _x1 , k' _y1 , k' _z1 are the waves propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z). Components of the wave vector of a spherical wave,
k' _x2 , k' _y2 , k' _z2 are the waves propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ). Components of the wave vector of a spherical wave,
k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
It is.

The fourth aspect of the present invention is
A data processing method for analyzing scattered waves of waves emitted to an object,
radiating the wave to the object from a plurality of transmission points p ₁ (x' ₁ , y' ₁ , z' ₁ ) arranged two-dimensionally on the xy plane,
The scattered waves reflected at the reflection point (x, y, z) on the object with a reflectance f(x, y, z) are transmitted to a plurality of receiving points p ₂ (x ' ₂ , y' ₂ , z' ₂ ) and receive the measured value s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k),
_{The measured value} s _{( x ' 1} , ' _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k),

Let k' _x1,nyq be the Nyquist wave number in the x' ₁ direction determined by the measurement interval in the x direction, and k' _{x2, nyq} be the Nyquist wave number in the x' ₂ direction, then -k _x1 ≦k' _x1 ≦k _x1 , and , -k _x2 ≦k' _x2 ≦k _x2 (however, k' _x1,nyq ≦k _x1 and k' _x2,nyq ≦k _x2 ) Perform variable substitution processing in the x direction,
If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 ), perform variable substitution processing in the y direction,
The reflectance f(x, y, z) is determined by triple inverse Fourier transform using equation (2).

This is a data processing method.
however,
z' ₁ = z' ₂ = 0
k is the wave number of the propagating wave,
k' _x1 , k' _y1 , k' _z1 are the waves propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z). Components of the wave vector of a spherical wave,
k' _x2 , k' _y2 , k' _z2 are the waves propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ). Components of the wave vector of a spherical wave,
k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
It is.

The fifth aspect of the present invention is
A measurement system that analyzes scattered waves of waves emitted to an object,
A transmitting/receiving section,
a transmitter that radiates the wave to the object from a plurality of transmission points p ₁ (x' ₁ , y' ₁ , z' ₁ ) arranged two-dimensionally on the xy plane;
The scattered waves reflected at the reflection point (x, y, z) on the object with a reflectance f(x, y, z) are transmitted to a plurality of receiving points p ₂ (x a receiving unit that receives a measured value s(x' ₁ , x' ₂ , y' ₁ , y' _{2 ,} z _{' 1 ,} z' ₂ , k) at ' 2 , y' 2 , z' 2 );
a transmitting/receiving unit having;
A processing device,
_{The measured value} s _{( x ' 1} , ' _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k);

Let k' _x1,nyq be the Nyquist wave number in the x' ₁ direction determined by the measurement interval in the x direction, and k' _{x2, nyq} be the Nyquist wave number in the x' ₂ direction, then -k _x1 ≦k' _x1 ≦k _x1 , and , -k _x2 ≦k' _x2 ≦k _x2 (however, k' _x1,nyq ≦k _x1 and k' _x2,nyq ≦k _x2 );
If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 );
A step of calculating the reflectance f(x, y, z) by performing triple inverse Fourier transform from equation (2);

a processing device that executes;
It is a measurement system with
however,
z' ₁ = z' ₂ = 0
k is the wave number of the propagating wave,
k' _x1 , k' _y1 , k' _z1 are the waves propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z). Components of the wave vector of a spherical wave,
k' _x2 , k' _y2 , k' _z2 are the waves propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ). Components of the wave vector of a spherical wave,
k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
It is.

The sixth aspect of the present invention is
A program that analyzes scattered waves of waves emitted to an object,
_{The measured value} s _{( x ' 1} , _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k);

Let k' _x1,nyq be the Nyquist wave number in the x' ₁ direction determined by the measurement interval in the x direction, and k' _{x2, nyq} be the Nyquist wave number in the x' ₂ direction, then -k _x1 ≦k' _x1 ≦k _x1 , and , -k _x2 ≦k' _x2 ≦k _x2 (however, k' _x1,nyq ≦k _x1 and k' _x2,nyq ≦k _x2 );
If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 );
Steps to calculate the reflectance f(x, y, z) by performing triple inverse Fourier transform from equation (2);

This is a program that causes a computer to execute.
however,
z' ₁ = z' ₂ = 0
k is the wave number of the propagating wave,
k' _x1 , k' _y1 , k' _z1 are the spherical waves of the wave propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z) The components of the wave vector of
k' _x2 , k' _y2 , k' _z2 is the spherical surface of the wave propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ) components of the wave vector of the wave,
k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
It is.

According to the data processing method, measurement system, and program of the present invention, it is possible to improve the spatial resolution in measurement while keeping the number of arranged antennas constant.

A diagram showing the configuration of a radar device according to the first embodiment Diagram showing the configuration of the array antenna shown in Fig. 1 A diagram explaining the positional relationship between the array antenna and the measurement target of the first embodiment Flowchart showing the data processing method of the first embodiment Conceptual diagram of aliasing An example of the measurement waveform of the first embodiment Fourier transform regarding the x-axis in the first embodiment Fourier transform regarding the y-axis in the first embodiment A diagram explaining the positional relationship between the array antenna and the measurement target according to the second embodiment Flowchart showing the data processing method of the second embodiment Fourier transform regarding x' ₁ axis and x' ₂ axis in second embodiment Fourier transform regarding y' ₁ axis and y' ₂ axis in second embodiment A diagram explaining the positional relationship between the array antenna and the measurement target according to the third embodiment Flowchart showing the data processing method of the third embodiment Fourier transform regarding the x-axis in the third embodiment Fourier transform regarding y' ₁ axis and y' ₂ axes in the third embodiment xz plane image of point target 1 simulated using the data processing method of comparative example 1 xz plane image of point target 1 simulated using the data processing method of Example 1 YZ plane image of point target 1 simulated using the data processing method of Comparative Example 1 YZ plane image of point target 1 simulated using the data processing method of Example 1 Waveforms in the x-axis direction of point target 1 in Comparative Example 1 and Example 1 Waveforms in the y-axis direction of point target 1 in Comparative Example 1 and Example 1 Waveforms in the z-axis direction of point target 1 in Comparative Example 1 and Example 1 xz plane image of point target 2 simulated using the data processing method of comparative example 1 xz plane image of point target 2 simulated using the data processing method of Example 1 YZ plane image of point target 2 simulated using the data processing method of Comparative Example 1 YZ plane image of point target 2 simulated using the data processing method of Example 1 Waveforms in the x-axis direction of point target 2 in Comparative Example 1 and Example 1 Waveforms in the y-axis direction of point target 2 in Comparative Example 1 and Example 1 Waveforms in the z-axis direction of point target 2 in Comparative Example 1 and Example 1 Conceptual diagram of cross-shaped k' _y1 -k' _y2 filter Simulation results for four point targets Point target recovery rate Transmit and receive antenna layout 3D image of the specimen without super-resolution processing 3D image of the specimen subjected to super-resolution processing Three-dimensional image of the specimen subjected to super-resolution processing by applying a cross-shaped k' _y1 - k' _y2 filter Waveforms in the x-axis direction of point target 1 in Comparative Example 2 and Example 2 Waveforms in the y-axis direction of point target 1 in Comparative Example 2 and Example 2 Waveforms in the z-axis direction of point target 1 in Comparative Example 2 and Example 2 Waveforms in the x-axis direction of point target 2 in Comparative Example 2 and Example 2 Waveforms in the y-axis direction of point target 2 in Comparative Example 2 and Example 2 Waveforms in the z-axis direction of point target 2 in Comparative Example 2 and Example 2 Conceptual diagram of cross-shaped k' _x1 - k' _x2 filter

<First embodiment>
The data processing method, measurement system, and program of the first embodiment will be described in detail below. FIG. 1 shows the configuration of a radar device according to this embodiment. FIG. 2 shows the configuration of the array antenna shown in FIG. FIG. 3 is a diagram illustrating the positional relationship between the array antenna of this embodiment and the object to be measured. In this embodiment, electromagnetic waves will be described as waves that are radiated into space, but instead of electromagnetic waves, waves that propagate in space, such as X-rays and ultrasonic waves, may be used.

The measurement system 1 of this embodiment includes a transmitting/receiving section and a processing device. The processing device may be provided integrally with the transmitting/receiving section, or may be provided at a separate location connected to the transmitting/receiving section via a network. In the following embodiments, an example will be described in which a processing device is provided integrally with a transmitting/receiving section.

The radar device 60 of this embodiment shown in FIG. 1 radiates electromagnetic waves from the transmitting antenna while sweeping the frequency of the electromagnetic waves using a transmitting array antenna and a receiving array antenna (transmitting/receiving section). Then, the radar device 60 receives the reflected wave from the object to be measured using a receiving antenna, and obtains measurement data s(x', y', z', k). The measurement data s(x', y', z', k) is data whose variables are an x-coordinate component, a y-coordinate component, a z-coordinate component, and the frequency of electromagnetic waves.

The radar device 60 includes a measurement unit 61, a data processing unit (processing device) 66, and an image display unit 68. The measurement unit 61 includes a transmitting array antenna 50, a receiving array antenna 52, high frequency switches 58 and 59, a high frequency circuit 62, and a system control circuit 64. The radar device 60 emits electromagnetic waves of 10 MHz or more, for example, 10 to 20 GHz, but the frequency of the electromagnetic waves is not particularly limited.

As shown in FIG. 2, the transmitting array antenna 50 has a plurality of transmitting antennas 10a arranged in one direction. Each transmitting antenna 10a radiates electromagnetic waves toward the object to be measured. The reception array antenna 52 has a plurality of reception antennas 10b arranged along the arrangement direction of the transmission antennas 10a. Each receiving antenna 10b receives electromagnetic waves reflected from the object to be measured.
The transmitting antenna 10a of the transmitting array antenna 50 and the receiving antenna 10b of the receiving array antenna 52 are arranged on one plane. A transmitting array antenna 50 and a receiving array antenna 52 are arranged so that the object to be measured faces this plane.

The data processing unit 66 processes a plurality of measurement data obtained by transmission toward the object to be measured by the plurality of transmitting antennas 10a and reception by the plurality of receiving antennas 10b, and calculates image data regarding the object to be measured. Although the transmitting antenna 10a and the receiving antenna 10b of this embodiment are planar antennas in which an antenna pattern is formed in a plane on a substrate, they are not limited to planar antennas.

The transmitting array antenna 50 and the receiving array antenna 52 move parallel to the surface of the object to be measured. That is, the transmitting array antenna 50 and the receiving array antenna 52 measure while scanning along the surface of the object to be measured. When the transmitting array antenna 50 and the receiving array antenna 52 move, the system control circuit 64 controls the operation of the high frequency circuit 62. Specifically, the system control circuit 64 radiates electromagnetic waves while switching the transmitting antenna 10a using the high frequency switch 58 for each unit length of the moving distance of the transmitting array antenna 50 and the receiving array antenna 52. Controls the operation of the high frequency circuit 62.

Radar device 60 has an encoder 69. Encoder 69 generates a pulse signal every fixed moving distance. The encoder 69 senses movement of the transmitting array antenna 50 and the receiving array antenna 52.
At this time, each time an electromagnetic wave is radiated from each transmitting antenna 10a, the high frequency switch 59 sequentially switches the plurality of receiving antennas 10b to cause each receiving antenna 10b to receive the electromagnetic wave.

Note that the frequency of the electromagnetic waves radiated from the transmitting array antenna 50 is swept at set frequency intervals within a range of 10 to 20 GHz at a certain time, and the electromagnetic waves are radiated. Therefore, the measurement data obtained from the high frequency circuit 62 is data whose value is determined by the position transmitted by the transmitting antenna 10a, the position received by the receiving antenna 10b, the frequency, and the position of the target.
At this time, the high frequency switch 59 is operated so that the receiving antenna 10b closest to the transmitting antenna 10a that radiated the electromagnetic wave receives the reflected wave of the electromagnetic wave when the electromagnetic wave radiated from the transmitting antenna 10a is reflected by the object to be measured. is controlled. The receiving microwave amplifier (RF amplifier) may be set to change the gain for each pair of the transmitting antenna 10a for transmitting and the receiving antenna 10b for receiving. At this time, the high frequency circuit 62 has a variable gain amplification function that switches the gain depending on the selection of the pair of transmitting antenna 10a and receiving antenna 10b. As a result, the depth at which defects and the like in the object to be measured can be inspected can be increased.

In this embodiment, the arrangement direction of the transmitting antenna 10a and the receiving antenna 10b is parallel, and as shown in FIG. 2, the arrangement direction is the y direction. On the other hand, the moving direction (scanning direction) of the transmitting array antenna 50 and the receiving array antenna 52 is assumed to be the x direction. When viewed from the transmitting array antenna 50 and the receiving array antenna 52, the direction in which the object to be measured is located (the direction in which electromagnetic waves are transmitted) is defined as the z direction.

Note that the moving direction (scanning direction) of the transmitting array antenna 50 and the receiving array antenna 52 may be the y direction. That is, it may be moved (scanned) in the same direction as the arrangement direction of the transmitting antenna 10a and the receiving antenna 10b.
Further, the transmitting array antenna 50 may have only one transmitting antenna 10a, and the receiving array antenna 52 may have a plurality of receiving antennas 10b. In this case as well, the moving direction (scanning direction) of the transmitting array antenna 50 and the receiving array antenna 52 may be set to the y direction. That is, it may be moved (scanned) in the same direction as the arrangement direction of the receiving antennas 10b.

The data processing unit 66 processes the measurement data s(x', y', z', k) obtained by transmitting and receiving electromagnetic waves by the transmitting array antenna 50 and the receiving array antenna 52 to understand the inside of the measurement target. Create image data to represent. The data processing unit 66 is configured by a computer, for example, and starts by calling a program stored in the storage section 66a. This allows the data processing unit 66 to perform its functions. That is, the data processing unit 66 is composed of software modules. The image display unit 68 displays an image of the inside of the object to be measured using the created image data.

FIG. 2 schematically shows a transmitting array antenna 50 and a receiving array antenna 52. The positions of the transmitting antenna 10a and the receiving antenna 10b in the x direction are shifted by ΔL, but in the following explanation, the positions of the transmitting antenna 10a and the receiving antenna 10b in the x direction are shifted from each other by ΔL between the transmitting antenna 10a and the receiving antenna 10b. Do what is at the middle circle point. This circled point is called the transmitting/receiving point.
Note that there may be cases where there is no deviation in the y direction between the transmitting antenna 10a and the receiving antenna 10b. That is, there are cases where Δy=0. Furthermore, the transmitting antenna 10a and the receiving antenna 10b may be shared. That is, Δy=0 and ΔL=0 in some cases.

Therefore, the positional relationship between the measurement object, the transmitting array antenna 50, and the receiving array antenna 52 can be expressed as shown in FIG.
Here, the coordinates of the transmitting/receiving point are p(x', y', z'). Let f(x, y, z) be the reflectance at the reflection point (x, y, z) of the object to be measured. Let s (x', y', z', k) be the measurement data at the transmission/reception point p (x', y', z'). Let the propagation wavelength of electromagnetic waves in vacuum be λ ₀ . Let the relative dielectric constant of the medium be ε _r . Let k be the wave number of the propagating electromagnetic wave.

At this time, the measurement data s(x', y', z', k) at the transmitting/receiving point p(x', y', z') can be expressed by the following formula.

however,

It is.

In equation (1-1), the electromagnetic wave is expressed as a spherical wave, and distance attenuation is omitted. This distance attenuation has been omitted because it has little effect on subsequent processing. When the exponent part of the integrand of the second stage equation in equation (1-1) is expressed in Fourier transform notation, it becomes the following equation. This is equivalent to decomposing the reciprocating spherical wave in equation (1-1) into a three-dimensional plane wave.

Here, (k _x , k _y , k _z ) is the wave number of the round trip spherical wave propagating between the transmitting/receiving point p (x', y', z') and the reflecting point (x, y, z). It is a component of a vector. however,

satisfy.

Hereinafter, the reflectance f(x, y, z) is derived from the measurement data s(x', y', z', k) based on equation (1-3). First, equation (1-3) is rearranged as follows.

Here, the integral inside { } is a triple Fourier transform with respect to (x, y, z). Also, the integral inside [ ] is a double inverse Fourier transform regarding (k _x , k _y ). Therefore, double Fourier transform is performed on both sides of equation (1-5) regarding (x', y'). Let F(k x , k y , k z ) be the function f( _x , _y , _z ) after triple Fourier transformation. Let S(k x , k y , z', k) be a function after double Fourier transform of measurement data s ( _x ', y _' , z', k). At this time, equation (1-5) is expressed by the following equation.

When both sides of the second line of equation (1-6) are subjected to triple inverse Fourier transform for (k _x , k _y , k _z ), the reflectance f(x, y, z) is obtained as follows. .

Here, in this embodiment, as shown in FIG. 3, the transmitting/receiving point p (x', y', z') is arranged on the xy plane, and z'=0, so equation (1-7) is It can be expressed as follows.

As described above, the data processing unit 66 determines the reflectance f(x, y, z) based on the measurement data s(x', y', z', k).

The data processing method and program of this embodiment will be described below with reference to FIG. FIG. 4 is a flowchart showing the data processing method of this embodiment.
First, the measurement unit 61 acquires measurement data s(x', y', 0, k) (step S1-1). Then, the data processing unit 66 performs Hilbert transformation on the measurement data s(x', y', 0, k) (step S1-2). Thereby, the imaginary component of the frequency data at each transmitting/receiving point is obtained.

Next, the data processing unit 66 performs a double Fourier transform on (x', y') on the measurement data s(x', y', 0, k) (step S1-3). As a result, S(k _x , k _y , 0, k) is obtained as shown in equation (1-6).
Next, the data processing unit 66 performs variable substitution on S(k _x , k _y , 0, k) (step S1-4). Specifically, using equation (1-4), the function of (k _x , k _y , k) is made into the function of (k _x , k _y , k _z ). This yields S(k _x , k _y , k _z ).
Next, the data processing unit 66 performs triple inverse Fourier transform on S(k _x , k _y , k _z ) and (k _x , k _y , k _z ) (step S1-5). . As a result, the reflectance f(x, y, z) is obtained as shown in equation (1-8).

The storage unit 66a stores a program for executing the data processing method of this embodiment. The program stored in the storage section 66a causes the data processing unit 66 to execute the data processing method of this embodiment.

Here, alias data, which is considered as aliasing noise in Fourier transform processing, will be explained. When measuring with an array antenna, aliasing will occur in the measured data unless the antenna size and scanning interval are below the Nyquist standard. Aliasing refers to the fact that continuous signals with different frequency components become indistinguishable after sampling.

FIG. 5 shows a conceptual diagram of aliasing. The Nyquist wave number determined from the spatial sampling interval is defined as k _nyq . The Nyquist wave number k _nyq is expressed by the following equation (1-10).

Here, Δx and Δy are measurement intervals in the x-axis direction and y-axis direction, respectively.
Signals exceeding the Nyquist wave number k _nyq are called aliases and are treated as aliasing noise during signal processing.

FIG. 6 shows an example of a measured waveform in the first embodiment. Here, the frequency of the electromagnetic wave is 15 GHz, the dielectric constant is 1, the measurement interval Δy = 7.5 mm, the measurement width is 562.5 mm, and the depth of the point target is 120 mm directly below the origin. Here, the wave number ky of the waveform measured at a distance y from the origin is expressed by the following equation (1-11).

Here, θ _ty is the angle between the line connecting the outermost transmitting/receiving antenna and the target and the z-axis. Further, the Nyquist wavenumber is expressed by the following equation (1-12).

From this, the Nyquist criterion is expressed by the following equation (1-13).

The y-coordinate waveform that does not satisfy this condition, that is, the data in the aliasing area in FIG. 6, all becomes alias data (aliasing noise).

From the above equation, it can be seen that the alias data increases as the wavelength of the electromagnetic waves becomes shorter, the measurement interval becomes wider, or the target gets closer to the measurement surface. Usually, alias data is also processed as aliasing noise during signal processing. Therefore, an increase in alias data causes a deterioration in image resolution, a decrease in image intensity of a shallow target, and a deterioration in the S/N ratio.

Hereinafter, the measurement interval Δx = 10 mm, the measurement interval Δy = 19.25 mm, and the wavelength of the highest frequency in the medium λ _fmax = 21.24 mm. At this time, the following equation (1-14) is obtained.

At this time, the Nyquist criterion is not satisfied, and aliasing occurs in both the x-axis direction and the y-axis direction. The Nyquist wavenumber in this case is expressed by the following equation (1-15).

The maximum wave number determined from the maximum frequency is expressed by the following equation (1-16).

Here, the Fourier transform regarding the x-axis in step S1-3 is shown in FIG. Furthermore, the Fourier transform regarding the y-axis in step S1-3 is shown in FIG.
From FIG. 7, in the Fourier transform regarding the x-axis, the wave number with the upper limit of ±k _{x, nyq} is output to range A. Data with a wave number exceeding ±k _x, nyq appears as alias data with a period of 2k _{x, nyq} . Therefore, in the variable replacement process in step S1-4, the variable replacement process is performed not in the range A but in the range A' having an upper limit of ±k _xmax . As a result, the alias data is not discarded as aliasing noise, but is actually added to the variable replacement process as meaningful data, improving the resolution in the x-axis direction and further increasing the image intensity of shallow targets. becomes possible.
Regarding the y-axis, variable replacement processing is similarly performed in range B' instead of range B, as shown in FIG.

<Second embodiment>
The data processing method, measurement system, and program of the second embodiment will be described in detail below. The transmitting array antenna 50 and the receiving array antenna 52 of the first embodiment are arranged in one direction (the y direction in FIG. 3), but in this embodiment, the transmitting array antenna 50 and the receiving array antenna 52 are arranged in one direction (the y direction in FIG. 3). The arrangement is different. The transmitting array antenna 50 and the receiving array antenna 52 of this embodiment are arranged in a plane.
Furthermore, in the first embodiment, the coordinates of the transmitting point and the receiving point are both p(x', y', z'), but in this embodiment, the coordinates of the transmitting point and the receiving point are different. _{In this embodiment ,} as _{shown in FIG .} Arranged on a plane.

Here, the reflectance at the reflection point (x, y, z) of the object to be measured is assumed to be f(x, y, z). The measurement data at the transmitting point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ) are expressed as s(x' ₁ , x' ₂ , y ' ₁ , y' ₂ , z' ₁ , z' ₂ , k). Let the propagation wavelength of electromagnetic waves in vacuum be λ ₀ . Let the relative dielectric constant of the medium be ε _r . Let k be the wave number of the propagating electromagnetic wave.

At this time, the measurement data s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) can be expressed by the following formula.

however,

It is.

In formula (2-1), the electromagnetic wave is expressed as a spherical wave, and distance attenuation is omitted. This distance attenuation has been omitted because it has little effect on subsequent processing. When the exponent part of the integrand in equation (2-1) is expressed in Fourier transform notation, it becomes the following equation. This is equivalent to decomposing the spherical wave in equation (2-1) into a three-dimensional plane wave.

Here, (k' _x1 , k' _y1 , k' _z1 ) is a component of the wave number vector of the spherical wave propagating from the transmission point to the reflection point. Furthermore, (k' _x2 , k' _y2 , k' _z2 ) is a component of the wave number vector of the spherical wave propagating between the reflection point and the reception point. however,

satisfy.

Below _, based on _equation ( _2-3 ), _the reflectance f ₍ x, y, _z ) is derived. First, quadruple Fourier transform is performed on both sides of equation (2-3) with respect to x' ₁ , x' ₂ , y' ₁ , and y' ₂ .

Rewrite the left side of equation (2-5) as shown in equation (2-6) below.

Then, equation (2-5) is expressed as equation (2-7).

Perform the following integration on both sides of equation (2-7).

Here, perform the following variable substitutions.

Here, the following equation can be obtained from equation (2-9) and equation (2-10).

The absolute value |J| of the Jacobian in this variable substitution is given by the following equations from equations (2-12) and (2-13), respectively.

By substituting formula (2-9), formula (2-10), formula (2-14), and formula (2-15) into formula (2-8) and performing variable substitution, the following formula is obtained. .

Here, the integral with respect to (u, v) on the right side of the second line of equation (2-16) is a constant, so it is omitted. When a triple inverse Fourier transform is performed on both sides of equation (2-16) for (k _x , k _y , k _z ), the reflectance f(x, y, z) is obtained as follows.

In order to solve equation (2-17), (k' _z1 , k' _z2 , k) is transformed into (k' _x1 , k' _x2 , k' _y1 , k' _y2 , k _z ) or (k _x , u, k _y , v, k _z ).
Organize using equations (2-4) and (2-11), and if both the transmitting point and receiving point are located on the xy plane passing through the origin, z' ₁ = z' ₂ = 0. Therefore, equation (2-17) can be expressed as follows.

As described _above , _the data processing unit ₆₆ calculates _the reflectance _{f (} x' , y, z).

The data processing method and program of this embodiment will be described below with reference to FIG. FIG. 10 is a flowchart showing the data processing method of this embodiment.
First, the measurement unit 61 acquires measurement data s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) (step S2-1). Then, the data processing unit 66 performs Hilbert transformation on the measurement data s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) (step S2- 2). Thereby, the imaginary component of the frequency data at each measurement point is obtained.

Next _{, the data} processing unit 66 calculates _{(x' 1 , x ' 2} , y' ₁ , y' ₂ ) is performed (step S2-3). As a result, S(k' _x1 , k' _x2 , k' _y1 , k' _y2 , 0, 0, k) is obtained as shown in equation (2-6).
Next, the data processing unit 66 performs variable substitution on S(k' _x1 , k' _x2 , k' _y1 , k' _y2 , 0, 0, k) (step S2-4). Specifically, (k' _x1 , k' _x2 , k' _y1 , k ' Make the function of _y2 , k) a function of (k _x , k _y , k _z ). This yields S(k _x , k _y , k _z ).
Next, the data processing unit 66 performs triple inverse Fourier transform on S(k _x , k _y , k _z ) and (k _x , k _y , k _z ) (step S2-5). . As a result, the reflectance f(x, y, z) is obtained as shown in equation (2-18).

The storage unit 66a stores a program for executing the data processing method of this embodiment. The program stored in the storage section 66a causes the data processing unit 66 to execute the data processing method of this embodiment.

Here, the frequency of the electromagnetic wave is 4.46484 GHz, the relative permittivity is 10, the measurement interval Δx' ₁ = 38.5 mm, the measurement interval Δx' ₂ = 38.5 mm, the measurement interval Δy' ₁ = 38.5 mm, the measurement interval Δy ' ₂ = 38.5 mm, the measurement width is 616 mm, and the wavelength of the highest frequency in the medium λ _fmax = 21.24 mm. At this time, the following equation (2-20) is obtained.

At this time, the Nyquist criterion is not satisfied, and aliasing occurs in each of the x' _1- axis direction, x' _2- axis direction, y' _1- axis direction, and y' _2- axis direction. The Nyquist wave number in this case is expressed by the following equation (2-21).

The maximum wave number determined from the highest frequency is expressed by the following equation (2-22).

Here, the Fourier transform regarding the x' ₁ axis and the x' ₂ axes in step S2-3 is shown in FIG. Further, the Fourier transform regarding the y' _1- axis and the y' _2- axis in step S2-3 is shown in FIG.
From FIG. 11, in the Fourier transform regarding the x' _1- axis and the x' _2- axis, wave numbers whose upper limits are ±k' _{x1, nyq} and ±k' _{x2, nyq} are output to range A, respectively. Data with wave numbers exceeding ±k' _{x1, nyq} and ±k' _x2, nyq appear as alias data in cycles of 2k' _{x1, nyq} and 2k' _{x2, nyq,} respectively. Therefore, in the variable replacement process in step S2-4, the variable replacement process is performed not in the range A but in the range A' whose upper limits are ±k' _x1max and ±k' _x2max , respectively. As a result, instead of discarding alias data as aliasing noise, it is added to the variable replacement process as actually meaningful data, improving the resolution in the x' _1- axis and x' _2- axis directions, and further improving the resolution of shallow targets. It becomes possible to further increase the image intensity.
Similarly, for the _y'1 axis _{and the y'2} axis, as shown in FIG. conduct.

<Third embodiment>
The data processing method, measurement system, and program of the third embodiment will be described in detail below. In the second embodiment, the transmitting array antenna 50 and the receiving array antenna 52 are arranged in a plane, but in this embodiment, the transmitting array antenna 50 and the receiving array antenna 52 are arranged differently. The transmitting array antenna 50 and the receiving array antenna 52 of this embodiment are arranged in a straight line.
Specifically, in this embodiment, the transmitting antenna 10a and the receiving antenna 10b are arranged in the y direction, as shown in FIG. 13. The moving direction (scanning direction) of the transmitting array antenna 50 and the receiving array antenna 52 is assumed to be the x direction. When viewed from the transmitting array antenna 50 and the receiving array antenna 52, the direction in which the object to be measured is located (the direction in which electromagnetic waves are transmitted) is defined as the z direction.
Note that the moving direction (scanning direction) of the transmitting array antenna 50 and the receiving array antenna 52 may be the y direction.

Therefore, the positional relationship between the measurement object, the transmitting array antenna 50, and the receiving array antenna 52 can be expressed as shown in FIG.
Here, the coordinates of the transmitting point are p ₁ (x' ₁ , y' ₁ , z' ₁ ), and the coordinates of the receiving point are p ₂ (x' ₂ , y' ₂ , z' ₂ ). Let f(x, y, z) be the reflectance at the reflection point (x, y, z) of the object to be measured. Let the measurement data at p ₂ (x' ₂ , y' ₂ , z' ₂ ) be s (x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k). Let the propagation wavelength of electromagnetic waves in vacuum be λ ₀ . Let the relative dielectric constant of the medium be ε _r . Let k be the wave number of the propagating electromagnetic wave.

At this time, the measurement data s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) can be expressed by the following formula.

however

It is.

In formula (3-1), the electromagnetic wave is expressed as a spherical wave, and distance attenuation is omitted. This distance attenuation has been omitted because it has little effect on subsequent processing. When the exponent part of the integrand of the second-stage equation in equation (3-1) is expressed in Fourier transform notation, it becomes the following equation. This is equivalent to decomposing the spherical wave in equation (3-1) into a three-dimensional plane wave.

Here, (k' _x1 , k' _y1 , k' _z1 ) is a component of the wave number vector of the spherical wave propagating from the transmission point to the reflection point. Furthermore, (k' _x2 , k' _y2 , k' _z2 ) is a component of the wave number vector of the spherical wave propagating between the reflection point and the reception point. however,

satisfy.

Here, since the x coordinates of the transmitting point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ) are equal, x'=x When ' ₁ =x' ₂ , equation (3-3) is expressed by the following equation.

Here, perform the following variable substitutions.

The following equation is obtained from equation (3-6).

By calculating the absolute value of the Jacobian from equation (3-7), the following equation is obtained.

By substituting equations (3-6) and (3-8) into equation (3-5) and performing variable substitution, the following equation is obtained.

Here, the integral with respect to u in the second line of equation (3-9) is omitted because it is a constant. When triple inverse Fourier transform is performed on both sides of equation (3-9) for (x, y ₁ , y ₂ ), the following equation is obtained.

Rewrite the left side of equation (3-10) as shown below.

Then, equation (3-10) is expressed by the following equation.

By performing the following integration on both sides of equation (3-12), the following equation is obtained.

Here, the following variable substitution is defined for (k' _y1 , k' _y2 ) and (k' _z1 , k' _z2 ).

Here, the following equation is obtained from equation (3-14).

The absolute value of the Jacobian in this variable substitution is given by the following formula.

By substituting equation (3-1), equation (3-15), and equation (3-17) into equation (3-13) and performing variable substitution, the following equation is obtained.

Here, the integral with respect to v on the right side of the second line of equation (3-18) is omitted because it is a constant. When triple inverse Fourier transform is performed on both sides of equation (3-18) for (k _x , k _y , k _z ), the reflectance f(x, y, z) is obtained as follows.

Here, in order to align the measurement plane with the x'-y' plane passing through the origin, if z'=0, then equation (3-19) can be expressed as follows.

In order to solve equation (3-20), it is necessary to express k as (k' _y1 , k' _y2 , k _z ) or (k _y , v, k _z ).
Solve the simultaneous equations of four equations: equation (3-4), equation (3-6), equation (3-15), and the following equation (3-21) obtained from the assumption.

From this, k is expressed by the following formula.

As described _above , _the data processing unit ₆₆ calculates _the reflectance _{f (} x' , y, z).

The data processing method and program of this embodiment will be described below with reference to FIG. FIG. 14 is a flowchart showing the data processing method of this embodiment.
First, the measurement unit 61 obtains measurement data s(x', y' ₁ , y' ₂ , 0, 0, k) (step S3-1). Then, the data processing unit 66 performs Hilbert transformation on the measurement data s(x', y' ₁ , y' ₂ , 0, 0, k) (step S3-2). Thereby, the imaginary component of the frequency data at each measurement point is obtained.

Next, the data processing unit 66 processes the measurement data s (x', y' ₁ , y' ₂ , 0, 0, k) by performing triple Fourier processing on (x', y' ₁ , y' ₂ ). Conversion is performed (step S3-3). As a result, S(k _x , k' _y1 , k' _y2 , 0, 0, k) is obtained as shown in equation (3-11).
Next, the data processing unit 66 performs variable substitution on S(k _x , k' _y1 , k' _y2 , 0, 0, k) (step S3-4). Specifically, using equations (3-14) and (3-15), the function of (k _x , k' _y1 , k' _y2 , k) is transformed into (k _x , k _y , v, k) Make it a function of This yields S(k _x , k _y , v, 0, 0, k).
Next, the data processing unit 66 performs a triple inverse Fourier transform on S(k _x , k _y , v, 0, 0, k) (k _x , k _y , k _z ) ( Step S3-5). As a result, the reflectance f(x, y, z) is obtained as shown in equation (3-20).

The storage unit 66a stores a program for executing the data processing method of this embodiment. The program stored in the storage section 66a causes the data processing unit 66 to execute the data processing method of this embodiment.

Here, the frequency of the electromagnetic wave is 4.46484 GHz, the relative dielectric constant is 10, the measurement interval Δx = 10 mm, the measurement interval Δy' ₁ = 38.5 mm, the measurement interval Δy' ₂ = 38.5 mm, the measurement width 616 mm, and the The wavelength of the highest frequency λ _fmax = 21.24 mm. At this time, the following equation (3-23) is obtained.

At this time, the Nyquist criterion is not satisfied, and aliasing occurs in each of the x-axis direction, y' _1- axis direction, and y' _2- axis direction. The Nyquist wave number in this case is expressed by the following equation (3-24).

The maximum wave number determined from the highest frequency is expressed by the following equation (3-25).

Here, the Fourier transform regarding the x-axis in step S3-3 is shown in FIG. Furthermore, the Fourier transform regarding the y' _1- axis and the y' _2- axis in step S3-3 is shown in FIG.
From FIG. 15, in the Fourier transform regarding the x-axis, wave numbers whose upper limits are ±k _{x, nyq} and ±k' _{x2, nyq} are output to range A, respectively. Data with a wave number exceeding ±k _{x, nyq} appears as alias data with a period of 2k _{x, nyq} . Therefore, in the variable replacement processing in step S3-4, the variable replacement processing is performed not in the range A but in the range A' having the upper limit of ±k _xmax . As a result, the alias data is not discarded as aliasing noise, but is actually added to the variable replacement process as meaningful data, improving the resolution in the x-axis direction and further increasing the image intensity of shallow targets. becomes possible.
Regarding the y' _1- axis and the y' _2- axis, similarly to the second embodiment, as shown in FIG. 16, instead of the range B, the range B' has the upper limit of ±k' _y1max and ±k' _y2max , respectively. Performs variable replacement processing.

(simulation result)
The results of a computer simulation of the data processing method of the third embodiment will be described below. The simulation conditions are as follows.

・Used frequency band f _min ~ f _max : DC ~ 4.46484 GHz
・Center frequency fc: 2.23242GHz
・Wavelength in vacuum of center frequency λ _0c : 134.38mm
・Relative dielectric constant ε _r of medium: 10
・Wavelength λc of center frequency in medium: 42.50 mm
・Wavelength of the highest frequency in the medium λ _fmax : 21.24 mm
・Number of measurements Nx: 128pt
・Measurement interval Δx: 10mm
・Measurement width: 1280mm
・Number of measurements N _y'1 : 16pt
・Measurement interval Δy'1: 38.5mm
・Measurement width: 616mm
・Number of measurements N _y'2 : 16pt
・Measurement interval Δy'2: 38.5mm
・Measurement width: 616mm
・Coordinates of point target 1 (unit: mm): (640, 307, 20)
・Coordinates of point target 2 (unit: mm): (640, 307, 200)

FIG. 17 shows an xz plane image of point target 1 simulated using the data processing method of comparative example 1. FIG. 18 shows an xz plane image of point target 1 simulated using the data processing method of Example 1. FIG. 19 shows a yz plane image of point target 1 simulated using the data processing method of comparative example 1. FIG. 20 shows a yz plane image of point target 1 simulated using the data processing method of Example 1.
FIG. 21 shows waveforms of the point target 1 in the x-axis direction of Comparative Example 1 and Example 1. FIG. 22 shows waveforms of the point target 1 in the y-axis direction of Comparative Example 1 and Example 1. FIG. 23 shows waveforms of the point target 1 in the z-axis direction of Comparative Example 1 and Example 1. In each of FIGS. 21 to 23, Example 1 is shown by a solid line, and Comparative Example 1 is shown by a broken line.

FIG. 24 shows an xz plane image of point target 2 simulated using the data processing method of comparative example 1. FIG. 25 shows an xz plane image of point target 2 simulated using the data processing method of Example 1. FIG. 26 shows a yz plane image of point target 2 simulated using the data processing method of comparative example 1. FIG. 27 shows a yz plane image of point target 2 simulated using the data processing method of Example 1.
FIG. 28 shows waveforms of the point target 2 in the x-axis direction of Comparative Example 1 and Example 1. FIG. 29 shows waveforms in the y-axis direction of the point target 2 of Comparative Example 1 and Example 1. FIG. 30 shows waveforms in the z-axis direction of the point target 2 of Comparative Example 1 and Example 1. In each of FIGS. 28 to 30, Example 1 is shown by a solid line, and Comparative Example 1 is shown by a broken line.

In Comparative Example 1, variable replacement processing was performed in range A in the variable replacement processing in step S3-4 in the present embodiment. In Example 1, in the variable replacement process of step S3-4 in this embodiment, the variable replacement process was performed not in the range A but in the range A' having the upper limit of ±k _xmax . Hereinafter, the variable replacement process in the range A' having ±k _xmax as the upper limit in Example 1 will be referred to as super-resolution process.

17 to 30, compared to Comparative Example 1, according to Example 1 which performs super-resolution processing, the resolution in the x-axis direction and the y-axis direction is improved for both point targets 1 and 2. confirmed.

<Fourth embodiment>
In this embodiment, noise that may occur in the third embodiment is reduced. The measurement data s in the y-axis direction is expressed by the following formula (4- 1).

At this time, if the amplitude components below the Nyquist frequency are t ₀ and r ₀ and the n folded amplitude components (alias data) are t _n and r _n , then t and r are respectively expressed by the following equations (4-2). It is expressed as

By substituting equation (4-2) into equation (4-1), the following equation (4-3) is obtained.

For example, in the super-resolution processing shown in the third embodiment, when replacing variables in step S3-4, in the (k' _y1 , k' _y2 ) space shown in FIG. As shown in the figure, if k' _y1 caused by switching the transmitting antenna and k' _y2 caused by the receiving antenna are both within the range of the Nyquist wave number or less, the component containing alias data (aliasing noise component) is , all become noise components in signal processing.

At this time, equation (4-3) is expressed as equation (4-5).

Here, r ₀ t ₀ surrounded by a solid line is a signal component, and the part surrounded by a broken line is a noise component.

Next, when only k' _y1 falls into the one-time aliasing noise region during variable substitution, the condition of equation (4-6) is met, and equation (4-3) is expressed as equation (4-7). .


Here, r ₀ t ₁ surrounded by a solid line is a signal component, and the part surrounded by a broken line is a noise component.

Next, when k' _y2 falls into the twice-folding noise region during variable substitution, the condition of equation (4-8) is met, and equation (4-3) is expressed by equation (4-9).


Here, r ₂ t ₁ surrounded by a solid line is a signal component, and the part surrounded by a broken line is a noise component.

Here, as shown in FIG. 5, when the y' ₁ coordinate of the transmitting antenna measurement point and the y' ₂ coordinate of the receiving antenna measurement point are close to the target, that is, when they are near directly above the target, y' The spatial frequencies of the measurement data in _{the 1} and y' _two- axis directions are respectively lower. Conversely, as the y' ₁ and y' ₂ coordinates become farther from the target, the spatial frequencies of the measurement data in the y' ₁ and y' ₂ axis directions become higher, respectively. Then, when the wave number at that point exceeds the Nyquist wave number, it is included in the measurement data as an aliasing noise (alias) component. In this way, alias data components (aliasing noise components) tend to occur at points relatively far from the target. Therefore, due to distance attenuation, the magnitude of the amplitude becomes smaller than that of data components below the Nyquist wave number. Therefore, the following equation (4-10) generally holds true.

From this, the magnitude of the quadratic term regarding r _n and t _n becomes smaller than the term including r ₀ or t ₀ . This is expressed by the following equation (4-11).

Furthermore, unlike in the air where the attenuation rate is small, when electromagnetic waves propagate inside concrete or underground where the attenuation rate is large, the magnitudes of r _n and t _n are determined by the attenuation rate due to the medium in addition to the distance attenuation. , becomes smaller. At this time, equation (4-11) is expressed as equation (4-12) below.

Therefore, from equation (4-5), equation (4-7), equation (4-9), and equation (4-12), the quadratic terms regarding t _n and r _n are relative to other terms. This results in data with a poor S/N ratio. Therefore, in order to improve the S/N ratio, the quadratic term regarding r _n and t _n may be selectively excluded in the variable replacement process of step S3-4. At this time, equation (4-3) is expressed by equation (4-13) below.

Here, only the components surrounded by solid lines are used in the variable replacement process, and the parts surrounded by broken lines are not used in the variable replacement process.

In this way, when electromagnetic waves propagate through a medium with a large attenuation rate, the deterioration of the S/N ratio of a three-dimensional image can be prevented by using the filter expressed by equation (4-13) during variable replacement processing. It can be suppressed. This is expressed as a "cruciform k' _y1 -k' _y2 filter" in the (k' _y1 , k' _y2 ) space. FIG. 31 is a conceptual diagram of a cross-shaped k' _y1 -k' _y2 filter.

k _xmax , k' _y1max , and k' _y2max shown in FIGS. 15 and 16 depend on the position of the target or the directivity of the antenna. Here, if targets near the measurement surface (z=0) are included in the measurement object, the maximum value may be determined by the antenna directivity θ _bx and θ _by as shown in equation (4-14).

however,

It is.

(Simulation result 1)
The results of a computer simulation of the data processing method of the fourth embodiment will be described below. First, the effect of improving the image intensity of a shallow target will be verified under the same simulation conditions as in the third embodiment.
Super-resolution processing uses the missing alias data. Therefore, the target image intensity (image amplitude) increases compared to the case without super-resolution processing. This effect is larger for shallow targets with more alias data than for deep targets. That is, the image intensity of a shallow target, which had been reduced by aliasing, is greatly restored by super-resolution processing.

Figure 32 shows simulation results for four point targets. The horizontal position of each point target is at the center of the measurement surface, and the depths are 20 mm, 200 mm, 400 mm, and 600 mm. In FIG. 32, Example 1 in which super-resolution processing was performed is shown by a solid line, and Comparative Example 1 in which super-resolution processing was not performed is shown in broken lines. The magnitude of the three-dimensional image amplitude is normalized by setting the amplitude of a target with a depth of 20 mm without super-resolution processing to 1.

The shallower the depth of the target, the stronger the image amplitude of the target becomes due to the effect of super-resolution processing. Because the simulation data takes into account distance attenuation, the amplitude of deep targets becomes very small, making it difficult to directly compare the recovery effects. Therefore, regarding the amplitude of a target at the same depth, the ratio of the three-dimensional image amplitude subjected to super-resolution processing to the three-dimensional image amplitude without super-resolution processing is defined as the recovery rate. FIG. 33 shows the recovery rate for each point target. From FIG. 33, it can be seen that the shallower the target, the more alias data there is, so the recovery rate by super-resolution processing is higher. At a depth of 20 mm, the recovery rate reaches 6 times. On the other hand, the deeper the target, the less alias data and the smaller the recovery rate. At a depth of 600 mm, the recovery rate is about 1.2 times.

From the above, it was confirmed that image amplitude intensity is improved by using alias data through super-resolution processing. Particularly noticeable improvements are seen in images of shallow targets with a lot of aliased data. This is mainly reflected in the three-dimensional image as an improvement in the S/N ratio of the shallow target image.

(Simulation result 2)
Hereinafter, the results of computer simulation of the data processing methods of the third and fourth embodiments will be described. The simulation conditions are as follows.
・Used frequency band f _min ~ f _max : DC ~ 4.46484 GHz
・Center frequency fc: 2.23242GHz
・Wavelength in vacuum of center frequency λ _0c : 134.38mm
・Relative dielectric constant ε _r of medium: 5
・Wavelength λc of center frequency in medium: 60.1 mm
・Wavelength of the highest frequency in the medium λ _fmax : 30.05 mm
・Measurement interval Δx in scanning direction: 10mm
・Number of measurements Nx: 115pt
・Measurement width: 1150mm
・Measurement interval Δy1 of transmitting antenna: 38.5mm
・Number of transmitting antennas N _y1 : 12
・Reception antenna measurement interval Δy2: 38.5mm
・Number of receiving antennas N _y2 : 12pt
・Measurement width: 481.25mm

The results of measuring the reinforced concrete test specimen are shown below. FIG. 34 shows the layout of the transmitting antenna and receiving antenna used for measurement. The test specimen has a plurality of reinforcing bars with different z-coordinates inside. FIG. 35 is a three-dimensional image of the test specimen without super-resolution processing. FIG. 36 is a three-dimensional image of the specimen subjected to super-resolution processing. In FIGS. 35 and 36, the z-coordinates of a plurality of reinforcing bars are shown. Comparing FIG. 35 and FIG. 36, in FIG. 36 in which super-resolution processing was performed, the image of the reinforcing bars was thinner overall, and it was confirmed that the image resolution was improved. In addition, the image intensity of reinforcing bars at depths of 35 mm and 45 mm became stronger, and the effect of restoring the image amplitude of shallow targets was also confirmed. On the other hand, it was also confirmed that unnecessary noise (aliasing noise) such as interference fringes was newly generated between reinforcing bars.

Next, the cross-shaped k' _y1 -k' _y2 filter described in the fourth embodiment and FIG. 31 was applied. FIG. 37 is a three-dimensional image of a test specimen subjected to super-resolution processing by applying a cross-shaped k' _y1 -k' _y2 filter. Figure 36 is a three-dimensional image of the specimen subjected to super-resolution processing without applying the cross-shaped k' _y1 - k' _y2 filter, and the super-resolution processing performed by applying the cross-shaped k' _y1 - k' _y2 filter. Comparing Fig. ₃₇ , _which was performed with It could be confirmed.

Next, the influence of the cross-shaped k' _y1 -k' _y2 filter on resolution was evaluated using simulation data. The simulation conditions are the same as in the third embodiment.
FIG. 38 shows the waveforms of the point target 1 in the x-axis direction of Comparative Example 2 and Example 2. FIG. 39 shows waveforms of the point target 1 in the y-axis direction of Comparative Example 2 and Example 2. FIG. 40 shows waveforms in the z-axis direction of the point target 1 of Comparative Example 2 and Example 2. FIG. 41 shows the waveforms of the point target 2 in the x-axis direction of Comparative Example 2 and Example 2. FIG. 42 shows waveforms of the point target 2 in the y-axis direction of Comparative Example 2 and Example 2. FIG. 43 shows waveforms in the z-axis direction of the point target 2 of Comparative Example 2 and Example 2.
In each of FIGS. 38 to 43, Example 2 is shown by a solid line, and Comparative Example 2 is shown by a broken line.

In Comparative Example 2, in the variable replacement process of step S3-4 in the third embodiment, the cross-shaped k' _y1 -k' _y2 filter is not applied, and variable replacement is performed in the range A' with the upper limit of ±k _xmax . processed. In Example 2, in the variable replacement process in step S3-4 in the third embodiment, a cross-shaped k' _y1 -k' _y2 filter is applied to perform variable replacement in a range A' whose upper limit is ±k _xmax . processed.
From Figures 38 to 43, even when super-resolution processing is performed by applying a cross-shaped k' _y1 - k' _y2 filter, there is almost no effect on the image resolution on any of the x-, y-, and z-axes. was confirmed.

<Fifth embodiment>
In this embodiment, noise that may occur in the second embodiment is reduced. In the fourth embodiment, super-resolution processing is performed by applying a cross-shaped k' _y1 -k' _y2 filter in the third embodiment. In this embodiment, super-resolution processing is performed by applying the cross-shaped k' _x1 -k' _x2 filter and the cross-shaped k' _y1 -k' _y2 filter in the second embodiment. FIG. 44 is a conceptual diagram of a cross-shaped k' _x1 -k' _x2 filter. That is, in this embodiment, super-resolution is achieved by applying both the cross-shaped k' _x1 -k' _x2 filter shown in FIG. 44 and the cross-shaped k' _y1 -k' _y2 filter shown in FIG. 31 in the second embodiment. Perform processing.

k' _x1max , k' _x2max , k' _y1max , and k' _y2max depend on the position of the target or the directivity of the antenna. Here, if targets near the measurement surface (z=0) are included in the measurement object, the maximum value may be determined by the antenna directivity θ _bx and θ _by as shown in equation (5-1).

however,

It is.

10a Transmitting antenna 10b Receiving antenna 50 Transmitting array antenna 52 Receiving array antenna 60 Radar device 61 Measurement unit 64 System control circuit 66 Data processing unit 68 Image display units 58, 59 High frequency switch 62 High frequency circuit 69 Encoder

Claims (17)

1. A data processing method for analyzing scattered waves of waves emitted to an object,
   radiating the wave to the object from a plurality of transmission points p ₁ (x' ₁ , y' ₁ , z' ₁ ) arranged on the y-axis;
   The scattered waves reflected at the reflection point (x, y, z) on the object with a reflectance f(x, y, z) are received at a plurality of receiving points p ₂ (x' ₂ , receive the measured value s(x _{' 1 ,} x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) at y' ₂ , z' 2 );
   The measured value s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) is subjected to triple Fourier transform using equation (1) to obtain S(k' _x1 , k ' _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k),
   
   If the Nyquist wave number in the x direction determined by the measurement interval in the x direction is k _x,nyq , then variable substitution processing in the x direction is performed in the range -k _x1 ≦k _x ≦k _x1 (however, k _x,nyq ≦k _x1 ). conduct,
   If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 ), perform variable substitution processing in the y direction,
   The reflectance f(x, y, z) is determined by triple inverse Fourier transform using equation (2).
   
   Data processing method.
   however,
   x' = x' ₁ = x' ₂
   z' ₁ = z' ₂ = 0
   k is the wave number of the propagating wave,
   k' _x1 , k' _y1 , k' _z1 are the waves propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z). Components of the wave vector of a spherical wave,
   k' _x2 , k' _y2 , k' _z2 are the waves propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ). Components of the wave vector of a spherical wave,
   k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
   It is.
2. If the maximum wave number in the x direction determined by equation (3) is k _xmax , perform variable substitution processing in the x direction in the range -k _xmax ≦k _x ≦k _xmax ,
   
   If the maximum wave number in the y' ₁ direction determined by equation (4) is k' _y1max , and the maximum wave number in the y' ₂ direction is k' _y2max , then -k' _y1max ≦k' _y1 ≦k' _y1max and -k ' _y2max ≦k' _y2 ≦k' Perform variable substitution processing in the y direction within the range of _y2max ,
   
   The data processing method according to claim 1.
3. In the variable substitution process in the y direction, the range of k' _y1,nyq ≦|k' _y1 |≦k' _y1max and k' _y2,nyq ≦|k' _y2 |≦k' _y2max is selectively excluded. and perform the variable replacement process in the y direction.
   The data processing method according to claim 2.
4. the plurality of transmitting points p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the plurality of receiving points p ₂ (x' ₂ , y' ₂ , z' ₂ ) move in the x direction;
   The data processing method according to any one of claims 1 to 3.
5. the plurality of transmission points p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the plurality of reception points p ₂ (x' ₂ , y' ₂ , z' ₂ ) move in the y direction;
   The data processing method according to any one of claims 1 to 3.
6. The wave number of the wave and the component of the wave number vector satisfy equation (5),
   The data processing method according to any one of claims 1 to 5.
   
7. A measurement system that analyzes scattered waves of waves emitted to an object,
   A transmitting/receiving section,
   a transmitter that radiates the wave to the object from a plurality of transmission points p ₁ (x' ₁ , y' ₁ , z' ₁ ) arranged on the y-axis;
   The scattered waves reflected at the reflection point (x, y, z) on the object with a reflectance f(x, y, z) are received at a plurality of receiving points p ₂ (x' ₂ , a receiving unit that receives the measured value s(x _{' 1 ,} x' ₂ , y' ₁ , y' ₂ , z _{' 1 ,} z' ₂ , k) at y' 2 , z' 2 );
   a transmitting/receiving unit having;
   A processing device,
   The measured value s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) is subjected to triple Fourier transform using equation (1) to obtain S(k' _x1 , k ' _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k);
   
   If the Nyquist wave number in the x direction determined by the measurement interval in the x direction is k _x,nyq , then variable substitution processing in the x direction is performed in the range -k _x1 ≦k _x ≦k _x1 (however, k _x,nyq ≦k _x1 ). The steps to take and
   If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 );
   A step of calculating the reflectance f(x, y, z) by performing triple inverse Fourier transform from equation (2);
   
   a processing device that executes;
   A measurement system with
   however,
   x' = x' ₁ = x' ₂
   z' ₁ = z' ₂ = 0
   k is the wave number of the propagating wave,
   k' _x1 , k' _y1 , k' _z1 are the waves propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z). Components of the wave vector of a spherical wave,
   k' _x2 , k' _y2 , k' _z2 are the waves propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ). Components of the wave vector of a spherical wave,
   k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
   It is.
8. The processing device includes:
   If the maximum wave number in the x direction determined by equation (3) is k _xmax , then the procedure for performing variable substitution processing in the x direction in the range -k _xmax ≦k _x ≦k _xmax ,
   
   If the maximum wave number in the y' ₁ direction determined by equation (4) is k' _y1max , and the maximum wave number in the y' ₂ direction is k' _y2max , then -k' _y1max ≦k' _y1 ≦k' _y1max and -k ' _y2max ≦k' _y2 ≦k' Procedure for performing variable substitution processing in the y direction within the range of _y2max ,
   
   execute,
   The measurement system according to claim 7.
9. In the variable substitution process in the y direction, the processing device calculates a range of k' _y1,nyq ≦|k' _y1 |≦k' _y1max and k' _y2,nyq ≦|k' _y2 |≦k' _y2max . selectively excluding and performing a procedure for performing the variable replacement process in the y direction;
   The measurement system according to claim 8.
10. the transmitter/receiver moves in the x direction;
    The measurement system according to any one of claims 7 to 9.
11. the transmitting and receiving unit moves in the y direction;
    The measurement system according to any one of claims 7 to 9.
12. The wave number of the wave and the component of the wave number vector satisfy equation (5),
    The measurement system according to any one of claims 7 to 11.
    
13. A program that analyzes scattered waves of waves emitted to an object,
    _{The measured value} s _{( x ' 1} , _x2 , k' _y1 , k' _y2 , z' ₁ , z' ₂ , k);
    
    If the Nyquist wave number in the x direction determined by the measurement interval in the x direction is k _x,nyq , then the variable substitution process in the x direction is performed in the range -k _x1 ≦k _x ≦k _x1 (k _x,nyq ≦k _x1 ). The steps to take and
    If the Nyquist wave number in the y' ₁ direction determined by the measurement interval in the y direction is k' _y1,nyq , and the Nyquist wave number in the y' ₂ direction is k' _y2,nyq , -k _y1 ≦k' _y1 ≦k _y1 , and , -k _y2 ≦k' _y2 ≦k _y2 (however, k' _y1,nyq ≦k _y1 and k' _y2,nyq ≦k _y2 );
    Steps to calculate the reflectance f(x, y, z) by performing triple inverse Fourier transform from equation (2);
    
    A program that causes a computer to execute.
    however,
    x' = x' ₁ = x' ₂
    z' ₁ = z' ₂ = 0
    k is the wave number of the propagating wave,
    k' _x1 , k' _y1 , k' _z1 are the waves propagating between the transmission point p ₁ (x' ₁ , y' ₁ , z' ₁ ) and the reflection point (x, y, z). Components of the wave vector of a spherical wave,
    k' _x2 , k' _y2 , k' _z2 are the waves propagating between the reflecting point (x, y, z) and the receiving point p ₂ (x' ₂ , y' ₂ , z' ₂ ). Components of the wave vector of a spherical wave,
    k _x = k' _x1 + k' _x2 , u = k' _x1 - k' _x2 , k _y = k' _y1 + k' _y2 , v = k' _y1 - k' _y2
    It is.
14. If the maximum wave number in the x direction determined by equation (3) is k _xmax , then the procedure for performing variable substitution processing in the x direction in the range -k _xmax ≦k _x ≦k _xmax ,
    
    If the maximum wave number in the y' ₁ direction determined by equation (4) is k' _y1max , and the maximum wave number in the y' ₂ direction is k' _y2max , then -k' _y1max ≦k' _y1 ≦k' _y1max and -k ' _y2max ≦k' _y2 ≦k' Procedure for performing variable substitution processing in the y direction within the range of _y2max ,
    
    The program according to claim 13, further causing a computer to execute.
15. In the variable substitution process in the y direction, the range of k' _y1,nyq ≦|k' _y1 |≦k' _y1max and k' _y2,nyq ≦|k' _y2 |≦k' _y2max is selectively excluded. and causing the computer to further execute a procedure for performing the variable substitution process in the y direction,
    The program according to claim 14.
16. The measurement value s(x' ₁ , x' ₂ , y' ₁ , y' ₂ , z' ₁ , z' ₂ , k) is obtained from multiple transmission points p ₁ (x' ₁ , y' ₁ , z' ₁ ) to the object, and the wave is reflected at the reflection point (x, y, z) on the object with the reflectance f(x, y, z). This is the value of the scattered waves received at multiple reception points p ₂ (x' ₂ , y' ₂ , z' ₂ ) arranged on the y-axis.
    The program according to any one of claims 13 to 15.
17. The wave number of the wave and the component of the wave number vector satisfy equation (5),
    The program according to any one of claims 13 to 16.
    

Images