This paper was converted on www.awesomepapers.org from LaTeX by an anonymous user.
Want to know more? Visit the Converter page.

RIS-Assisted Joint Uplink Communication and Imaging: Phase Optimization and Bayesian
Echo Decoupling

Shengyu Zhu, Zehua Yu, Qinghua Guo, , Jinshan Ding, ,
Qiang Cheng, , and Tie Jun Cui
This work was supported by the Basic Scientific Center of Information Metamaterials of the National Natural Science Foundation of China under Grant 6228810001 and also supported by the National Natural Science Foundation of China under Grant 62171358.Shengyu Zhu and Jinshan Ding are with the National Laboratory of Radar Signal Processing, Xidian University, Xi’an 710071, China (e-mail: [email protected]). Zehua Yu is with the School of Electronic Engineering, Xidian University, Xi’an 710071, China (e-mail: [email protected]). Qinghua Guo is with the School of Electrical, Computer, and Telecommunications Engineering, University of Wollongong, Wollongong, Australia (e-mail: [email protected]). Qiang Cheng and Tie Jun Cui are with the State Key Laboratory of Millimeter Waves, Southeast University, Nanjing 210096, China (e-mail: [email protected]).
Abstract

Achieving integrated sensing and communication (ISAC) via uplink transmission is challenging due to the unknown waveform and the coupling of communication and sensing echoes. In this paper, a joint uplink communication and imaging system is proposed for the first time, where a reconfigurable intelligent surface (RIS) is used to manipulate the electromagnetic signals for echo decoupling at the base station (BS). Aiming to enhance the transmission gain in desired directions and generate required radiation pattern in the region of interest (RoI), a phase optimization problem for RIS is formulated, which is high dimensional and nonconvex with discrete constraints. To tackle this problem, a back propagation based phase design scheme for both continuous and discrete phase models is developed. Moreover, the echo decoupling problem is tackled using the Bayesian method with the factor graph technique, where the problem is represented by a graph model which consists of difficult local functions. Based on the graph model, a message-passing algorithm is derived, which can efficiently cooperate with the adaptive sparse Bayesian learning (SBL) to achieve joint communication and imaging. Numerical results show that the proposed method approaches the relevant lower bound asymptotically, and the communication performance can be enhanced with the utilization of imaging echoes.

Index Terms:
integrated sensing and communication (ISAC), joint communication and imaging, reconfigurable intelligent surface (RIS), factor graph, information metamaterials.

I Introduction

Integrated sensing and communication (ISAC) has been recognized as a promising technology for the next-generation wireless networks[1]. ISAC not only allows communication and radar systems to share spectrum resources, but also enables communication and radar sensing functionalities simultaneously, which significantly reduces hardware resource consumption. However, this new framework introduces great challenges to hardware design and signal processing, thus bringing a research upsurge in system design 2, 3, 4, 5, 6, optimal joint waveform design [7, 8, 9, 10, 11], exploration of joint processing algorithms [12, 13, 14], and resource allocation[15, 16, 17].

In the past years, information metamaterials have opened a new era of real-time digital regulation of electromagnetic waves [18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29]. As one representative of information metamaterials, reconfigurable intelligent surface (RIS) is expected to break the dependence of traditional sensing and communication on the channel environment. RIS bridges the gap between ISAC and information metamaterials by configuring the electromagnetic environment more flexibly, quickly and intelligently, which is highly expected to better serve the next generation of wireless communication and senisng. Essentially, it is a low-cost passive system that achieves control of electromagnetic beam energy by tuning the phase of array elements. It is shown in [18] that the highly controllable reflection of RIS can be practically achieved by leveraging the existing digitally reconfigurable or programmable metasurface. RIS was first applied as a wireless relay to redirect electromagnetic signals in the field of wireless communication [30, 31, 32], thereby overcoming the adverse effects in natural environments. In addition, it can also be used as a new communication architecture transmitter to implement various modulation [33, 34, 35, 36, 37, 38, 39, 40], e.g., FSK [33], PSK [35], and MIMO [38]. Recently, RIS has been employed to achieve environment awareness and parameters estimation. A new self-sensing RIS architecture was proposed in [41], where the performance of different benchmark sensing systems in the cases of with and without RIS was compared. The authors in [42] and [43] shed light on the interplay among the system parameters, including the radar-RIS distance, the RIS size and the location of the prospective target. The results show that the radar system can achieve the optimal performance when the RIS is deployed in the near field of the radar arrays on both the transmitter and receiver sides. The authors in [44] and [45] proposed to sense humans, recognize their gestures and physiological state simultaneously by utilizing the programmable metasurface and Wi-Fi signals.

Recently, RIS-assisted ISAC system has attracted significant research interest. A Dual-Functional Radar-Communication (DFRC) system with RIS deployed near the communication devices has been proposed in [46, 47]. By optimizing the Cramer-Rao bound (CRB) of DOA estimation, the constant-waveform and RIS phase shifts are designed jointly to mitigate the mitigating multi-user interference (MUI). The RIS-assisted ISAC system in the cases of congested and obstructed channels has been investigated in [48] and [49], respectively, where the PSM of the RIS and the precoding of the base station (BS) are optimized jointly to improve the SNR at receiver. To explore the potential of multiple RISs in assisting ISAC, the authors in [50] propose a double-RIS-assisted ISAC system, where two RISs are deployed to enhance the communication performance while suppressing mutual interference. The beam patterns of RISs and radar are optimized jointly based on penalized dual decomposition (PDD) to enhane the system performance. In [51], the authors propose an RIS-aided localization and communication system, where the theoretical performance limits of localization and communication is derived for both near-field and far-field scenarios. Numerical results show that with the assistance of multiple RISs, both spectral efficiency and localization accuracy can be improved significantly.

Most of the existing literature in ISAC focus on the utilization of downlink transmission. On the contrary, ISAC in uplink transmission is underexplored, technically because it is much more challenging compared with downlink counterparts due to unknown signal waveform and coupling of communication and sensing echoes. In this article, we propose a new RIS-assisted ISAC system, where the uplink transmission is exploited to achieve joint communication and imaging. In contrast to existing works [46, 47, 48, 49, 50] that achieve communication and sensing at different receivers, the proposed uplink ISAC system performs communication and imaging at one receiver, where the communication performance can be enhanced with the utilization of sensing echo. We show that the radiation pattern can be controlled expectedly with the proposed phase design scheme and the factor graph technique is very suitable for echo decoupling for uplink ISAC system. The contribution of this work can be summarized as follows

  • A novel RIS-assisted uplink joint communication and imaging system is presented for the first time, to the best of our knowledge. It is shown that by properly designing the phase shift of the RIS, the radiation pattern can be modulated as desired, thereby allowing joint communication and imaging at the BS with only one RF chain.

  • A phase optimization problem based on the requirements of the system is formulated. A back propagation based phase design scheme with the combination of temperature parameter is proposed to achieve efficient phase optimization for both continuous and discrete phase models.

  • A factor graph representation of the joint communication and imaging is established, and then, an efficient message passing algorithm is successfully developed to decoupling echoes at the BS. It is also demonstrated that the communication performance can be enhanced by making full use of imaging echoes.

The article is organized as follows. Section II briefly introduces the system and the signal model used. In Section III, the phase design problem is formulated and a back propagation based optimization method is proposed. In Section IV, a joint maximum a posteriori estimation problem is established. Then the problem is represented by a graph model, based on which we propose an efficient message passing algorithm to achieve joint communication and imaging. Numerical results are provided in Section V. Section VI concludes this article.

Notations: Throughout this paper, column vectors and matrices are denoted by bold lower-case and bold upper-case letters, respectively. The notation ()H(\cdot)^{\mathrm{H}} denotes the conjugate operator and ()1(\cdot)^{-1} denotes inverse operator. |||\cdot|, \|\cdot\| and F\|\cdot\|_{\mathrm{F}} are the l2l_{2} norm, the l2l_{2} norm and the Frobenius norm, respectively. The transpose operation is denoted by ()T(\cdot)^{\mathrm{T}}. We use 𝔼{}\mathbb{E}\{\cdot\} to denote the expectation operator. The notation 𝒩(x;𝐦,𝐕)\mathcal{N}(\textbf{x};\mathbf{m},\mathbf{V}) denotes a Gaussian distribution of x with mean vector 𝐦\mathbf{m} and covariance matrix 𝐕\mathbf{V}. In many cases, we also use the inverse of the covariance matrix 𝐕\mathbf{V}, which is denoted by 𝐖\mathbf{W} and called weight matrix as in [52, 53]. The notation 𝐈a\mathbf{I}_{a} denotes an identity matrix of a size a×aa\times a. The notation \propto denotes equality of functions up to a scale factor. We use \odot to represent the inner product, and a()\int_{\sim a}\mathcal{F}(\cdot) denotes integral over all variables in ()\mathcal{F}(\cdot) except aa. The notation diag()\mathrm{diag(\cdot)} denotes the diagonal operation and δ()\delta(\cdot) is the Dirac delta function.

II System Model

Refer to caption
Figure 1: Diagram of a RIS-assisted joint uplink imaging and communication system.

We consider an uplink RIS-assisted ISAC system in 2D plane as shown in Fig.1, which consists of a user equipment (UE), a BS, a RIS and a region of interest (RoI). Assume there is an obstacle between the UE and the BS. The RIS is deployed between the UE and BS to customize the radio environment for the signal from the UE side and reflect to the BS, and its center is considered as the origin of the coordinate system. The UE sends communication signals to the BS with the assistant of the RIS. Due to the high path loss, it is also assumed that the signals with more than two reflections are negligible. We aim to perform communication at the BS while realizing reconstruction of the RoI.

The UE and the BS are located at the far field of the RIS, in directions of θU{\theta_{U}} and θB{\theta_{B}} from the origin, respectively, and both of which are assumed to adopt the single-antenna structure. The RIS is a uniform linear array (ULA) with NN elements that are half a wavelength apart. The RoI is divided into an equi-spaced grid with MM pixel units. We denote the angle vector of the RoI from the origin as 𝜽S=[θS,1,θS,2,,θS,M]T\bm{\theta}_{S}=\left[\theta_{S,1},\theta_{S,2},\ldots,\theta_{S,M}\right]^{T}.

It is noted that at the transmitted signals reach the BS through two paths, i.e., UE-RIS-BS, namely communication link and UE-RoI-RIS-BS, namely imagine link. Supposing that the UE transmits a communication frame of length LL, due to the difference of propagation delay, the signal received at the BS falls into 3 categories: 1) non-overlapped communication echo, 2) overlapped communication and imaging echo and 3) non-overlapped imaging echo. A generic model of received signal is shown in Fig.2. Specifically, the sampled base band received signal at time index tt can be expressed as

y(t)={αc𝐠T𝚯(t)𝐡𝐜x(t)UERISBS+w(t),tkαc𝐠T𝚯(t)𝐡𝐜x(t)UERISBS+αI𝐠T𝚯(t)𝐇𝐈𝝈x(tk)UERoIRISBS+w(t),k<t<LαI𝐠T𝚯(t)𝐇𝐈𝝈x(tk)UERoIRISBS+w(t),L<tL+k\displaystyle y(t)\!=\!\left\{\!\!{\begin{array}[]{l}{\underbrace{{\alpha_{c}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}x(t)}_{UE\to RIS\to BS}+w(t),\quad t\leq k}\\ \\ {\underbrace{{\alpha_{c}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}x(t)}_{UE\to RIS\to BS}+\underbrace{{\alpha_{I}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}{\bm{\sigma}}x(t-k)}_{UE\to RoI\to RIS\to BS}+w(t)},\\ \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad{k<t<L}\\ {\underbrace{{\alpha_{I}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}{\bm{\sigma}}x(t-k)}_{UE\to RoI\to RIS\to BS}+w(t),\quad L<t\leq L+k}\end{array}}\right. (1)

where we define the variables and symbols list below

  • αc\alpha_{c} and αI\alpha_{I} denote the attenuation coefficients of communication and imaging links, respectively, which are assumed to have been obtained by the pilot signals.

  • The steering vectors are 𝐠=𝐚(θB)N×1{\bf{g}}=\mathbf{a}({\theta_{B}})\in{\mathbb{C}^{N\times 1}}, 𝐡𝐜=𝐚(θU)N×1{{\bf{h}}_{\bf{c}}}=\mathbf{a}({\theta_{U}})\in{\mathbb{C}^{N\times 1}} and 𝐇𝐈=𝐚(𝜽S)N×M{{\mathbf{H}}_{\bf{I}}}=\mathbf{a}(\bm{\theta}_{S})\in{\mathbb{C}^{N\times M}}, where 𝐚(θ)=[1,ejπsin(θ),,ejπ(N1)sin(θ)]T{\mathbf{a}}(\theta)={\left[{1,{e^{j{{\pi}}\sin(\theta)}},\ldots,{e^{j{{\pi}}\left({N-1}\right)\sin(\theta)}}}\right]^{T}}.

  • 𝚯(t)=diag(ejθt,1,,ejθt,N){{\bf{\Theta}}(t)={\rm{diag}}\left({{e^{j{\theta_{t,1}}}},\ldots,{e^{j{\theta_{t,N}}}}}\right)} denotes the phase shift matrix of RIS at the time index tt, where θt,n{0,2π2nbit,,2π(2nbit1)2nbit}{\theta_{t,n}}\in\{0,\frac{{2\pi}}{{2^{n_{bit}}}},...,\frac{{2\pi(2^{n_{bit}}-1)}}{{2^{n_{bit}}}}\}, nbitn_{bit} denotes the number of discrete values [54, 45]. When nbit+n_{bit}\to+\infty, it means that the phase shift of RIS can be controlled continuous.

  • kk indicates the time delay between communication and imaging links, which is assumed to has been estimated by the pilot signals [55].

  • x(t)x(t) denotes the communication symbol sent by the UE. In this paper, we consider QPSK codes for communication, i.e., x(t){exp(jπ4),exp(j3π4),exp(j5π4),exp(j7π4)}x(t)\in\{\exp(j\frac{\pi}{4}),\exp(j\frac{3\pi}{4}),\exp(j\frac{5\pi}{4}),\exp(j\frac{7\pi}{4})\}.

  • 𝝈M×1\bm{\sigma}\in{\mathbb{C}^{M\times 1}} denotes the vector of scattering coefficients of the RoI, which is assumed to be sparse.

  • w(t)w(t) is the i.i.d. complex Gaussian noise with zero mean and variance ξ2\xi^{2}.

Refer to caption
Figure 2: Diagram of overlapped receive signal model.

III Phase Shift Design Scheme

In this section, a phase optimization problem is established based on the requirement of the joint communication and imaging system. We then present a back propagation based phase design scheme in the both cases of continuous and discrete phase models.

III-A Phase Shifts Optimization Formulation

The phase shift {𝚯(t)}\left\{\mathbf{\Theta}(t)\right\} is critical since it determines the gain of UE-RIS-BE and RoI-RIS-BS links simultaneously, and both communication and imaging performance should be considered in optimizing {𝚯(t)}\left\{\mathbf{\Theta}(t)\right\}. For communication functionality, we hope to reflect more radiated energy from the UE to the BS after modulation by the RIS, i.e., maximizing 𝐠T𝚯(t)𝐡𝐜{\bf{g}}^{\rm{T}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}. For imaging functionality, on the one hand, we hope to gather more radiated energy reflected by the RoI and reflect it to the BS, i.e., maximizing 𝐠T𝚯(t)𝐇𝐈{{{\bf{g}}^{T}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}}. On the other hand, according to the compressed sensing theory [56, 57], when reconstructing 𝝈\bm{\sigma}, we need to modify {𝚯(t),k<tL+k}\left\{\mathbf{\Theta}(t),k<t\leq L+k\right\} to minimize the correlation between columns of the sensing matrix which given by

𝐆=[𝐠T𝚯(k)𝐇𝐈𝐠T𝚯(L+k)𝐇𝐈].\begin{array}[]{*{20}{l}}{\bf{G}}=\left[{\begin{array}[]{*{20}{c}}{{{\bf{g}}^{T}}{\bf{\Theta}}(k){{\bf{H}}_{\bf{I}}}}\\ {\vdots}\\ {{{\bf{g}}^{T}}{\bf{\Theta}}(L+k){{\bf{H}}_{\bf{I}}}}\end{array}}\right].\end{array} (2)

The orthogonality of 𝐆\bf{G} can be evaluated by the mean of non-diagonal elements of (𝐆)\mathcal{R}(\bf{G}), i.e., (𝐆)𝐈MF/(M2M)\left\|{{\mathcal{R}}({\bf{G}})-{\bf{I}}_{M}}\right\|_{\mathrm{F}}/(M^{2}-M), where (𝐆)\mathcal{R}(\bf{G}) denotes the correlation coefficient matrix of columns of 𝐆\bf{G}. Based on the above consideration, we can formulate the following phase shift optimization problem

max𝚯(t)t=1L+kρ𝐠T𝚯(t)𝐡𝐜2+(1ρ)𝐠T𝚯(t)𝐇𝐈2\displaystyle\!\!\!\!\!\!\!\!\!\!\!\max_{\bm{\Theta}(t)}\sum\limits_{t=1}^{L+{\rm{k}}}{\rho{{{\left\|{{{\bf{g}}^{T}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}}\right\|}}^{2}}+(1-\rho){{{\left\|{{{\bf{g}}^{T}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}}\right\|}}^{2}}} (3a)
s.t. (𝐆)𝐈MFM2Mη0,\displaystyle\frac{\left\|{{\mathcal{R}}({\bf{G}})-{\bf{I}}_{M}}\right\|_{\mathrm{F}}}{M^{2}-M}\leq{\eta_{0}}, (3b)
θt,n{0,2π2nbit,,2π(2nbit1)2nbit}.\displaystyle{{\theta_{t,n}}\in\left\{0,\frac{{2\pi}}{{{2^{{n_{bit}}}}}},...,\frac{{2\pi({2^{{n_{bit}}}}-1)}}{{{2^{{n_{bit}}}}}}\right\}}. (3c)

The cost function (3a) indicates that we expect to concentrate energy in the RoI and the direction of communication as much as possible, where ρ[0,1]\rho\in[0,1] is used to realize a trade-off between communication and imaging performance. η0{\eta_{0}} is a threshold between 0 and 1 to constrain the orthogonality of 𝐆\bf{G}. The constraint (3c) denotes the feasible set of phase shift. We noted that (3) is non-convex, highly nonlinear and there are numerous coupled variables to be optimized, making it highly challenging to solve the problem directly.

III-B Back Propagation Based Phase Shift Optimization

Refer to caption
Figure 3: The flowchart of the proposed phase optimization method for continuous and discrete phase model.

Note that (3b) indicates a strongly coupled constraint on {𝚯(t),1<tL+k}\left\{\mathbf{\Theta}(t),1<t\leq L+k\right\}, while the constraint (3c) makes the problem even more challenging. To solve this problem, we relax the original problem (3) into two subproblems and turn to find the sub-optimal solution. Specifically, we propose to first initialize {𝚯(t),1tL+k}\left\{\mathbf{\Theta}(t),1\leq t\leq L+k\right\} by minimizing (𝐆)𝐈M/(M2M)\left\|{{\mathcal{R}}({\bf{G}})-{\bf{I}}_{M}}\right\|/(M^{2}-M) subject to the constraint (3c). With an appropriate start point of {𝚯(t),1tL+k}\left\{\mathbf{\Theta}(t),1\leq t\leq L+k\right\}, we can further modify it to maximize (3a) until (3b) and (3c) can not be satisfied. With this in mind, we first need to solve the following problem

min𝚯(t)\displaystyle\min_{\bm{\Theta}(t)} (𝐆)𝐈MFM2M\displaystyle\ \frac{\left\|{{\mathcal{R}}({\bf{G}})-{\bf{I}}_{M}}\right\|_{\mathrm{F}}}{M^{2}-M} (4)
s.t. (3c).\displaystyle\eqref{const_1b}.

However, the above problem is still challenging due to the coupled cost function. Since Gaussian matrix is recognized to be suit for a sensing matrix [58], we convert (4) into the following problem

min𝚯(t)\displaystyle\min_{\bm{\Theta}(t)} 𝐑L×N𝐆F\displaystyle\ \ {\left\|{\mathbf{R}_{L\times N}-{\bf{G}}}\right\|_{\rm{F}}} (5)
s.t. (3c).\displaystyle\eqref{const_1b}.

where 𝐑L×N\mathbf{R}_{L\times N} is a Gaussian matrix in the size of L×NL\times N.

Inspired by the color section method in [59], we introduce an increasing temperature parameter α\alpha and the softmax function into back propagation framework to solve (5), where the temperature parameter in the llth iteration is calculated by α=1+(rl)2\alpha=1+{(rl)^{2}} and rr is a factor to adjust the increasing rate, and the softmax operation of α|𝐰nl|\alpha|\mathbf{w}_{n}^{l}| is defined as

Fsoftmax(α|𝐰nl|)=exp(α|𝐰nl|)s=12nbitexp(α|wn,sl|)F_{soft\max}(\alpha|\mathbf{w}_{n}^{l}|)=\frac{{\exp(\alpha|\mathbf{w}_{n}^{l}|)}}{\sum_{s=1}^{{2^{{n_{bit}}}}}{\exp(\alpha\left|{w_{n,s}^{l}}\right|)}} (6)

The corresponding flowchart of the proposed method is shown in Fig. 3, where {𝚯(t)}\left\{\mathbf{\Theta}(t)\right\} can be updated by alternating iterative forward and backward propagation procedures. In the forward propagation, assume that 𝐰n\mathbf{w}_{n} in the (l1)(l-1)th iteration has been obtained, which denoted by 𝐰nl1\mathbf{w}_{n}^{l-1}, then 𝐰~nl\mathbf{\tilde{w}}_{n}^{l} is calculated by making softmax operation for α|𝐰nl1|\alpha|\mathbf{w}_{n}^{l-1}|. It is noteworthy that α\alpha increases with the iteration, and the softmax operation essentially plays a role of selector. Specifically, the softmax operation will push the larger element in 𝐰𝐧\mathbf{w_{n}} approaching 1 while the smaller elements approaching 0 with the increase of α\alpha. Then θn\theta_{n} can be obtained by making inner product of 𝐰~nl\mathbf{\tilde{w}}_{n}^{l} and 𝐬\mathbf{s}. After applying the above procedure for all {θ}\{\theta\}, we can update {𝚯(t)}\left\{\mathbf{\Theta}(t)\right\} in parallel thereby calculate the loss in the llth iteration. In the backward propagation, the chain rule is utilized to update the gradients of all variables except 𝐬\mathbf{s} and thereby update 𝐰\mathbf{w} by stochastic gradient descent (SGD) [60]. We can finally obtain {𝚯(t)}\left\{\mathbf{\Theta}(t)\right\} along with corresponding 𝐰\mathbf{w} after the forward/backward propagation iterated alternately until convergence. The steps of the proposed method for solving (5) is summarized in Algorithm 1.

Based the initialized {𝚯(t)}\left\{\mathbf{\Theta}(t)\right\} , we need to further modify {𝚯(t)}\left\{\mathbf{\Theta}(t)\right\} by solving the following problem

min𝚯(t)\displaystyle\min_{\bm{\Theta}(t)} t=1L+k1ρ𝐠T𝚯(t)𝐡𝐜2+(1ρ)𝐠T𝚯(t)𝐇𝐈2\displaystyle\sum\limits_{t=1}^{L+{\rm{k}}}\frac{1}{{\rho{{{\left\|{{{\bf{g}}^{T}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}}\right\|}}^{2}}+(1-\rho){{{\left\|{{{\bf{g}}^{T}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}}\right\|}}^{2}}}} (7)
s.t. (3c).\displaystyle\eqref{const_1b}.

Note that (3a) is rewritten to be a minimized form to facilitate the calculation of loss. Fortunately, it is also worth noting that (7) can be solved in the same way as (4). All we need to do is execute Algorithm 1, where the loss is calculated by the cost function in (7) and the iteration condition Γ\Gamma is given by (3b).

In practical, when the phase shift can be controlled continuously, i.e., nbit+n_{bit}\to+\infty, the constraint (3c) can be removed from (3), which makes (3) greatly simplified. In this case, we can first initialize 𝚯\mathbf{\Theta} by minimizing (𝐆)𝐈M/(M2M)\left\|{{\mathcal{R}}({\bf{G}})-{\bf{I}}_{M}}\right\|/(M^{2}-M). Then {𝚯(t)}\left\{\mathbf{\Theta}(t)\right\} can be further optimized by minimize t=1L+k1ρ𝐠T𝚯(t)𝐡𝐜2+(1ρ)𝐠T𝚯(t)𝐇𝐈2\sum\limits_{t=1}^{L+{\rm{k}}}\frac{1}{{\rho{{{\left\|{{{\bf{g}}^{T}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}}\right\|}}^{2}}+(1-\rho){{{\left\|{{{\bf{g}}^{T}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}}\right\|}}^{2}}}}, which can be solved efficient by back propagation procedure in the dashed box of Fig.3.

Algorithm 1 Discrete Phase Optimization Algorithm
1:  Initialize the parameters: weights 𝐰t,n1\mathbf{w}_{t,n}^{1}, learning rate η\eta and scale factor rr
2:  while Γ\Gamma  do
3:     α=1+(rl)2\alpha=1+{(rl)^{2}}
4:     n,𝐰~t,nl=softmax(α|𝐰t,nl|){\forall n},{{\bf{\tilde{w}}}_{t,n}^{l}}=soft\max(\alpha|{{{\bf{w}}_{t,n}^{l}}}|)
5:     n,θt,nl=𝐰~t,nl[0,2π2nbit,,2π(2nbit1)2nbit]{\forall n},{{\theta_{t,n}^{l}}}={{\bf{\tilde{w}}}_{t,n}^{l}}\odot[0,\frac{{2\pi}}{{{2^{{n_{bit}}}}}},...,\frac{{2\pi({2^{{n_{bit}}}}-1)}}{{{2^{{n_{bit}}}}}}]
6:     𝚯l(t)=diag(ejθt,1l,,ejθt,Nl){{\bf{\Theta}}^{l}(t)={\rm{diag}}\left({{e^{j{\theta_{t,1}^{l}}}},\ldots,{e^{j{\theta_{t,N}^{l}}}}}\right)}
7:     Ctl=Costfunction(𝚯l(t))C_{t}^{l}=Cost\ function({\bf{\Theta}}^{l}(t))
8:     𝐖tl+1=SGD(𝐖tl,η,Ctl𝐖tl)\mathbf{W}_{t}^{l+1}=SGD\left(\mathbf{W}_{t}^{l},\eta,\frac{\partial C_{t}^{l}}{\partial\mathbf{W}_{t}^{l}}\right)
9:     l=l+1l=l+1
10:  end while

IV Message Passing Based Joint Communication and Imaging

In this part, we investigate joint communication and imaging at the BS. The problem is first formulated into a maximum posteriori probability problem, which is nonconvex and strongly coupled. To solve it, a factor graph representation is then established, based on which an efficient iterative message passing algorithm is derived to estimate 𝝈\bm{\sigma} and 𝐱\mathbf{x} jointly .

IV-A Joint Communication and Imaging Formulation

According to the Bayes’ theorem, the joint posterior distribution of 𝐱\mathbf{x} and 𝝈\bm{\sigma} conditioned on 𝐲\mathbf{y} is given by

p(𝐱,𝝈|𝐲)\displaystyle p(\mathbf{x},\bm{\sigma}|\mathbf{y}) p(𝐲|𝐱,𝝈)p(𝝈)p(𝐱)\displaystyle\propto p(\mathbf{y}|\mathbf{x},\bm{\sigma})p(\bm{\sigma})p(\mathbf{x}) (8)

where 𝐲=[y(1),,y(L+k)]T(L+k)×1\mathbf{y}=[y(1),...,y(L+k)]^{\mathrm{T}}\in\mathbb{C}^{(L+k)\times 1} and 𝐱=[x(1),,x(L)]TL×1\mathbf{x}=[x(1),...,x(L)]^{\mathrm{T}}\in\mathbb{C}^{L\times 1}; p(𝐱)p(\mathbf{x}) and p(𝝈)p(\bm{\sigma}) denote the prior for 𝐱\mathbf{x} and 𝝈\bm{\sigma}, respectively. For QPSK modulation, we have

p(𝐱)=t=1Lp(x(t))=t=1L14i=14δ[x(t)ej(π2iπ4)].p(\mathbf{x})=\prod_{t=1}^{L}p(x(t))=\prod_{t=1}^{L}\frac{1}{4}\sum_{i=1}^{4}\delta[x(t)-e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}]. (9)

Following [61], it is assumed that the elements in 𝝈\bm{\sigma} are independent and following the two-layer sparsity-promoting prior, i.e.,

p(𝝈)\displaystyle p(\bm{\sigma}) =p(𝝈|𝜸)p(𝜸|ϵ,η)\displaystyle=p(\bm{\sigma}|\bm{\gamma})p(\bm{\gamma}|\epsilon,\eta) (10)
=m=1M𝒩(σm|0,γm1)m=1MGa(γm1|ϵ,η),\displaystyle=\prod_{m=1}^{M}\mathcal{N}\left(\sigma_{m}|0,\gamma_{m}^{-1}\right)\cdot\prod_{m=1}^{M}Ga(\gamma_{m}^{-1}|\epsilon,\eta),

where the precision vector 𝜸=[γ1,,γM]H\bm{\gamma}=[\gamma_{1},...,\gamma_{M}]^{\mathrm{H}}.

Based on (1), p(𝐲|𝐱,𝝈)p(\mathbf{y}|\mathbf{x},\bm{\sigma}) can be decomposed into the product of three categories of probability, i.e.,

p(𝐲|𝐱,𝝈)\displaystyle p(\mathbf{y}|\mathbf{x},\bm{\sigma}) =t=1kp(y(t)|x(t))t=k+1Lp((y(t)|x(t),x(tk),𝝈)\displaystyle\!=\!\prod_{t=1}^{k}p\left(y(t)|x(t)\right)\!\cdot\!\prod_{t=k+1}^{L}p\left((y(t)|x(t),x(t-k),\bm{\sigma}\right) (11)
t=L+1L+kp(y(t)|x(tk),𝝈),\displaystyle\cdot\!\prod_{t=L+1}^{L+k}p\left(y(t)|x(t-k),\bm{\sigma}\right),

where

p(y(t)|x(t))exp{[y(t)αc𝐠T𝚯(t)𝐡𝐜x(t)]2ξ2},tk\displaystyle\small p\left(y(t)|x(t)\right)\propto\exp\left\{-\frac{\left[y(t)-{\alpha_{c}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}x(t)\right]^{2}}{\xi^{2}}\right\},t\leq k (12)
p(y(t)|x(t),x(tk),𝝈)\displaystyle p(y(t)|x(t),x(t-k),\bm{\sigma})\propto (13)
exp{[y(t)αc𝐠T𝚯(t)𝐡𝐜x(t)αI𝐠T𝚯(t)𝐇𝐈𝝈x(tk)]2ξ2},\displaystyle\exp\left\{-\frac{\left[y(t)-{\alpha_{c}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}x(t)-{\alpha_{I}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}{\bm{\sigma}}x(t-k)\right]^{2}}{\xi^{2}}\right\},
k<tL\displaystyle\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad k<t\leq L
p(\displaystyle p( y(t)|x(tk),𝝈)exp{[y(t)αI𝐠T𝚯(t)𝐇𝐈𝝈x(tk)]2ξ2},\displaystyle y(t)|x(t-k),\bm{\sigma})\propto\exp\left\{-\frac{\left[y(t)-{\alpha_{I}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}{\bm{\sigma}}x(t-k)\right]^{2}}{\xi^{2}}\right\}, (14)
L<tL+k\displaystyle\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad L<t\leq L+k

Based on the decomposition in (9) - (14), we aim to find 𝐱\mathbf{x} and 𝝈\bm{\sigma} that maximize the joint posterior distribution p(𝐱,𝝈|𝐲)p(\mathbf{x},\bm{\sigma}|\mathbf{y}) in (8). However, it is a challenging problem as 𝝈\bm{\sigma} and 𝐱\mathbf{x} is strongly coupled in different time segments.

IV-B Factor Graph Representation for Joint Imaging and Communication

TABLE I: Definitions of Variables and Functions
at=Δαc𝐠T𝚯(t)𝐡𝐜xta_{t}\buildrel\Delta\over{=}{\alpha_{c}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}x_{t}
bt=ΔαI𝐠T𝚯(t)𝐇𝐈xtk{b_{t}}\buildrel\Delta\over{=}{{\alpha_{I}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t)}{{\bf{H}}_{\bf{I}}}x_{t-k}
ct=Δat+bt{c_{t}}\buildrel\Delta\over{=}a_{t}+{b_{t}}
ht=ΔαI𝐠T𝚯(t+k)𝐇𝐈𝝈{h_{t}}\buildrel\Delta\over{=}{\alpha_{I}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t+k){{\bf{H}}_{\bf{I}}}\bm{\sigma}
fA,t1(at)=Δ𝒩(yt;at,ξ2)f_{A,t}^{1}(a_{t})\buildrel\Delta\over{=}{\cal N}(y_{t};a_{t},{\xi^{2}})
fA,t2(ct)=Δ𝒩(yt;ct,ξ2)f_{A,t}^{2}(c_{t})\buildrel\Delta\over{=}{\cal N}(y_{t};c_{t},{\xi^{2}})
fA,t3(bt)=Δ𝒩(yt;bt,ξ2)f_{A,t}^{3}(b_{t})\buildrel\Delta\over{=}{\cal N}(y_{t};b_{t},{\xi^{2}})
fB,t(at,xt)=Δδ(αc𝐠T𝚯(t)𝐡𝐜xtat)f_{B,t}\left(a_{t},x_{t}\right)\!\buildrel\!\Delta\over{=}\delta({\alpha_{c}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}x_{t}-a_{t})
fD,t(at,bt,ct)=Δδ(at+btct){f_{D,t}}(a_{t},b_{t},c_{t})\buildrel\Delta\over{=}\delta(a_{t}+{b_{t}}-{c_{t}})
fE,t(xt)=Δ14i=14δ[xtej(π2iπ4)]{f_{E,t}}\left(x_{t}\right)\buildrel\Delta\over{=}\frac{1}{4}\sum_{i=1}^{4}\delta[x_{t}-e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}]
fK,t(ht,xt,bt+k)=Δδ(htxtbt+k){f_{K,t}}\left(h_{t},x_{t},b_{t+k}\right)\buildrel\Delta\over{=}\delta({h_{t}}x_{t}-{b_{t+k}})
fH,t(𝝈,ht)=Δδ(αI𝐠T𝚯(t+k)𝐇𝐈𝝈ht){f_{H,t}}(\bm{\sigma},h_{t})\buildrel\Delta\over{=}\delta({\alpha_{I}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t+k){{\bf{H}}_{\bf{I}}}{\bm{\sigma}}-{h_{t}})
fG(𝝈,𝜸)=ΔM𝒩(σm;0,γm1){f_{G}}(\bm{\sigma},\bm{\gamma})\buildrel\Delta\over{=}\prod_{M}{{\cal N}({\sigma_{m}};{0,{\gamma}_{m}^{-1}})}
fF(𝜸)MGa(γm;ϵ,η)f_{F}(\bm{\gamma})\triangleq\prod_{M}Ga\left(\gamma_{m};\epsilon,\eta\right)

Factor graphs and message passing are a powerful tool for inference and estimation. In this subsection, a factor graph model for (8) is established, based on which we can then jointly estimate 𝐱\mathbf{x} and 𝝈{\bm{\sigma}}. To simplify the notations, we use xtx_{t}, yty_{t} and 𝚯t\bm{\Theta}_{t} to represent x(t)x(t), y(t)y(t) and 𝚯(t)\bm{\Theta}(t) respectively and we also define some variables and functions in Table I. With these abbreviations and definitions, (9), (10) and (12) - (14) can be respectively rewritten as

p(𝐱)=t=1LfE,t(xt),\displaystyle p(\mathbf{x})=\prod_{t=1}^{L}{f_{E,t}}\left(x_{t}\right), (15)
p(𝝈)\displaystyle p(\bm{\sigma}) =fG(𝝈,𝜸)fF(𝜸),\displaystyle=f_{G}(\bm{\sigma},\bm{\gamma})f_{F}(\bm{\gamma}), (16)
p(yt|xt)\displaystyle p(y_{t}|x_{t}) =p(yt|at)p(at|xt)\displaystyle=p(y_{t}|a_{t})p(a_{t}|x_{t}) (17)
=fA,t1(at)fB,t(at,xt),tk\displaystyle=f_{A,t}^{1}({a_{t}})f_{B,t}({a_{t}},x_{t}),t\leq k
p(yt|xt,xtk,𝝈)\displaystyle p\left(y_{t}|x_{t},x_{t-k},\bm{\sigma}\right) (18)
=\displaystyle= p(yt|ct)p(ct|at,bt)p(at|xt)p(bt|xtk,htk)p(htk|𝝈)\displaystyle\ p(y_{t}|c_{t})p(c_{t}|a_{t},b_{t})p(a_{t}|x_{t})p(b_{t}|x_{t-k},h_{t-k})p(h_{t-k}|\bm{\sigma})
=\displaystyle= fA,t2(ct)fD,t(at,bt,ct)fB,t(at,xt)\displaystyle\ f_{A,t}^{2}({c_{t}})f_{D,t}({a_{t}},{b_{t}},{c_{t}})f_{B,t}\left({a_{t}},x_{t}\right)
fK,tk(htk,xtk,bt)fH,tk(𝝈,htk),k<tL\displaystyle\!\cdot f_{K,t-k}\left({h_{t-k}},x_{t-k},{b_{t}}\right)f_{H,t-k}(\bm{\sigma},h_{t-k}),k<t\leq L

and

p(yt|xtk,𝝈)=p(yt|bt)p(bt|xtk,htk)p(htk|𝝈)\displaystyle p\left(y_{t}|x_{t-k},\bm{\sigma}\right)=p\left(y_{t}|b_{t}\right)p\left(b_{t}|x_{t-k},h_{t-k}\right)p(h_{t-k}|\bm{\sigma}) (19)
=fA,t3(bt)fK,tk(htk,xtk,bt)fH,tk(𝝈,htk),\displaystyle=f_{A,t}^{3}({b_{t}})f_{K,t-k}\left({h_{t-k}},x_{t-k},{b_{t}}\right)f_{H,t-k}(\bm{\sigma},h_{t-k}),
L<kL+k.\displaystyle L<k\leq L+k.

Then the joint distribution of 𝐱,𝝈,at,bt,ct,ht\mathbf{x},\bm{\sigma},a_{t},b_{t},c_{t},h_{t} and 𝜸\bm{\gamma} conditioned on 𝐲\mathbf{y} can be given by

p(𝐱,𝝈,at,bt,ct,ht,𝜸|𝐲)\displaystyle p(\mathbf{x},\bm{\sigma},a_{t},b_{t},c_{t},h_{t},\bm{\gamma}|\mathbf{y}) (20)
\displaystyle\propto fE(𝐱)fG(𝝈,𝜸)fF(𝜸)\displaystyle f_{E}(\mathbf{x})f_{G}(\bm{\sigma},\bm{\gamma})f_{F}(\bm{\gamma})\cdot
t=1kfA,t1(at)fB,t(at,xt)\displaystyle\cdot\prod_{t=1}^{k}f_{A,t}^{1}({a_{t}})f_{B,t}({a_{t}},x_{t})
t=k+1LfA,t2(ct)fD,t(at,bt,ct)fB,t(at,xt)\displaystyle\cdot\prod_{t=k+1}^{L}f_{A,t}^{2}({c_{t}})f_{D,t}({a_{t}},{b_{t}},{c_{t}})f_{B,t}({a_{t}},x_{t})
fK,tk(htk,xtk,bt)fH,tk(𝝈,htk)\displaystyle\qquad\quad\ \cdot f_{K,t-k}({h_{t-k}},x_{t-k},{b_{t}})\cdot f_{H,t-k}(\bm{\sigma},h_{t-k})
t=L+1L+kfA,t3(bt)fK,tk(htk,xtk,bt)fH,tk(𝝈,htk).\displaystyle\cdot\prod_{t=L+1}^{L+k}f_{A,t}^{3}({b_{t}})f_{K,t-k}({h_{t-k}},x_{t-k},{b_{t}})f_{H,t-k}(\bm{\sigma},h_{t-k}).

Accordingly, the marginals of 𝐱\mathbf{x} and 𝝈\bm{\sigma} can be expressed by

(𝐱)={𝐱}p(𝐱,𝝈,at,bt,ct,ht,𝜸𝐲)\mathcal{M}(\mathbf{x})=\int_{\sim\{\mathbf{x}\}}p\left(\mathbf{x},\bm{\sigma},a_{t},b_{t},c_{t},h_{t},\bm{\gamma}\mid\mathbf{y}\right) (21)

and

(𝝈)={𝝈}p(𝐱,𝝈,at,bt,ct,ht,𝜸𝐲).\mathcal{M}(\bm{\sigma})=\int_{\sim\{\bm{\sigma}\}}p\left(\mathbf{x},\bm{\sigma},a_{t},b_{t},c_{t},h_{t},\bm{\gamma}\mid\mathbf{y}\right). (22)

Based on the factorization in (20), we take k=1k=1 as an example and establish the factor graph for joint communication and imaging as shown in Fig. 4, where we drop the augments of the functions for simplicity of notations.

It is noteworthy that though (𝐱)\mathcal{M}(\mathbf{x}) and (𝝈)\mathcal{M}(\bm{\sigma}) are given in (21) and (22), it is difficult to maximize them due to the high-dimensional integral. However, the factor graph in Fig.4 provides us a possible way to find approximation marginal functions of (𝐱)\mathcal{M}(\mathbf{x}) and (𝝈)\mathcal{M}(\bm{\sigma}), which are expected to allow easy maximization.

Refer to caption
Figure 4: Factor graph representation for joint imaging and communication problem as formulated by (20).

IV-C Message Passing Based Echoes Decoupling Algorithm

Considering that the factor graph contains loops, so iterative message passing is needed. We divide the graph into upper and lower parts by the dotted line, where the lower part represents the message update for 𝐱\mathbf{x} and the upper part represents the message update for 𝝈\bm{\sigma}. To estimate 𝐱\mathbf{x} and 𝝈\bm{\sigma} jointly, the upper and lower parts of Fig.4 are updated alternatively, and the adaptive Sparse Bayesian learning (SBL) algorithm in [61] is adopted to recover 𝝈\bm{\sigma}.

Based on the sum-product algorithm (SPA) [52, 53], the belief (or marginal function) of xtx_{t} can be calculated by the product of all incoming messages to xtx_{t}, i.e.,

B(xt)=μfE,txt(xt)μfB,txt(xt)μfK,txt(xt).B(x_{t})=\mu_{{f_{E,t}}\to{x_{t}}}(x_{t})\cdot{\mu_{{f_{B,t}}\to{x_{t}}}}(x_{t})\cdot\mu_{{f_{K,t}}\to{x_{t}}}(x_{t}). (23)

where

μfE,txt=14i=14δ[xtej(π2iπ4)].\mu_{{f_{E,t}}\to{x_{t}}}=\frac{1}{4}\sum_{i=1}^{4}\delta[x_{t}-e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}]. (24)

To obtain B(xt)B(x_{t}), we first need to derive messages in Fig. 4 from time index 1 to L+kL+k successively to update μfB,txt(xt){\mu_{{f_{B,t}}\to{x_{t}}}}(x_{t}). After that, we need to derive messages in reverse order of time index, i.e., form time index L+kL+k to 11, to update μfK,txt(xt)\mu_{{f_{K,t}}\to{x_{t}}}(x_{t}).

1) When tkt\leq k: Based on (1) and (12), we have

μfA,t1at(at)\displaystyle{\mu_{{f^{1}_{A,t}}\to{a_{t}}}}(a_{t}) =𝒩(at;mfA,t1at,WfA,t1at1)\displaystyle=\mathcal{N}\left({{a_{t}};{m_{{f^{1}_{A,t}}\to{a_{t}}}},W_{{f^{1}_{A,t}}\to{a_{t}}}^{-1}}\right) (25)
=𝒩(at;yt,ξ2).\displaystyle=\mathcal{N}\left({a_{t}};y_{t},\xi^{2}\right).

According to the structure of the factor graph, the message μatfB,t(at)\mu_{{a_{t}}\to{f_{B,t}}}(a_{t}) is forward to fB,tf_{B,t}, i.e.,

μatfB,t(at)=μfA,t1at(at)=𝒩(at;yt,ξ2).\mu_{{a_{t}}\to{f_{B,t}}}(a_{t})={\mu_{{f^{1}_{A,t}}\to{a_{t}}}}(a_{t})=\mathcal{N}\left({a_{t}};y_{t},\xi^{2}\right). (26)

Based on SPA, we have

μfB,txt(xt)\displaystyle{\mu_{{f_{B,t}}\to{x_{t}}}}(x_{t}) =fB,t(at,xt)μatfB,t(at)dat.\displaystyle=\int{{f_{B,t}}}\left({{a_{t}},{x_{t}}}\right)\cdot{\mu_{{a_{t}}\to{f_{B,t}}}}(a_{t}){\rm{d}}{a_{t}}. (27)

Since μatfB,t(at){\mu_{{a_{t}}\to{f_{B,t}}}}(a_{t}) is Gaussian message and fB,t(at,xt)f_{B,t}\left({{a_{t}},{x_{t}}}\right) indicates a linear relationship between ata_{t} and xtx_{t}, μfB,txt(xt)\mu_{{f_{B,t}}\to{x_{t}}}(x_{t}) also has a Gaussian form, i.e.,

μfB,txt(xt)=𝒩(xt;mfB,txt,WfB,txt1),\displaystyle{\mu_{{f_{B,t}}\to{x_{t}}}}(x_{t})=\mathcal{N}\left({{x_{t}};{m_{{f_{B,t}}\to{x_{t}}}},W_{{f_{B,t}}\to{x_{t}}}^{-1}}\right), (28)

where mfB,txt{m_{{f_{B,t}}\to{x_{t}}}} and WfB,txtW_{{f_{B,t}}\to{x_{t}}} are calculated by

WfB,txt=(αc𝐠T𝚯t𝐡𝐜)HWfA,t1at(αc𝐠T𝚯t𝐡𝐜),{W_{{f_{B,t}}\to{x_{t}}}}={({{\alpha_{c}\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}})^{\mathrm{H}}}{W_{{f^{1}_{A,t}}\to{a_{t}}}}({\alpha_{c}{\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}}), (29)

and

mfB,txt=WfB,txt1(αc𝐠T𝚯t𝐡𝐜)H(WfA,t1atmfA,t1at).{m_{{f_{B,t}}\to{x_{t}}}}=W_{{f_{B,t}}\to{x_{t}}}^{-1}{({\alpha_{c}{\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}})^{\mathrm{H}}}({W_{{f^{1}_{A,t}}\to{a_{t}}}}{m_{{f^{1}_{A,t}}\to{a_{t}}}}). (30)

The message outgoing from xtx_{t} can be calculated by

μxtfK,t(xt)=μfE,txt(xt)μfB,txt(xt)\displaystyle\mu_{x_{t}\to f_{K,t}}(x_{t})=\mu_{f_{E,t}\to x_{t}}(x_{t})\cdot\mu_{{f_{B,t}}\to{x_{t}}}(x_{t})\propto (31)
i=14δ(xtej(π2iπ4))𝒩(xt;mfB,txt,WfB,txt1)\displaystyle\sum_{i=1}^{4}\delta(x_{t}-e^{j(\frac{\pi}{2}i-\frac{\pi}{4})})\cdot\mathcal{N}\left({{x_{t}};{m_{{f_{B,t}}\to{x_{t}}}},W_{{f_{B,t}}\to{x_{t}}}^{-1}}\right)

Unfortunately, the message in (31) is no longer Gaussian, which makes the subsequent message calculation intractable. To this end, we approximate μxtfK,t(xt)\mu_{x_{t}\to f_{K,t}}(x_{t}) to a Gaussian message μ~xtfK,t=𝒩(xt;mxtfK,t,VxtfK,t)\tilde{\mu}_{x_{t}\to f_{K,t}}=\mathcal{N}\left({{x_{t}};{m_{x_{t}\to f_{K,t}}},V_{x_{t}\to f_{K,t}}}\right) by minimizing the Kullback-Leibler (KL) divergence between them. The optimal solution is given by moment matching, i.e.,

mxtfK,t=𝔼μxtfK,t(xt)\displaystyle m_{x_{t}\to f_{K,t}}=\mathbb{E}_{\mu_{x_{t}\to f_{K,t}}}(x_{t}) (32)
=14i=14exp(|ej(π2iπ4)mfB,txt|2WfB,txt),\displaystyle=\frac{1}{4}\sum\limits_{i=1}^{4}{\exp\left(-{\left|{e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}-{m_{{f_{B,t}}\to{x_{t}}}}}\right|^{2}}{W_{{f_{B,t}}\to{x_{t}}}}\right)},

and

VxtfK,t=𝔼μxtfK,t(|xt𝔼μxtfK,t(xt)|2)=\displaystyle V_{x_{t}\to f_{K,t}}=\mathbb{E}_{\mu_{x_{t}\rightarrow f_{K,t}}}\left(\left|x_{t}-\mathbb{E}_{\mu_{x_{t}\rightarrow f_{K,t}}}\left(x_{t}\right)\right|^{2}\right)= (33)
14i=14[exp(|ej(π2iπ4)mfB,txt|2WfB,txt)mxtfK,t]2\displaystyle\frac{1}{4}\sum_{i=1}^{4}\!\left[\exp\!\left(\!\!-\!\left|e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}-m_{f_{B,t}\rightarrow x_{t}}\right|^{2}\!W_{f_{B,t}\rightarrow x_{t}}\!\!\right)\!-\!m_{x_{t}\rightarrow f_{K,t}}\!\right]^{2}

2) When k<tLk<t\leq L : Suppose that the message from the variable node htkh_{t-k} to the factor node fK,tkf_{K,{t-k}}, i.e., μhtkfK,tk(htk)\mu_{h_{t-k}\to f_{K,{t-k}}}(h_{t-k}) is given in the last round of iteration, with its belief denoted by B(htk)B(h_{t-k}). According to the BP rule, the message outgoing from fK,tkf_{K,{t-k}} to btb_{t} can be expressed by

μfK,tkbt(bt)=fK,tk(htk,xtk,bt)\displaystyle\mu_{f_{K,{t-k}}\to b_{t}}(b_{t})=\int\!\!\!\int{{f_{K,{t-k}}}}\left({{h_{t-k}},{x_{t-k}},b_{t}}\right) (34)
μhtkfK,tk(htk)μxtkfK,tk(xtk)dhtkdxtk.\displaystyle\cdot{\mu_{{h_{t-k}}\to{f_{K,{t-k}}}}}(h_{t-k})\cdot{\mu_{{x_{t-k}}\to{f_{K,{t-k}}}}}(x_{t-k}){\rm{d}}{h_{t-k}}{\rm{d}}x_{t-k}.

It can be seen that the computation in (34) is difficult due to the involved complex integral. Although numerical calculation seems feasible, it will suffer from high computation complexity and make the subsequent message update intractable due to the lack of closed-form expression. To simplify the message computation, we approximate the incoming message μhtkfK,tk(htk)\mu_{{h_{t-k}}\to{f_{K,{t-k}}}}(h_{t-k}) to the Dirac delta function with proper location, i.e.,

μhtkfK,tk(htk)δ(htkmhtk),\mu_{{h_{t-k}}\to{f_{K,{t-k}}}}(h_{t-k})\approx\delta(h_{t-k}-m_{h_{t-k}}), (35)

where B(htk)B(h_{t-k}) is maximized at mhtkm_{h_{t-k}}. Substituting (35) into (34) yields

μfK,tkbt(bt)μ~fK,tkbt(bt)\displaystyle\mu_{f_{K,{t-k}}\to b_{t}}(b_{t})\approx\tilde{\mu}_{f_{K,{t-k}}\to b_{t}}(b_{t}) (36)
=\displaystyle= fK,tk(mhtk,xtk,bt)μxtkfK,tk(xtk)dxtk.\displaystyle\int{{f_{K,{t-k}}}}\left({m_{h_{t-k}},{x_{t-k}},b_{t}}\right)\cdot{\mu_{{x_{t-k}}\to{f_{K,{t-k}}}}}(x_{t-k}){\rm{d}}x_{t-k}.

It is not hard to show that μ~fK,tkbt(bt)\tilde{\mu}_{f_{K,{t-k}}\to b_{t}}(b_{t}) still preserves Gaussianity, i.e.,

μ~fK,tkbt(bt)=𝒩(bt;mfK,tkbt,VfK,tkbt),\tilde{\mu}_{f_{K,{t-k}}\to b_{t}}(b_{t})=\mathcal{N}\left({{b_{t}};{m_{{f_{K,{t-k}}}\to{b_{t}}}},V_{{f_{K,{t-k}}}\to{b_{t}}}}\right), (37)

where

VfK,tkbt=mhtkVxtkfK,tk(mhtk)HV_{{f_{K,{t-k}}}\to{b_{t}}}=m_{h_{t-k}}V_{x_{t-k}\to f_{K,{t-k}}}(m_{h_{t-k}})^{\text{H}} (38)

and

mfK,tkbt=mhtkmxtkfK,tk.m_{{f_{K,{t-k}}}\to{b_{t}}}=m_{h_{t-k}}\cdot m_{x_{t-k}\to f_{K,{t-k}}}. (39)

It is clear that

μbtfD,t(bt)=μfK,tkbt(bt).\displaystyle\mu_{b_{t}\to f_{D,t}}(b_{t})=\mu_{f_{K,{t-k}}\to b_{t}}(b_{t}). (40)

The message outgoing from the function node fD,tf_{D,t} to the variable node ata_{t} can be given by

μfD,tat(at)=\displaystyle\mu_{f_{D,t}\to a_{t}}(a_{t})= (41)
fD,t(at,bt,ct)μbtfD,t(bt)μctfD,t(ct)dbtdct\displaystyle\int\!\!\int{{f_{D,t}}}\left({a_{t},b_{t},c_{t}}\right)\cdot{\mu_{{b_{t}}\to{f_{D,t}}}}(b_{t})\cdot{\mu_{{c_{t}}\to{f_{D,t}}}}(c_{t}){\rm{d}}{b_{t}}{\rm{d}}c_{t}

where

μctfD,t(ct)=μfA,t2ct(ct)\displaystyle{\mu_{{c_{t}}\to{f_{D,t}}}}(c_{t})={\mu_{{f^{2}_{A,t}}\to{c_{t}}}}(c_{t}) (42)
=\displaystyle= 𝒩(ct;mfA,t2ct,VfA,t2ct)=𝒩(ct;yt,ξ2)\displaystyle\mathcal{N}\left({{c_{t}};{m_{{f^{2}_{A,t}}\to{c_{t}}}},V_{{f^{2}_{A,t}}\to{c_{t}}}}\right)=\mathcal{N}\left({c_{t}};y_{t},\xi^{2}\right)

With the combination of (40) and (42), we have

μfD,tat(at)=𝒩(at;mfD,tat,VfD,tat),\mu_{f_{D,t}\to a_{t}}(a_{t})=\mathcal{N}\left({a_{t}};{m_{f_{D,t}\to a_{t}},V_{f_{D,t}\to a_{t}}}\right), (43)

where

mfD,tat=mctfD,tmbtfD,t,{m_{{f_{D,t}}\to{a_{t}}}}={m_{{c_{t}}\to{f_{D,t}}}}-{m_{{b_{t}}\to{f_{D,t}}}}, (44)
VfD,tat=VctfD,t+VbtfD,t.{V_{{f_{D,t}}\to{a_{t}}}}={V_{{c_{t}}\to{f_{D,t}}}}+{V_{{b_{t}}\to{f_{D,t}}}}. (45)

and

Based on the structure of the factor graph, we have

μatfB,t(at)=μfD,tat(at),\mu_{a_{t}\to f_{B,t}}(a_{t})=\mu_{f_{D,t}\to a_{t}}(a_{t}), (46)

and

μfB,txt(xt)\displaystyle\mu_{f_{B,t}\to x_{t}}(x_{t}) =fB,t(at,xt)μatfB,t(at)𝑑at\displaystyle=\int f_{B,t}(a_{t},x_{t})\cdot\mu_{a_{t}\to f_{B,t}}(a_{t})da_{t} (47)
=𝒩(xt;mfB,txt,WfB,txt1)\displaystyle=\mathcal{N}\left({{x_{t}};{m_{{f_{B,t}}\to{x_{t}}}},W_{{f_{B,t}}\to{x_{t}}}^{-1}}\right)

where

WfBtxt=(αc𝐠T𝚯t𝐡𝐜)HVfD,tat1(αc𝐠T𝚯t𝐡𝐜),{W_{{f_{Bt}}\to{x_{t}}}}={({{\alpha_{c}\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}})^{\rm{H}}}{V_{{f_{D,t}}\to{a_{t}}}^{-1}}({\alpha_{c}{\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}}), (48)

and

mfB,txt=WfB,txt1(αc𝐠T𝚯t𝐡𝐜)H(VfD,tat1mfD,tat).{m_{{f_{B,t}}\to{x_{t}}}}=W_{{f_{B,t}}-{x_{t}}}^{-1}{({\alpha_{c}{\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}})^{\text{H}}}({V_{{f_{D,t}}\to{a_{t}}}^{-1}}{m_{{f_{D,t}}\to{a_{t}}}}). (49)

Note that the non-Gaussian incoming message to xtx_{t} will make the outgoing message μxtfK,t(xt)\mu_{x_{t}\to f_{K,t}}(x_{t}) no longer keeps the Gaussian form. We here approximate μxtfK,t(xt)\mu_{x_{t}\to f_{K,t}}(x_{t}) to a Gaussian form based on moment matching with the same way in (31) - (33), i.e.,

μxtfK,t(xt)=μfE,txt(xt)μfB,txt(xt)\displaystyle\mu_{x_{t}\to f_{K,t}}(x_{t})=\mu_{f_{E,t}\to x_{t}}(x_{t})\cdot\mu_{{f_{B,t}}\to{x_{t}}}(x_{t}) (50)
i=14δ[xtej(π2iπ4)]𝒩(xt;mfB,txt,WfB,txt1)\displaystyle\propto\sum_{i=1}^{4}\delta\left[x_{t}-e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}\right]\cdot\mathcal{N}\left({{x_{t}};{m_{{f_{B,t}}\to{x_{t}}}},W_{{f_{B,t}}\to{x_{t}}}^{-1}}\right)
μ~xtfK,t(xt)=𝒩(xt;mxtfK,t,VxtfK,t),\displaystyle\approx\tilde{\mu}_{x_{t}\to f_{K,t}}(x_{t})=\mathcal{N}\left({{x_{t}};{m_{x_{t}\to f_{K,t}}},V_{x_{t}\to f_{K,t}}}\right),

where mxtfK,tm_{x_{t}\to f_{K,t}} and VxtfK,tV_{x_{t}\to f_{K,t}} have same expressions as (32) and (33).

We have derived μxtfK,t(xt)\mu_{x_{t}\to f_{K,t}}(x_{t}) for tkt\leq k and k<tLk<t\leq L. We then need to derive μfK,txt(xt)\mu_{{f_{K,t}}\to{x_{t}}}(x_{t}) in different time segments. When L<tL+kL<t\leq L+k, based on (1) and (14), we have

μbtfK,tk(bt)=μfA,t3bt(bt)\displaystyle\mu_{b_{t}\to f_{K,t-k}}(b_{t})=\mu_{f^{3}_{A,t}\to b_{t}}(b_{t}) (51)
=𝒩(bt;mfA,t3bt,WfA,t3bt+k1)\displaystyle=\mathcal{N}\left({{b_{t}};{m_{{f^{3}_{A,t}}\to b_{t}}},W_{{f^{3}_{A,t}}\to b_{t+k}}^{-1}}\right)
=𝒩(bt;yt,ξ2).\displaystyle=\mathcal{N}\left(b_{t};y_{t},\xi^{2}\right).

When Lk<tLL-k<t\leq L, (51) can be recast as μbt+kfK,t(bt+k)=𝒩(bt+k;yt+k,ξ2)\mu_{b_{t+k}\to f_{K,t}}(b_{t+k})=\mathcal{N}\left(b_{t+k};y_{t+k},\xi^{2}\right). With the combination of the same approximation in (35), we have

μfK,txt(xt)=\displaystyle\mu_{f_{K,t}\to x_{t}}(x_{t})= fK,t(ht,xt,bt+k)μhtfK,t(ht)\displaystyle\int\!\!\int{{f_{K,t}}}\left({{h_{t}},{x_{t}},b_{t+k}}\right)\cdot{\mu_{{h_{t}}\to{f_{K,t}}}}(h_{t}) (52)
μbt+kfK,t(bt+k)dhtdxt\displaystyle\cdot\mu_{b_{t+k}\to f_{K,t}}(b_{t+k}){\rm{d}}{h_{t}}{\rm{d}}x_{t}
fK,t\displaystyle\approx\!\int{{f_{K,t}}} (mht,xt,bt+k)μbt+kfK,t(bt+k)dxt\displaystyle\left({{m_{h_{t}}},{x_{t}},b_{t+k}}\right)\cdot\mu_{b_{t+k}\to f_{K,t}}(b_{t+k}){\rm{d}}x_{t}
=𝒩\displaystyle=\mathcal{N} (xt;mfK,txt,WfK,txt1).\displaystyle\left({{x_{t}};{m_{{f_{K,t}}\to{x_{t}}}},W_{{f_{K,t}}\to{x_{t}}}^{-1}}\right).

where

mfK,txt=yt+k/mht,m_{{f_{K,t}}\to{x_{t}}}=y_{t+k}/m_{h_{t}}, (53)
WfK,txt=mfK,txtHmfK,txtξ2.W_{{f_{K,t}}\to{x_{t}}}=\frac{m^{\rm{H}}_{{f_{K,t}}\to{x_{t}}}m_{{f_{K,t}}\to{x_{t}}}}{\xi^{2}}. (54)

When k<tLkk<t\leq L-k, the message μxtfB,t(xt)\mu_{x_{t}\to f_{B,t}}(x_{t}) can be calculated by

μxtfB,t(xt)=μfE,txt(xt)μfK,txt(xt)\displaystyle\mu_{x_{t}\to f_{B,t}}(x_{t})=\mu_{f_{E,t}\to x_{t}}(x_{t})\cdot\mu_{{f_{K,t}}\to{x_{t}}}(x_{t}) (55)
i=14δ(xtej(π2iπ4))𝒩(xt;mfK,txt,WfK,txt1)\displaystyle\propto\sum_{i=1}^{4}\delta\left(x_{t}-e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}\right)\cdot\mathcal{N}\left({{x_{t}};{m_{{f_{K,t}}\to{x_{t}}}},W_{{f_{K,t}}\to{x_{t}}}^{-1}}\right)
𝒩(xt;mxtfB,t,VxtfB,t).\displaystyle\approx\mathcal{N}\left({{x_{t}};{m_{x_{t}\to f_{B,t}}},V_{x_{t}\to f_{B,t}}}\right).

where

mxtfB,t=𝔼μxtfB,t(xt)\displaystyle m_{x_{t}\to f_{B,t}}=\mathbb{E}_{\mu_{x_{t}\to f_{B,t}}}(x_{t}) (56)
=14i=14exp(|ej(π2iπ4)mfK,txt|2WfK,txt),\displaystyle=\frac{1}{4}\sum\limits_{i=1}^{4}{\exp\left(-{\left|{e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}-{m_{{f_{K,t}}\to{x_{t}}}}}\right|^{2}}{W_{{f_{K,t}}\to{x_{t}}}}\right)},

and

VxtfB,t=𝔼μxtfB,t(|xt𝔼μxtfB,t(xt)|2)=\displaystyle V_{x_{t}\to f_{B,t}}=\mathbb{E}_{\mu_{x_{t}\rightarrow f_{B,t}}}\left(\left|x_{t}-\mathbb{E}_{\mu_{x_{t}\rightarrow f_{B,t}}}\left(x_{t}\right)\right|^{2}\right)= (57)
14i=14[exp(|ej(π2iπ4)mfK,txt|2WfK,txt)mxtfB,t]2\displaystyle\frac{1}{4}\sum_{i=1}^{4}\!\left[\exp\!\left(\!\!-\left|e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}\!-\!m_{f_{K,t}\rightarrow x_{t}}\right|^{2}\!W_{f_{K,t}\rightarrow x_{t}}\!\right)\!-\!m_{x_{t}\rightarrow f_{B,t}}\!\right]^{2}

Note that the same approximation method in (31) - (33) is adopted here to approximate μxtfB,t(xt)\mu_{x_{t}\to f_{B,t}}(x_{t}) to a Gaussian form. Then μfB,tat(at)\mu_{f_{B,t}\to a_{t}}(a_{t}) can be calculated by BP, i.e.,

μfB,tat(at)=fB,t(at,xt)μxtfB,t(xt)𝑑xt\displaystyle\mu_{f_{B,t}\to a_{t}}(a_{t})=\int f_{B,t}\left(a_{t},x_{t}\right)\mu_{x_{t}\to f_{B,t}}(x_{t})dx_{t} (58)
=𝒩(at;mfB,tat,VfB,tat).\displaystyle=\mathcal{N}\left({{a_{t}};{m_{f_{B,t}\to a_{t}}},V_{f_{B,t}\to a_{t}}}\right).

where

mfB,tat=(αc𝐠T𝚯t𝐡𝐜)mxtfB,t{m_{{f_{B,t}}\to{a_{t}}}}=({\alpha_{c}}{{\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}}){m_{{x_{t}}\to{f_{B,t}}}} (59)

and

VfB,tat=(αc𝐠T𝚯t𝐡𝐜)VxtfB,t(αc𝐠T𝚯t𝐡𝐜)H{V_{{f_{B,t}}\to{a_{t}}}}=({\alpha_{c}}{{\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}}){V_{{x_{t}}\to{f_{B,t}}}}{({\alpha_{c}}{{\bf{g}}^{\rm{T}}}{{\bf{\Theta}}_{t}}{{\bf{h}}_{\bf{c}}})^{\rm{H}}} (60)

Based on SPA, we have μatfD,t(at)=μfB,tat(at)\mu_{a_{t}\to f_{D,t}}(a_{t})=\mu_{f_{B,t}\to a_{t}}(a_{t}) and

μfD,tbt(bt)=\displaystyle\mu_{f_{D,t}\to b_{t}}(b_{t})= (61)
fD,t(at,bt,ct)μatfD,t(at)μctfD,t(ct)dbtdct\displaystyle\int\!\!\int{{f_{D,t}}}\left({a_{t},b_{t},c_{t}}\right)\cdot{\mu_{{a_{t}}\to{f_{D,t}}}}(a_{t})\cdot{\mu_{{c_{t}}\to{f_{D,t}}}}(c_{t}){\rm{d}}{b_{t}}{\rm{d}}c_{t}
=𝒩(bt;mfD,tbt,VfD,tbt)\displaystyle=\mathcal{N}\left({{b_{t}};{m_{f_{D,t}\to b_{t}}},V_{f_{D,t}\to b_{t}}}\right)

and

μbtfK,tk(bt)=μfD,tbt(bt).\mu_{b_{t}\to f_{K,t-k}}(b_{t})=\mu_{f_{D,t}\to b_{t}}(b_{t}). (62)

It can be verified that

μfK,tkxtk(xtk)=fK,tk(htk,xtk,bt)\displaystyle\mu_{f_{K,t-k}\to x_{t-k}}(x_{t-k})=\int\!\!\!\int{{f_{K,{t-k}}}}\left({{h_{t-k}},{x_{t-k}},b_{t}}\right) (63)
μhtkfK,tk(htk)μbtfK,tk(bt)dhtkdxtk.\displaystyle\cdot{\mu_{{h_{t-k}}\to{f_{K,{t-k}}}}}(h_{t-k})\cdot\mu_{b_{t}\to f_{K,t-k}}(b_{t}){\rm{d}}{h_{t-k}}{\rm{d}}x_{t-k}.
Algorithm 2 Message Passing Based Joint Communication and Imaging Algorithm.
1:  Initialize the maximum number of iterations KmaxK_{max} , and threshold ϵ0\epsilon_{0}
2:  while Δϵ0\Delta\leq\epsilon_{0} and p<Kmaxp<K_{max}  do
3:     For tkt\leq k, compute μfB,txt(xt){\mu_{{f_{B,t}}\to{x_{t}}}}(x_{t}) with (28).
4:     For k<tLk<t\leq L, compute μfB,txt(xt)\mu_{f_{B,t}\to x_{t}}(x_{t}) with (47).
5:     For Lk<tLL-k<t\leq L, compute μfK,txt(xt)\mu_{f_{K,t}\to x_{t}}(x_{t}) with (52).
6:     For k<tLkk<t\leq L-k, compute μfK,txt(xt)\mu_{f_{K,t}\to x_{t}}(x_{t}) with (64).
7:     Compute B(xt)B(x_{t}) with (67).
8:     The estimate of xtx_{t} in the ppth iteration is given by xtp=x^tx_{t}^{p}=\hat{x}_{t}, where x^tp\hat{x}_{t}^{p} maximizes B(xt)B(x_{t}).
9:     Compute 𝝈p\bm{\sigma}^{p} by iterating (69)- (73) until convergence.
10:     Compute Δ=𝐱p𝐱p1+𝝈p𝝈p1\Delta=\left\|\mathbf{x}^{p}-\mathbf{x}^{p-1}\right\|+\left\|\bm{\sigma}^{p}-\bm{\sigma}^{p-1}\right\|.
11:  end while

Substituting (35) into (63), we have

μfK,tkxtk(xtk)fK,tk(htk,xtk,bt)\displaystyle\mu_{f_{K,t-k}\to x_{t-k}}(x_{t-k})\approx\int\!\!\!\int{{f_{K,{t-k}}}}\left({{h_{t-k}},{x_{t-k}},b_{t}}\right) (64)
μhtkfK,tk(htk)δ(htkmhtk)dhtkdxtk\displaystyle\cdot{\mu_{{h_{t-k}}\to{f_{K,{t-k}}}}}(h_{t-k})\cdot\delta(h_{t-k}-m_{h_{t-k}}){\rm{d}}{h_{t-k}}{\rm{d}}x_{t-k}
=fK,tk(mhtk,xtk,bt)μhtkfK,tk(htk)dxtk\displaystyle=\int{{f_{K,{t-k}}}}\left({m_{h_{t-k}},{x_{t-k}},b_{t}}\right)\cdot{\mu_{{h_{t-k}}\to{f_{K,{t-k}}}}}(h_{t-k})\cdot{\rm{d}}x_{t-k}
=𝒩(xtk;mfK,tkxtk,WfK,tkxtk1)\displaystyle=\mathcal{N}\left({x_{t-k}};m_{{f_{K,{t-k}}}\to{x_{t-k}}},W_{{f_{K,{t-k}}}\to x_{t-k}}^{-1}\right)

where

mfK,tkxtk=mfD,tbt/mhtkm_{{f_{K,{t-k}}}\to{x_{t-k}}}={m_{f_{D,t}\to b_{t}}}/m_{h_{t-k}} (65)
WfK,tkxtk=mfK,tkxtkHWfD,tbtmfK,tkxtkW_{{f_{K,{t-k}}}\to x_{t-k}}=m^{\rm{H}}_{{f_{K,{t-k}}}\to{x_{t-k}}}{W_{f_{D,t}\to b_{t}}}m_{{f_{K,{t-k}}}\to{x_{t-k}}} (66)

Based on (28), (47), (52), and (64), B(xt)B(x_{t}) in (23) can be calculated by

B(xt)\displaystyle B(x_{t}) =14i=14δ[xtej(π2iπ4)]\displaystyle=\frac{1}{4}\sum_{i=1}^{4}\delta[x_{t}-e^{j(\frac{\pi}{2}i-\frac{\pi}{4})}] (67)
𝒩(xt;mfB,txt,WfB,txt1)\displaystyle\cdot\mathcal{N}\left({{x_{t}};{m_{{f_{B,t}}\to{x_{t}}}},W_{{f_{B,t}}\to{x_{t}}}^{-1}}\right)
𝒩(xt;mfK,txt,WfK,txt1)\displaystyle\cdot\mathcal{N}\left({{x_{t}};{m_{{f_{K,t}}\to{x_{t}}}},W_{{f_{K,t}}\to{x_{t}}}^{-1}}\right)

The estimate of xtx_{t} is given by x^t\hat{x}_{t}, which maximize B(xt)B(x_{t}). With the estimate of 𝐱\mathbf{x}, the adaptive SBL algorithm [61] can be adopted to update 𝝈{\bm{\sigma}}. We define

𝐡=[yk+1a^k+1,,yLa^L,yL+1,,yL+k]T.\small\mathbf{h}=\left[y_{k+1}\!-\!\hat{a}_{{k}\!+\!1},\cdots,y_{L}\!-\!\hat{a}_{{L}},\ y_{L+1},\cdots,y_{L+k}\right]^{\rm{T}}. (68)

where at^=αc𝐠T𝚯(t)𝐡𝐜x^t\hat{a_{t}}={\alpha_{c}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}\hat{x}_{t}. Then 𝝈{\bm{\sigma}} can be updated by iterating the following steps until 𝝈{{\bm{\sigma}}} converges.

𝐖𝝈=𝐆H𝐆ξ^+diag(𝜸^),{{\bf{W}}_{\bm{\sigma}}}={\frac{{\bf{G}}^{\text{H}}{\bf{G}}}{\hat{\xi}}+{\mathop{\rm diag}\nolimits}(\hat{\bm{\gamma}})}, (69)
𝝈=𝐖𝝈𝐆H𝐡ξ^,{{\bm{\sigma}}}=\frac{{{\bf{W}}_{\bm{\sigma}}}{{\bf{G}}^{\text{H}}}{\bf{h}}}{\hat{\xi}}, (70)
γ^m=(2ϵ+1)/(2ζ+|𝝈(m)|2+𝐖𝝈1(m,m)),m=1,,M{\hat{\gamma}_{m}}=(2\epsilon+1)/\left({2\zeta+{{\left|{{{\bm{\sigma}}}(m)}\right|}^{2}}+{{\bf{W}}^{-1}_{\bm{\sigma}}}(m,m)}\right),m=1,\ldots,M (71)
ϵ=12log(1Mγ^m)1Mlog(γ^m)\epsilon=\frac{1}{2}\sqrt{\text{log}(\frac{1}{M}\sum{\hat{\gamma}_{m}})-\frac{1}{M}\sum{\text{log}(\hat{\gamma}_{m})}} (72)
ξ^=𝐡𝐆𝝈2Lm=1Mγ^m\hat{\xi}=\frac{{{{\left\|{{\bf{h}}-{\bf{G\bm{\sigma}}}}\right\|}^{2}}}}{{L-\sum\limits_{m=1}^{M}{{{\hat{\gamma}}_{m}}}}} (73)
Refer to caption
Figure 5: The receive beam pattern of the BS, i.e., 𝐠T𝚯(t)a(θp){{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t)a(\theta_{p}), in the case of different time index tt.
Refer to caption
Figure 6: The receive beam pattern of the BS, i.e., 𝐠T𝚯(t)a(θp){{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t)a(\theta_{p}), in the case of different weight factor ρ\rho.
Refer to caption
Figure 7: The receive beam pattern of the BS, i.e., 𝐠T𝚯(t)a(θp){{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t)a(\theta_{p}), in the case of different nbitn_{bit}.
Refer to caption
(a) The ground truth of 𝝈\bm{\sigma} in scenario 1.
Refer to caption
(b) NMSE = -15.16 dB.
Refer to caption
(c) NMSE = -20.49 dB
Refer to caption
(d) NMSE = -25.34 dB.
Refer to caption
(e) The ground truth of 𝝈\bm{\sigma} in scenario 2.
Refer to caption
(f) NMSE = -15.14 dB.
Refer to caption
(g) NMSE = -20.42 dB.
Refer to caption
(h) NMSE = -25.03 dB.
Refer to caption
(i) The ground truth of 𝝈\bm{\sigma} in scenario 3.
Refer to caption
(j) NMSE = -15.24 dB.
Refer to caption
(k) NMSE = -20.11 dB.
Refer to caption
(l) NMSE = -25.10 dB.
Figure 8: The ground truth of rearranged 𝝈\bm{\sigma} and corresponding imaging results under different NMSE levels.

The proposed message passing based joint communication and imaging algorithm is summarized in Algorithm 2.

V Numerical Results

In this section, numerical experiments are conducted to evaluate the uplink ISAC system and the performance of proposed method. We consider a RIS with 150 elements and a BS located at θB=45{\theta_{B}}={-45^{\circ}} from the RIS. A single UE is located at θU=70{\theta_{U}}={70^{\circ}}, and the transmitted signal is QPSK modulated with the length L=1024L=1024. The RoI is equidistantly divided into M=64M=64 grids, spacing from 1515^{\circ} to 50{50^{\circ}} from the RIS. In Algorithm 1, we set η=0.001\eta=0.001 and r=0.005r=0.005.

We define the communication-to-noise ratio (CNR) and the imaging-to-noise ratio (INR) as

CNR=10log10𝔼t[|αc𝐠T𝚯(t)𝐡𝐜x(t)|2]ξ2CNR=10\text{log}_{10}\frac{\mathbb{E}_{t}\left[\left|{\alpha_{c}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{h}}_{\bf{c}}}x(t)\right|^{2}\right]}{\xi^{2}} (74)

and

INR=10log10𝔼t[|αI𝐠T𝚯(t)𝐇𝐈𝝈x(tk)|2]ξ2.INR=10\text{log}_{10}\frac{\mathbb{E}_{t}\left[\left|{{\alpha_{I}}{{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t){{\bf{H}}_{\bf{I}}}{\bm{\sigma}}x(t-k)}\right|^{2}\right]}{\xi^{2}}. (75)

The performance of communication and imaging is assessed by the symbol-error-rate (SER) and normalized mean squared error (NMSE), which are respectively defined as

SER=q=1QEqQLSER=\sum_{q=1}^{Q}\frac{E_{q}}{QL} (76)

and

NMSE=10log10𝝈𝒒𝝈2𝝈2,NMSE=10\text{log}_{10}\frac{\|\bm{\sigma_{q}}-\bm{\sigma}\|^{2}}{\|\bm{\sigma}\|^{2}}, (77)

where QQ is the number of trials, EqE_{q} and 𝝈q\bm{\sigma}_{q} denote the number of error symbol and the estimate of 𝝈\bm{\sigma} for the qqth trial, respectively.

V-A Phase Optimization Result

We first illustrate the performance of the proposed phase design scheme. Firstly, we set ρ=0.5\rho=0.5 and investigate the receive beam pattern of the BS at different time index tt, i.e., 𝐠T𝚯(t)a(θp),θp(90,90){{\bf{g}}^{\rm{T}}}{\bf{\Theta}}(t)a(\theta_{p}),\theta_{p}\in(-90^{\circ},90^{\circ}). In Fig.5, we show the receive beam pattern at time indexes t=1,t=200,t=500,t=800t=1,t=200,t=500,t=800 and t=1000t=1000. We observe that the pattern is relatively flat in the directions of the RoI, and the maximum gain is obtained in the direction of the UE, i.e., θU=70\theta_{U}=70^{\circ}. Besides, it is also noted that to meet the orthogonality requirement of the sensing matrix, the beam pattern is disorderly in the directions of RoI at different time indexes.

Then we investigate the receive beam pattern under different weight factor ρ\rho. As shown in Fig. 6, the gain in the direction of the UE increases with ρ\rho, while the gain in the direction of RoI is decreases accordingly. The results illustrate that the choice of ρ\rho is a trade-off between communication and imaging performance.

In Fig. 7, we fix ρ=0.5\rho=0.5 and investigate the effect of nbitn_{bit} on the receive pattern. Compared with discrete phase model, a higher gain can be obtained in the direction of the UE. When nbit=1n_{bit}=1, it is shown that there is a gain loss about 0.5 dB compared with the continuous phase. We also observe that the gain loss at the main lobe is only about 0.05dB when the phase shift is 2-bit quantized.

V-B Performance of Uplink Joint Communication and Imaging

In this subsection, we investigate the performance of the proposed message passing based joint communication and imaging method. We rearrange the scattering coefficients 𝝈\bm{\sigma} into a 8×88\times 8 matrix to make it intuitive. As shown in Fig. 8 (a), (e) and (i), the letters X, D and U are respectively set to be the grand truth of 𝝈\bm{\sigma} in different scenarios of our simulations.

Refer to caption
Figure 9: The SER performance of the proposed method versus CNR and INR when the CNR is 5 dB higher than the INR.
Refer to caption
Figure 10: The NMSE performance of the proposed method versus CNR and INR when the CNR is 5 dB higher than the INR.
Refer to caption
Figure 11: The SER performance of the proposed method versus INR when CNR = 12 dB.
Refer to caption
Figure 12: The NMSE performance of the proposed method versus INR when CNR = 12 dB.

Scenario 1: In the first scenario, we investigate communication and imaging performance versus the noise level. Note that the CNR and INR will change synchronous with the noise level, we fix CNR-INR = 5dB. The ground truth of 𝝈\bm{\sigma} is shown in Fig. 8 (a). In Fig. 9, the SER performance versus CNR is plotted, along with three benchmarks to better illustrate the performance of the proposed method. It is not surprising that we observe the worst SER performance when we ignore the imaging echo directly, i.e., the echoes are not decouple. Besides, we pretend that there is no imaging echo, i.e., the received signals are simple uncoupled QPSK modulated signals, and we take the SER performance in this case as another benchmark. It is shown that under low CNR cases, the SER of the proposed method is higher than the QPSK signal. However, with the increase of CNR, the proposed method achieves a better SER performance than the QPSK signals. This result shows that under high SNR levels, the decoupled echoes are conducive to improving communication performance. Moreover, we also plot the SER performance under given 𝝈\bm{\sigma} as a low bound. It is shown that, with the increase of CNR, the SER of the proposed method can approaches this low bound gradually.

To illustrate the imaging performance of the proposed method, the NMSE versus CNR is investigated in Fig.10, where the NMSE performance under given 𝐱\mathbf{x} is also included as a low bound. It can be seen that the proposed method can attain the low bound with the increase of CNR. To visually view the imaging results, we also show the 𝝈q\bm{\sigma}_{q} in different NMSE levels in Fig. 8 (b)-(d).

Scenario 2: In this scenario, we fix CNR = 12 dB and vary INR from 4 dB to 8 dB to evaluate the performance under different CNR levels. The ground truth of rearranged 𝝈\bm{\sigma} is given in Fig. 8 (e). As shown in Fig. 11, we observe that the communication performance improves with the increase of INR and gradually approaches the low bound in the case of given 𝝈\bm{\sigma}. It also can be seen that, with the increase of INR, the SER performance transcends the performance boundary of uncoupled QPSK signal, which also indicates that the imaging echo is utilized to enhance the communication performance. In addition, the NMSE performance is shown in Fig. 12, along with the NMSE in the case of given 𝐱\mathbf{x} as a benchmark. In Fig. 8 (f)-(h), we show the 𝝈q\bm{\sigma}_{q} in different NMSE levels to view the imaging results visually. The ground truth of rearranged 𝝈\bm{\sigma} is given in Fig. 8 (e). We can draw the conclusion that the imaging performance of the proposed method improves with the increase of CNR and it is able to reach the low bound asymptotically.

Scenario 3: Lastly, we fix INR = 6 dB and investigate the performance of the proposed method under different CNRs. The ground truth of the rearranged 𝝈\bm{\sigma} is given in Fig. 8 (i), and the SER performance of the proposed method is shown in Fig. 13. Besides, we consider the NMSE performance in the case of given 𝝈\bm{\sigma} as a benchmark. As depicted in Fig. 13, the performance of the proposed method approaches the benchmark in a wide range of CNR. In Fig. 14, we evaluate the imaging performance under different CNR levels. It can be noted that with the increase of CNR, the NMSE of the proposed method decreases and gradually approaches the NMSE performance under given 𝐱\mathbf{x}.

Refer to caption
Figure 13: The SER performance of the proposed method versus CNR when INR = 6 dB.
Refer to caption
Figure 14: The NMSE performance of the proposed method versus CNR when INR = 6 dB.

VI Conclusion

We have presented a novel RIS-assisted uplink joint communication and imaging system in this article. A phase optimization problem is established with the consideration of both communication and imaging performance. Since the problem is non-convex and strongly coupled, we propose a back propagation based method, where temperature parameters and discrete prior are jointly utilized to solve it. After that, an efficient iterative message passing based method is proposed, cooperating the adaptive SBL to decouple echoes and achieve joint communication and imaging. Simulation results have demonstrated that the proposed phase design scheme can make a trade-off between communication and imaging functionality. It is also shown that the proposed method approaches the low bound asymptotically, and the communication performance can be enhanced by making full use of imaging echoes.

References

  • [1] Y. Cui, F. Liu, X. Jing, and J. Mu, “Integrating sensing and communications for ubiquitous iot: Applications, trends, and challenges,” IEEE Network, vol. 35, no. 5, pp. 158–167, 2021.
  • [2] F. Peng, Z. Jiang, S. Zhou, Z. Niu, and S. Zhang, “Sensing and Communication Co-Design for Status Update in Multiaccess Wireless Networks,” IEEE Transactions on Mobile Computing, pp. 1–1, 2021.
  • [3] S. Biswas, K. Singh, O. Taghizadeh, and T. Ratnarajah, “Coexistence of MIMO Radar and FD MIMO Cellular Systems With QoS Considerations,” IEEE Transactions on Wireless Communications, vol. 17, no. 11, pp. 7281–7294, 2018.
  • [4] P. Kumari, J. Choi, N. González-Prelcic, and R. W. Heath, “IEEE 802.11ad-Based Radar: An Approach to Joint Vehicular Communication-Radar System,” IEEE Transactions on Vehicular Technology, vol. 67, no. 4, pp. 3012–3027, Apr. 2018.
  • [5] J. Qian, M. Lops, L. Zheng, X. Wang, and Z. He, “Joint System Design for Coexistence of MIMO Radar and MIMO Communication,” IEEE Transactions on Signal Processing, vol. 66, no. 13, pp. 3504–3519, Jul. 2018.
  • [6] J. Qian, Z. He, N. Huang, and B. Li, “Transmit Designs for Spectral Coexistence of MIMO Radar and MIMO Communication Systems,” IEEE Transactions on Circuits and Systems II: Express Briefs, vol. 65, no. 12, pp. 2072–2076, 2018.
  • [7] N. Cao, Y. Chen, X. Gu, and W. Feng, “Joint Radar-Communication Waveform Designs Using Signals From Multiplexed Users,” IEEE Transactions on Communications, vol. 68, no. 8, pp. 5216–5227, 2020.
  • [8] F. Dong, W. Wang, Z. Hu, and T. Hui, “Low-Complexity Beamformer Design for Joint Radar and Communications Systems,” IEEE Communications Letters, vol. 25, no. 1, pp. 259–263, 2021.
  • [9] F. Liu, L. Zhou, C. Masouros, A. Li, W. Luo, and A. Petropulu, “Toward Dual-functional Radar-Communication Systems: Optimal Waveform Design,” IEEE Transactions on Signal Processing, vol. 66, no. 16, pp. 4264–4279, 2018.
  • [10] Y. Luo, J. A. Zhang, X. Huang, W. Ni, and J. Pan, “Optimization and Quantization of Multibeam Beamforming Vector for Joint Communication and Radio Sensing,” IEEE Transactions on Communications, vol. 67, no. 9, pp. 6468–6482, Sep. 2019.
  • [11] S. Shi, Z. Wang, Z. He, and Z. Cheng, “Constrained waveform design for dual-functional MIMO radar-Communication system,” Signal Processing, vol. 171, pp. 107 530–107 541, Jun. 2020.
  • [12] A. Hassanien, M. G. Amin, Y. D. Zhang, and F. Ahmad, “Signaling strategies for dual-function radar communications: An overview,” IEEE Aerospace and Electronic Systems Magazine, vol. 31, no. 10, pp. 36–45, Oct. 2016.
  • [13] J. A. Zhang, F. Liu, C. Masouros, R. W. Heath, Z. Feng, L. Zheng, and A. Petropulu, “An Overview of Signal Processing Techniques for Joint Communication and Radar Sensing,” IEEE Journal of Selected Topics in Signal Processing, vol. 15, no. 6, pp. 1295–1315, 2021.
  • [14] L. Zheng, M. Lops, Y. C. Eldar, and X. Wang, “Radar and Communication Coexistence: An Overview: A Review of Recent Methods,” IEEE Signal Processing Magazine, vol. 36, no. 5, pp. 85–99, Sep. 2019.
  • [15] B. K. Chalise, M. G. Amin, and B. Himed, “Performance Tradeoff in a Unified Passive Radar and Communications System,” IEEE Signal Processing Letters, vol. 24, no. 9, pp. 1275–1279, Sep. 2017.
  • [16] A. Liu, Z. Huang, M. Li, Y. Wan, W. Li, T. X. Han, C. Liu, R. Du, D. T. K. Pin, J. Lu, Y. Shen, F. Colone, and K. Chetty, “A Survey on Fundamental Limits of Integrated Sensing and Communication,” IEEE Communications Surveys Tutorials, vol. 24, no. 2, pp. 994–1034, 2022.
  • [17] P. Ren, A. Munari, and M. Petrova, “Performance Tradeoffs of Joint Radar-Communication Networks,” IEEE Wireless Communications Letters, vol. 8, no. 1, pp. 165–168, 2019.
  • [18] T. J. Cui, M. Q. Qi, X. Wan, J. Zhao, and Q. Cheng, “Coding metamaterials, digital metamaterials and programmable metamaterials,” Light: science & applications, vol. 3, no. 10, pp. e218–e218, 2014.
  • [19] S. Liu and T. J. Cui, “Concepts, working principles, and applications of coding and programmable metamaterials,” Advanced optical materials, vol. 5, no. 22, p. 1700624, 2017.
  • [20] T. J. Cui, L. Li, S. Liu, Q. Ma, L. Zhang, X. Wan, W. X. Jiang, and Q. Cheng, “Information metamaterial systems,” iScience, vol. 23, no. 8, p. 101403, 2020. [Online]. Available: https://www.sciencedirect.com/science/article/pii/S2589004220305939
  • [21] X. Q. Chen, L. Zhang, S. Liu, and T. J. Cui, “Artificial Neural Network for Direction-of-Arrival Estimation and Secure Wireless Communications Via Space-Time-Coding Digital Metasurfaces,” Advanced Optical Materials, vol. 10, p. 2201900, 2022.
  • [22] J. Y. Dai, W. Tang, M. Wang, M. Z. Chen, Q. Cheng, S. Jin, T. J. Cui, and C. H. Chan, “Simultaneous in situ Direction Finding and Field Manipulation Based on Space-Time-Coding Digital Metasurface,” IEEE Transactions on Antennas and Propagation, vol. 70, no. 6, pp. 4774–4783, Jun. 2022.
  • [23] Z. A. Huang, J. W. Wang, Q. Xiao, B. Y. Li, W. H. Wang, X. Wan, Q. Cheng, and T. J. Cui, “A Single Metasurface Can Perform Range-Velocity Detection and Target Imaging Simultaneously at Single Frequency,” Advanced Optical Materials, vol. 10, no. 21, p. 2201382, 2022.
  • [24] H. Li, Y. B. Li, G. Chen, S. Y. Dong, J. L. Shen, C. Y. Gong, S. Y. Wang, H. P. Wang, and T. J. Cui, “High-Resolution Near-Field Imaging and Far-Field Sensing Using a Transmissive Programmable Metasurface,” Advanced Materials Technologies, vol. 7, no. 5, p. 2101067, 2022.
  • [25] Z. Wang, H. Zhang, H. Zhao, T. J. Cui, and L. Li, “Intelligent electromagnetic metasurface camera: System design and experimental results,” Nanophotonics, vol. 11, no. 9, pp. 2011–2024, Apr. 2022.
  • [26] X. Wan, Z. A. Huang, J. W. Wang, B. Y. Li, W. H. Wang, Q. Cheng, and T. J. Cui, “Joint Radar and Communication Empowered by Digital Programmable Metasurface,” Advanced Intelligent Systems, vol. 4, no. 9, pp. 2 200 083–2 200 092, 2022.
  • [27] X. Wan, C. K. Xiao, H. Huang, Q. Xiao, W. Xu, J. W. Wang, Z. A. Huang, Q. Cheng, S. Jin, and T. J. Cui, “User Tracking and Wireless Digital Transmission through a Programmable Metasurface,” Advanced Materials Technologies, vol. 6, no. 7, p. 2001254, 2021.
  • [28] J. Wei Wang, Z. Ai Huang, Q. Xiao, W. Han Li, B. Yang Li, X. Wan, and T. Jun Cui, “High-Precision Direction-of-Arrival Estimations Using Digital Programmable Metasurface,” Advanced Intelligent Systems, vol. 4, no. 4, pp. 2 100 164–2 100 171, 2022.
  • [29] Q. Y. Zhou, J. W. Wu, S. R. Wang, Z. Q. Fang, L. J. Wu, J. C. Ke, J. Y. Dai, T. J. Cui, and Q. Cheng, “Two-dimensional direction-of-arrival estimation based on time-domain-coding digital metasurface,” Applied Physics Letters, vol. 121, no. 18, pp. 181 702–181 708, Oct. 2022.
  • [30] E. Basar, M. Di Renzo, J. De Rosny, M. Debbah, M.-S. Alouini, and R. Zhang, “Wireless Communications Through Reconfigurable Intelligent Surfaces,” IEEE Access, vol. 7, pp. 116 753–116 773, 2019.
  • [31] H. Lu, Y. Zeng, S. Jin, and R. Zhang, “Enabling panoramic full-angle reflection via aerial intelligent reflecting surface,” in 2020 IEEE International Conference on Communications Workshops (ICC Workshops), 2020, pp. 1–6.
  • [32] W. Tang, M. Z. Chen, X. Chen, J. Y. Dai, Y. Han, M. Di Renzo, Y. Zeng, S. Jin, Q. Cheng, and T. J. Cui, “Wireless Communications With Reconfigurable Intelligent Surface: Path Loss Modeling and Experimental Measurement,” IEEE Transactions on Wireless Communications, vol. 20, no. 1, pp. 421–439, Jan. 2021.
  • [33] J. Zhao, X. Yang, J. Y. Dai, Q. Cheng, X. Li, N. H. Qi, J. C. Ke, G. D. Bai, S. Liu, S. Jin, A. Alù, and T. J. Cui, “Programmable time-domain digital-coding metasurface for non-linear harmonic manipulation and new wireless communication systems,” National Science Review, vol. 6, no. 2, pp. 231–238, Mar. 2019.
  • [34] J. Y. Dai, W. Tang, M. Z. Chen, C. H. Chan, Q. Cheng, S. Jin, and T. J. Cui, “Wireless Communication Based on Information Metasurfaces,” IEEE Transactions on Microwave Theory and Techniques, vol. 69, no. 3, pp. 1493–1510, Mar. 2021.
  • [35] J. Y. Dai, W. K. Tang, J. Zhao, X. Li, Q. Cheng, J. C. Ke, M. Z. Chen, S. Jin, and T. J. Cui, “Wireless Communications through a Simplified Architecture Based on Time-Domain Digital Coding Metasurface,” Advanced Materials Technologies, vol. 4, no. 7, pp. 1 900 044–1 900 055, Jul. 2019.
  • [36] J. Y. Dai, W. Tang, L. X. Yang, X. Li, M. Z. Chen, J. C. Ke, Q. Cheng, S. Jin, and T. J. Cui, “Realization of Multi-Modulation Schemes for Wireless Communication by Time-Domain Digital Coding Metasurface,” IEEE Transactions on Antennas and Propagation, vol. 68, no. 3, pp. 1618–1627, Mar. 2020.
  • [37] Y. Li, X. Wan, J. Eisenbeis, X. Long, M. Jozwicka, T. Cui, and T. Zwick, “A novel wireless communication system using programmable metasurface operating at 28 ghz,” in 2020 International Conference on Microwave and Millimeter Wave Technology (ICMMT), 2020, pp. 1–3.
  • [38] W. Tang, J. Y. Dai, M. Z. Chen, K.-K. Wong, X. Li, X. Zhao, S. Jin, Q. Cheng, and T. J. Cui, “MIMO Transmission Through Reconfigurable Intelligent Surface: System Design, Analysis, and Implementation,” IEEE Journal on Selected Areas in Communications, vol. 38, no. 11, pp. 2683–2699, Nov. 2020.
  • [39] X. Wan, Q. Zhang, T. Yi Chen, L. Zhang, W. Xu, H. Huang, C. Kun Xiao, Q. Xiao, and T. Jun Cui, “Multichannel direct transmissions of near-field information,” Light: Science & Applications, vol. 8, no. 1, pp. 60–68, Dec. 2019.
  • [40] L. Zhang, M. Z. Chen, W. Tang, J. Y. Dai, L. Miao, X. Y. Zhou, S. Jin, Q. Cheng, and T. J. Cui, “A wireless communication scheme based on space- and frequency-division multiplexing using digital surfaces,” Nature Electronics, vol. 4, no. 3, pp. 218–227, Mar. 2021.
  • [41] X. Shao, C. You, W. Ma, X. Chen, and R. Zhang, “Target sensing with intelligent reflecting surface: Architecture and performance,” IEEE Journal on Selected Areas in Communications, vol. 40, no. 7, pp. 2070–2084, 2022.
  • [42] S. Buzzi, E. Grossi, M. Lops, and L. Venturino, “Radar target detection aided by reconfigurable intelligent surfaces,” IEEE Signal Processing Letters, vol. 28, pp. 1315–1319, 2021.
  • [43] S. Buzzi, E. Grossi, M. Lops, and L. Venturino, “Foundations of MIMO Radar Detection Aided by Reconfigurable Intelligent Surfaces,” IEEE Transactions on Signal Processing, vol. 70, pp. 1749–1763, 2022.
  • [44] L. Li, Y. Shuang, Q. Ma, H. Li, H. Zhao, M. Wei, C. Liu, C. Hao, C.-W. Qiu, and T. J. Cui, “Intelligent metasurface imager and recognizer,” Light: Science & Applications, vol. 8, no. 1, p. 97, Dec. 2019.
  • [45] L. Dai, B. Wang, M. Wang, X. Yang, J. Tan, S. Bi, S. Xu, F. Yang, Z. Chen, M. Di Renzo et al., “Reconfigurable intelligent surface-based wireless communications: Antenna design, prototyping, and experimental results,” IEEE Access, vol. 8, pp. 45 913–45 923, 2020.
  • [46] X. Wang, Z. Fei, J. Huang, and H. Yu, “Joint Waveform and Discrete Phase Shift Design for RIS-Assisted Integrated Sensing and Communication System Under Cramer-Rao Bound Constraint,” IEEE Transactions on Vehicular Technology, vol. 71, no. 1, pp. 1004–1009, Jan. 2022.
  • [47] X. Wang, Z. Fei, Z. Zheng, and J. Guo, “Joint Waveform Design and Passive Beamforming for RIS-Assisted Dual-Functional Radar-Communication System,” IEEE Transactions on Vehicular Technology, vol. 70, no. 5, pp. 5131–5136, 2021.
  • [48] Z.-M. Jiang, M. Rihan, P. Zhang, L. Huang, Q. Deng, J. Zhang, and E. M. Mohamed, “Intelligent reflecting surface aided dual-function radar and communication system,” IEEE Systems Journal, vol. 16, no. 1, pp. 475–486, 2022.
  • [49] X. Song, D. Zhao, H. Hua, T. X. Han, X. Yang, and J. Xu, “Joint transmit and reflective beamforming for irs-assisted integrated sensing and communication,” in 2022 IEEE Wireless Communications and Networking Conference (WCNC), 2022, pp. 189–194.
  • [50] Y. He, Y. Cai, H. Mao, and G. Yu, “Ris-assisted communication radar coexistence: Joint beamforming design and analysis,” IEEE Journal on Selected Areas in Communications, vol. 40, no. 7, pp. 2131–2145, 2022.
  • [51] Z. Wang, Z. Liu, Y. Shen, A. Conti, and M. Z. Win, “Location awareness in beyond 5g networks via reconfigurable intelligent surfaces,” IEEE Journal on Selected Areas in Communications, vol. 40, no. 7, pp. 2011–2025, 2022.
  • [52] H. Loeliger, “An introduction to factor graphs,” IEEE Signal Processing Magazine, vol. 21, no. 1, pp. 28–41, 2004.
  • [53] H. Loeliger, J. Dauwels, J. Hu, S. Korl, L. Ping, and F. R. Kschischang, “The factor graph approach to model-based signal processing,” Proceedings of the IEEE, vol. 95, no. 6, pp. 1295–1322, 2007.
  • [54] M. Di Renzo, K. Ntontin, J. Song, F. H. Danufane, X. Qian, F. Lazarakis, J. De Rosny, D.-T. Phan-Huy, O. Simeone, R. Zhang et al., “Reconfigurable intelligent surfaces vs. relaying: Differences, similarities, and performance comparison,” IEEE Open Journal of the Communications Society, vol. 1, pp. 798–807, 2020.
  • [55] F. Liu, A. Garcia-Rodriguez, C. Masouros, and G. Geraci, “Interfering channel estimation in radar-cellular coexistence: How much information do we need?” IEEE Transactions on Wireless Communications, vol. 18, no. 9, pp. 4238–4253, 2019.
  • [56] R. G. Baraniuk, “Compressive sensing [lecture notes],” IEEE signal processing magazine, vol. 24, no. 4, pp. 118–121, 2007.
  • [57] S. Ji, Y. Xue, and L. Carin, “Bayesian compressive sensing,” IEEE Transactions on signal processing, vol. 56, no. 6, pp. 2346–2356, 2008.
  • [58] E. J. Candes and T. Tao, “Near-optimal signal recovery from random projections: Universal encoding strategies?” IEEE Transactions on Information Theory, vol. 52, no. 12, pp. 5406–5425, 2006.
  • [59] A. Chakrabarti, “Learning sensor multiplexing design through back-propagation,” CoRR, vol. abs/1605.07078, 2016. [Online]. Available: http://arxiv.org/abs/1605.07078
  • [60] L. Bottou, F. E. Curtis, and J. Nocedal, “Optimization methods for large-scale machine learning,” Siam Review, vol. 60, no. 2, pp. 223–311, 2018.
  • [61] M. Luo, Q. Guo, M. Jin, Y. C. Eldar, D. Huang, and X. Meng, “Unitary approximate message passing for sparse bayesian learning,” IEEE transactions on signal processing, vol. 69, pp. 6023–6039, 2021.