Multi-frequency calibration for DOA estimation with distributed sensors

2020·Arxiv

Abstract

Abstract

In this work, we investigate direction finding in the presence of sensor gain uncertainties and directional perturbations for sensor array processing in a multi-frequency scenario. Specifically, we adopt a distributed optimization scheme in which coherence models are incorporated and local agents exchange information only between connected nodes in the network, i.e., without a fusion center. Numerical simulations highlight the advantages of the proposed parallel iterative technique in terms of statistical and computational efficiency.

Keywords: Calibration, source localization, multi-frequency, sensor array processing, distributed optimization.

1 Introduction

Calibration and Direction-of-Arrival (DoA) estimation is a major issue in array processing [2, 3]. The latter has been studied in several applications, e.g., radar, sonar, satellite, wireless communication and radio interferometric systems [4, 5], where we commonly use largely distributed sensors elements aiming to achieve high resolution. In all these sensor network applications, calibration is required as some parameters are not exactly known due to imperfect instrumentation or propagation conditions [6]. Let us note that calibration algorithms are distinguished by the presence [7] or absence [8] of one or more cooperative sources, named calibrator sources. Indeed, prior source information can be available [6] and consists mainly in the true/nominal directions and powers of calibrator sources (i.e., without any perturbation effects or antenna imperfections). Furthermore, most calibration algorithms are based on the least squares approach, with a sequential procedure updating each parameter alternatively [4]. The least squares estimator is indeed equivalent to the Maximum Likelihood (ML) method under a (unrealistic) Gaussian noise model.

The aim of the proposed methodology here is to estimate successively the unknown sensor gains and phase errors, along with the calibrator and noise parameters, through minimization of a proper weighting cost function. In this work, uncertainties are estimated from the array covariance matrix, since dealing directly with time series data and operating on the signal domain quickly becomes computationally unfeasible for a large number of samples [9]. The scenario under study is general but could be adapted to any practical application as in the radio astronomy context, where the number of parameters to estimate is tremendous and frequency bands are wide.

In the multi-frequency scenario, a suboptimal way to perform calibration is to consider one wavalength bin at a time, with only one centralized processor, which has access to data in the whole available range of wavelengths. In this work, we study an accelerated version based on the scalable form of the Alternating Direction Method of Multipliers (ADMM) [10, 11] with a specific network topology: there is no fusion center and agents exchange information only among themselves. The goal being to reduce the complexity in operation flow and signaling exchanging [12, 13, 14, 15, 16, 17]. For estimation of the directional gains, the compressive sensing framework, especially the sparse representation method, is well-adapted and has already been applied for source localization in fully and partially calibrated arrays [18, 19, 20, 21].

The notation used through this paper is the following: denote, respectively, the complex conjugate, transpose, Hermitian operator, element-wise raising to , real part and the n-th element of a vector. The expectation operator is respectively, the Kronecker, the Khatri-Rao and the Hadamard product. The operator diag(.) converts a vector into a diagonal matrix, blkdiag(.) is the block-diagonal operator, whereas vecdiag(.) produces a vector from the main diagonal of a matrix and vec(.) stacks the columns of a matrix on top of one another. The operators refer to the and Frobenius norms, respectively. Finally, identity matrix and refers to the cardinality of a set.

2 Model setup

Let us consider Q emitting signal sources and P sensor elements in the array. Each source direction is defined by a 2-dimensional vector , s.t., all nominal/true known directions, without any

disturbances, are stacked in Propagation conditions induce wavelength dependent distortions, leading to apparent source directions different from the true ones. Under the narrowband assumption, the array response matrix reads

includes the known Cartesian coordinates describing each sensor location in the array, s.t., for . Therefore, the narrowband signals measured by all antennas is written as follows, for the n-th time sample and wavelength

where the undirectional antenna gains are collected in the complex diagonal matrix and the directional gain responses, assumed identical for all antennas, are modeled by the diagonal matrix Finally, are the i.i.d. calibrator source signal and additive Gaussian thermal noise vectors with their corresponding diagonal covariance matrices , respectively. From (1), we deduce the following covariance matrix

where In this context, the calibration problem consists in estimating the parameter vector of interest the total number of available wavelengths and . To this end, we exploit sample co- variance matrices , defined as for wavelength

In estimation theory, the ML estimator is well-known for its statistical efficiency but not always easy to implement in practice. The Weighting Least Squares approach is an appropriate alternative as it is asymptotically equivalent to the ML for a large number of samples N. Therefore, we wish to minimize the following local cost function, associated to wavelength

where . Most sources are assumed buried beneath the noise and antennas are identical in the array with negligible mutual coupling. The aim of the designed calibration algorithm is to minimize the global cost function in a parallel and step-wise approach, with the total set of available wavelengths. Usually, minimization is conducted w.r.t. one specific parameter while fixing the others in [22].Here, our approach is different: we propose an accelerated version where estimation is performed directly w.r.t. the consensus (hidden) variables, as described in Algorithm 1 and detailed in the following.

3 Description of the proposed estimator

To achieve multi-frequency calibration in the sensor array, coherence is imposed along wavelength subbands for both directional and undirectional gains, by imposing available constraints or enforcing smooth variation. The choice of the basis functions is motivated by the application under analysis and can be adapted accordingly.

3.1 Coherence model for the undirectional antenna gains

To impose coherence along subbands, we introduce a set of smooth wavelength dependent basis functions and express the gains as linear combinations. Let us define , the consensus vector for the p-th sensor with unknown linear coefficients. Therefore, for

stands for the polynomial terms, describing the variation of the undirectional gains w.r.t. wavelength. For instance, we

can consider the typical basis function is the studied frequency of interest with c the speed of light and reference frequency [22, 23]. By stacking all vectors , we obtain the global consensus vector leading to

with

3.2 Coherence model for the directional gains

Similarly as for the undirectional gains, the coherence model is defined as follows: let us consider , such that for

in which is the vector of hidden variables for the q-th calibrator source, associated to directional gains is the corresponding basis vector. As in section 3.1, all are stacked in

, finally leading to

with . We assume identical behavior for all sources but the process can be straightforwardly adapted to different behavior. In [22], the directional gains in were assumed inversely proportional to the algorithm can be adjusted to any general existing models.

3.3 Distributed network with a fusion center

Dealing with large data volumes delivered by advanced sensor array systems requires computationally efficient calibration algorithms, with a huge number of unknowns to solve. To improve both computational cost and estimation accuracy, distributed calibration has been proposed by exploiting data parallelism across frequency. Contrary to a centralized hardware architecture which processes all frequency bands at a single location and is therefore computationally challenging, distributed optimization introduces more than one compute agents and analyzes the data simultaneously across smaller frequency intervals [10]. By distributing the total computations across the network, we gain a significant reduction in operational and energy cost and each agent receives information indirectly across the whole frequency range, thus improving the calibration accuracy. To handle this, let us consider Z computational agents disposed on a network. Each agent has access to some wavelengths . The corresponding unknown pa- rameters in p are estimated locally and consensus is enforced among agents by imposing constraints in (4) and (6).

To start with, let us focus on estimation of the undirectionnal sensor gains in section 3.1. We define as the local copy of the common optimization variable -th agent and we note the set of all in the network. Calibration is reformulated as the following constrained problem

where is the cost function for the z-th agent, i.e., for depends on the local variable and is associated to data To solve this problem, we use the augmented Lagrangian, given by [24]

where Lagrange multipliers and is the regularization term. We resort to the consensus ADMM in the scaled form by introducing the scaled dual variable [10]. The three updates of the iterative

algorithm are therefore given by

where t is the iteration counter. Minimization (9) leads to the following average, computed at the fusion center and sent to all agents in the network,

from which the undirectional gains can be directly deduced with (4). The local minimization step in (8) is the computationally most expensive one. To this end, we adopt an iterative approach and notice that the problem is separable w.r.t. each , i.e., w.r.t. each agent. Let us assume as two independent variables [25]. We then minimize w.r.t. , considering as fixed and neglecting the diagonal elements in the cost function. In this case, the local cost function becomes separable w.r.t. the sub-vectors of local consensus vector for the p-th sensor at the z-th agent. The following decompositions w.r.t. the sensor elements are also possible

and ρ

corresponds to the cost function for the p-th row of, which only depends on since the remaining parameters are considered as fixed in this step. Let us define the operator converts to a vector the p-th row of a matrix and removes the p-th element of this selected vector. We also introduce the quantity (reference source model) and the following vectors

in which In addition, let us consider the

in (12) as and finally obtain the following estimate

3.4 Distributed network with no fusion center

We consider a specific formulation of the ADMM where every node in the network performs calibration locally and consensus is only reached with clearly identified neighbours without fusion center [12]. We note the index set that corresponds to the neighbours of the z-th agent. The considered network architecture is exposed in Figure 1 where for example, define the quantity as the copy available at the z-th agent, transferred to the y-th agent. In such context, the minimization problem becomes

where the auxiliary variables impose consensus contraints on two neighboring agents and are meant to be local copies of . The decentralized strategy enables to cooperatively minimize a sum of local objective functions, the final aim being to converge to a common value, with fast convergence speed and good estimation performance [26]. To obtain a more compact form of the problem in (15), we define leading to

with . As in section 3.3, the scaled version of the ADMM leads to

and through decomposition of the problem in (17) w.r.t. sensor dependence, we obtain

with . The selected variables are obtained from via an appropriate selection matrix. After considering the projection onto B and denoting the messages passed between the agent as

we solve (18) thanks to

The steps of the proposed distributed method for calibration of sensor gains are exposed in Algorithm 1.2.

3.5 Estimation of directional gains

In this section, we describe the part of the algorithm dedicated to the estimation of DoA and directional gains , for fixed sensor gains, with a sparse and distributed implementation. Assuming a sparse observed scene, we define dictionaries of steering matrices for

notes the total number of directions on the grid. The sparse vectors in contain the corresponding squared direc-

tion dependent gains. The covariance model is rewritten as

estimation and satisfy both sparsity and positivity requirements, we use the Distributed Iterative Hard Thresholding (IHT) [27, 28]. But contrary to [22],

the following hard-thresholding operator

ered to provide access to the DoA of the q-th source, and a first estimate of the directional gain . The quantity refers to the q-th column of a ma- trix, the expression discards the elements corresponding to the diagonal of and the hard thresholding operator -largest components and sets the remaining entries equal to zero. Finally, thanks to (5) and dealing with the consensus variables as in section 3.4, the minimization problem becomes

where we benefit from the previous hard-thresholding estimate to define

. As previously, we impose consensus between neighbours thanks to some auxiliary variables but due to lack of space, we only present here the resulting local update for

where , we obtain an estimate of and process the next source, as shown in Algorithm 1.3.

4 Numerical simulations

In order to evaluate the method, we consider realistic simulations for the radio astronomy context where the new generation of phased array systems such as the Low Frequency Array (LOFAR) and the Square Kilometre Array (SKA) requires the development of new advanced signal processing techniques for calibration purpose [4, 29]. Indeed, lack of calibration leads to dramatic effects and distortions in the reconstructed images. We consider P = 60 antennas spread over a five-armed spiral [30, 31], which corresponds to the LOFAR’s Initial Test Station. Let us assume a sky model with Q = 3 strong calibrator sources and weak unknown sources in the background. The reference frequency and we consider frequencies ranging from 29.6 MHz to 30.4 MHz, with Z = 3 agents in the network and . The polynomial orders are chosen as The consensus variables are initialized as zeros and the squared directional gains are generated thanks to power law functions

4.1 Influence of the number of frequency channels

First of all, we investigate the statistical performance of the proposed distributed algorithm as a function of the number of samples N or the Signal-to-Noise Ratio (SNR). The SNR is defined as the ratio between the sum of apparent powers for all Q sources and the noise power. Results are averaged for 100 Monte-Carlo runs. In Figure 2, we plot the three following cases: F = 3 and each agent handles one frequency, i.e., (green curve), (blue curve) and (red curve). In Figure 2 (a), we plot the Root Mean Square Error (RMSE) as a function of N for the undirectional gains , defined as for fixed SNR = dB. A similar figure is presented in Figure 2 (b), for the source directions , as a function of the SNR and fixed illustrate the performance by comparing with the mono-calibration scenario where each agent handles one single frequency, independently. We notice that mono-calibration is clearly improved, by using a distributed procedure where the whole information is flowing through the entire network.

4.2 Influence of the network architecture

We aim to show the advantages of the proposed distributed network with no fusion center and only exchange of local information between neighboring agents, in terms of complexity. With similar number of iterations in all loops of the algorithm, different estimation performance are attained in Figure 2 (a) while similar RMSE is reachable in Figure 2 (b) but with an additional computational cost if there is a fusion center (an increase of at least a factor 5 in computing time).

4.3 Convergence analysis

We illustrate the convergence behavior of the proposed algorithm by analyzing the following residuals as function of the iteration number. Depending on the iteration in Algorithm 1, we plot the primal residual as a function of the iteration number of Algorithm 1.2, defined as

Likewise, we also study the different estimates between agents through

Similar statistical behavior can be ob- tained for corresponding residuals in Algorithm 1.3

5 Conclusion

In this work, we proposed an iterative algorithm for parallel calibration, applied in a general context of sensor array processing: complex electronic gains are imprecisely known and propagation disturbances lead to deviations in the source locations. In order to reduce the communication overhead, the spe-cific variation of parameters across wavelength is exploited in a distributed network with no fusion center and local exchange of information between adjacent connected nodes. The two main steps of the algorithm are based on the scalable form of the ADMM and distributed IHT procedures. We highlighted the effectiveness and time efficiency of the proposed method using simulated data, even in the presence of non-calibrator sources at unknown directions.

Figure 1: Example of distributed network with no fusion center.

References

[1] M. Brossard, V. Ollier, M. N. El Korso, R. Boyer, P. Larzabal, Fast decentralized multi-frequency calibration for doa estimation with distributed sensors, in: DAT 2020, Proceeding Conf. Algiers, Algeria, 2020.

[2] P. Stoica, K. C. Sharman, Maximum likelihood methods for direction-of- arrival estimation, IEEE Transactions on Acoustics, Speech and Signal Processing 38 (7) (1990) 1132–1143.

[3] S. Vorobyov, A. B. Gershman, K. M. Wong, Maximum likelihood direction-of-arrival estimation in unknown noise fields using sparse sensor arrays, IEEE Transactions on Signal Processing 53 (1) (2005) 34–43.

[4] A.-J. van der Veen, S. J. Wijnholds, Signal processing tools for radio astronomy, in: Handbook of Signal Processing Systems, Springer, 2013, pp. 421–463.

[5] L. C. Godara, Application of antenna arrays to mobile communications, part II: Beam-forming and direction-of-arrival considerations, Proceedings of the IEEE 85 (8) (1997) 1195–1245.

[6] B. C. Ng, C. M. S. See, Sensor-array calibration using a maximum- likelihood approach, IEEE Transactions on Antennas and Propagation 44 (6) (1996) 827–835.

[7] B. C. Ng, A. Nehorai, Active array sensor localization, Elsevier Signal processing 44 (3) (1995) 309–327.

[8] A. J. Weiss, B. Friedlander, Array shape calibration using sources in unknown locations-a maximum likelihood approach, IEEE Transactions on Acoustics, Speech, and Signal Processing 37 (12) (1989) 1958–1966.

[9] S. J. Wijnholds, S. Chiarucci, Blind calibration of phased arrays us- ing sparsity constraints on the signal model, in: 24th European Signal Processing Conference (EUSIPCO), 2016, pp. 270–274.

[10] S. Boyd, N. Parikh, E. Chu, B. Peleato, J. Eckstein, Distributed opti- mization and statistical learning via the alternating direction method of multipliers, Foundations and Trends in Machine Learning 3 (1) (2011) 1–122.

[11] S. Kazemi, S. Yatawatta, S. Zaroubi, Clustered calibration: an im- provement to radio interferometric direction-dependent self-calibration, Monthly Notices of the Royal Astronomical Society 430 (2) (2013) 1457– 1472.

[12] T. Erseghe, A distributed and scalable processing method based upon ADMM, IEEE Signal Processing Letters 19 (9) (2012) 563–566.

[13] W. Shi, Q. Ling, K. Yuan, G. Wu, W. Yin, On the linear convergence of the ADMM in decentralized consensus optimization., IEEE Transactions on Signal Processing 62 (7) (2014) 1750–1761.

[14] J. F. C. Mota, J. M. F. Xavier, P. M. Q. Aguiar, M. Puschel, D-ADMM: a communication-efficient distributed algorithm for separable optimization, IEEE Transactions on Signal Processing 61 (10) (2013) 2718–2723.

[15] V. Ollier, M. N. El Korso, A. Ferrari, R. Boyer, P. Larzabal, Robust distributed calibration of radio interferometers with direction dependent distortions, accepted for publication by Elsevier Signal Processing.

[16] S. Yatawatta, Fine tuning consensus optimization for distributed radio interferometric calibration, in: 24th European Signal Processing Conference (EUSIPCO), Budapest, Hungary, 2016.

Figure 2: (a) RMSE on the undirectional gains as function of the number of samples N, (b) RMSE on the apparent source directions as function of the SNR.

[17] S. Yatawatta, S. Kazemi, S. Zaroubi, GPU accelerated nonlinear op- timization in radio interferometric calibration, in: Innovative Parallel Computing (InPar), IEEE, San Jose, CA, 2012, pp. 1–6.

[18] D. Malioutov, M. Cetin, A. S. Willsky, A sparse signal reconstruction perspective for source localization with sensor arrays, IEEE Transac-

Figure 3: (a) Statistical comparison between different network topologies for same computational cost, (b) Statistical comparison between different network topologies for different computational cost.

tions on Signal Processing 53 (8) (2005) 3010–3022.

[19] C. Steffens, P. Parvazi, M. Pesavento, Direction finding and array cal- ibration based on sparse reconstruction in partly calibrated arrays, in: 8th Sensor Array and Multichannel Signal Processing Workshop (IEEE SAM), 2014, pp. 21–24.

Figure 4: (a) Primal residual and (b) estimates difference of the local consensus variables among agents as function of the iteration t in Algorithm 2, for different values of the iteration i in Algorithm 1.

[20] M. Haardt, M. Pesavento, F. Roemer, M. N. El Korso, Subspace meth- ods and exploitation of special array structures, in: Electronic Reference in Signal Processing:Array and Statistical Signal Processing (M. Virberg, ed.), Vol. 3, Academic Press Library in Signal Processing, Elsevier Ltd., 2014.

[21] V. Ollier, M. N. El Korso, R. Boyer, P. Larzabal, M. Pesavento, Joint ML calibration and DOA estimation with separated arrays, in: International Conference on Acoustics, Speech and Signal Processing (IEEE ICASSP), Shanghai, China, 2016.

[22] M. Brossard, M. N. El Korso, M. Pesavento, R. Boyer, P. Larzabal, S. J. Wijnholds, Parallel multi-wavelength calibration algorithm for radio astronomical arrays, Elsevier Signal Processing Journal 145 (2018) 258–271.

[23] S. Yatawatta, Distributed radio interferometric calibration, Monthly Notices of the Royal Astronomical Society 449 (4) (2015) 4506–4514.

[24] L. Li, X. Wang, G. Wang, Alternating Direction Method of Multipliers for separable convex optimization of real functions in complex variables, Mathematical Problems in Engineering 2015 (Article ID 104531) (2015) 1–14.

[25] S. Salvini, S. J. Wijnholds, Fast gain calibration in radio astronomy using alternating direction implicit methods: Analysis and applications, Astronomy & Astrophysics 571 (A97).

[26] T. Erseghe, D. Zennaro, E. Dall’Anese, L. Vangelista, Fast consensus by the alternating direction multipliers method, IEEE Transactions on Signal Processing 59 (11) (2011) 5523–5537.

[27] T. Blumensath, M. E. Davies, Normalized iterative hard thresholding: Guaranteed stability and performance, IEEE Journal of selected topics in signal processing 4 (2) (2010) 298–309.

[28] S. Patterson, Y. C. Eldar, I. Keidar, Distributed compressed sensing for static and time-varying networks, IEEE Transactions on Signal Processing 62 (19) (2014) 4931–4946.

[29] V. Ollier, M. N. El Korso, R. Boyer, P. Larzabal, M. Pesavento, Relaxed concentrated MLE for robust calibration of radio interferometers, in: 24th European Signal Processing Conference (EUSIPCO), Budapest, Hungary, 2016, pp. 280–284.

[30] S. J. Wijnholds, J. D. Bregman, A. J. Boonstra, Sky noise limited snap- shot imaging in the presence of RFI with Lofar’s Initial Test Station, Experimental Astronomy 17 (1-3) (2004) 35–42.

[31] V. Ollier, M. N. El Korso, R. Boyer, P. Larzabal, M. Pesavento, Robust calibration of radio interferometers in non-Gaussian environment, IEEE Transactions on Signal Processing 65 (21) (2017) 5649–5660.

[32] E. Ghadimi, A. Teixeira, I. Shames, M. Johansson, Optimal parameter selection for the alternating direction method of multipliers (ADMM): quadratic problems, IEEE Transactions on Automatic Control 60 (3) (2015) 644–658.