CN112147608A - Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method - Google Patents

Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method Download PDF

Info

Publication number
CN112147608A
CN112147608A CN202011000874.9A CN202011000874A CN112147608A CN 112147608 A CN112147608 A CN 112147608A CN 202011000874 A CN202011000874 A CN 202011000874A CN 112147608 A CN112147608 A CN 112147608A
Authority
CN
China
Prior art keywords
uniform
wall
imaging
radar
gaussian
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202011000874.9A
Other languages
Chinese (zh)
Inventor
李家强
陈焱博
刘然
陈金立
史佳惠
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing University of Information Science and Technology
Original Assignee
Nanjing University of Information Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing University of Information Science and Technology filed Critical Nanjing University of Information Science and Technology
Priority to CN202011000874.9A priority Critical patent/CN112147608A/en
Publication of CN112147608A publication Critical patent/CN112147608A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/887—Radar or analogous systems specially adapted for specific applications for detection of concealed objects, e.g. contraband or weapons
    • G01S13/888—Radar or analogous systems specially adapted for specific applications for detection of concealed objects, e.g. contraband or weapons through wall detection
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging

Abstract

The invention discloses a rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method, which comprises the following steps: establishing a through-wall radar detection scene model, arranging antenna array elements parallel to a wall body, transmitting electromagnetic waves, and uniformly sampling an obtained target echo signal to obtain echo signal data e; gridding pixel points in the imaging area of the through-wall radar, projecting echo signal data onto each pixel point grid in the imaging area to obtain an amplitude value I (x) of a BP imaging pixel point of the through-wall radarp,yp) (ii) a Constructing an expression of the non-uniform FFT, and pre-calculating and storing each coefficient in the expression of the non-uniform FFT; transforming uniformly sampled data to non-uniformly sampled data, pair I (x)p,yp) And (4) performing fast Gaussian gridding non-uniform FFT calculation and accumulating to obtain a radar imaging graph. The invention can effectively reduce the calculation complexity while ensuring the imaging quality, and solves the problem that the calculation amount of the BP method is along with the working frequencyThe number of points, antennas and pixel points are increased.

Description

Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method
Technical Field
The invention belongs to the technical field of through-wall radar imaging, and particularly relates to a rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method.
Background
The through-wall imaging radar can detect targets behind hidden objects such as buildings, forts, grasses and the like by utilizing the good penetration characteristic of electromagnetic waves, has low hardware cost, light weight, small volume and good portability, can perform real-time target positioning, classification and identification, and is widely applied to military and civil fields such as law enforcement action, anti-terrorism and fight, criminal investigation, emergency rescue and the like.
Common imaging methods for the through-wall radar include a compressed sensing imaging method, a backscattering tomography method, a boundary estimation method, a Back Projection (BP) method, and the like. The compressed sensing imaging method has the advantages that sparse or compressible echo signals are subjected to compressed sensing imaging, the norm optimization problem with constraint is solved, the original echo signal imaging is reconstructed by using the sampling number far lower than that required by the Nyquist sampling theorem, and the problems of additive noise, model error interference solving process, image result and the like still exist although the compressed sensing imaging application prospect is wide. The backscattering chromatography method can accurately reconstruct the image of an imaging scene, but the backscattering chromatography method needs multiple iterative operations, has huge calculated amount and is difficult to apply in engineering. The boundary estimation method that the inverse boundary scattering change exists between the shape of the target boundary and the receiving pulse delay is utilized, the target boundary can be clearly imaged, but the propagation path and the time delay of the electromagnetic wave can be changed when an obstacle exists, and the method does not meet the inverse boundary scattering change any more, so that the target cannot be accurately imaged.
In recent years, the BP method projects radar echo data to each pixel point of an imaging area, calculates the time delay of the distance between a radar antenna and an image pixel of the radar echo, and performs coherent accumulation on a time domain to realize high-resolution imaging. The BP method is widely applied to through-wall imaging radars, but has large calculation data amount and redundancy phenomenon, so that the calculation complexity is high, and the memory requirement of a computer is high, so that domestic and foreign scholars propose some related BP acceleration methods.
Document [1] proposes a fast factorization back projection method, which projects a sub-aperture signal onto an imaging grid under local polar coordinates, corrects to obtain a coarse sub-image through deviation of distance dimension and rotation of angle dimension, then fuses the sub-images to obtain a final imaging image, and subsequently Kyra Moon et al proposes a new method of decomposing BP, namely document [2], which divides the image into separately processed columns and coherently adds the sub-images. The method is easy to realize parallelization, and can reduce the computational complexity because each column can be formed independently of other columns. The above methods are all performed by fusing after dividing the sub-images, and although the computational complexity is reduced, the imaging accuracy is also reduced. Document [3] proposes a BP imaging method for iterative sub-images, which reduces complexity by performing image reconstruction on iterative sub-images, but is not suitable for large-scale radar imaging. The document [4] proposes a fast BP method, which is based on the principle that an imaging region is divided into blocks, and the calculation amount of the BP method is reduced in a hierarchical coherent accumulation mode. The BP imaging method based on the non-uniform fast Fourier transform, namely, document [5] utilizes a pixel point amplitude expression in the BP method to meet the non-uniform Fourier transform expression, adopts the fast Fourier transform to calculate, reduces the calculation complexity and ensures the imaging quality, is suitable for large-scale real-time imaging, but the method needs to repeatedly calculate coefficients in the amplitude expression, and has large requirements on a computer memory. Therefore, it is a research focus of the present invention to reduce the computational complexity and memory requirement of the BP method while ensuring the imaging quality.
[1]Zuo S,Sun G,Xing M,et al.A modified fast factorized back projection algorithm for the spotlight SAR imaging[C]//2015IEEE 5th Asia-Pacific Conference on Synthetic Aperture Radar (APSAR).IEEE,2015.
[2]K.Moon and D.G.Long,A new factorized backprojection algorithm for stripmap synthetic aperture radar[J],Positioning,2013,4(1):42-56.
[3]Kelly S I,Rilling G,Davies M,et al.Iterative image formation using fast (Re/Back)-projection for spotlight-mode SAR[C]//IEEE Radarcon.IEEE,2011.
[4]Boag A,Bresler Y,Michielssen E.A multilevel domain decomposition algorithm for fast O(N/sup 2/logN)reprojection of tomographic images[J].IEEE Transactions on Image Processing, 2000,9(9):1573-1582.
[5]Qu LL,Yin Y Q,Sun Y P,et al.Efficient back projection imaging approach for airborne GPR using NUFFT technique[C]//2016 16th International Conference on Ground Penetrating Radar(GPR).IEEE,2016.
Disclosure of Invention
The technical problem to be solved by the invention is to provide a fast Gaussian gridding non-uniform FFT through-wall imaging radar BP method aiming at the defects of the prior art, firstly deconvolving the pixel point amplitude matrix coefficient and the Gaussian kernel function in the BP method to eliminate the influence of Gaussian smoothness, then carrying out fast Fourier transform on uniform data, and finally carrying out convolution operation on the obtained data to realize uniform and smooth output of the data. Through processing the through-the-wall radar data obtained by simulation software GprMax2D/3D based on a Finite-Time difference method (FDTD), simulation experiments prove that the method effectively reduces the computational complexity and the memory requirement under the condition of ensuring the imaging quality, and can effectively solve the problems of higher computational complexity and large memory requirement of a computer due to the fact that the BP method in the through-the-wall imaging radar has larger computation data amount and redundancy phenomenon.
In order to achieve the technical purpose, the technical scheme adopted by the invention is as follows:
a rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method comprises the following steps:
step 1, establishing a through-wall radar detection scene model, arranging antenna array elements parallel to a wall body, transmitting electromagnetic waves, and uniformly sampling an obtained target echo signal to obtain echo signal data e;
step 2, gridding and dividing the pixels in the imaging area of the through-wall radar, projecting the echo signal data to each pixel grid in the imaging area to obtain an amplitude value I (x) of the BP imaging pixel of the through-wall radarp,yp)；
And step 3: constructing an expression of non-uniform FFT, converting the uniform sampling data into the non-uniform sampling data, and pre-calculating and storing each coefficient in the expression of the non-uniform FFT;
step 4, converting the uniform sampling data into non-uniform sampling data, and comparing I (x)p,yp) And performing fast Gaussian gridding non-uniform FFT calculation, and finally accumulating to obtain a radar imaging graph.
In order to optimize the technical scheme, the specific measures adopted further comprise:
the through-wall radar detection scene model established in the step 1 is as follows:
the wall thickness of the front wall and the rear wall is d, and the relative dielectric constant iswThe detection target is an ideal electric conductor with the radius r, and is arranged between the two walls, and the vertical distance between the circle center and the front wall is equal to the vertical distance;
is parallel to the wall body and has a distance h from the wall body1The N antenna array elements transmit electromagnetic waves, and M times of sampling is carried out on target echo signals of the N antenna array elements to obtain an M multiplied by N dimensional matrix e:
e＝[z1,z2,…zN] (1)
in the formula, zn＝[zn(1),zn(2),…,zn(M)]T，n＝1,2,…N，zn(M), M is 1,2, … M represents the echo signal amplitude of the same antenna at different sampling points.
In the step 2, the pixels in the through-wall radar imaging area are divided into grids, the whole imaging area is divided into P pixel grids, and z is dividedn(m) projecting the data to each pixel point grid of the imaging area to obtain a through-wall radar BP imaging pixel point amplitude value:
Figure 7
in the formula, I (x)p,yp) Complex amplitude values for the P (P ═ 1,2, …, P-1) th pixel grid point, fm＝f0+ M Δ f is the M (M ═ 1,2, … M) th operating frequency point, f0For the starting frequency of the transmitted signal, Δ f is the frequency interval, τpFor the two-way time delay between the imaging pixel grid point p and the nth (N ═ 1,2, … N) antenna:
τp＝2(l1+l3)/c+2l2/v (3)
in the formula (I), the compound is shown in the specification,
Figure BDA0002694274030000041
is the propagation velocity of the electromagnetic wave in the wall body, c is the propagation velocity of the electromagnetic wave in the air,wis the relative dielectric constant of the wall;
l1，l2，l3respectively showing the slant distance between the antenna array element for transmitting the electromagnetic wave and the front wall, the distance of the electromagnetic wave transmitted on the wall and the distance of the electromagnetic wave reaching a target after passing through the wall.
The constructing of the non-uniform FFT expression in step 3 above, which is used to transform the uniformly sampled data into the non-uniformly sampled data, and pre-calculate and store each coefficient in the non-uniformly FFT expression, includes the following steps:
step 3.1: for uniform sampling sequence
Figure BDA0002694274030000042
Defining a non-uniform sampling sequence xjDiscrete fourier transform of (d):
Figure BDA0002694274030000043
wherein xj∈[0,2π] (4)
gτ(x) Is [0,2 π]The above one-dimensional periodic gaussian kernel function has the following expression:
Figure BDA0002694274030000044
in the formula, tau is a Gaussian kernel function parameter and determines the exponential decay rate of the Gaussian kernel function;
step 3.2: constructing non-uniform FFT expressions
Figure BDA0002694274030000045
Step 3.3: precomputing and storing f (x)j) To obtain the expected value f (x)j)。
The step 3.2 comprises the following steps:
step 3.2.1: for Gaussian kernel function gτ(x) Fourier transform to obtain Gτ(k) Let Gτ(k) Deconvoluting with coefficients in the expression, i.e. with uniform data F (k), to remove the effect of Gaussian smoothing and obtain an auxiliary function F-τ(k)：
Figure BDA0002694274030000046
Step 3.2.2: to F-τ(k) F is obtained by fast discrete Fourier inverse transformation-τ(x)：
Figure BDA0002694274030000047
In the formula, Kr＝R×K，KrIs the oversampled grid number, σ is the oversampled grid point, and R is the oversampling coefficient. Due to f-τ(x) The uniformity is satisfied and is distributed on uniform grid points, so that the calculation can be carried out by using fast Fourier transform;
step 3.2.3: to f-τ(x) Convolution is carried out to realize smooth output, and an expected value f (x) is obtainedj)：
Figure BDA0002694274030000051
In step 3.3 above, only x is calculatedjMost recent KspA plurality of uniform grid points, f (x) is calculated and stored in advancej) The obtained expected value f (x) of each coefficient in (b)j) Comprises the following steps:
Figure BDA0002694274030000052
in the formula, -Ksp＜σ'＜Ksp；
Figure BDA0002694274030000053
KspThe number of grid points covered by the unilateral extendable coverage of the Gaussian kernel function is shown, and sigma' is the grid points of the bilateral extendable coverage of the Gaussian kernel function; let ζ 2 π σ/Kr(ζ≤xj) And xi is expressed in xjVery close uniform grid points.
In the above step 4, the amplitude value I (x) of the pixel point is calculatedp,yp) Is deformed into
Figure BDA0002694274030000054
Wherein
Figure 6
Amplitude I of each pixeln(xp,yp) Then corresponds to non-uniform sample point data f (x)j) Uniform target echo data Zn(m) corresponding to the uniform sampling points F (k), the product of the frequency interval and the two-way delay 2 π Δ f τpCorresponding to non-uniform data xj。
The invention has the following beneficial effects:
1. the invention can effectively reduce the calculation complexity while ensuring the imaging quality, solves the problem that the calculation amount of the BP method is increased along with the increase of working frequency points, the number of antennas and pixel points, and provides a new thinking and solving way for large-scale imaging and real-time imaging of the through-wall radar.
2. The invention can pre-calculate and store the coefficient of the pixel point amplitude matrix in the BP method, avoids repeated calculation and effectively reduces the memory requirement on a computer.
Drawings
FIG. 1 is a flow chart of an implementation of the present invention;
FIG. 2 is a through-wall radar detection scene model;
FIG. 3 is a signal time domain echo plot;
FIG. 4 is a grid of pixel points of an imaging region;
FIG. 5 is a result of imaging based on the BP method;
FIG. 6 is an imaging result based on the non-uniform FFT BP method;
FIG. 7 is an imaging result based on a fast gridding non-uniform FFT BP method;
FIG. 8 is a graph of normalized amplitude of imaging pixels for the BP method, the non-uniform FFT BP method, and the fast gridding non-uniform FFT BP method;
FIG. 9 is a graph of the computational complexity and the change of the grid point of a pixel point for the BP method, the non-uniform FFT BP method, and the fast gridding non-uniform FFT BP method;
FIG. 10 is a graph of memory requirements and pixel point grid point changes for the BP method, the non-uniform FFT BP method, and the fast gridding non-uniform FFT BP method
Detailed Description
Embodiments of the present invention are described in further detail below with reference to the accompanying drawings.
Referring to fig. 1, the fast gaussian gridding non-uniform FFT through-wall imaging radar BP method of the present invention includes:
step 1, establishing a through-wall imaging radar model as shown in fig. 2: the wall thickness of the front wall and the rear wall is d, and the relative dielectric constant iswThe ideal electric conductor with radius r as the detection target is arranged between the two walls, and the vertical distance between the circle center and the front wall is equal to the vertical distance between the circle center and the front wall. A group of antenna array elements is arranged, which are parallel to the wall body and have a distance h from the wall body1And carrying out omnidirectional scanning to obtain echo signals.
Respectively sampling the echo signals of the N antenna array elements, recording the sampling times as M, and respectively forming an M multiplied by N dimensional matrix e by two groups of echo signal data:
e＝[z1,z2,…zN] (1)
in the formula, zn＝[zn(1),zn(2),…,zn(M)]T，n＝1,2,…N，zn(M), M is 1,2, … M represents the echo signal amplitude of the same antenna at different sampling points.
Step 2, fig. 3 shows the time domain echo of the signal, which needs to be imaged by using a BP method. FIG. 4 is a graph of gridding pixel points in a through-wall radar imaging area, dividing a whole detection scene into P pixel grid points l1，l2，l3Respectively representing the slant distance between the antenna array element for transmitting the electromagnetic wave and the front wall body, the distance of the electromagnetic wave transmitted on the wall body and the distance of the electromagnetic wave reaching a target after passing through the wall body, the through-wall radar BP imaging can be represented as follows:
in the formula, I (x)p,yp) Complex amplitude values for the P (P ═ 1,2, …, P-1) th pixel grid point, fm＝f0+ M Δ f is the M (M ═ 1,2, … M) th operating frequency point, f0For the starting frequency of the transmitted signal, Δ f is the frequency interval, τpFor the two-way time delay between the imaging pixel grid point p and the nth (N ═ 1,2, … N) antenna:
τp＝2(l1+l3)/c+2l2/v (3)
in the formula (I), the compound is shown in the specification,
Figure BDA0002694274030000072
is the propagation velocity of the electromagnetic wave in the wall body, c is the propagation velocity of the electromagnetic wave in the air,wis the relative dielectric constant of the wall.
And 3, the non-uniform FFT can directly and quickly process the non-uniform data, namely the uniformly sampled data can be converted into the non-uniformly sampled data.
For uniform sampling sequence
Figure BDA0002694274030000073
The discrete fourier transform of the non-uniform sampling sequence is defined as:
Figure BDA0002694274030000074
wherein xj∈[0,2π] (4)
gτ(x) Is [0,2 π]The above one-dimensional periodic gaussian kernel function has the following expression:
Figure BDA0002694274030000075
in the formula, tau is a Gaussian kernel function parameter and determines the exponential decay rate of the Gaussian kernel function.
For Gaussian kernel function gτ(x) Fourier transform to obtain Gτ(k) Let Gτ(k) Deconvoluting with the uniform data F (k) to eliminate the influence of Gaussian smoothing to obtain an auxiliary function F-τ(k)：
Figure BDA0002694274030000081
Further performing inverse discrete Fourier transform on the formula (6) to obtain f-τ(x)：
Figure BDA0002694274030000082
In the formula, Kr＝R×K，KrIs the oversampled grid number, σ is the oversampled grid point, and R is the oversampling coefficient. Due to f-τ(x) Homogeneity is satisfied and is distributed over uniform grid points, so it can be calculated using the fast fourier transform.
To f-τ(x) Convolution is carried out to realize smooth output, and an expected value f (x) is obtainedj)：
Figure BDA0002694274030000083
When solving equation (8), non-uniform sampling point x is calculated each timejAll uniform grid points are traversed in time, so that the calculation amount is huge, and due to the exponential attenuation characteristic of a Gaussian function, the distance x is far awayjThe grid points of (a) are ignored, so the grid spread range can be set, considering only xjNearby KspDot, KspNumber of grid points covered by one-sided extension of Gaussian kernel function (K)sp6 stands for single precision, Ksp12 for double precision).
In the formula (8) < f >-τ、gτCan be expressed as:
Figure BDA0002694274030000084
Figure BDA0002694274030000085
in the formulae (9) and (10),
Figure 1
and
Figure 2
all contain xjTo do so
Figure BDA0002694274030000088
And xjIs irrelevant. Thus, for non-uniform sampling point xjOnly two exponents need to be calculated and stored
Figure 3
And
Figure BDA00026942740300000810
while
Figure BDA00026942740300000811
Figure BDA00026942740300000812
The calculation is performed only once in the whole calculation process, so that repeated calculation is avoided, and the storage space is reduced. Let ζ 2 π σ/Kr(ζ≤xj) And xi is expressed in xjVery close uniform grid points. Calculating a storage index
Figure BDA00026942740300000813
E1＝eikζ，
Figure BDA00026942740300000814
The expected value f (x)j) Can be expressed as:
Figure BDA0002694274030000091
in the formula, -Ksp＜σ'＜Ksp。
Step 4, echo data z due to the targetn(m) is uniform data, time delay taupAmplitude I (x) for a single pixelp,yp) The method is non-uniform and cannot be directly calculated by fast Fourier transform, and the fast Gaussian gridding non-uniform FFT method provided by the invention is adaptive to the non-uniformity and can be used for carrying out the fast Fourier transform on non-uniform data.
The method has the main idea that the uniform sample and the Gaussian function are deconvoluted to eliminate the influence of Gaussian smoothing, then fast Fourier transform is carried out, and finally smooth output is realized by convolving the sample with the Gaussian function. Therefore, the amplitude I (x) of the pixel point is processed by using the fast Gaussian gridding non-uniform FFT methodp,yp) And performing quick calculation.
The BP imaging method in equation (2) can be deformed as follows:
Figure BDA0002694274030000092
wherein the content of the first and second substances,
Figure 4
equation (13) can be calculated by fast Gaussian gridding non-uniform FFT, where the amplitude I of each pixel pointn(xp,yp) Non-uniform sampling point f (x) in corresponding equation (4)j) Uniform target echo data Zn(m) corresponds to the uniform sampling point F (k) in the formula (4), the product of the frequency interval and the two-way time delay is 2 pi delta f taupNon-uniform data x in corresponding equation (4)j。
Simulation experiment 1:
a BP imaging method, a non-uniform FFT BP imaging method and the method provided by the invention are respectively used for imaging the through-wall radar detection area. The frequency range of the electromagnetic wave is selected to be 1 GHz-2 GHz, the frequency interval delta f is 0.49MHz, and the frequency point M is 2036. The imaging area is set to be at a transverse distance of 2.2m and a longitudinal distance of 2.1m, and 279 × 293 pixels are divided according to the transverse direction and the longitudinal direction. The imaging results of the three methods are shown in fig. 5, 6 and 7. In order to compare the imaging resolutions of the three methods, the pixel point matrix is normalized in amplitude, and fig. 8 is a graph of normalized amplitude of the imaging pixel points of the three methods. The width of a pixel point at the position of-3 dB reduction of the maximum imaging amplitude is taken, BP imaging is 3.8 pixel points, and nonuniform FFT BP imaging is 5.3 pixel points, the imaging of the method is 4.9 pixel points, the width difference of the three methods is small, accurate imaging can be carried out on a target area, and the method keeps good effect and is better than the nonuniform FFT BP imaging.
Simulation experiment 2:
the calculation complexity of BP imaging, non-uniform FFT BP imaging and imaging of the method provided by the invention are briefly analyzed. Supposing to neglect the influence of the number of the antennas, taking N as 1, and because the sampling point of the antenna is far less than the number of the pixel points, namely M & lt P. The calculation complexity of the BP method is C0MP, the computational complexity of the non-uniform FFT BP is C1M log M + | log | P ≈ | log | P, where is the ideal calculation accuracy. The invention sets the Gaussian kernel function single side to extend and cover 6 grid points, namely the calculation precision Ksp6; the oversampling coefficient R is 2. In the process of the invention E, E1，E2， E3，E4Can be calculated and stored in advance without repeated calculation, so the calculation complexity of the method is C2＝(Ksp+1)M+(Krlog Kr) + P ≈ P. Fig. 9 is a graph of the computational complexity of the three methods and the change of the grid point of the pixel point, and comparison shows that the computational complexity of the method of the present invention is much lower than that of the other two methods for the same pixel point under the same condition.
Simulation experiment 3:
the memory requirements of BP imaging, non-uniform FFT BP imaging and imaging of the method provided by the invention are briefly analyzed. The real part and the imaginary part are both expressed by floating point numbers, the memory requirement of each real part is 4 bytes, and the memory requirement of each imaginary part is 8 bytes. BP pixel point amplitude matrix
Figure BDA0002694274030000101
Thus the memory requirement for direct computation of BP imaging is S08 MP. The memory requirements for non-uniform FFT BP imaging are: s18(M log M + | log | P) ≈ 8| log | P. The method only needs to store the coefficient for imaging
Figure BDA0002694274030000102
Then S2＝4P(1+4Ksp)+8Ksp(1+2M) ≈ 100P. Fig. 10 is a graph of memory requirements of three methods and changes of pixel point grid points, and comparison shows that under the same condition, the memory requirements of the method of the present invention are lower than those of the other two methods for the same pixel point.
The above is only a preferred embodiment of the present invention, and the protection scope of the present invention is not limited to the above-mentioned embodiments, and all technical solutions belonging to the idea of the present invention belong to the protection scope of the present invention. It should be noted that modifications and embellishments within the scope of the invention may be made by those skilled in the art without departing from the principle of the invention.

Claims (7)

1. A rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method is characterized by comprising the following steps:
step 1, establishing a through-wall radar detection scene model, arranging antenna array elements parallel to a wall body, transmitting electromagnetic waves, and uniformly sampling an obtained target echo signal to obtain echo signal data e;
step 2, gridding and dividing the pixels in the imaging area of the through-wall radar, projecting the echo signal data to each pixel grid in the imaging area to obtain an amplitude value I (x) of the BP imaging pixel of the through-wall radarp,yp)；
And step 3: constructing an expression of non-uniform FFT, converting the uniform sampling data into the non-uniform sampling data, and pre-calculating and storing each coefficient in the expression of the non-uniform FFT;
step 4, converting the uniform sampling data into non-uniform sampling data, and comparing I (x)p,yp) And performing fast Gaussian gridding non-uniform FFT calculation, and finally accumulating to obtain a radar imaging graph.
2. The fast Gaussian gridding non-uniform FFT through-wall imaging radar BP method according to claim 1, wherein the through-wall radar detection scene model established in the step 1 is as follows:
the wall thickness of the front wall and the rear wall is d, and the relative dielectric constant iswThe detection target is an ideal electric conductor with the radius r, and is arranged between the two walls, and the vertical distance between the circle center and the front wall is equal to the vertical distance;
is parallel to the wall body and has a distance h from the wall body1The N antenna array elements transmit electromagnetic waves, and M times of sampling is carried out on target echo signals of the N antenna array elements to obtain an M multiplied by N dimensional matrix e:
e＝[z1,z2,…zN] (1)
in the formula, zn＝[zn(1),zn(2),…,zn(M)]T，n＝1,2,…N，zn(M), M is 1,2, … M represents the echo signal amplitude of the same antenna at different sampling points.
3. According to claim 2The fast Gaussian gridding non-uniform FFT through-wall imaging radar BP method is characterized in that in the step 2, pixel points in an imaging area of the through-wall radar are gridded and divided, the whole imaging area is divided into P pixel point grids, and z is divided into P pixel point gridsn(m) projecting the data to each pixel point grid of the imaging area to obtain a through-wall radar BP imaging pixel point amplitude value:
Figure FDA0002694274020000011
in the formula, I (x)p,yp) Complex amplitude values for the P (P ═ 1,2, …, P-1) th pixel grid point, fm＝f0+ M Δ f is the M (M ═ 1,2, … M) th operating frequency point, f0For the starting frequency of the transmitted signal, Δ f is the frequency interval, τpFor the two-way time delay between the imaging pixel grid point p and the nth (N ═ 1,2, … N) antenna:
τp＝2(l1+l3)/c+2l2/v (3)
in the formula (I), the compound is shown in the specification,
Figure FDA0002694274020000021
is the propagation velocity of the electromagnetic wave in the wall body, c is the propagation velocity of the electromagnetic wave in the air,wis the relative dielectric constant of the wall;
l1，l2，l3respectively showing the slant distance between the antenna array element for transmitting the electromagnetic wave and the front wall, the distance of the electromagnetic wave transmitted on the wall and the distance of the electromagnetic wave reaching a target after passing through the wall.
4. The fast Gaussian gridding non-uniform FFT through-wall imaging radar BP method according to claim 3, wherein the step 3 of constructing the expression of the non-uniform FFT, pre-calculating and storing each coefficient in the expression of the non-uniform FFT comprises the following steps:
step 3.1: for uniform sampling sequence
Figure FDA0002694274020000022
Defining a non-uniform sampling sequence xjDiscrete fourier transform of (d):
Figure FDA0002694274020000023
gτ(x) Is [0,2 π]The above one-dimensional periodic gaussian kernel function has the following expression:
Figure FDA0002694274020000024
in the formula, tau is a Gaussian kernel function parameter and determines the exponential decay rate of the Gaussian kernel function;
step 3.2: constructing non-uniform FFT expressions
Figure FDA0002694274020000025
Step 3.3: precomputing and storing f (x)j) To obtain the expected value f (x)j)。
5. The fast Gaussian gridding non-uniform FFT through-wall imaging radar BP method according to claim 4, wherein the step 3.2 comprises the following steps:
step 3.2.1: for Gaussian kernel function gτ(x) Fourier transform to obtain Gτ(k) Let Gτ(k) Deconvoluting with coefficients in the expression, i.e. with uniform data F (k), to remove the effect of Gaussian smoothing and obtain an auxiliary function F-τ(k)：
Figure FDA0002694274020000031
Step 3.2.2: to F-τ(k) F is obtained by fast discrete Fourier inverse transformation-τ(x)：
Figure FDA0002694274020000032
In the formula, Kr＝R×K，KrIs the oversampling grid number, R is the oversampling coefficient;
step 3.2.3: to f-τ(x) Convolution is carried out to realize smooth output, and an expected value f (x) is obtainedj)：
Figure FDA0002694274020000033
6. The fast Gaussian gridding non-uniform FFT through-wall imaging radar BP method of claim 5, wherein in step 3.3, only x is calculatedjMost recent KspA plurality of uniform grid points, f (x) is calculated and stored in advancej) The obtained expected value f (x) of each coefficient in (b)j) Comprises the following steps:
Figure FDA0002694274020000034
in the formula, -Ksp＜σ'＜Ksp；
E1＝eikζ，
Figure FDA0002694274020000035
KspThe number of grid points covered by the unilateral extendable coverage of the Gaussian kernel function is shown, and sigma' is the grid points of the bilateral extendable coverage of the Gaussian kernel function; let ζ 2 π σ/Kr(ζ≤xj) And xi is expressed in xjVery close uniform grid points.
7. The fast Gaussian gridding non-uniform FFT through-wall imaging radar BP method of claim 1, wherein in step 4, the pixel point amplitude value I (x) is obtainedp,yp) Is deformed into
Figure FDA0002694274020000036
Wherein
Figure FDA0002694274020000037
Figure FDA0002694274020000038
Figure FDA0002694274020000041
Amplitude I of each pixeln(xp,yp) Then corresponds to non-uniform sample point data f (x)j) Uniform target echo data Zn(m) corresponding to the uniform sampling points F (k), the product of the frequency interval and the two-way delay 2 π Δ f τpCorresponding to non-uniform data xj。
CN202011000874.9A 2020-09-22 2020-09-22 Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method Pending CN112147608A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011000874.9A CN112147608A (en) 2020-09-22 2020-09-22 Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011000874.9A CN112147608A (en) 2020-09-22 2020-09-22 Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method

Publications (1)

Publication Number Publication Date
CN112147608A true CN112147608A (en) 2020-12-29

Family

ID=73892682

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011000874.9A Pending CN112147608A (en) 2020-09-22 2020-09-22 Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method

Country Status (1)

Country Link
CN (1) CN112147608A (en)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113238221A (en) * 2021-05-11 2021-08-10 沈阳航空航天大学 MIMO through-wall radar imaging method based on two-dimensional minimum phase coherence factor
CN114063176A (en) * 2021-12-09 2022-02-18 中国石油大学(华东) Radar imaging method, device and computer readable storage medium
CN114463638A (en) * 2022-02-23 2022-05-10 中国科学院空天信息创新研究院 Geometric correction method for airborne interferometric synthetic aperture radar image

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103364768A (en) * 2012-03-31 2013-10-23 中国科学院电子学研究所 Compressive sensing radar reconstruction method
CN107462887A (en) * 2017-07-07 2017-12-12 清华大学 Wide cut satellite-borne synthetic aperture radar imaging method based on compressed sensing
CN109884627A (en) * 2019-02-22 2019-06-14 中国人民解放军军事科学院国防科技创新研究院 The short range millimeter wave rapid three dimensional imaging process of any linear array configuration

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103364768A (en) * 2012-03-31 2013-10-23 中国科学院电子学研究所 Compressive sensing radar reconstruction method
CN107462887A (en) * 2017-07-07 2017-12-12 清华大学 Wide cut satellite-borne synthetic aperture radar imaging method based on compressed sensing
CN109884627A (en) * 2019-02-22 2019-06-14 中国人民解放军军事科学院国防科技创新研究院 The short range millimeter wave rapid three dimensional imaging process of any linear array configuration

Non-Patent Citations (6)

* Cited by examiner, † Cited by third party
Title
LESLIE: "Accelerating the Nonuniform Fast Fourier Transform", 《SIAM REVIEW》, pages 1 - 3 *
文智波: "FPGA与DSP在小波变换应用中的对比研究", 《仪表技术》 *
曹倩: "一种基于特征选择优化的图像检索方法", 《济南大学学报(自然科学版)》 *
狄培: "便携式穿墙雷达成像技术研究", 《中国优秀硕士学位论文全文数据库 信息科技辑》, pages 1 - 5 *
钱峰;胡光岷;: "基于滑动时窗的小波变换实时算法", 信号处理, no. 03 *
韩裕生;张延厚;王硕;姜兆祯;: "基于谱分析的分布式无源雷达非均匀成像", 雷达科学与技术, no. 02 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113238221A (en) * 2021-05-11 2021-08-10 沈阳航空航天大学 MIMO through-wall radar imaging method based on two-dimensional minimum phase coherence factor
CN113238221B (en) * 2021-05-11 2023-09-22 沈阳航空航天大学 MIMO through-wall radar imaging method based on two-dimensional minimum phase coherence factor
CN114063176A (en) * 2021-12-09 2022-02-18 中国石油大学(华东) Radar imaging method, device and computer readable storage medium
CN114063176B (en) * 2021-12-09 2023-09-15 中国石油大学(华东) Radar imaging method, radar imaging device and computer readable storage medium
CN114463638A (en) * 2022-02-23 2022-05-10 中国科学院空天信息创新研究院 Geometric correction method for airborne interferometric synthetic aperture radar image

Similar Documents

Publication Publication Date Title
CN110275166B (en) ADMM-based rapid sparse aperture ISAR self-focusing and imaging method
CN107193003B (en) Sparse singular value decomposition scanning radar foresight imaging method
CN110244303B (en) SBL-ADMM-based sparse aperture ISAR imaging method
CN112147608A (en) Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method
CN110068805B (en) High-speed target HRRP reconstruction method based on variational Bayesian inference
CN106680815B (en) MIMO radar imaging method based on tensor sparse representation
CN111551928B (en) Wall-penetrating radar imaging method based on low-rank sparse constraint of wall
WO2023045431A1 (en) Inverse synthetic aperture radar imaging method and apparatus, electronic device and storage medium
CN109507666B (en) ISAR sparse band imaging method based on off-network variational Bayesian algorithm
Moses et al. An autoregressive formulation for SAR backprojection imaging
CN107390215A (en) A kind of high speed super-resolution MIMO array imaging method
CN112099007B (en) Azimuth multi-channel SAR fuzzy suppression method suitable for non-ideal antenna directional diagram
CN111856465A (en) Forward-looking sea surface target angle super-resolution method based on sparse constraint
CN115453528A (en) Method and device for realizing segmented observation ISAR high-resolution imaging based on rapid SBL algorithm
US20130216151A1 (en) Systems and methods for image sharpening
CN113608218A (en) Frequency domain interference phase sparse reconstruction method based on back projection principle
CN111948652B (en) SAR intelligent parameterized super-resolution imaging method based on deep learning
CN111722227B (en) Bunching SAR compressed sensing imaging method based on approximate observation matrix
CN114779247B (en) Strabismus imaging method of spaceborne sparse synthetic aperture radar based on frequency scaling
CN112946644B (en) Based on minimizing the convolution weight l1Norm sparse aperture ISAR imaging method
CN116027293A (en) Rapid sparse angle super-resolution method for scanning radar
CN115453523A (en) Scanning radar sparse target batch processing super-resolution method
CN115902886A (en) ITR-DTV radar forward-looking super-resolution imaging method based on Renyi entropy
CN114966687A (en) Sparse ISAR imaging method and system based on low rank and non-local self-similarity
CN114910905A (en) GEO satellite-machine bistatic SAR moving target intelligent imaging method under similarity constraint

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination