CN113359095A - Coherent passive MIMO radar Clarithrome boundary calculation method - Google Patents

Coherent passive MIMO radar Clarithrome boundary calculation method Download PDF

Info

Publication number
CN113359095A
CN113359095A CN202110543516.0A CN202110543516A CN113359095A CN 113359095 A CN113359095 A CN 113359095A CN 202110543516 A CN202110543516 A CN 202110543516A CN 113359095 A CN113359095 A CN 113359095A
Authority
CN
China
Prior art keywords
target
signal
matrix
representing
crb
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.)
Granted
Application number
CN202110543516.0A
Other languages
Chinese (zh)
Other versions
CN113359095B (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.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
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 University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Publication of CN113359095A publication Critical patent/CN113359095A/en
Application granted granted Critical
Publication of CN113359095B publication Critical patent/CN113359095B/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
    • 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/40—Means for monitoring or calibrating
    • 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/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G01S7/418—Theoretical aspects

Abstract

The invention discloses a coherent passive MIMO radar Clarmero bound calculation method, and belongs to the technical field of radars. The Cramer-Rao bound calculated by the method can be used for evaluating the joint estimation performance of the target speed and the position parameters of the coherent passive MIMO radar, the parameter estimation performance under the influence of direct wave signals under the conditions that the signals are unknown and the signals have constant modulus constraints is considered, and the method is closer to the practical working process of parameter estimation of the coherent passive MIMO radar, and in the process, the availability of the direct waves to the parameter estimation and the unknown property of the signals are considered.

Description

Coherent passive MIMO radar Clarithrome boundary calculation method
Technical Field
The invention belongs to the technical field of radar, and particularly relates to calculation of a parameter estimation performance evaluation index Cramer-Lo boundary (CRB) in radar signal processing, which is suitable for passive radar signal processing.
Background
The passive radar uses the existing opportunistic signals as transmitting signals to complete tasks such as target detection, classification and parameter estimation. The advantages of Multiple Input Multiple Output (MIMO) techniques may be exploited to improve performance when multiple opportunistic signal sources and receiving stations are used. According to document 1 (e.fisher, a.haimovich, r.s.blum, l.j.cimini, d.chizhik and r.a.valenczuela, "Spatial diversity in radar-modules and detection performance," IEEE Transactions on Signal Processing, vol.54, No.3, pp.823-838, 2006.), in a passive MIMO radar system, coherent or non-coherent Processing may be employed depending on the relative relationship between the antenna position and the target beam width. There is a large body of literature available that demonstrates the advantages of coherent processing in parameter estimation.
The cramer-perot lower bound (CRB) is an important indicator for evaluating the performance of a parameter estimation. Document 2(C.Shi, F.Wang and J.Zhou, "frame-random basis for joint target location and velocity estimation in frequency modulated base radar networks," IET Signal Processing, vol.10, No.7, pp.780-790, 2016.) shows a CRB that is a joint estimate of the coherent and passive radar target position and velocity for known transmitted signals. Reference 3(c.shi, f.wang, m.sellatthurai and j.zhou, "Transmitter subset selection in FM-based passive radio networks for joint target estimation," IEEE Sensors Journal, vol.16, No.15, pp.6043-6052, 2016.) considers the CRB-based Transmitter selection problem for coherent passive MIMO radar with unknown transmit signals. Document 4(l.wang, q.he, r.s.blum and h.li, "Joint parameter estimation applying coherent MIMO radar," The Journal of Engineering, vol.2019, No.20, pp.6859-6862, 2019.) analyzes The Joint estimation of target position and velocity assuming that The transmit signal of The coherent passive MIMO radar is known. When the statistics of the transmitted codewords of coherent passive multi-base radar are known, reference 5(b.d. rising and a.nehiri, "Cramer-Rao bases for UMTS-based passive multi-static radar," IEEE Transactions on Signal Processing, vol.62, No.1, pp.95-106, 2014.) calculates the CRB for target parameter estimation. To our knowledge, current work on coherent passive MIMO radar parameter estimation assumes known transmitted signals or statistical information of known transmitted signals.
In passive radar, it is generally assumed that a target path (transmit-target-receive) signal and a direct path (transmit-receive) signal at a receiving end can be separated, and thus a direct path can be used to observe a transmission signal. However, in some cases, it is not possible to separate the target path and direct path signals, which makes it difficult to obtain the transmit signal information directly from the direct signal. Therefore, in a coherent passive MIMO radar system, it is necessary to consider joint estimation of the position and velocity of an object whose transmitted signal is completely unknown and which is affected by a direct wave signal. Since the constant modulus constraint is a common practical constraint for the transmitted signal, it is also important to further consider the joint parameter estimation that the coherent passive radar knows that the transmitted signal has a constant modulus value.
Disclosure of Invention
The invention aims at solving the technical problem of the prior art that the coherent passive MIMO radar joint target speed and position parameter estimation with unknown transmitting signals and direct wave influence considered is obtained, maximum likelihood estimation is carried out, and the Clarmetro bound is calculated.
The technical scheme of the invention is a coherent passive MIMO radar Clarmero bound calculation method, which comprises the following steps:
step 1: arranging signal sampling values from the mth transmitter received by the nth receiver of the coherent passive MIMO radar into a line in sequence to form a received signal rnm；
rnm＝[rnm[1]，…，rnm[K]]T
＝Γd，nmsd，nm+Γt，nmst，nm⊙a(ft，nm)+wnm
Wherein N is 1, …, N, M is 1, …, M, rnm[k](K-1, …, K) is kTsThe received signal at the time of day is,
Figure BDA0003072692150000021
wherein the mth transmitter is at kTsSampled value of time being
Figure BDA0003072692150000022
sm(kTs) Representing the transmission signal of the mth transmitter, EmIs the transmitted energy of the mth transmitter, TsFor the sampling interval, K (K ═ 1, …, K) is the sampling number, K is the total number of samples;
Figure BDA0003072692150000023
is the target reflection coefficient; f. ofcIs the carrier frequency; dd，nmRepresenting a distance corresponding to a direct signal path from the mth transmitter to the nth receiver; dt，mIndicating a distance to the target corresponding to the mth transmitter; dr，nIndicating a distance to the target corresponding to the nth receiver; tau isd，nmRepresenting a time delay corresponding to a direct signal path from the mth transmitter to the nth receiver; tau ist，nmRepresenting a time delay corresponding to a target path from the mth transmitter to the target and back to the nth receiver; f. oft，nmIndicating the doppler frequencies corresponding to the mth transmitter and the nth receiver; w is anm[k]The nth receiver receives the clutter plus noise of the signal path of the mth transmitter; m represents the total number of transmitters, and N represents the total number of receivers; p0Indicating a direct path at dd，nmThe ratio of the received energy to the transmitted energy is 1; p1Indicates the target path is at dt，m＝dr，nThe ratio of the received energy to the transmitted energy is 1;
Figure BDA0003072692150000031
sd，nm＝[sm(Ts-τd，nm)，…，sm(KTs-τd，nm)]T，st，nm＝[sm(Ts-τt，nm)，…，sm(KTs-τt，nm)]T，
Figure BDA0003072692150000032
(·)Trepresenting a transpose;
step 2: the time domain signal r is obtained according tonmConversion to frequency domain signal Rnm；
Rnm＝Γd，nmZ(τd，nm)Sm+Γt，nmA(ft，nm)Z(τt，nm)Sm+Wnm
Wherein s ism＝[sm[1]，sm[2]，…，sm[K]]T，Sm＝Tsm，Wnm＝Twnm，A(ft，nm) Is composed of Ta (f)t，nm) The cyclic matrix of the construction is constructed,
Figure BDA0003072692150000033
Figure BDA0003072692150000034
t is a discrete time transformation matrix having elements of
Figure BDA0003072692150000035
fΔ＝1/(TsK) Diag {. } represents a diagonal matrix;
and step 3: arranging all the frequency domain signals of M transmitters received by N receivers into a line in sequence to form a total received signal vector
Figure BDA0003072692150000036
Wherein the content of the first and second substances,
Figure BDA0003072692150000037
Γd＝Diag{Γd，11IK×K，Γd， 12IK×K，…，Γd，NMIK×K}Γt＝Diag{Γt，11IK×K，Γt，12IK×K，…，Γt，NMIK×K}，Zd，n＝Diag{Zd，n1，Zd，n2，…，Zd，nM}Zt，n＝Diag{Zt，n1，Zt，n2，…，Zt，nM}，
Figure BDA0003072692150000038
Af＝Diag{Af，11，Af，12，…，Af，NM}
Figure BDA0003072692150000039
diag {. is said to constitute a block diagonal matrix, IK×KAn identity matrix representing a K dimension;
and 4, step 4: determining a covariance matrix Q of a frequency domain received signal for maximum likelihood estimation
Figure BDA00030726921500000310
Wherein
Figure BDA00030726921500000311
Figure BDA00030726921500000312
Represents the kronecker product (·)HDenotes conjugate transpose, QrIs a zero mean value Gaussian distribution matrix
Figure BDA0003072692150000041
The covariance matrix of (a);
and 5: according to the formula
Figure BDA0003072692150000042
Obtaining an estimated value of theta
Figure BDA0003072692150000043
Where θ is the target parameter we want to estimate, including: target position x, y, target velocity vx，vyAnd reflection coefficient
Figure BDA0003072692150000044
Real part of
Figure BDA0003072692150000045
And imaginary part
Figure BDA0003072692150000046
Expressed as:
Figure BDA0003072692150000047
B＝ΓdZd+ΓtAfZt；
step 6: according to the formula
Figure BDA0003072692150000048
Determining an estimated value S of the transmission signal SMLAccording to
Figure BDA0003072692150000049
Obtaining an estimate of the time domain transmit signal
Figure BDA00030726921500000410
The real and imaginary parts are arranged as vectors in the following manner
Figure BDA00030726921500000411
And 7: repeating steps 1 to 6 according to the estimated
Figure BDA00030726921500000412
The root mean square error was found to be:
Figure BDA00030726921500000413
where num is the number of repetitions,
Figure BDA00030726921500000414
Figure BDA00030726921500000415
wherein
Figure BDA00030726921500000416
Are respectively a transmitting signal sm[k]The real and imaginary parts of (c);
and 8: let τ be ═ τt，11，τt，12，…，τt，NM]T，f＝[ft，11，ft，12，…，ft，NM]T，dt＝[dt，1，dt，2，…，dt，M]T，dr＝[dr，1，dr，2，…，dr，N]TObtaining a matrix:
Figure BDA00030726921500000417
wherein
Figure BDA00030726921500000418
Respectively, the derivative of the time delay tau with respect to the target position x, y,
Figure BDA00030726921500000419
respectively representing the Doppler frequency f vs. x, y, vx，vyThe derivative of (a) of (b),
Figure BDA00030726921500000420
respectively, the distance d from the transmitter to the targettThe derivative of x, y,
Figure BDA00030726921500000421
respectively representing the distance d of the receiver to the targetrDerivatives of x, y;
and step 9: the ij (i, j ═ 1,2, …, (2NM + M + N +2+2MK)) th element of matrix B is obtained as:
Figure BDA0003072692150000051
wherein
Figure BDA0003072692150000052
The representation takes the imaginary part of a complex number,
Figure BDA0003072692150000053
step 10: according to the formula:
CRB(β)＝(ABAH)-1
calculating the corresponding x, y, vx，vyCRB of (b), the first four of the diagonal elements of CRB (β) being the target position x, y and target velocity v, respectivelyx，vyLower cramer-mello boundary of (c);
step 11: according to the formula
Figure BDA0003072692150000054
s.t.fi(β)＝|sm[k]|2-1＝0，i＝K(m-1)+k，m＝1，2，...，M
Finding an estimate of beta under constant modulus constraints
Figure BDA0003072692150000055
Step 12: repeat step 11 based on the estimated betaCMAnd solving the root mean square error under the constant modulus constraint as:
Figure BDA0003072692150000056
step 13: calculate the matrix U as
Figure BDA0003072692150000057
Wherein U is1＝[A1，A2]T，
Figure BDA0003072692150000058
Figure BDA0003072692150000059
Figure BDA0003072692150000061
Step 14: according to the formula:
CRBCM(β)＝U(UTJU)-1UT
calculating the x, y, v under the constant modulus constraintx，vyWherein J ═ { CRB (β) }-1，CRBCMThe first four diagonal elements of (beta) are respectively target position x, y and target speed v under constant modulus constraintx，vyLower boundary of cramer.
The Cramer-Rao bound calculated by the method can be used for evaluating the joint estimation performance of the target speed and the position parameters of the coherent passive MIMO radar, the parameter estimation performance under the influence of direct wave signals under the conditions that the signals are unknown and the signals have constant modulus constraints is considered, and the method is closer to the practical working process of parameter estimation of the coherent passive MIMO radar, and in the process, the availability of the direct waves to the parameter estimation and the unknown property of the signals are considered.
Drawings
FIG. 1 is a graph of x, y, v computed for the unconstrained case and constant modulus constrained scenario under different SCNRs when DSR is 0dBx，vySchematic of RMSE and RCRB of (a).
FIG. 2 is a plot of x, y, v computed for the unconstrained case and constant modulus constrained scenario under different DSRs when SCNR is 20dBx，vySchematic representation of the RCRB of (a).
Detailed Description
For convenience of description, the following definitions are first made:
()Tis a transposition ofHIs a conjugate transpose of the original image,
Figure BDA0003072692150000067
representing a mathematical expectation, Diag {. denotes the diagonal matrix, Diag {. denotes the building Block diagonal matrix, IK×KAn identity matrix of dimension K is expressed, det (-) expresses determinant of matrix,
Figure BDA0003072692150000062
the representation takes the real part of a complex number,
Figure BDA0003072692150000063
representing the kronecker product.
Considering a coherent passive MIMO radar with M single-antenna transmitters and N single-antenna receivers, the M (M-1, …, M) th transmitting station and the N (N-1, …, N) th receiving station are located at the respective positions in a cartesian coordinate system
Figure BDA0003072692150000064
And
Figure BDA0003072692150000065
m transmitter at kTsSampled value of time being
Figure BDA0003072692150000066
sm[k]For the transmitted signal of the m-th transmitting antenna, TsFor the sampling interval, K (K ═ 1, …, K) is the number of samples, K being the total number of samples, assuming that the transmitted signals of different transmitters are orthogonal. EmIs the transmitted signal energy of the mth transmitter. Assume a target position of (x, y) and a velocity of (v)x，vy) Target reflection coefficient of
Figure BDA0003072692150000071
Then kTsThe nth receiver receives the received signal from the mth transmitter at the time of
Figure BDA0003072692150000072
Wherein f iscIs the carrier frequency; dd，nmRepresenting a distance corresponding to a direct signal path from the mth transmitter to the nth receiver; dt，mIndicating a distance to the target corresponding to the mth transmitter; dr，nIndicating a distance to the target corresponding to the nth receiver; tau isd，nmRepresenting a time delay corresponding to a direct signal path from the mth transmitter to the nth receiver; tau ist，nmRepresenting a time delay corresponding to a target path from the mth transmitter to the target and back to the nth receiver; f. oft，nmIndicating the doppler frequencies corresponding to the mth transmitter and the nth receiver; w is anm[k]The nth receiver receives the clutter plus noise of the signal path of the mth transmitter; m represents the total number of transmitters, and N represents the total number of receivers; p0Indicating a direct path at dd，nmThe ratio of the received energy to the transmitted energy is 1; p1Indicates the target path is at dt，m＝dr，nThe ratio of received energy to transmitted energy is 1. Suppose the parameters to be estimated are position (x, y), velocity (v)x，vy) And a transmission signal sm[k]Where M is 1, …, M, K is 1, …, K, which may be written as
Figure BDA0003072692150000073
ft，nmFor the position (x, y) and velocity (v) of the object to be estimatedx，vy) Function of (2)
Figure BDA0003072692150000074
Where λ represents the carrier wavelength.
Distance dt，mAnd dr，nIs a function of the position (x, y) of the object to be estimated
Figure BDA0003072692150000075
Figure BDA0003072692150000076
Signal path delay τt，nmAnd direct path delay τd，nmAs a function of the target position (x, y) to be estimated:
Figure BDA0003072692150000081
Figure BDA0003072692150000082
where c represents the speed of light.
The received signal from the M (M1, …, M) th transmitter received by the N (N1, …, N) th receiver is rnm
Figure BDA0003072692150000083
Wherein
Figure BDA0003072692150000084
Figure BDA0003072692150000085
sd，nm＝[sm(Ts-τd，nm)，…，sm(KTs-τd，nm)]T (10)
st，nm＝[sm(Ts-τt，nm)，…，sm(KTs-τt，nm)]T (11)
Figure BDA0003072692150000086
wnm＝[wnm[1]，...，wnm[K]]T (13)
Will r isnmConversion to frequency domain signal RnmWhich is a
Rnm＝Γd，nmZ(τd，nm)Sm+Γt，nmA(ft，nm)Z(τt，nm)Sm+Wnm (14)
Wherein
sm＝[sm[1]，sm[2]，…，sm[K]]T (15)
Sm＝Tsm (16)
Wnm＝Twnm (17)
Figure BDA0003072692150000087
Figure BDA0003072692150000088
A(ft，nm) Is composed of Ta (f)t，nm) Constructed circulant matrix, T being a discrete time transformation matrix whose elements are
Figure BDA0003072692150000089
Wherein f isΔ＝1/(TsK)。
The total frequency domain signals of M transmitters received by N receivers are
Figure BDA00030726921500000810
Wherein
Figure BDA00030726921500000811
Figure BDA0003072692150000091
Γd＝Diag{Γd，11IK×K，Γd，12IK×K，…，ΓdNMIK×K} (24)
Γt＝Diag{Γt，11IK×K，Γt，12IK×K，…，Γt，NMIK×K} (25)
Zd，n＝Diag{Zd，n1，Zd，n2，…，Zd，nM} (26)
Zt，n＝Diag{Zt，n1，Zt，n2，…，Zt，nM} (27)
Figure BDA0003072692150000092
Af＝Diag{Af，11，Af，12，…，Af，NM} (29)
Figure BDA0003072692150000093
Assuming time domain noise
Figure BDA0003072692150000094
Obeying a 0-mean covariance matrix of QrThe complex Gaussian random variable of (1) is a covariance matrix with zero mean of frequency domain noise W of
Figure BDA0003072692150000095
Of complex Gaussian random variables, wherein
Figure BDA0003072692150000096
The invention adopts the following steps to calculate the maximum likelihood estimation and CRB of the external radiation source MIMO radar:
step 1, arranging the signal sampling values received by N receivers into a line in sequence by the signal model (21) to determine a received signal R,
R＝ΓdZdS+ΓtAfZtS+W (31)
step 2: according to the formula
Figure BDA0003072692150000097
Obtaining an estimated value of the target parameter theta
Figure BDA0003072692150000098
Where θ is the target parameter we want to estimate, including: target position x, y, target velocity vx，vyAnd reflection coefficient
Figure BDA0003072692150000099
Real part of
Figure BDA00030726921500000910
And imaginary part
Figure BDA00030726921500000911
Expressed as:
Figure BDA00030726921500000912
B＝ΓdZd+ΓtAfZt；
and step 3: according to the formula
Figure BDA00030726921500000913
Obtaining an estimated value S of the transmitted frequency domain signal SMLFurther according to
Figure BDA00030726921500000914
Obtaining an estimate of the time domain transmit signal
Figure BDA00030726921500000915
The real and imaginary parts are arranged as vectors in the following manner
Figure BDA00030726921500000916
And 4, step 4: repeating steps 1 to 3 according to the estimated
Figure BDA00030726921500000917
The root mean square error was found to be:
Figure BDA00030726921500000918
where num is the number of repetitions,
Figure BDA00030726921500000919
Figure BDA00030726921500000920
wherein
Figure BDA00030726921500000921
Are respectively a transmitting signal sm[k]The real and imaginary parts of (c);
and 5: let τ be ═ τt，11，τt，12，…，τt，NM]T，f＝[ft，11，ft，12，…，ft，NM]T，dt＝[dt，1，dt，2，…，dt，M]T，dr＝[dr，1，dr，2，…，dr，N]TObtaining a matrix
Figure BDA0003072692150000101
Wherein
Figure BDA0003072692150000102
Respectively, the derivative of the time delay tau with respect to the target position x, y,
Figure BDA0003072692150000103
respectively representing the Doppler frequency f vs. x, y, vx，vyThe derivative of (a) of (b),
Figure BDA0003072692150000104
respectively, the distance d from the transmitter to the targettThe derivative of x, y,
Figure BDA0003072692150000105
respectively representing the distance d of the receiver to the targetrDerivative to x, y.
Step 6: the ij (i, j ═ 1,2, …, (2NM + M + N +2+2MK)) th element of matrix B is obtained as:
Figure BDA0003072692150000106
wherein
Figure BDA0003072692150000107
The representation takes the imaginary part of a complex number,
Figure BDA0003072692150000108
and 7: according to the formula:
CRB(β)＝(ABAH)-1 (37)
calculating the corresponding x, y, vx，vyCRB of (b), the first four of the diagonal elements of CRB (β) being the target position x, y and target velocity v, respectivelyx，vyLower boundary of cramer.
And 8: according to the formula
Figure BDA0003072692150000109
s.t.fi(β)＝|sm[k]|2-1＝0，i＝K(m-1)+k，m＝1，2，...，M
Finding an estimate of beta under constant modulus constraints
Figure BDA00030726921500001010
And step 9: repeating step 8 according to the estimated betaCMAnd solving the root mean square error under the constant modulus constraint as:
Figure BDA0003072692150000111
step 10: calculate the matrix U as
Figure BDA0003072692150000112
Wherein U is1＝[A1，A2]T，
Figure BDA0003072692150000113
Figure BDA0003072692150000114
Figure BDA0003072692150000115
Step 11: according to the formula:
CRBCM(β)＝U(UTJU)-1UT (43)
calculating constant modulusCorresponding to x, y, v under constraintx，vyWherein J ═ { CRB (β) }-1，CRBCMThe first four diagonal elements of (beta) are respectively target position x, y and target speed v under constant modulus constraintx，vyLower boundary of cramer.
Working principle of the invention
The likelihood function is expressed as a function of the frequency domain received signal model (21)
Figure BDA0003072692150000116
Further obtaining a log likelihood function of
Figure BDA0003072692150000117
Wherein C is0Is a term that is independent of the parameter to be estimated. L (R | β) is derived from S and the reciprocal is made 0, giving an ML estimate of S as
Figure BDA0003072692150000118
Wherein B ═ ΓdZd+ΓtAfZt. Will be provided with
Figure BDA0003072692150000119
Bringing in (46) available target parameters
Figure BDA00030726921500001110
Is estimated as
Figure BDA0003072692150000121
Order to
Figure BDA0003072692150000122
Wherein τ ═ τ [ τ ]t，11，τt，12，…，τt，NM]T，f＝[ft，11，ft，12，…，ft，NM]T，dt＝[dt，1，dt，2，…，dt，M]T，dr＝[dr，1，dr，2，…，dr，N]TAccording to the chain rule
CRB(β)＝(ABAH)-1 (49)
Wherein
Figure BDA0003072692150000123
According to the literature (S.Kay, "Fundamentals of Statistical Signal Processing: Estimation Theory," Prentice-Hall.Englewood Cli _ s, NJ, 1993), the compounds of formula I are available
Figure BDA0003072692150000124
Under constant modulus constraint, write the constant modulus constraint set as
fi(β)＝|sm[k]|2-1＝0，i＝K(m-1)+k (52)
Wherein M1, 2, 1., M, K1, 2. According to the literature (Terference J.Moore, Brian M.Sadler and Richard J.Kozick, "Maximum-similarity-probability estimation, the Cramer-Rao bound, and the method of estimating with parameter constraints," IEEE Transactions on Signal Processing, vol.54, No.3, pp.895-908, 2008.) ML estimates for the parameters to be estimated are obtained
Figure BDA0003072692150000125
s.t.fi(β)＝|sm[k]|2-1＝0，i＝K(m-1)+k，m＝1，2，...，M
Defining the constraint set vector as f (β) ═ f1(β), …, fKM(β)]TThen a gradient matrix is obtained
Figure BDA0003072692150000131
According to the literature (Brian M.Sadler, Richard J.Kozick and Terference J.Moore, "On the performance of source with constraint modules," IEEE International Conference On Acoustics, Speech and Signal Processing (ICASSP), pp.2373-2376, 2002.), if a gradient matrix Y (β) is a row full-rank matrix for a certain β satisfying a constraint f (β), then there is a (2KM +6) x (KM +6) dimensional matrix U, whose columns are the orthogonal bases of the null space of the matrix Y (β), i.e., the columns of U are the columns of the matrix Y (β)
YU＝0，UTU＝I (55)
Wherein U is1＝[A1，A2]T，
Figure BDA0003072692150000132
Figure BDA0003072692150000133
Figure BDA0003072692150000134
Further, a CRB under constant modulus constraint can be written as
CRBCM(β)＝U(UTJU)-1UT (58)
Wherein J ═ { CRB (β) }-1Is the FIM matrix in the unconstrained case.
The maximum likelihood estimation and the CRB are calculated in an unconstrained scene and a constant modulus constrained scene based on a coherent passive MIMO radar signal model, wherein the maximum likelihood estimation is obtained by adopting 1000 Monte Carlo experiments, simulation results are shown in figures 1 and 2, and simulation parameters are set as follows:
assuming that the target position is (200, -1500) m, it is moved at a speed of (-70, 20) m/s. A coherent passive MIMO radar has 2 transmitters with positions (-340, -290) km and (-330, 80) km, and 3 receivers with positions (270, 220) km, (260, 200) km and (100, -350) km.
The transmitted constant modulus signal is sm[k]＝s′m[k]/|s′m[k]1,2,. M, K1, 2,. K, where s'm[k]Is a signal that is an OFDM signal and,
Figure BDA0003072692150000141
NcΔ f is the subcarrier spacing, α, for the total number of subcarriersmiE { -1, 1} is the ith (i ═ 1, 2.. Nc) A data symbol, fbIs the frequency offset between different transmitters. Assume carrier frequency fc＝1GHz，Nc＝64，Δf＝125Hz，fb1600Hz, and the sampling frequency is set to fs7800Hz and K46 total samples.
Defining the reflected signal to clutter plus noise ratio as
Figure BDA0003072692150000142
The ratio of direct wave signal to clutter plus noise is defined as
Figure BDA0003072692150000143
In fig. 1, the RMSE and RCRB of ML estimation as a function of SCNR are plotted at DSR 0dB, taking into account two scenarios, radar vs. constant modulus constraint, no-a-priori (unconditioned) and known (constrained). As can be seen from the figure, all RMSEs in the liangz scenario decrease with increasing SCNR, and all RMSE curves have a threshold, and above which the RMSE starts to approach the RCRB, proving the correctness of the CRB. It can also be seen that the case of known constant modulus constraints has a lower threshold, and that the RCRB of the known constant modulus constraints is smaller than the RCRB of the unconstrained. Therefore, the performance under the known constant modulus constraint condition is superior to the performance without prior to the constant modulus constraint, mainly because the known constant modulus constraint condition utilizes the prior information of the constant modulus constraint of the unknown signal.
Figure 2 SCNR-20 dB, gives the curves of RCRB as a function of DSR for both cases. We see that for a sufficiently small DSR (DSR ≦ -20dB), the RCRB remains almost unchanged because the DSR of the direct-path signal is too small in this case to obtain useful information from the direct-path signal. As DSR becomes larger, RCRB becomes smaller for both cases, indicating that the direct path signal helps to improve estimation performance. When DSR increases to some extent (DSR ≧ 20dB in this example), RCRB hardly decreases with increasing DSR, since we have obtained enough transmitted signal information from the direct-path signal. Thus, when DSR is large enough, with the direct path signal, the estimation performance can be greatly improved in both cases. At the same time, similar to the conclusions of FIG. 1, under different DSRs, the performance under constrained conditions is better than the performance under unconstrained conditions.

Claims (1)

1. A method for computing a coherent passive MIMO radar clalmelo boundary, the method comprising:
step 1: arranging signal sampling values from the mth transmitter received by the nth receiver of the coherent passive MIMO radar into a line in sequence to form a received signal rnm；
Figure FDA0003072692140000011
Wherein N is 1, …, N, M is 1, …, M, rnm[k](K-1, …, K) is kTsThe received signal at the time of day is,
Figure FDA0003072692140000012
wherein the mth transmitter is at kTsSampled value of time being
Figure FDA0003072692140000013
sm(kTs) Representing the transmission signal of the mth transmitter, EmIs the transmitted energy of the mth transmitter, TsFor the sampling interval, K (K ═ 1, …, K) is the sampling number, K is the total number of samples;
Figure FDA0003072692140000015
is the target reflection coefficient; f. ofcIs the carrier frequency; dd,nmRepresenting a distance corresponding to a direct signal path from the mth transmitter to the nth receiver; dt,mIndicating a distance to the target corresponding to the mth transmitter; dr,nIndicating a distance to the target corresponding to the nth receiver; tau isd,nmRepresenting a time delay corresponding to a direct signal path from the mth transmitter to the nth receiver; tau ist,nmRepresenting a time delay corresponding to a target path from the mth transmitter to the target and back to the nth receiver; f. oft,nmIndicating the doppler frequencies corresponding to the mth transmitter and the nth receiver; w is anm[k]The nth receiver receives the clutter plus noise of the signal path of the mth transmitter; m represents the total number of transmitters, and N represents the total number of receivers; p0Indicating a direct path at dd,nmThe ratio of the received energy to the transmitted energy is 1; p1Indicates the target path is at dt,m＝dr,nThe ratio of the received energy to the transmitted energy is 1;
Figure FDA0003072692140000014
sd,nm＝[sm(Ts-τd,nm),…,sm(KTs-τd,nm)]T，st,nm＝[sm(Ts-τt,nm),…,sm(KTs-τt,nm)]T，
Figure FDA0003072692140000016
wnm＝[wnm[1],...,wnm[K]]T，(·)Trepresenting a transpose;
step 2: the time domain signal r is obtained according tonmConversion to frequency domain signal Rnm；
Rnm＝Γd,nmZ(τd,nm)Sm+Γt,nmA(ft,nm)Z(τt,nm)Sm+Wnm
Wherein s ism＝[sm[1],sm[2],…,sm[K]]T，Sm＝Tsm，Wnm＝Twnm，A(ft,nm) Is composed of Ta (f)t,nm) The cyclic matrix of the construction is constructed,
Figure FDA0003072692140000021
Figure FDA0003072692140000022
t is a discrete time transformation matrix having elements of
Figure FDA0003072692140000023
fΔ＝1/(TsK) Diag {. } represents a diagonal matrix;
and step 3: arranging all the frequency domain signals of M transmitters received by N receivers into a line in sequence to form a total received signal vector
Figure FDA0003072692140000024
Wherein the content of the first and second substances,
Figure FDA0003072692140000025
Γt＝Diag{Γt,11IK×K,Γt,12IK×K,…,Γt,NMIK×K}，Zd,n＝Diag{Zd,n1,Zd,n2,…,Zd,nM}
Figure FDA0003072692140000026
Af＝Diag{Af,11,Af,12,…,Af,NM}
Figure FDA0003072692140000027
diag {. is said to constitute a block diagonal matrix, IK×KAn identity matrix representing a K dimension;
and 4, step 4: determining a covariance matrix Q of a frequency domain received signal for maximum likelihood estimation
Figure FDA00030726921400000216
Wherein
Figure FDA0003072692140000028
Figure FDA0003072692140000029
Represents the kronecker product (·)HDenotes conjugate transpose, QrIs a zero mean value Gaussian distribution matrix
Figure FDA00030726921400000210
The covariance matrix of (a);
and 5: according to the formula
Figure FDA00030726921400000211
Obtaining an estimated value of theta
Figure FDA00030726921400000212
Where θ is the target parameter we want to estimate, including: target position x, y, target velocity vx,vyAnd vice versaCoefficient of radiation
Figure FDA00030726921400000217
Real part of
Figure FDA00030726921400000218
And imaginary part
Figure FDA00030726921400000219
Expressed as:
Figure FDA00030726921400000220
B＝ΓdZd+ΓtAfZt；
step 6: according to the formula
Figure FDA00030726921400000213
Determining an estimated value S of the transmission signal SMLAccording to
Figure FDA00030726921400000214
Obtaining an estimate of the time domain transmit signal
Figure FDA00030726921400000215
The real and imaginary parts are arranged as vectors in the following manner
Figure FDA0003072692140000031
And 7: repeating steps 1 to 6 according to the estimated
Figure FDA0003072692140000032
The root mean square error was found to be:
Figure FDA0003072692140000033
where num is the number of repetitions,
Figure FDA0003072692140000034
Figure FDA0003072692140000035
wherein
Figure FDA0003072692140000036
Are respectively a transmitting signal sm[k]The real and imaginary parts of (c);
and 8: let τ be ═ τt,11,τt,12,…,τt,NM]T，f＝[ft,11,ft,12,…,ft,NM]T，dt＝[dt,1,dt,2,…,dt,M]T，dr＝[dr,1,dr,2,…,dr,N]TObtaining a matrix:
Figure FDA0003072692140000037
wherein
Figure FDA0003072692140000038
Respectively, the derivative of the time delay tau with respect to the target position x, y,
Figure FDA0003072692140000039
respectively representing the Doppler frequency f vs. x, y, vx,vyThe derivative of (a) of (b),
Figure FDA00030726921400000310
respectively, the distance d from the transmitter to the targettThe derivative of x, y,
Figure FDA00030726921400000311
respectively representing the distance d of the receiver to the targetrDerivatives of x, y;
and step 9: the ij (i, j ═ 1,2, …, (2NM + M + N +2+2MK)) th element of matrix B is obtained as:
Figure FDA00030726921400000312
wherein
Figure FDA00030726921400000313
The representation takes the imaginary part of a complex number,
Figure FDA00030726921400000314
step 10: according to the formula:
CRB(β)＝(ABAH)-1
calculating the corresponding x, y, vx,vyCRB of (b), the first four of the diagonal elements of CRB (β) being the target position x, y and target velocity v, respectivelyx,vyLower cramer-mello boundary of (c);
step 11: according to the formula
Figure FDA0003072692140000041
s.t.fi(β)＝|sm[k]|2-1＝0,i＝K(m-1)+k,m＝1,2,...,M
Finding an estimate of beta under constant modulus constraints
Figure FDA0003072692140000042
Step 12: repeat step 11 based on the estimated betaCMAnd solving the root mean square error under the constant modulus constraint as:
Figure FDA0003072692140000043
step 13: calculate the matrix U as
Figure FDA0003072692140000044
Wherein U is1＝[A1,A2]T，
Figure FDA0003072692140000045
Figure FDA0003072692140000046
Figure FDA0003072692140000047
Step 14: according to the formula:
CRBCM(β)＝U(UTJU)-1UT
calculating the x, y, v under the constant modulus constraintx,vyWherein J ═ { CRB (β) }-1，CRBCMThe first four diagonal elements of (beta) are respectively target position x, y and target speed v under constant modulus constraintx,vyLower boundary of cramer.
CN202110543516.0A 2021-04-27 2021-05-19 Coherent passive MIMO radar Clarithrome boundary calculation method Active CN113359095B (en)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN202110460402 2021-04-27
CN202110460402X 2021-04-27

Publications (2)

Publication Number Publication Date
CN113359095A true CN113359095A (en) 2021-09-07
CN113359095B CN113359095B (en) 2022-10-14

Family

ID=77526948

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110543516.0A Active CN113359095B (en) 2021-04-27 2021-05-19 Coherent passive MIMO radar Clarithrome boundary calculation method

Country Status (1)

Country Link
CN (1) CN113359095B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113777575A (en) * 2021-09-14 2021-12-10 电子科技大学长三角研究院(衢州) MIMO radar multi-target parameter estimation and evaluation method applied to complex environment

Citations (17)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2008070671A2 (en) * 2006-12-04 2008-06-12 Ntt Docomo Inc. Method for optimum threshold selection of time-of-arrival estimators
CN103852750A (en) * 2014-01-28 2014-06-11 大连大学 Robust MIMO radar waveform optimization method for improving worst estimated performance
US20140347211A1 (en) * 2013-05-24 2014-11-27 Robert Bosch Gmbh Method for operating a mimo radar
CN104833959A (en) * 2015-04-09 2015-08-12 大连大学 MIMO radar waveform optimization method based on target prior information
CN105527610A (en) * 2015-12-05 2016-04-27 中国人民解放军信息工程大学 Multi-antenna joint optimization clutter suppression method based on fractional time delay estimation
CN106646410A (en) * 2017-01-06 2017-05-10 天津大学 Learning-perception-decision making-responding method under broadband cognitive passive radar architecture
CN106886011A (en) * 2017-01-19 2017-06-23 电子科技大学 A kind of MIMO radar Cramér-Rao lower bound computational methods for reflecting through wave action
CN106909779A (en) * 2017-02-17 2017-06-30 电子科技大学 MIMO radar Cramér-Rao lower bound computational methods based on distributed treatment
CN107144827A (en) * 2017-04-12 2017-09-08 北京理工大学 A kind of distributed radar preferred configuration construction method based on Combined estimator CRLB
CN108469601A (en) * 2018-01-17 2018-08-31 大连大学 Improve the transmit-receive combination Robust Optimization of the worst detection performances of MIMO-STAP
CN108957422A (en) * 2018-06-01 2018-12-07 电子科技大学 A kind of root carat Metro lower bound calculation method of the cloud MIMO radar based on quantized data
CN109239683A (en) * 2018-08-28 2019-01-18 天津大学 The Cramér-Rao lower bound analysis method of the passive MIMO radar in broadband
CN110133635A (en) * 2019-04-03 2019-08-16 电子科技大学 A kind of method of cooperation MIMO radar and communication system calculating target positioning and mutual information
EP3690483A1 (en) * 2019-02-04 2020-08-05 Fraunhofer-Gesellschaft zur Förderung der angewandten Forschung e.V. A method for synthesis of antenna array layouts or selection of waveform in a set of mutually incoherent apertures for radar and radio-frequency applications
CN111650563A (en) * 2020-06-15 2020-09-11 桂林电子科技大学 System and method for quickly estimating co-channel interference time delay and energy of external radiation source radar
CN111866709A (en) * 2020-06-29 2020-10-30 重庆邮电大学 Indoor Wi-Fi positioning error bound estimation method for moving target
CN111948619A (en) * 2020-06-30 2020-11-17 西安电子科技大学 Passive cooperative detection method and system for target under irradiation of multi-type external radiation sources

Patent Citations (17)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2008070671A2 (en) * 2006-12-04 2008-06-12 Ntt Docomo Inc. Method for optimum threshold selection of time-of-arrival estimators
US20140347211A1 (en) * 2013-05-24 2014-11-27 Robert Bosch Gmbh Method for operating a mimo radar
CN103852750A (en) * 2014-01-28 2014-06-11 大连大学 Robust MIMO radar waveform optimization method for improving worst estimated performance
CN104833959A (en) * 2015-04-09 2015-08-12 大连大学 MIMO radar waveform optimization method based on target prior information
CN105527610A (en) * 2015-12-05 2016-04-27 中国人民解放军信息工程大学 Multi-antenna joint optimization clutter suppression method based on fractional time delay estimation
CN106646410A (en) * 2017-01-06 2017-05-10 天津大学 Learning-perception-decision making-responding method under broadband cognitive passive radar architecture
CN106886011A (en) * 2017-01-19 2017-06-23 电子科技大学 A kind of MIMO radar Cramér-Rao lower bound computational methods for reflecting through wave action
CN106909779A (en) * 2017-02-17 2017-06-30 电子科技大学 MIMO radar Cramér-Rao lower bound computational methods based on distributed treatment
CN107144827A (en) * 2017-04-12 2017-09-08 北京理工大学 A kind of distributed radar preferred configuration construction method based on Combined estimator CRLB
CN108469601A (en) * 2018-01-17 2018-08-31 大连大学 Improve the transmit-receive combination Robust Optimization of the worst detection performances of MIMO-STAP
CN108957422A (en) * 2018-06-01 2018-12-07 电子科技大学 A kind of root carat Metro lower bound calculation method of the cloud MIMO radar based on quantized data
CN109239683A (en) * 2018-08-28 2019-01-18 天津大学 The Cramér-Rao lower bound analysis method of the passive MIMO radar in broadband
EP3690483A1 (en) * 2019-02-04 2020-08-05 Fraunhofer-Gesellschaft zur Förderung der angewandten Forschung e.V. A method for synthesis of antenna array layouts or selection of waveform in a set of mutually incoherent apertures for radar and radio-frequency applications
CN110133635A (en) * 2019-04-03 2019-08-16 电子科技大学 A kind of method of cooperation MIMO radar and communication system calculating target positioning and mutual information
CN111650563A (en) * 2020-06-15 2020-09-11 桂林电子科技大学 System and method for quickly estimating co-channel interference time delay and energy of external radiation source radar
CN111866709A (en) * 2020-06-29 2020-10-30 重庆邮电大学 Indoor Wi-Fi positioning error bound estimation method for moving target
CN111948619A (en) * 2020-06-30 2020-11-17 西安电子科技大学 Passive cooperative detection method and system for target under irradiation of multi-type external radiation sources

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
QIAN HE;JIANBIN HU;RICK S. BLUM;YONGGANG WU: "Generalized Cramér–Rao Bound for Joint Estimation of Target Position and Velocity for Active and Passive Radar Networks", 《IEEE TRANSACTIONS ON SIGNAL PROCESSING》 *
SANDEEP GOGINENI;MURALIDHAR RANGASWAMY;BRIAN D. RIGLING;ARYE NEH: "Cramér-Rao Bounds for UMTS-Based Passive Multistatic Radar", 《IEEE TRANSACTIONS ON SIGNAL PROCESSING》 *
XIN ZHANG;HONGBIN LI;JUN LIU;BRAHAM HIMED: "Joint Delay and Doppler Estimation for Passive Sensing With Direct-Path Interference", 《IEEE》 *
吴永刚: "外辐射源MIMO雷达的参数估计问题研究", 《中国优秀硕士学位论文全文数据库 信息科技辑》 *
孙文: "多接收站外辐射源雷达目标检测与跟踪", 《中国博士学位论文全文数据库 信息科技辑》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113777575A (en) * 2021-09-14 2021-12-10 电子科技大学长三角研究院(衢州) MIMO radar multi-target parameter estimation and evaluation method applied to complex environment
CN113777575B (en) * 2021-09-14 2023-08-25 电子科技大学长三角研究院(衢州) MIMO radar multi-target parameter estimation evaluation method applied to complex environment

Also Published As

Publication number Publication date
CN113359095B (en) 2022-10-14

Similar Documents

Publication Publication Date Title
CN104181517B (en) Method for running MIMO radar
Yardibi et al. Source localization and sensing: A nonparametric iterative adaptive approach based on weighted least squares
CN105068049B (en) A kind of Cramér-Rao lower bound computational methods for splitting antenna MIMO radar
CN108957422B (en) Quantitative data-based root-caramello lower bound calculation method for cloud MIMO radar
CN109507641B (en) Cooperative MIMO radar and MIMO communication system performance calculation method
CN113777575B (en) MIMO radar multi-target parameter estimation evaluation method applied to complex environment
CN104678372A (en) Joint estimation method for super-resolution distance value and angle value by using orthogonal frequency division multiplexing radar
US20090204362A1 (en) Indoor location determination system and method
CN108828504B (en) MIMO radar target direction fast estimation method based on partial correlation waveform
CN114286307B (en) Channel state information parameter estimation method based on matrix beams
CN113359095B (en) Coherent passive MIMO radar Clarithrome boundary calculation method
CN108416105A (en) Steady adaptive beam-forming algorithm under pulse and Gaussian noise
CN113360841B (en) Distributed MIMO radar target positioning performance calculation method based on supervised learning
CN113406562B (en) TOA and DOA combined estimation dimension reduction method in Beidou and ultra-wideband system
Bayat et al. Joint TOA/DOA wireless position location using matrix pencil
Hyder et al. A new target localization method for bistatic FDA radar
CN110391820B (en) Communication receiving method for avoiding same frequency interference based on DFT
CN110231590B (en) Array target angle measurement method based on DFT (discrete Fourier transform)
Gogineni et al. Adaptive design for distributed MIMO radar using sparse modeling
Velickovic et al. The performance of the modified GCC technique for differential time delay estimation in the cooperative sensor network
Rawat et al. RLS-Based adaptive time-varying RCS estimation and imaging in MIMO radar systems
Javed et al. Characterizing the impacts of rice fading on a WiFi-based passive multistatic radar using Cramer-Ráo bound
WO2022168294A1 (en) Radar device, radar system, and radar method
KR102308983B1 (en) Oscillatory pulse sequence generation, detection method and apparatus therefor
CN115514436B (en) Inter-satellite link channel prediction method for Massive MIMO LEO communication

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