CN111665500B - Pulse through-wall radar imaging method based on single-bit compressed sensing - Google Patents

Pulse through-wall radar imaging method based on single-bit compressed sensing Download PDF

Info

Publication number
CN111665500B
CN111665500B CN202010533605.2A CN202010533605A CN111665500B CN 111665500 B CN111665500 B CN 111665500B CN 202010533605 A CN202010533605 A CN 202010533605A CN 111665500 B CN111665500 B CN 111665500B
Authority
CN
China
Prior art keywords
pulse
matrix
channel
layer iteration
imaging
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.)
Active
Application number
CN202010533605.2A
Other languages
Chinese (zh)
Other versions
CN111665500A (en
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.)
Shenyang Aerospace University
Original Assignee
Shenyang Aerospace University
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 Shenyang Aerospace University filed Critical Shenyang Aerospace University
Priority to CN202010533605.2A priority Critical patent/CN111665500B/en
Publication of CN111665500A publication Critical patent/CN111665500A/en
Application granted granted Critical
Publication of CN111665500B publication Critical patent/CN111665500B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

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/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • 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
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/28—Details of pulse systems
    • G01S7/285—Receivers
    • G01S7/292—Extracting wanted echo-signals
    • G01S7/2923—Extracting wanted echo-signals based on data belonging to a number of consecutive radar periods
    • G01S7/2927—Extracting wanted echo-signals based on data belonging to a number of consecutive radar periods by deriving and controlling a threshold value
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • General Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

The invention provides a pulse through-wall radar imaging method based on single bit compressed sensing, which belongs to the technical field of radar imaging and comprises the following steps: step 1: establishing a full-sampling echo signal representation of an IQ channel of the pulse through-wall radar; step 2: space grid division is carried out on the target imaging area, and a dictionary matrix is generated; step 3: generating a measurement matrix to obtain down-sampled IQ channel single-bit compressed sensing measurement signals; step 4: generating an observation matrix according to the dictionary matrix and the measurement matrix; step 5: image reconstruction is performed using an enhanced binary iterative hard threshold algorithm. The imaging method provided by the invention can randomly downsample the time-domain echo signals and quantize the time-domain echo signals into single-bit data, so that the data is more convenient to transmit and process, the pixels of a target area are more concentrated, the background clutter is less, the imaging result quality is higher, the performance is more excellent, and the imaging method is more suitable for pulse through-wall radar imaging in an actual scene.

Description

Pulse through-wall radar imaging method based on single-bit compressed sensing
Technical Field
The invention belongs to the technical field of radar imaging, and particularly relates to a pulse through-wall radar imaging method based on single-bit compressed sensing.
Background
The pulse through-wall radar is one radar system for directly collecting time domain data, and consists of pulse source, transceiver antenna and data collecting part. The pulse through-wall radar detects a target area through a wall body by time domain pulses emitted by an antenna. The pulse system-based through-wall radar has the characteristics of low hardware implementation cost, strong penetrating capacity and the like, but the pulse through-wall radar is extremely large in collected data volume while pursuing high resolution and high measurement accuracy, so that the data storage and transmission volume are increased, the data quantization is difficult, and a heavy burden is caused on hardware. However, with the rise of the compressed sensing theory, the above problems have been solved properly.
Compressed sensing theory is an emerging signal processing technology, and development in the last decade is fully applied to various fields. Compressed sensing can recover a signal with a much smaller amount of data than the nyquist sampling law. The single bit compressed sensing is further used for quantizing the measurement data into single bits on the basis of compressed sensing, so that the transmission and processing of the data become convenient, and meanwhile, the single bit compressed sensing has better robustness to noise. The hardware aspect is degraded by a complex quantizer into a simple comparator, so that the quantization process becomes simple, and the burden of the hardware is greatly reduced. The single-bit compressed sensing pulse through-wall radar imaging method directly utilizes downsampled single-bit data to realize accurate imaging and positioning of a target.
The binary iterative hard threshold algorithm is widely applied to the field of single-bit compressed sensing radar imaging because of the advantages of simple calculation, high reconstruction quality and the like. However, the traditional single-bit compressed sensing radar imaging method only considers the sparsity of the scene and ignores the structural sparsity of the radar signal.
Disclosure of Invention
Aiming at the problems existing in the prior art, the invention provides a pulse through-wall radar imaging method based on single-bit compressed sensing, which utilizes an enhanced binary iterative hard threshold algorithm to realize accurate imaging of the pulse through-wall radar on a block target in a detection scene. The enhanced binary iterative hard threshold algorithm is a derivative algorithm of the binary iterative hard threshold algorithm, fully applies double-layer sparsity in radar imaging, and specifically comprises the following steps: (1) When a target is detected, the I channel and the Q channel of the pulse through-wall radar system receiver correspond to the same scene reflection coefficient, namely the IQ channel has joint sparsity. (2) In an actual scene, the detection targets are generally clustered instead of punctiform, i.e. the image to be reconstructed has block sparsity.
In order to achieve the above object, the present invention adopts the following technical scheme:
a pulse through-wall radar imaging method based on single bit compressed sensing comprises the following steps:
step 1: and establishing a full-sampling echo signal expression of the pulse through-wall radar IQ channel.
Assume pulse through-wall radarM antenna positions are provided, and sampling is carried out on the IQ channel for N times for each antenna position; in the pulse duration period, the sampling start time is t 1 The sampling interval is Δt, and the N (n=1, 2,) th IQ channel sampling signal of the M (m=1, 2,) th antenna position is shown in formula (1):
Figure BDA0002536282980000021
wherein P is the target number, sigma p Reflection coefficient for P (p=1, 2,., P) target, t n Is the nth sampling time and t n ＝t 1 +(n-1)Δt，τ mp Delay omega of the p-th target corresponding to the m-th antenna position c ＝2πf c For the carrier angular frequency, f c Is the carrier frequency, S (t n ) A gaussian pulse that is first-order differential, as shown in equation (2):
Figure BDA0002536282980000022
wherein,,
Figure BDA0002536282980000023
χ＝1/f 0 ，f 0 is the center frequency of the first order differential gaussian pulse.
The m-th antenna position collects N multiplied by 1 dimension time domain echo signal beta m ＝[β m (t 1 ),β m (t 2 ),...,β m (t N )] T Stacking the time domain echo signals acquired by all antenna positions, and then the MN multiplied by 1 dimension full sampling echo signal beta of the pulse through-wall radar IQ channel is shown as formula (3):
Figure BDA0002536282980000024
step 2: and carrying out space grid division on the target imaging area to generate a dictionary matrix.
The target imaging region is divided into H×L space grids, and the HL×1-dimensional reflection coefficient vector sigma is obtained through column stacking operation.
Dictionary matrix ψ corresponding to mth antenna position m N. row i (i=1, 2,) column elements are expressed as:
Figure BDA0002536282980000025
wherein τ mi Representing the two-way transmission delay of the ith spatial grid corresponding to the mth antenna position.
The n×1-dimensional time-domain echo signal collected at the mth antenna position can be expressed as:
β m ＝Ψ m σ (5)
as M antenna positions are used, MN X HL-dimensional dictionary matrix can be obtained
Figure BDA0002536282980000026
The relationship between the pulse through-wall radar IQ channel full-sampling echo signal β and the dictionary matrix ψ is:
β＝Ψσ (6)
step 3: and generating a measurement matrix to obtain a down-sampled IQ channel single-bit compressed sensing measurement signal.
N×1-dimensional time-domain echo signal β collected for mth antenna position m Generating a corresponding measurement matrix Φ m ，Φ m Is a J×N Gaussian random matrix, where J is a measurement matrix Φ m Measured number of phi m Is subject to a gaussian distribution with a mean value of 0; for the fully sampled echo signal beta of MN multiplied by 1, a MJ multiplied by MN dimension measurement matrix phi is generated, wherein phi is represented by phi 1 ,Φ 2 ,...,Φ M Is a diagonal array of diagonal elements, as shown in formula (7):
Φ＝diag[Φ 1 ,Φ 2 ,…,Φ M ] (7)
performing single-bit compressed sensing measurement on an IQ channel of a pulse through-wall system to obtain downsamplingIQ channel single bit compressed sensing measurement signal
Figure BDA0002536282980000031
The down-sampled IQ channel single bit compressed sensing measurement signal>
Figure BDA0002536282980000032
As shown in formula (8):
Figure BDA0002536282980000033
wherein beta is I As the real part of beta, beta Q An imaginary part of β; sign is a sign function, and if the given data is positive, the sign output is 1, and if the given data is negative, the sign output is-1.
Step 4: and generating an observation matrix according to the dictionary matrix and the measurement matrix.
The observation matrix Θ is shown as formula (9):
Θ＝ΦΨ (9)
finally obtaining the observation matrix theta with MJ×HL dimensions.
After the observation matrix Θ is obtained, the equation (8) is further deformed as shown in the equation (10):
Figure BDA0002536282980000034
wherein,,
Figure BDA0002536282980000035
Θ I ＝ΦΨ I ，Θ Q ＝ΦΨ Q ，Ψ I 、Ψ Q the real part and the imaginary part of the dictionary matrix ψ are respectively, and the reflection coefficient vectors corresponding to the I channel and the Q channel are the same, namely the reflection coefficient vectors corresponding to the I channel and the Q channel
Figure BDA0002536282980000036
Step 5: image reconstruction is performed using an enhanced binary iterative hard threshold algorithm. The method comprises the following specific steps:
and (3) outer layer iteration input:
Figure BDA0002536282980000037
sparsity K, gradient descent step μ, maximum number of iterations t max Iterative accuracy epsilon.
Outer layer iteration initialization:
Figure BDA0002536282980000038
the number of iterations t=0.
Outer layer iteration:
(1) Updating the outer layer iteration times: t=t+1.
(2) And (3) performing outer layer iteration gradient descent calculation according to the formula (11) to obtain an outer layer iteration intermediate vector a.
Figure BDA0002536282980000041
(3) And calculating inner layer iteration parameters gamma and lambda.
Let a be ii For the ii-th element in vector a, pair sets
Figure BDA0002536282980000042
(ii=1, 2, HL) are arranged in descending order, γ being the K-th element in the set, λ=0.4γ 2 。
(4) Performing inner layer iterative computation
Inner layer iteration input: a, gamma, lambda, sparsity K, gradient descent step size
Figure BDA0002536282980000043
Maximum number of iterations->
Figure BDA0002536282980000044
Iterative accuracy->
Figure BDA0002536282980000045
Inner layer iterationInitializing: b 0 Number of iterations =a
Figure BDA0002536282980000046
Inner layer iteration:
(1) updating the number of inner layer iterations
Figure BDA0002536282980000047
Introduction of intermediate variable +.>
Figure BDA0002536282980000048
(2) Performing an inner layer iterative gradient descent calculation according to equation (12):
Figure BDA0002536282980000049
wherein,,
Figure BDA00025362829800000410
represents->
Figure BDA00025362829800000411
The j-th element of (2)>
Figure BDA00025362829800000412
The descending size of the iteration gradient of the inner layer is controlled,
Figure BDA00025362829800000413
as shown in formulas (13), (14):
if j is less than or equal to HL,
Figure BDA00025362829800000414
if 2HL is more than or equal to j and is more than HL,
Figure BDA00025362829800000415
wherein a is j (j=1, 2,., 2 HL) represents the j-th element in vector a, N j Represents the set of adjacent element points corresponding to the jth element, j' represents N j Any element in the collection; the function ζ (x) is shown in equation (15):
Figure BDA0002536282980000051
wherein ( * Representing a conjugation operation.
(3) The greedy selection is adopted to ensure joint sparsity, and the calculation method is as shown in formula (16):
Figure BDA0002536282980000052
if j is less than or equal to HL,
Figure BDA0002536282980000053
as shown in the formulas (17), (18):
Figure BDA0002536282980000054
Figure BDA0002536282980000055
if 2HL is more than or equal to j and is more than HL,
Figure BDA0002536282980000056
as shown in the formulas (19), (20):
Figure BDA0002536282980000057
Figure BDA0002536282980000058
wherein ρ is a set
Figure BDA0002536282980000059
(ii = 1,2,) HL, K-th element after descending order.
(4) And judging whether the inner layer iteration is continued or not.
If it is
Figure BDA00025362829800000510
Or->
Figure BDA00025362829800000511
The inner layer iteration is terminated and output +.>
Figure BDA00025362829800000512
Otherwise, returning to (1) and continuing the inner layer iteration.
(5) Will be
Figure BDA00025362829800000513
Assign->
Figure BDA00025362829800000514
Judging whether the outer layer iteration termination condition is met.
If t is greater than or equal to t max Or (b)
Figure BDA00025362829800000515
The outer layer iteration is terminated and the reconstructed reflection coefficient vector +.>
Figure BDA00025362829800000516
Namely, a reflection coefficient vector corresponding to the region to be imaged; otherwise, returning to (1) and continuing the outer layer iteration.
(6) Reflection coefficient vector to be reconstructed
Figure BDA00025362829800000517
The former HL elements are rearranged into an H multiplied by L dimensional matrix through column stacking operation, and then a single bit compressed sensing pulse through-wall radar imaging result can be obtained.
The invention has the beneficial effects that:
according to the technical scheme, the invention provides a pulse through-wall radar imaging method based on single-bit compressed sensing. And reconstructing the reflection coefficient of the detection scene by using an enhanced binary iterative hard threshold algorithm according to the sparse characteristics of the double-layer structure of the radar imaging target. The imaging method provided by the invention can randomly downsample the time-domain echo signals and quantize the time-domain echo signals into single-bit data, so that the data is more convenient to transmit and process, the pixels of a target area are more concentrated, the background clutter is less, the imaging result quality is higher, the performance is more excellent, and the imaging method is more suitable for pulse through-wall radar imaging in an actual scene by adopting the double-layer structure sparse feature to be applied to the imaging reconstruction process.
Drawings
Fig. 1 is a schematic flow chart of a pulse through-wall radar imaging method based on single bit compressed sensing according to an embodiment of the present invention;
FIG. 2 is a schematic flow chart of an outer layer iteration of an enhanced binary iterative hard threshold algorithm according to an embodiment of the present invention;
FIG. 3 is a schematic flow chart of an inner layer iteration of an enhanced binary iterative hard threshold algorithm according to an embodiment of the present invention;
fig. 4 is a schematic diagram of a pulse through-wall radar simulation scene based on single bit compressed sensing according to an embodiment of the present invention;
FIG. 5 is a graph of the real reflection coefficient vector imaging result of a pulse through-wall radar imaging method based on single bit compressed sensing according to an embodiment of the present invention;
FIG. 6 is a view of the result of reconstructing reflectance vector imaging of a pulse through-the-wall radar imaging method based on single bit compressed sensing according to an embodiment of the present invention;
fig. 7 is a graph showing the comparison of the reconstructed SNR of the method of the present invention with a conventional single-layer sparse imaging method at different data volumes in an embodiment of the present invention.
Fig. 8 is a graph comparing reconstructed SNR of the method of the present invention with a common single-layer sparse imaging method under different echo signal SNRs in an embodiment of the present invention.
Detailed Description
The following detailed description of the invention proceeds with reference to the accompanying drawings and examples. The following examples are illustrative of the invention and are not intended to limit the scope of the invention. As shown in fig. 1, the method of this embodiment is as follows. A pulse through-wall radar imaging method based on single bit compressed sensing comprises the following steps:
step 1: and establishing a full-sampling echo signal expression of the pulse through-wall radar IQ channel.
Assuming that the pulse through-wall radar has M antenna positions, sampling the IQ channel for N times at each antenna position; in the pulse duration period, the sampling start time is t 1 The sampling interval is Δt, and the N (n=1, 2,) th IQ channel sampling signal of the M (m=1, 2,) th antenna position is shown in formula (1):
Figure BDA0002536282980000061
wherein P is the target number, sigma p Reflection coefficient for P (p=1, 2,., P) target, t n Is the nth sampling time and t n ＝t 1 +(n-1)Δt，τ mp Delay omega of the p-th target corresponding to the m-th antenna position c ＝2πf c For the carrier angular frequency, f c Is the carrier frequency, S (t n ) A gaussian pulse that is first-order differential, as shown in equation (2):
Figure BDA0002536282980000062
wherein,,
Figure BDA0002536282980000063
χ＝1/f 0 ，f 0 is the center frequency of the first order differential gaussian pulse.
The m-th antenna position collects N multiplied by 1 dimension time domain echo signal beta m ＝[β m (t 1 ),β m (t 2 ),...,β m (t N )] T Time domain echo signals acquired by all antenna positionsStacking the numbers, and then, the mn×1 dimension full-sampling echo signal β of the pulse through-wall radar IQ channel is shown in formula (3):
Figure BDA0002536282980000071
step 2: and carrying out space grid division on the target imaging area to generate a dictionary matrix.
The target imaging region is divided into H×L space grids, and the HL×1-dimensional reflection coefficient vector sigma is obtained through column stacking operation.
Dictionary matrix ψ corresponding to mth antenna position m N. row i (i=1, 2,) column elements are expressed as:
Figure BDA0002536282980000072
wherein τ mi Representing the two-way transmission delay of the ith spatial grid corresponding to the mth antenna position.
The n×1-dimensional time-domain echo signal collected at the mth antenna position can be expressed as:
β m ＝Ψ m σ (5)
as M antenna positions are used, MN X HL-dimensional dictionary matrix can be obtained
Figure BDA0002536282980000073
The relationship between the pulse through-wall radar IQ channel full-sampling echo signal β and the dictionary matrix ψ is:
β＝Ψσ (6)
step 3: and generating a measurement matrix to obtain a down-sampled IQ channel single-bit compressed sensing measurement signal.
N×1-dimensional time-domain echo signal β collected for mth antenna position m Generating a corresponding measurement matrix Φ m ，Φ m Is a J×N Gaussian random matrix, where J is a measurement matrix Φ m Measured number of phi m Each element of (2)The diathesis is distributed from a Gaussian with the mean value of 0; for the fully sampled echo signal beta of MN multiplied by 1, a MJ multiplied by MN dimension measurement matrix phi is generated, wherein phi is represented by phi 1 ,Φ 2 ,...,Φ M Is a diagonal array of diagonal elements, as shown in formula (7):
Φ＝diag[Φ 1 ,Φ 2 ,…,Φ M ] (7)
performing single-bit compressed sensing measurement on an IQ channel of the pulse through-wall system to obtain a down-sampled IQ channel single-bit compressed sensing measurement signal
Figure BDA0002536282980000074
The down-sampled IQ channel single bit compressed sensing measurement signal>
Figure BDA0002536282980000075
As shown in formula (8):
Figure BDA0002536282980000076
wherein beta is I As the real part of beta, beta Q An imaginary part of β; sign is a sign function, and if the given data is positive, the sign output is 1, and if the given data is negative, the sign output is-1.
Step 4: and generating an observation matrix according to the dictionary matrix and the measurement matrix.
The observation matrix Θ is shown as formula (9):
Θ＝ΦΨ (9)
finally obtaining the observation matrix theta with MJ×HL dimensions.
After the observation matrix Θ is obtained, the equation (8) is further deformed as shown in the equation (10):
Figure BDA0002536282980000081
wherein,,
Figure BDA0002536282980000082
Θ I ＝ΦΨ I ，Θ Q ＝ΦΨ Q ，Ψ I 、Ψ Q the real part and the imaginary part of the dictionary matrix ψ are respectively, and the reflection coefficient vectors corresponding to the I channel and the Q channel are the same, namely the reflection coefficient vectors corresponding to the I channel and the Q channel
Figure BDA0002536282980000083
Step 5: image reconstruction is performed using an enhanced binary iterative hard threshold algorithm. The method comprises the following specific steps:
the outer layer iterative calculation, the flow is as shown in figure 2:
and (3) outer layer iteration input:
Figure BDA0002536282980000084
sparsity K, gradient descent step μ, maximum number of iterations t max Iterative accuracy epsilon.
Outer layer iteration initialization:
Figure BDA0002536282980000085
the number of iterations t=0.
Outer layer iteration:
(1) Updating the outer layer iteration times: t=t+1.
(2) And (3) performing outer layer iteration gradient descent calculation according to the formula (11) to obtain an outer layer iteration intermediate vector a.
Figure BDA0002536282980000086
(3) And calculating inner layer iteration parameters gamma and lambda.
Let a be ii For the ii-th element in vector a, pair sets
Figure BDA0002536282980000087
(ii=1, 2, HL) are arranged in descending order, γ being the K-th element in the set, λ=0.4γ 2 。
(4) And (3) performing inner layer iterative computation, wherein the flow is as shown in fig. 3:
inner layer iteration input: a, gamma, lambda, sparsity K, gradient descent step size
Figure BDA0002536282980000088
Maximum number of iterations->
Figure BDA0002536282980000089
Iterative accuracy->
Figure BDA00025362829800000810
Inner layer iteration initialization: b 0 Number of iterations =a
Figure BDA00025362829800000811
Inner layer iteration:
(1) updating the number of inner layer iterations
Figure BDA0002536282980000091
Introduction of intermediate variable +.>
Figure BDA0002536282980000092
(2) Performing an inner layer iterative gradient descent calculation according to equation (12):
Figure BDA0002536282980000093
wherein,,
Figure BDA0002536282980000094
represents->
Figure BDA0002536282980000095
The j-th element of (2)>
Figure BDA0002536282980000096
The descending size of the iteration gradient of the inner layer is controlled,
Figure BDA0002536282980000097
as shown in formulas (13), (14):
if j is less than or equal to HL,
Figure BDA0002536282980000098
if 2HL is more than or equal to j and is more than HL,
Figure BDA0002536282980000099
wherein a is j (j=1, 2,., 2 HL) represents the j-th element in vector a, N j Represents the set of adjacent element points corresponding to the jth element, j' represents N j Any element in the collection; the function ζ (x) is shown in equation (15):
Figure BDA00025362829800000910
wherein ( * Representing a conjugation operation.
(3) The greedy selection is adopted to ensure joint sparsity, and the calculation method is as shown in formula (16):
Figure BDA00025362829800000911
if j is less than or equal to HL,
Figure BDA00025362829800000912
as shown in (17), (18)
Figure BDA00025362829800000913
Figure BDA00025362829800000914
If 2HL is more than or equal to j and is more than HL,
Figure BDA00025362829800000915
as shown in (19), (20)
Figure BDA00025362829800000916
Figure BDA00025362829800000917
Wherein ρ is a set
Figure BDA0002536282980000101
(ii = 1,2,) HL, K-th element after descending order.
(4) And judging whether the inner layer iteration is continued or not.
If it is
Figure BDA0002536282980000102
Or->
Figure BDA0002536282980000103
The inner layer iteration is terminated and output +.>
Figure BDA0002536282980000104
Otherwise, returning to (1) and continuing the inner layer iteration.
(5) Will be
Figure BDA0002536282980000105
Assign->
Figure BDA0002536282980000106
Judging whether the outer layer iteration termination condition is met.
If t is greater than or equal to t max Or (b)
Figure BDA0002536282980000107
The outer layer iteration is terminated and the reconstructed reflection coefficient vector +.>
Figure BDA0002536282980000108
Namely, a reflection coefficient vector corresponding to the region to be imaged; otherwise, returning to (1) and continuing the outer layer iteration.
(6) Reflection coefficient vector to be reconstructed
Figure BDA0002536282980000109
The former HL elements are rearranged into an H multiplied by L dimensional matrix through column stacking operation, and then a single bit compressed sensing pulse through-wall radar imaging result can be obtained.
The simulation computing environment is configured to be an internal memory of Intel (R) Core (TM) i5-4210 CPU@2.60GHz,4GB, a 64-bit Microsoft Windows 7 system, and simulation software adopts MATLAB (R2017 a) to carry out simulation experiments.
As shown in fig. 4, a coordinate system is established by taking the upper left edge point of the wall body as a coordinate origin, the coordinate system extends from the origin to the wall body direction, and an extension line is drawn to be an x-axis, namely an azimuth direction, by clinging to the surface of the wall body facing the receiving and transmitting antenna; a straight line extending perpendicular to the x-axis imaging region is drawn from the origin as the y-axis, i.e., distance-wise.
The modulated first-order differential Gaussian pulse is adopted as pulse excitation, the center frequency is 1200MHz, and the carrier frequency is 2400MHz. The thickness of the wall body is 0.2m, the relative dielectric constant is 4, the system adopts a single-station measurement mode, the receiving and transmitting antenna has 41 antenna positions in total, the y-direction (distance direction) position of the receiving and transmitting antenna is same as-2 m, the x-direction (azimuth direction) position of the receiving and transmitting antenna takes 2.3m as the initial position, and the position interval is 0.06m. The imaging area size was 4.125m×4.125m, divided into 33×33 grids. Assume that there are two block targets in the detection scene, one of which consists of 6 point targets and the other consists of 5 point targets, and the 11 point target positions are (3.5,2.2) m, (3.5,2.075) m, (3.375,2.2) m, (3.5,2.325) m, (3.625,2.2) m, (2.5,0.7) m, (2.625,0.7) m, (2.75,0.7) m, (2.5,0.825) m, (2.625,0.825) m and (2.75,0.825) m, respectively. The reflection coefficient of the target is assumed to follow a gaussian distribution. An additive gaussian white noise is added to the received echo signal before single bit quantization, and the signal-to-noise ratio (SNR) of the echo signal is 30dB. In the imaging process, the number of time domain sampling points of each antenna position is n=101, the data of each antenna position is downsampled by using a measurement matrix with the measurement number of 20, and finally, 41×20=820 bits of data is used for imaging. In this embodiment, fig. 5 shows an imaging result of the real reflection coefficient of the detection scene, and fig. 6 shows an imaging result of the reconstructed reflection coefficient vector of the imaging method provided by the invention, from which it can be seen that the target position is accurately reconstructed. In order to further prove the superiority of the imaging method provided by the invention, the reconstruction SNR is taken as an index, and the binary iterative hard threshold algorithm (traditional single bit compressed sensing imaging) is compared with the enhanced binary iterative hard threshold algorithm (single bit compressed sensing imaging with double-layer sparsity) provided by the invention in a quantization way, wherein the reconstruction SNR is shown as a formula (21):
Figure BDA0002536282980000111
wherein,,
Figure BDA0002536282980000112
representing the true reflection coefficient vector corresponding to the I channel and the Q channel in the implementation process of the embodiment +.>
Figure BDA0002536282980000113
Representing the reconstructed reflection coefficient vector.
Fig. 7 shows a graph of reconstructed SNR versus data volume for two imaging methods when the echo signal SNR is 20 dB. The reconstructed SNR for the two imaging methods when the data are 615bit, 820bit, 1025bit and 1230bit, respectively, is shown. When the data volume is increased, the imaging quality of both imaging methods is improved, but it can be obviously seen that the imaging quality of the imaging method provided by the invention is better.
Fig. 8 shows a graph of reconstructed SNR as a function of echo signal SNR for two imaging methods. As shown, the reconstructed SNR for the two imaging methods is given when the echo signal SNR is 10dB, 15dB, 20dB, 25dB, 30dB, respectively. According to the image, the performance of both imaging methods is improved along with the increase of the SNR of the echo signals, but the reconstruction SNR of the imaging result of the enhanced binary iterative hard threshold algorithm provided by the invention is higher, and the reconstruction image quality is more excellent.
Experimental results show that the imaging reconstruction time of the binary iterative hard threshold algorithm is 10 seconds, the imaging reconstruction time of the enhanced binary iterative hard threshold algorithm provided by the invention is 14 seconds, the reconstruction time is not greatly different, but the imaging quality of the imaging method provided by the invention is higher. According to the results, the conclusion can be drawn that the pulse through-wall radar imaging method based on single bit compressed sensing has higher reconstruction quality and more excellent performance in the imaging reconstruction time with little difference.

Claims (1)

1. A pulse through-wall radar imaging method based on single bit compressed sensing is characterized by comprising the following steps:
step 1: establishing a full-sampling echo signal representation of an IQ channel of the pulse through-wall radar;
step 2: space grid division is carried out on the target imaging area, and a dictionary matrix is generated;
step 3: generating a measurement matrix to obtain down-sampled IQ channel single-bit compressed sensing measurement signals;
step 4: generating an observation matrix according to the dictionary matrix and the measurement matrix;
step 5: using an enhanced binary iterative hard threshold algorithm to reconstruct imaging;
the specific method of the step 1 is as follows:
assuming that the pulse through-wall radar has M antenna positions, sampling the IQ channel for N times at each antenna position; in the pulse duration period, the sampling start time is t 1 The sampling interval is Δt, and the sample signal of the nth IQ channel at the mth antenna position is shown in formula (1):
Figure FDA0004278830910000011
wherein P is the target number, sigma p The reflection coefficient of the p-th target, t n Is the nth sampling time and t n ＝t 1 +(n-1)Δt，τ mp Delay omega of the p-th target corresponding to the m-th antenna position c ＝2πf c For the carrier angular frequency, f c For carrier frequencies, m=1, 2, M, n=1, 2, N, p=1, 2, P, S (t n ) A gaussian pulse that is first-order differential, as shown in equation (2):
Figure FDA0004278830910000012
wherein,,
Figure FDA0004278830910000013
χ＝1/f 0 ，f 0 the center frequency of the first-order differential Gaussian pulse;
the m-th antenna position collects N multiplied by 1 dimension time domain echo signal beta m ＝[β m (t 1 ),β m (t 2 ),,β m (t N )] T Stacking the time domain echo signals acquired by all antenna positions, and then the MN multiplied by 1 dimension full sampling echo signal beta of the pulse through-wall radar IQ channel is shown as formula (3):
Figure FDA0004278830910000014
the specific method of the step 2 is as follows:
dividing a target imaging region into H multiplied by L space grids, and obtaining an HL multiplied by 1-dimensional reflection coefficient vector sigma through column stacking operation;
dictionary matrix ψ corresponding to mth antenna position m The nth row and ith column elements of (a) are expressed as:
Figure FDA0004278830910000015
where n=1, 2, N, i=1, 2, HL, τ mi Double pass transmission representing the ith spatial grid corresponding to the mth antenna positionTime delay is carried out;
the n×1-dimensional time-domain echo signal collected at the mth antenna position can be expressed as:
β m ＝Ψ m σ (5)
as M antenna positions are used, MN X HL-dimensional dictionary matrix can be obtained
Figure FDA0004278830910000021
The relationship between the pulse through-wall radar IQ channel full-sampling echo signal β and the dictionary matrix ψ is:
β＝Ψσ (6)
the specific method of the step 3 is as follows:
n×1-dimensional time-domain echo signal β collected for mth antenna position m Generating a corresponding measurement matrix Φ m ，Φ m Is a J×N Gaussian random matrix, where J is a measurement matrix Φ m Measured number of phi m Is subject to a gaussian distribution with a mean value of 0; for the fully sampled echo signal beta of MN multiplied by 1, a MJ multiplied by MN dimension measurement matrix phi is generated, wherein phi is represented by phi 1 ,Φ 2 ,...,Φ M Is a diagonal array of diagonal elements, as shown in formula (7):
Φ＝diag[Φ 1 ，Φ 2 ，…，Φ M ] (7)
performing single-bit compressed sensing measurement on an IQ channel of the pulse through-wall system to obtain a down-sampled IQ channel single-bit compressed sensing measurement signal
Figure FDA0004278830910000022
The down-sampled IQ channel single bit compressed sensing measurement signal>
Figure FDA0004278830910000023
As shown in formula (8):
Figure FDA0004278830910000024
wherein beta is I As the real part of beta, beta Q An imaginary part of β; sign is a sign function, if the given data is positive, the sign output is 1, and if the given data is negative, the sign output is-1;
the specific method of the step 4 is as follows:
the observation matrix Θ is shown as formula (9):
Θ＝ΦΨ (9)
finally obtaining an observation matrix theta with MJ×HL dimensions;
after the observation matrix Θ is obtained, the equation (8) is further deformed as shown in the equation (10):
Figure FDA0004278830910000025
wherein,,
Figure FDA0004278830910000026
Θ I ＝ΦΨ I ，Θ Q ＝ΦΨ Q ，Ψ I 、Ψ Q the real part and the imaginary part of the dictionary matrix ψ are respectively, and the reflection coefficient vectors corresponding to the I channel and the Q channel are the same, namely the reflection coefficient vectors corresponding to the I channel and the Q channel +.>
Figure FDA0004278830910000031
The specific method in the step 5 is as follows:
and (3) outer layer iteration input:
Figure FDA0004278830910000032
sparsity K, gradient descent step μ, maximum number of iterations t max Iteration precision epsilon;
outer layer iteration initialization:
Figure FDA0004278830910000033
the number of iterations t=0;
outer layer iteration:
(1) Updating the outer layer iteration times: t=t+1;
(2) Performing outer layer iteration gradient descent calculation according to the formula (11) to obtain an outer layer iteration intermediate vector a;
Figure FDA0004278830910000034
(3) Calculating inner layer iteration parameters gamma and lambda;
let a be ii For the ii-th element in vector a, pair sets
Figure FDA0004278830910000035
Arranged in descending order, γ is the K-th element in the set, ii=1, 2,..hl, λ=0.4γ 2 ；
(4) Performing inner layer iterative computation;
inner layer iteration input: a, gamma, lambda, sparsity K, gradient descent step size
Figure FDA0004278830910000036
Maximum number of iterations->
Figure FDA0004278830910000037
Iterative accuracy->
Figure FDA0004278830910000038
Inner layer iteration initialization: b 0 Number of iterations =a
Figure FDA0004278830910000039
Inner layer iteration:
(1) updating the number of inner layer iterations
Figure FDA00042788309100000310
Introduction of intermediate variable +.>
Figure FDA00042788309100000311
(2) Performing an inner layer iterative gradient descent calculation according to equation (12):
Figure FDA00042788309100000312
wherein,,
Figure FDA00042788309100000313
represents->
Figure FDA00042788309100000314
J=1, 2,..2 HL,/-j>
Figure FDA00042788309100000315
Controlling the gradient decrease of the inner layer iteration>
Figure FDA00042788309100000316
As shown in formulas (13), (14):
if j is less than or equal to HL,
Figure FDA00042788309100000317
if 2HL is more than or equal to j and is more than HL,
Figure FDA0004278830910000041
wherein a is j Represents the j-th element, N in the vector a j Represents a set of adjacent element points corresponding to the j-th element, j=1, 2,..2 hl, j' represents N j Any element in the collection; the function ζ (x) is shown in equation (15):
Figure FDA0004278830910000042
wherein ( * Represents a conjugation operation;
(3) the greedy selection is adopted to ensure joint sparsity, and the calculation method is as shown in formula (16):
Figure FDA0004278830910000043
if j is less than or equal to HL,
Figure FDA0004278830910000044
as shown in the formulas (17), (18):
Figure FDA0004278830910000045
Figure FDA0004278830910000046
if 2HL is more than or equal to j and is more than HL,
Figure FDA0004278830910000047
as shown in the formulas (19), (20):
Figure FDA0004278830910000048
Figure FDA0004278830910000049
wherein ρ is a set
Figure FDA00042788309100000410
The K-th element after descending order, ii=1, 2, HL;
(4) judging whether the inner layer iteration is continued or not;
if it is
Figure FDA00042788309100000411
Or->
Figure FDA00042788309100000412
The inner layer iteration is terminated and output +.>
Figure FDA00042788309100000413
Otherwise, returning to (1) continuing the inner layer iteration;
(5) Will be
Figure FDA00042788309100000414
Assign->
Figure FDA00042788309100000415
Judging whether an outer layer iteration termination condition is met;
if t is greater than or equal to t max Or (b)
Figure FDA00042788309100000416
The outer layer iteration is terminated and the reconstructed reflection coefficient vector +.>
Figure FDA00042788309100000417
Namely, a reflection coefficient vector corresponding to the region to be imaged; otherwise, returning to (1) continuing the outer layer iteration;
(6) The reconstructed reflection coefficient vector
Figure FDA00042788309100000418
The former HL elements are rearranged into an H multiplied by L dimensional matrix through column stacking operation, and then a single bit compressed sensing pulse through-wall radar imaging result can be obtained.
CN202010533605.2A 2020-06-12 2020-06-12 Pulse through-wall radar imaging method based on single-bit compressed sensing Active CN111665500B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010533605.2A CN111665500B (en) 2020-06-12 2020-06-12 Pulse through-wall radar imaging method based on single-bit compressed sensing

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010533605.2A CN111665500B (en) 2020-06-12 2020-06-12 Pulse through-wall radar imaging method based on single-bit compressed sensing

Publications (2)

Publication Number Publication Date
CN111665500A CN111665500A (en) 2020-09-15
CN111665500B true CN111665500B (en) 2023-07-07

Family

ID=72387321

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010533605.2A Active CN111665500B (en) 2020-06-12 2020-06-12 Pulse through-wall radar imaging method based on single-bit compressed sensing

Country Status (1)

Country Link
CN (1) CN111665500B (en)

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
IL245366A0 (en) * 2016-05-01 2016-08-31 Technion Res & Dev Foundation Mimo radar and method of using thereof
CN106772365B (en) * 2016-11-25 2019-07-12 南京理工大学 A kind of multipath based on Bayes's compressed sensing utilizes through-wall radar imaging method
CN108896990B (en) * 2018-05-10 2022-06-03 桂林电子科技大学 Building wall imaging method and device by using coupled mode dictionary learning

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
超宽带穿墙雷达未知墙体参数估计研究;朱江;中国优秀硕士学位论文全文数据库（信息科技辑）;第16页 *

Also Published As

Publication number Publication date
CN111665500A (en) 2020-09-15

Similar Documents

Publication Publication Date Title
CN110275166B (en) ADMM-based rapid sparse aperture ISAR self-focusing and imaging method
CN105676168B (en) A kind of acoustic vector sensor array direction estimation method
CN109298383B (en) Mutual-prime array direction-of-arrival estimation method based on variational Bayes inference
CN106646344B (en) A kind of Wave arrival direction estimating method using relatively prime battle array
CN111736131B (en) Method for eliminating false targets of one-bit signal harmonic waves and related components
CN109375154B (en) Coherent signal parameter estimation method based on uniform circular array in impact noise environment
CN110244303B (en) SBL-ADMM-based sparse aperture ISAR imaging method
CN109633538B (en) Maximum likelihood time difference estimation method of non-uniform sampling system
CN108802705B (en) Space-time adaptive processing method and system based on sparsity
CN111562545B (en) PD-ALM algorithm-based sparse array DOA estimation method
CN110687528B (en) Adaptive beam former generation method and system
CN107064896B (en) MIMO radar parameter estimation method based on truncation correction SL0 algorithm
CN112147608A (en) Rapid Gaussian gridding non-uniform FFT through-wall imaging radar BP method
CN110954860B (en) DOA and polarization parameter estimation method
CN112305537A (en) Single-bit random frequency control array radar target distance-angle joint estimation method
CN112285647A (en) Signal orientation high-resolution estimation method based on sparse representation and reconstruction
CN111913155A (en) Two-dimensional DOA estimation method based on array radar
CN115453528A (en) Method and device for realizing segmented observation ISAR high-resolution imaging based on rapid SBL algorithm
CN117092585B (en) Single-bit quantized DoA estimation method, system and intelligent terminal
CN111665500B (en) Pulse through-wall radar imaging method based on single-bit compressed sensing
CN112800599A (en) Non-grid DOA estimation method based on ADMM under array element mismatch condition
CN116660883A (en) Single-bit compressed sensing radar target delay Doppler DOA estimation method
CN111044996A (en) LFMCW radar target detection method based on dimension reduction approximate message transfer
CN113589265B (en) Block near-end gradient dual-sparse dictionary learning beam forming method and system
CN115453527A (en) Periodic sectional observation ISAR high-resolution imaging method

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
GR01 Patent grant
GR01 Patent grant