同济大学：《传播信道特征估计和建模》课程教学资源（教案讲义）Chapter 07 Statistical channel parameter estimation.pdf_大学文库

7 Statistical channel parameter estimation In Chapter 6, different estimators of the deterministic parameters of propagation paths have been introduced. These methods are widely used to process measurement channel sounding data for characterization and modeling of channels in many scenarios. Based on their applications, channel models can be divided into two categories, i.e. the deterministic channel models, and the stochastic channel models. The former, e.g. the time-delay-line (TDL) models in 3GPP channel standards Spatial channel model for Multiple Input Multiple Output (MIMO) simulations (Release 7) (2007), can be used as calibration models for conformity testing of communication technologies and systems. The latter is more suitable for generating random channel realizations for the performance testing in either link level or system level. The stochastic channel models are usually obtained by first estimating the parameters of channel, such as the composite spreading parameters of a channel, e.g. delay spread, Doppler frequency spread, and the parameters of individual propagation paths from individual channel snapshots, and then secondly, extracting the statistics of the channel parameters for modeling. In this chapter, we describe some recently developed parameter estimation algorithms, which can be used to estimate directly the statistical channel parameters from the measurement data. Some of results presented in this chapter have been published in Yin et al. (2007a,b, 2008, 2006a, 2011). 7.1 Dispersive component estimation algorithms In wireless environments, the individual contributions of distributed scatterers to the received signal can be spread in delay, direction of departure, direction of arrival, Doppler frequency and polarization. In this section, the 1 st-order generalized array manifold (GAM) model proposed in Asztély et al. (1997) is used to approximate the individual contributions of slightly distributed scatterers (SDSs) to the signal received with a multiple-element antenna. A new definition of SDS is proposed based on the two largest eigenvalues of these signal components. With this definition, a distributed scatterer is an SDS when its signal component is closely approximated by the 1 st-order GAM model. A measure is described that quantitatively assesses the degree of the approximation. Azimuth dispersion caused by an SDS is characterized by means of two parameters, i.e. the nominal azimuth (NA) and the azimuth spread (AS). Based on the 1 st-order GAM model, we derive NA and AS estimators using the standard deterministic and stochastic maximum likelihood (ML) methods, as well as a new extension of the conventional MUSIC algorithm. Due to the discrepancy between the 1 st-order GAM model and the effective signal model, the AS estimators are biased. Instead of reducing the bias by using a multi-dimensional look-up table Bengtsson and Ottersten (2000), which is usually cumbersome to generate when the used antennas are non-isotropic, an empirical technique is proposed that adaptively selects the size of the array aperture in such a way to guarantee a good agreement between the two models. Simulation results demonstrate the performance improvement achieved with these methods compared to previously published related techniques. In this section, an approximate probability distribution is erived for the error of the ML azimuth estimator derived based on the conventional specular-scatterer model, when this estimator is used in a scenario with a single SDS. This result is an extension to the incoherent-SDS case of a previous result derived for the coherent-SDS case Asztély and Ottersten (1998). It is shown that in the incoherent-SDS case, the probability of occurrence of large estimation errors is significantly reduced by considering more observation samples. 7.1.1 A brief review of dispersive parameter estimators In a radio propagation environment, situations frequently occur where scatterers have a certain geometrical extent which is small in the view of the receiver (Rx), or local scattering exists around a transmitter (Tx) located far away This is a Book Title Name of the Author/Editor c XXXX John Wiley & Sons, Ltd

176 Statistical channel parameter estimation from the Rx. In both cases, the received signal contributed by each of these scatterers or clusters of local scatterers can be conceived as the sum of contributions originating from multiple sub-paths with slightly different directions of arrival. Retaining the terminology introduced in Asztély and Ottersten (1998); Bengtsson and Ottersten (2000); Fleury (2000) we refer to such scatterers or clusters of local scatterers as slightly distributed scatterers (SDSs). A quantitative definition of a SDS matching the theoretical framework and purpose of this paper is provided. Notice that this definition differs from the one given in these references. Direction dispersion in the channel induced by an SDS can be characterized by means of the nominal direction and the spread of the direction power spectrum of the signal contributed by this SDS Fleury (2000). For simplicity, we consider horizontal-only propagation in this contribution. In this case, dispersion is characterized by a nominal azimuth (NA) and an azimuth spread (AS) Bengtsson and Ottersten (2000); Shahbazpanahi et al. (2001); Trump and Ottersten (1996). By abuse of terminology, we refer to the azimuth power spectrum, NA and AS as the azimuth power spectrum, the NA and the AS of the SDS. In this contribution, we are interested in the design and performance investigation of NA and AS estimators for channel sounding. Practical limitations and constraints in the sounding equipment, such as non-isotropic antenna arrays, limited processing time, calibration errors, etc., put the following stringent constraints on these estimators: i) applicability with arbitrary transmit and receive arrays, e.g. arbitrary array configuration, non-isotropic element responses, ii) low computational complexity, iii) robust performance against calibration errors in the measured array response, and iv) applicability in the case where the components in the channel response contributed by an SDS are incoherent (i.e. totally uncorrelated). The situation considered in Item iv) is referred to as the incoherent distributed (ID) case in Shahbazpanahi et al. (2001). The coherent distributed (CD) case where the samples of the component contributed by an SDS are coherent (i.e. fully correlated) is not considered in this contribution. In the sequel, we briefly review the NA and AS estimators available in open literature, and comment on their applicability in the sought application scenario, i.e. MIMO channel sounding. We will see that most hitherto published estimators do not fulfill all the above four constraints and therefore cannot be used in this scenario. Conventional high-resolution channel parameter estimation algorithms Bresler (1986); Fleury et al. (1999); Roy and Kailath (1989); Schmidt (1986); Viberg (1991) are derived based on the specular-scatterer (SS) model, which assumes that the propagation environment consists of point scatterers located in the far field of the Tx and Rx arrays Lee et al. (1997). As a result, the contribution of each specular scatterer to the received signal is modeled as a (specular) plane wave. These algorithms are not suitable for the estimation of the NAs and the ASs of SDSs, due to two reasons. First, they do not provide with an estimate of the AS. Secondly, the probability of occurrence of large estimation errors is significant Bengtsson and Völcker (2001). In recent years, various model-based methods have been proposed for the estimation of the NA and the AS of SDSs. These estimators can be categorized into two groups: the estimators based on a distribution model and the estimators based on a linear approximation model. The estimators based on a distribution model rely on a probability distribution characterizing the shape of the azimuth power spectrum of individual SDSs Besson and Stoica (1999); Han et al. (2006); Meng et al. (1996); Qiang and Zhishun (2007); Ribeiro et al. (2004); Trump and Ottersten (1996); Valaee et al. (1995); Wang and Zoubir (2007); Zoubir et al. (2007). The uniform distribution Besson and Stoica (1999); Meng et al. (1996); Sieskul (2006); Valaee et al. (1995); Zoubir et al. (2007), the (truncated) Gaussian distribution Besson and Stoica (1999); Meng et al. (1996); Qiang and Zhishun (2007); Sieskul (2006); Trump and Ottersten (1996); Wang and Zoubir (2007), the Laplacian distribution Sieskul (2006) and the von-Mises distribution Fleury (2000); Ribeiro et al. (2004) have been proposed for horizontal-only propagation. These estimators are not appropriate for applications in channel sounding due to two reasons. First, no prior information of the shape of the azimuth power spectrum of the SDS is available. Second, these algorithms exhibit high computational complexity. Their computation requires to solve at least a onedimensional optimization problem for a certain objective function: the objective function is usually expressed in form of an integral that needs to be calculated numerically when the antenna array response is non-isotropic. Some recently proposed algorithms within this category achieve low computational complexity by exploiting the specific structure or applying a series expansion of the sample covariance matrix Han et al. (2006); Qiang and Zhishun (2007); Wang and Zoubir (2007); Zoubir et al. (2007). However these properties of the sample covariance matrix do not hold when the array is nonlinear and the array elements are non-isotropic. Thus, these estimators are not suitable for the scenarios of interest

Statistical channel parameter estimation 177 The estimators based on a linear approximation model include the estimators derived using the two-specularscatterer model Bengtsson and Ottersten (2000); Souden et al. (2008), the two-SDS model (Shahbazpanahi et al. 2001, Section VI), and the generalized array manifold (GAM) model Asztély et al. (1997); Jeong et al. (2001); Shahbazpanahi et al. (2001); Shahbazpanahi (2004); Tan et al. (2003). The two-specular-scatterer model is motivated from the observation that in the ID case, the covariance matrix of the signal contributed by a single SDS is well approximated by the covariance matrix of the signals resulting from two incoherent specular scatterers. The Spread-F algorithm, with F denoting a standard azimuth estimator, like Root-MUSIC, ESPRIT, MODE, etc., is derived based on this model Bengtsson and Ottersten (2000). This algorithm calculates the NA of an individual SDS to be the mean of two azimuth estimates, while the AS estimate is defined to be half the distance between these two azimuth estimates. Pairing of the two azimuth estimates corresponding to the same SDS in a multi-SDS scenario is a major implementation issue, especially in the case where the SDSs are closely spaced in azimuth. Furthermore, a look-up-table solution has to be used to reduce the bias of the AS estimator that results due to the discrepancy between the two-specular-scatterer model and the effective model Bengtsson and Ottersten (2000). In the case where channel sounding is performed with non-isotropic array responses, this look-up table needs to be generated for multiple entries. For instance, in the case of horizontal-only propagation, the NA, AS and signal-to-noise ratio (SNR) are entries of the look-up table. The high complexity involved prohibits the implementation of this approach in the scenario of interest. The GAM model has been used to derive NA and AS estimators in Asztély et al. (1997); Shahbazpanahi et al. (2001); Shahbazpanahi (2004). This model utilizes the first-order Taylor series expansion of the array responses to approximate the effective impact of each individual SDS on the received signal. For brevity, we refer to a GAM model using the n th-order Taylor series expansion to be a n th-order GAM model. We now briefly review these methods. The estimators proposed in Asztély et al. (1997) rely on the MUSIC and the noise subspace fitting algorithms. They cannot be employed in the ID case. The estimators proposed in Shahbazpanahi et al. (2001) make use of the ESPRIT algorithm. They are only applicable when the antenna arrays have two identical antennas in each matched pair Roy and Kailath (1989) and, additionally the distance between the two antennas in each pair is much less than the wavelength Shahbazpanahi et al. (2001). This condition is usually not fulfilled in the scenario of interest. The NA and AS estimators derived in Shahbazpanahi (2004) apply the covariance fitting principle. A major implementation issue of these estimators is the knowledge of the second-order derivative of the array response required for the computation of the AS estimate. Calibration errors in the measurement of the array response generate significant discrepancies in the high-order derivatives of the response. As a consequence, the performance of the estimator is degraded dramatically. To summarize, among all existing estimation methods, only the standard ML method applied in combination with the 1 st-order GAM model lead to NA and AS estimators that satisfy our four requirements. In this section, we present these estimators as well as a new MUSIC algorithm utilizing the 1 st-order GAM model, and propose a new technique to improve the performance of the AS estimators. The main contributions in this work are i) Derivation of an approximate probability density function (pdf) of the error of the ML azimuth estimator derived based on the conventional specularscatterer model when this estimator is used in a scenario with a single SDS. The approximate pdf reported Bengtsson and Völcker (2001) is a special case with a single observation of the new pdf; ii) Extension of the standard MUSIC algorithm for estimation of the parameters characterizing SDSs in the IS case; iii) Definition of a new measure that quantitatively assesses the degree of the agreement between the 1 st-order GAM model and the effective signal model. This measure underlies the derivation of an array size adaptation (ASA) technique that selectively adjusts the array size in such a way to guarantee a good agreement. The AS Estimators augmented with this ASA technique exhibit low bias. It is worth mentioning that this work has a definite practical significance, since the scenario of interest here is frequently encountered in mobile radio channel sounding. Due to various reasons, situations where uniform linear arrays with isotropic array elements are used, that make it possible to use some of the already proposed estimation methods are seldom. 7.1.2 Effective Signal Model Without loss of generality we restrict attention to azimuth dispersion generated by SDS at the Rx site, when the Rx is equipped with an antenna array. We introduce a model describing the effective contribution of SDSs to the signal received with the Rx array. This model is referred to as the effective signal model throughput the paper. Note that the model and the estimation methods derived based on the model immediately apply to assess direction dispersion by SDS at the Tx site as well, provided the Tx is equipped with an antenna array. In a propagation scenario with a single SDS, the signal at the output of an M-element Rx array is viewed as composed of the contributions of signals originating from waves propagating along multiple sub-paths distributed with respect to (w.r.t.) the direction of arrival. For simplicity consideration, we assume horizontal-only propagation

178 Statistical channel parameter estimation and neglect all dispersion effects but dispersion in azimuth of arrival. We consider a system where the Rx is equipped with M correlators. In the narrow-band scenario considered here the vector y(t) ∈ CM containing the outputs of the correlators at time t is expressed as y(t) = X L ℓ=1 aℓ(t)c(φ¯ + φ˜ ℓ) + w(t). (7.1) In (7.1), the total number L of sub-paths originating from the SDS is assumed to be large, aℓ(t) denotes the complex weight of the ℓth sub-path, and c(φ) ∈ CM represents the array response in azimuth φ. The azimuth of the ℓth subpath is decomposed as the sum of a NA φ¯ and a deviation φ˜ ℓ from φ¯. Throughout the paper, angles are expressed in radian with range [−π, π) in the theoretical investigations, while they are given in degree with range [−180◦ , +180◦ ) in the simulation studies. Addition of angles is defined in such a way that the resultant angle lies in these ranges. The noise vector w(t) in (7.1) is a circularly symmetric, spatially and temporally white M-dimensional Gaussian process with component spectral height σ 2 w. We assume that totally N observation samples are considered. These samples are collected at the time instances t1, . . . , tN , i.e. t ∈ {t1, . . . , tN } ⊂ R. The following assumptions are made 1. The azimuth deviations φ˜ ℓ, ℓ = 1, . . . , L, are independent and identically distributed with zero-mean. The azimuths of the sub-paths are concentrated with high probability around the NA φ¯, i.e. the azimuth deviations φ˜ ℓ, ℓ = 1, . . . , L, are small with high probability. 2. The sub-path weight processes a1(t), . . . , aL(t) are uncorrelated zero-mean complex circularly-symmetric widesense-stationary processes with autocorrelation function Raℓ (τ) .= E aℓ(t)a ∗ ℓ (t + τ) . Here, (·) ∗ denotes complex conjugation. Moreover, the sub-path weights have equal variance, i.e. Ra1 (0) = · · · = RaL (0). 3. Any two random elements in the set consisting of the azimuth deviations and the sub-path weights are uncorrelated. 4. Temporal samples of the sub-path weights are uncorrelated, i.e. Raℓ (tn′ − tn) = 0, n6=n ′ . In a scenario with D SDSs, (7.1) becomes y(t) = X D d=1 X Ld ℓ=1 ad,ℓ(t)c(φ¯ d+φ˜ d,ℓ)+ w(t), (7.2) where Ld denotes the number of sub-paths originating from the dth SDS. We make the additional assumption: 5) The random elements characterizing the signal contributed by each of the D SDSs satisfy Assumptions 1)–4). The probability distributions of the azimuth deviations of distinct SDSs may be different. In addition, any two distinct SDSs are uncorrelated, i.e. more specifically, any two random elements related each to distinct SDSs are uncorrelated. Assumptions 1)–5) corresponds to the case with ID sources described in Shahbazpanahi et al. (2001). 7.1.3 Specular-Scatterer (SS) Model Estimation The SS model is widely used in conventional high-resolution parameter estimation algorithms Bresler (1986); Fleury et al. (1999); Roy and Kailath (1989); Schmidt (1986). In this section, we briefly introduce the SS model and the ML azimuth estimator derived based on this model. Then we derive an approximation of the error distribution of this estimator when this estimator is applied in a scenario with one SDS. We consider the ID case with N (N ≥ 1) independent observation samples. The case with N = 1, corresponding to the CD case, was discussed previously in Asztély and Ottersten (1998).

Statistical channel parameter estimation 179 The SS Model and the ML Azimuth Estimator The SS model assumes that point scatterers re-radiate impinging waves and that the Rx is in the far field of these scatterers. Thus, specular plane waves are incident at the Rx site. Assuming a scenario with D point scatterers, the output signal vector (7.1) is approximated according to y(t) ≈ ySS(t) .= X D d=1 αd(t)c(φ¯ d) + w(t) = C(φ¯)α(t) + w(t), (7.3) where φ¯ .= [φ¯ 1, . . . , φ¯D], C(φ¯) .= [c(φ¯ 1), . . . , c(φ¯D)] and α(t) .= [α1(t), . . . , αD(t)]T. The subscript SS emphasizes that the SS assumption is used to approximate the signals contributed by scatterers. The deterministic ML estimator of φ¯ derived based on the SS model (SS-ML) reads Böhme (1984); Krim and Viberg (1996) φ ˆ¯ SS−ML=arg max φ¯ {tr[ΠC(φ¯)Σˆ y]}, (7.4) where tr[·] denotes the trace of the matrix given as an argument, ΠC(φ¯) .=C(φ¯)C(φ¯) † is the projection operator onto the column space of C(φ¯), C(φ¯) † .= C(φ¯) HC(φ¯) −1C(φ¯) H is the pseudo-inverse of C(φ¯) and Σˆ y = 1 N PtN t=t1 y(t)y(t) H is the sample covariance matrix. Here, [·] H denotes Hermitian transposition. Distribution of the Estimation Error The pdf of the azimuth estimate φ ˆ¯ SS−ML is affected by the discrepancy between the effective signal model (7.2) and the SS model (7.3). Appendix ?? shows that in a noiseless single-SDS scenario, the NA estimation error φˇ .= φ ˆ¯ SS−ML − φ¯ can be approximated as φˇ ≈ PtN t=t1 |α(t)| 2Re{ β(t) α(t) } PtN t=t1 |α(t)| 2 , (7.5) where α(t) .= PL ℓ=1 aℓ(t), β(t) .= PL ℓ=1 φ˜ ℓaℓ(t) are computed using the true values of aℓ(t) and φ˜ ℓ, and Re{·} denotes the real part of the argument. Making use of Assumptions 2) and 3) in Subsection 7.1.2 and invoking the Central Limit Theorem (Shanmugan and Breipohl 1988, Section 2.8.2), α(t) and β(t) can be approximated as uncorrelated complex circularly-symmetric Gaussian random processes. Based on this assumption, the pdf of the right-hand side (r.h.s.) in (7.5) can be calculated. This pdf provides an approximation of the pdf of φˇ: fφˇ(φ) ≈ Γ(N + 1 2 ) √ π Γ(N) · 1 σφ˜ · 1 1 + φ 2 σ 2 φ˜ (N+ 1 2 ) , (7.6) where Γ(·) is the gamma function. The pdf in the r.h.s. in (7.6) can be used to compute approximations of the moments of φˇ. For instance, Var[φˇ] ≈ Γ(N − 1) 2Γ(N) · σ 2 φ˜ (7.7) E[|φˇ|] ≈ Γ(N − 1 2 ) √ π Γ(N) · σφ˜. (7.8) The derivations of (7.6)–(7.8) are given in Appendix ??. For N = 1, the r.h.s. term in (7.7) is infinite due to the heavy tails of the r.h.s. pdf in (7.6), as reported in Asztély and Ottersten (1998). Actually, the variance of φˇ is always finite as the NA estimation error φˇ is confined within [−π, π). However, the fact that the r.h.s. of (7.7) is infinite indicates that the variance of φˇ is large. As a consequence, large estimation errors occur with significant probability. When N > 1, the r.h.s. of (7.7) is finite and decreases as N increases. Thus, in the ID case, the probability of large estimation error can be reduced by considering more observation samples.

180 Statistical channel parameter estimation 7.1.4 Parameter Estimation using the 1 st-order GAM Model In this section, we introduce the 1 st-order GAM model and derive estimators of the model parameters using standard deterministic and stochastic ML methods, as well as a novel MUSIC algorithm. An AS estimator based on these parameter estimators is also proposed. The 1 st-order GAM Model The GAM model Asztély et al. (1997) makes use of the fact that the deviations φ˜ d,ℓ are small with high probability. In the 1 st-order GAM model, the first-order Taylor series expansion of the array response is used to approximate the effective impact of the SDSs on the received signal. We regard a distributed scatterer as an SDS when its contribution to the output signal of the Rx array is closely approximated using the 1 st-order GAM model. A measure of the fit between the approximation model and the effective signal model is provided in Section 7.1.5. We first consider a single-SDS scenario. The function c(φ¯ + φ˜ ℓ) in (7.1) can be approximated by its first-order Taylor series expansion at φ¯. Inserting this approximation for each c(φ¯ + φ˜ ℓ) in (7.1) yields the 1 st-order GAM model Asztély et al. (1997) y(t) ≈ yGAM(t) .= X L ℓ=1 aℓ(t) c(φ¯)+φ˜ ℓc ′ (φ¯) + w(t) =α(t)c(φ¯) + β(t)c ′ (φ¯) + w(t), (7.9) where c ′ (φ¯) .= dc(φ) dφ φ=φ¯ . In matrix notation, (7.9) reads yGAM(t) = F(φ¯)ξ(t) + w(t) (7.10) with F(φ¯) .=[c(φ¯) c ′ (φ¯)] and ξ(t) .=[α(t), β(t)]T. The autocorrelation functions of α(t) and β(t) is calculated to be respectively, Rα(τ) .= E[α(t)α ∗ (t + τ)] = X L ℓ=1 Raℓ (τ) and Rβ(τ) .= E[β(t)β ∗ (t + τ)] = σ 2 φ˜ · Rα(τ), (7.11) where σ 2 φ˜ .= E[φ˜2 ℓ ]. Note that E[φ˜ ℓ] = 0 according to Assumption 1). The parameter σ 2 φ˜ is the second-central moment of the azimuth deviation. By denoting the variances of α(t) and β(t) with σ 2 α and σ 2 β respectively, we conclude from (7.11) that σ 2 β = σ 2 φ˜ · σ 2 α. (7.12) This equality can also be obtained using the results given in (Shahbazpanahi et al. 2001, (49)–(51)). We refer to the parameter σφ˜ as the AS of the SDS. Note that as shown in Fleury (2000), the natural figure for characterizing direction dispersion is the direction spread. However, in a scenario with horizontal-only propagation and small azimuth deviations, the direction spread can be approximated by σφ˜ expressed in radian Fleury (2000). For example, in the case where the azimuth power spectrum of an SDS is proportional to the von-Mises probability density function Mardia (1975) fφ˜ ℓ (φ) = 1 2πI0(κ) exp{κ cos(φ − φ¯)}, (7.13) where κ denotes the concentration parameter and I0(·) represents the modified Bessel function of the first kind and order 0, the approximation is close provided κ ≥ 7, i.e. σφ˜ ≤ 10◦ Fleury (2000). In a scenario with D SDSs, (7.9) extends to y(t) ≈ yGAM(t) .= X D d=1 αd(t)c(φ¯ d) + βd(t)c ′ (φ¯ d) + w(t) = B(φ¯)γ(t) + w(t), (7.14)

Statistical channel parameter estimation 181 where B(φ¯) .= [c(φ¯ 1), c ′ (φ¯ 1), . . . , c(φ¯D), c ′ (φ¯D)] and γ(t) .= [α1(t), β1(t), . . . , αD(t), βD(t)]T. Under Assumptions 3)–5) in Subsection 7.1.2 the elements in the vector γ(t) are uncorrelated. NA Estimators In this subsection, the standard deterministic and stochastic ML estimation methods as well as a novel MUSIC algorithm are applied using the 1 st-order GAM model to derive estimators of the NAs of SDSs. Deterministic ML (DML) NA Estimator The DML NA estimator based on the 1 st-order GAM model can be derived similarly to the SS-ML azimuth estimator (7.4). Assuming that the weight samples αd(t) and βd(t), t = t1, . . . , tN , d = 1, . . . , D in (7.14) are deterministic, the ML estimator of φ¯ is calculated as Krim and Viberg (1996) φ ˆ¯ DML=arg max φ¯ {tr[ΠB(φ¯)Σˆ y]}. (7.15) The parameters γ(t), t = t1, . . . , tN are estimated as (\γ(t))DML=B(φ ˆ¯) †y(t), t = t1, . . . , tN . (7.16) Stochastic ML (SML) NA Estimator The SML azimuth estimator derived based on the SS model was introduced in Jaffer (1988). We obtain the SML NA estimator based on the 1 st-order GAM model in a similar manner. Making use of the assumptions 1)–5) in Section 7.1.2 and invoking the central limit theorem, the weight samples αd(t) and βd(t), t = t1, . . . , tN , d = 1, . . . , D are uncorrelated complex circularly-symmetric Gaussian random processes with variances σ 2 αd and σ 2 βd respectively. Let Ω be the vector containing the parameters to be estimated: Ω .= [σ 2 w, φ¯ d, σ2 αd , σ2 βd ; d = 1, . . . , D]. (7.17) The ML estimator of Ω is a solution to the maximization problem Krim and Viberg (1996) Ωb SML = arg max Ω {−ln[|ΣyGAM |] − tr (ΣyGAM ) −1Σˆ y }, (7.18) where the covariance matrix ΣyGAM of yGAM(t) in (7.14) reads ΣyGAM = B(φ¯)RγB(φ¯) H + σ 2 wIM . (7.19) Here, IM denotes the M × M identity matrix and Rγ = diag(σ 2 α1 , σ2 β1 , . . . , σ2 αD , σ2 βD ) is the covariance matrix of γ(t). Here, diag(·) denotes a diagonal matrix with diagonal elements listed as argument. The maximization operations in (7.15) and (7.18) require respectively a D-dimensional and a (3D + 1)-dimensional search. The high computational complexity of these search procedures prohibits the implementation of φ ˆ¯ DML and Ωb SML in real applications. As an alternative, the SAGE algorithm Fleury et al. (1999); Yin and Fleury (2005) provides with a low-complexity approximation of these ML estimators. MUSIC NA Estimator The standard MUSIC algorithm Schmidt (1986) derived based on the SS model (7.3) uses the pseudo-spectrum fMUSIC(φ)= kc(φ)k 2 F kc(φ)HEwk 2 F . (7.20) Here, k · kF denotes the Frobenius norm and Ew is an orthonormal basis of the estimated noise subspace calculated from Σˆ y. The azimuths of the D scatterers are estimated to be the arguments of the pseudo-spectrum corresponding to its D highest peaks.

182 Statistical channel parameter estimation We propose a natural extension of the standard MUSIC algorithm for the estimation of the NAs of SDSs based on the 1 st-order GAM model. The extension considers the following generalization of the pseudo-spectrum in (7.20): fMUSIC(φ) = 1 kF˜(φ)HEwk 2 F . (7.21) In the r.h.s. term in (7.21), F˜(φ) is an orthonormal basis of the space spanned by the columns of F(φ). The NAs of the D SDSs are estimated to be the arguments of the pseudo-spectrum corresponding to its D highest peaks. Both the standard MUSIC algorithm and the proposed extension rely on the same principle, i.e. parameter estimates are obtained by minimizing the distance between the subspace spanned by the signal originating from single scatterer and an estimate of this subspace computed from the sample covariance matrix. In the SS case, the signal subspace induced by an SS is spanned by the steering vector c(φ), while in the SDS scenario the subspace induced by an SDS is spanned by the columns of F(φ). In that sense, the latter algorithm is a natural extension of the former one. Following (Edelman et al. 1998, p. 337), the distance between the subspace spanned by the columns of F(φ) and the estimated signal subspace coincides with the Frobenius norm of the difference between the projection matrices of the two subspaces. It can be shown that this distance is proportional to the Frobenius norm of the projection of one subspace onto the null space of the other subspace, i.e. kF˜(φ) HEwk 2 F in our case. Thus, the inverse of the pseudospectrum (7.21) provides with a measure of the distance between the signal subspace spanned by the columns of F(φ) and the estimated signal subspace. A thorough discussion of the relationships between this extended MUSIC algorithm and other previously published extensions of the standard MUSIC algorithm is given in Subsection 7.1.4. Generalization of the Proposed MUSIC Algorithm and its Relation to other Extensions of the Standard MUSIC Algorithm The proposed MUSIC algorithm, which makes use of the pseudo-spectrum (7.21), can be generalized to the scenario where the signals contributed by one component (e.g. a scatterer) span a subspace of any arbitrary dimension. In this case, F˜(φ) is an orthonormal basis of the signal subspace. The argument φ of F˜(φ) may be also multi-dimensional and it is not required that a closed–form expression exists which relates F˜(φ) to φ. For instance, in the case where azimuth dispersion of an SDS is characterized using a pdf, F˜(φ) can be obtained by the eigenvalue decomposition of the covariance matrix calculated using this pdfs. We propose now an alternative interpretation of the proposed MUSIC algorithm using the concept of principal angles between subspaces Golub and Loan (1996), which allows for a comparison with the variant of the MUSIC algorithm published in Christensen et al. (2004). As shown in (Edelman et al. 1998, p. 337), minimizing the distance between two subspaces is equivalent to minimizing the norm of sin(θ), where θ represents the vector containing all principal angles between these two subspaces and sin(·) is the operator computing the element–wise sin of θ. Thus, in our case the NA estimates obtained by maximizing the pseudo-spectrum (7.21) in fact minimize k sin(θ)k where the components of θ are the principal angles between the subspace spanned by the columns of F(φ) and the signal subspace estimated from the sample covariance matrix. This is a reasonable approach in the ID case where the dimension of the signal subspace induced by an SDS is larger than 1. By contrast, the variant of the MUSIC algorithm proposed in Asztély et al. (1997) computes the NA estimates by maximizing the smallest principal angle between the 2-dimensional subspace spanned by the columns of F(φ) and the estimated signal subspace. This maximization is indeed equivalent to the maximization of the objective function λ −1 min(F(φ) HEwEH w F(φ)), with λmin(·) denoting the smallest eigenvalue of the matrix given as argument, described in Asztély et al. (1997) to compute the NA estimates Drmac (2000). The resulting algorithm is applicable when the dimension of the subspace effectively induced by an SDS is equal to one, e.g. in the CD case for which the algorithm was initially designed. The pseudo-spectrum (7.21) can be recast as: fMUSIC(φ)= 1 tr{EH w F(φ)W(φ)F(φ)H Ew} , (7.22) where W(φ) is an azimuth-dependent weighting matrix defined as W(φ) .= F(φ) †F˜(φ)F˜(φ) H(F(φ) † ) H. (7.23) At first glance the representation in (7.22) seems to be similar to the pseudo-spectrum (Krim and Viberg 1996, Eq. (37)) of the weighted MUSIC algorithm. However, the proposed MUSIC algorithm and the standard weighted MUSIC

Statistical channel parameter estimation 183 algorithm have fundamental differences. First of all, it is impossible to recast the pseudo-spectrum (7.22) in exactly the same form as the pseudo-spectrum of the weight MUSIC algorithm. More specifically, the weighting matrix in the weighted MUSIC algorithm is inserted between Ew and EH w , while the weighting matrix is placed between F(φ) and F(φ) H in the proposed MUSIC algorithm. Furthermore, the criteria for the selection of the weighting matrices are fundamentally different. In the standard weighted MUSIC algorithm Krim and Viberg (1996), the weighting matrix is computed from the eigenvalues and eigenvectors of the sample covariance matrix Σˆ y and is constant. By contrast the weighting matrix in (7.22) is explicitly computed as a function of F˜(φ) and as a consequence, depends on the parameter to be estimated. The pseudo-spectrum (7.21) looks similar to the objective functions maximized in the pseudo-subspace fitting (PSF) method (Bengtsson 1999, Subsection 4.5.1). However, an essential difference between this method and the proposed MUSIC algorithm is that the latter computes the NA estimates by “scanning” a measure of the distance between a multi-dimensional subspace (induced by a single SDS in our case) and the estimated signal subspace, while in the PSF method, the NA estimates are the values providing the best “fit” between the estimated signal subspace and the subspace spanned by all signals. Thus, a one-dimensional search is required in the proposed MUSIC algorithm, while the PSF method requires a multi-dimensional search. Only in a single-SDS scenario is the objective function maximized in the PSF method identical to the pseudo-spectrum (7.21) calculated in the proposed MUSIC algorithm. In Christensen et al. (2004) another extension of the standard MUSIC algorithm is proposed which considers the projection of the columns of F(φ) on to Ew. This method relies on the pseudo-spectrum kF(φ)k 2 F kF(φ)HEwk 2 F . (7.24) The inverse of the pseudo-spectrum (7.24) corresponds to the distance between the multi-dimensional subspace spanned by the columns of F(φ) and the signal subspace, if and only if, the columns of F(φ) are orthonormal for any values of φ. This condition is usually not satisfied in real applications. Simulation results also show that the NA estimator derived from (7.21) outperform the estimator obtained from (7.24) in terms of lower root mean square estimation error. To our best knowledge, the proposed extension of the MUSIC algorithm according to (7.21) has not been reported in any published work yet. This algorithm is indeed the natural extension of the standard MUSIC algorithm to the case where the subspace induced by each individual signal component is multi-dimensional. One application example is the IC case in the SDS scenario considered in this contribution. Another example is fundamental frequency estimation for signals with a harmonic structure Christensen et al. (2004). AS Estimator Identity (7.12) inspires the following estimator of the AS of SDS: σcφ˜ = q σc2 β σc2 α. (7.25) The estimates σc2 β and σc2 α can be directly obtained for each of the D SDSs from (7.18) when using the SML estimators, or computed as σc2 β = 1 N XtN t=t1 |βˆ(t)− | 2 and σc2 α = 1 N XtN t=t1 |αˆ(t)− | 2 (7.26) when DML estimation is used. In (7.26), βˆ(t) and αˆ(t), t = t1, . . . , tN are calculated from (7.16) and denotes averaging. In the case where the proposed MUSIC algorithm (7.21) is applied, σc2 β and σc2 α can be obtained by applying the least-square covariance matrix fitting method (Johnson and Dudgeon n.d., Section 7.1.2) that we shortly describe below. First we rewrite the covariance matrix ΣyGAM in (7.19) according to vec(ΣyGAM) = D(φ¯)e (7.27)

184 Statistical channel parameter estimation with vec(·) denoting the vectorization Minka (2000) D(φ¯) .= [c(φ¯ 1) ⊗ c(φ¯ 1) ∗ , c ′ (φ¯ 1) ⊗ c ′ (φ¯ 1) ∗ , . . . , c(φ¯D) ⊗ c(φ¯D) ∗ , c ′ (φ¯D) ⊗ c ′ (φ¯D) ∗ , vec(IM )], where ⊗ is the Kronecker product, and e .= [σ 2 α1 , σ2 β1 , σ2 α2 , σ2 β2 , . . . , σ2 w] T. In the covariance matrix fitting method, the estimate Ωˆ of Ω in (7.17) minimizes the Euclidean distance between Σˆ y and ΣyGAM. Thus, the identity ∂kΣˆ y − ΣyGAMk 2 F ∂eH e=eˆ = 0 (7.28) holds for Ω = Ωˆ . Solving (7.28) yields the close-form expression for the solution eˆ: eˆ = D(φ ˆ¯) † vec(Σˆ y). (7.29) It is worth mentioning that the AS estimator (7.25) does not require knowledge of the pdf of the azimuth deviation. In the case where some assumption is made on the pdf in form of a parametric model, the AS estimate can be used to calculate the model parameters. In (Shahbazpanahi 2004, (39)–(43)), the AS σφ˜ is related to the parameters controlling the spread of the truncated Gaussian, the Laplacian and the confined uniform distributions. In addition, when the von-Mises pdf Ribeiro et al. (2005) is used, the relation between the AS and the concentration parameter κ of this pdf is approximated according to σφ˜ ≈ p 1 − |I1(κ)/I0(κ)| 2. In the sequel, the notation “AS(F) estimator” is used to denote the AS estimator calculated using (7.25), where σd2 βd and σd2 αd are computed from the estimates obtained using method “F”. More specifically, F can be “DML”, “SML” or “MUSIC”. 7.1.5 A New Definition of SDS and the Array Size Adaptation Technique In this section, a new definition of SDS is first provided, which makes use of the ratio of the largest eigenvalue by the second largest eigenvalue of the covariance matrix of the individual signal components contributed by scatterers. Then, based on this ratio we introduce a measure that quantitatively assesses the degree of the fit between the effective signal model (7.1) and the 1 st-order GAM model (7.9). Finally, we present a technique, called array size adaptation (ASA), which ensures a good fit between the two models by appropriately selecting the array size. The array size here is referred to as the size of the array aperture. Except when explicitly mentioned the single-SDS scenario is considered in this section. The Traditional Definition and a New Definition of SDS In this subsection, we present two definitions of SDS, the traditional definition relying on the effective rank of the signal subspace and a new definition based on the ratio of the largest signal eigenvalue by the second largest signal eigenvalue. We then discuss some issues arising when we apply the criteria induced by these definitions to decide based on measurement data whether a distributed scatterer is an SDS or not. The discussion reveals the advantages of the new definition. We first briefly review the conventional definition of SDS. In a scenario with a single distributed scatterer, the signal subspace is spanned by the L vectors in the sum in (7.1). For L > M, the signal covariance matrix has rank M with probability one. However, for small to moderate values of the AS the signal energy is concentrated in a few eigenvalues of the signal covariance matrix Meng et al. (1996); Shahbazpanahi et al. (2001); Xu (2003). The effective dimension of the signal subspace is determined by the so-called effective rank of the covariance matrix, which is defined to be the number of eigenvalues larger than twice the noise spectral height Zatman (1998). As shown in Xu (2003), the effective rank of the signal subspace induced by a distributed scatterer increases along with its AS. An SDS is a distributed scatterer which results in a signal subspace with small effective rank Bengtsson and Ottersten (2000); Bengtsson and Völcker (2001); Meng et al. (1996); Shahbazpanahi (2004); Trump and Ottersten (1996). The new definition of SDS is based on the following experimental evidence. For distributed scatterers with small AS the largest eigenvalue of the signal covariance matrix significantly dominates the other eigenvalues. A quantitative