1. Introduction

Information and communication technology (ICT) is rapidly developed. It has been widely discussed that the traditional all-electronic solutions may be insufficient for the ever stringent I/O interconnection requirement. Inspired by the success of the integrated circuits (IC) industry and the practical demand for faster interconnection, photonic IC (PIC) technologies are actively studied and developed in the last decades. With the growing momentum of the PIC industry, tremendous trials/tests of applying PIC to various applications have also been extensively performed.

Free-space optics is crucial for applications like optical wireless communication (OWC), light detection and ranging (LiDAR), and image scanning. For their advantages of programmability and flexibility, optical phased arrays (OPAs) become a popular option for manipulating the free-space optical beams. A simple OPA architecture built on silicon photonics is introduced in 2009 [1]. Since then, many PIC based OPAs have been implemented in various platforms, and their performance and scales are quickly improved [2–10]. More detailed illustration of the background and comparisons for different state-of-the-art works can be found in [5,9,11]. Due to the potential fabrication imperfection, the optical signals passing through different waveguides/paths typically result in random phase errors. These phase errors can seriously deviate the interference pattern of an OPA, and radically change the spatial distribution of the optical field. The resultant beams are typically broadened and have poor directivity; hence greatly deteriorate the scanning and transmission performance. The problems are commonly solved by applying thermally or electrically tunable phase shifters to calibrate the phase errors after the fabrication of the OPAs [2–11].

A feedback sub-system is typically applied for the calibration processes. Specific beam parameters are monitored and the phase shifters are tuned according to the corresponding responses until the best beam performance is detected. Pre-calibration can be performed for the factory OPAs by using the images detected by remote cameras as the feedback signals. Real-time calibration can be performed by monitoring the feedback signals of additional on-chip detectors as Ref. [3,12]. As the solution set of all possible phase states exponentially increases with the number of the OPA’s elements [12], it is surely inefficient to perform the calibration through nested sweeps over all the possible phase states. Many optimization algorithms have been proposed to reduce the time consumption of the calibration processes. Hill climber [2,4] and gradient descent (GD) [8] are two common algorithms for the optimization problems. However, as the phase error calibration is not a convex problem, local extreme uncertainty may exist. While the number of antenna elements increases, extra processing time for sequentially querying the responses from different antenna elements (channels) may also be significant. Stochastic parallel gradient descent (SPGD) is a variant algorithm which can get the stochastic gradient though a single-step calculation [12]. Hence, the time expended on the sequential querying processes can be saved. Based on SPGD, a two-stage algorithm named deterministic-stochastic gradient descent (DSGD) is proposed in [12]. DSGD first coarsely adjusts the OPA with a smaller possible solution set. This tunes the OPA to be at phase states close to the optimal solution. Then SPGD is applied to derive the final optimal solution. According to the published results, DSGD can efficiently mitigate the local extreme uncertainty problems, and get accurate results of the optimal points. Particle swarm optimization (PSO) [3] and genetic algorithm (GA) [10] have also been common options for the phase error calibration. These algorithms are more resistant to local extreme traps. In [10] and [12], it is shown that these two algorithms typically can more quickly and accurately get the optimal solutions compared to the GD and hill-climber based options. Though these algorithms have shown tremendous efficiency improvement compared to the nested sweeps, the increasing computation time for a large OPA can still be a cost concern. On the other hand, as the OPAs may work at different environment conditions, calibration needs to be constantly performed to dynamically adjust its best driving parameters. Hence, it is highly desirable to further push the efficiency of the algorithms. In [13], a modified rotating element electric field vector (mREV) method is proposed. This algorithm shows great potential in enhancing the computation efficiency. However, the in-phase condition cannot be exactly guaranteed (though could eventually reach after sufficient iterations in a theoretical viewpoint of the gradually reduced average errors). Moreover, according to the histograms of [13], the final phase deviation as high as >30° can still exist after the calibration processes.

In this paper, novel and efficient radix-p calibration processing is proposed. The proposed scheme can be directly combined with most of the already reported/developed optimization algorithms as a speed-up turbo. It is shown that the computation time of the proposed radix-p processing grows linearly with the number of the OPA’s elements. In comparison to the conventional particle swarm optimization (PSO) algorithm exhibiting exponential computation time (according to the numerical results), the radix-p processing shows great efficiency improvement. The basic concept of the proposed radix-p processing is to transform the optimization problems for a large antenna array into much simplified optimization problems for smaller sub-arrays composed of less antenna elements. Compared to other schemes (e.g. mREV) building on similar concepts, the results derived by the proposed radix-p method theoretically match with the exact optimal point. Slight final errors are observed, which may mainly come from the numerical limitation. The errors majorly influence the field pattern at regions several tens of dB lower than the main peak. Considering the potential noise levels, this may only cause negligible degradation of the systems. In the numerical results, ∼13.256 dB side-mode suppression ratio (SMSR) is obtained for the uniformly spaced OPA and more than 12 times efficiency improvement is shown for the proposed radix-p processing compared to the traditional PSO algorithms for an OPA with 64 elements.

2. Algorithm

An imaging system composed of an antenna array with N elements is considered as shown in Fig. 1. The antenna array is located on the object plane, which is described by the coordinate system of ${\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\rightharpoonup$}} \over r} _{obj}}$ = (x’, y’). The emission field of the n-th antenna element is denoted as ${\Phi _n}\left( {{{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\rightharpoonup$}} \over r} }_{obj}}} \right){e^{i \cdot {\phi _n}}}$. Considering general cases allowing arbitrary emission patterns, ${\Phi _n}\left( {{{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\rightharpoonup$}} \over r} }_{obj}}} \right)$ can be a complex number. However, ${\phi _n}$ is purposely separated from ${\Phi _n}\left( {{{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\rightharpoonup$}} \over r} }_{obj}}} \right)$ to describe the phase errors. ${\phi _n}$ is the variable that will be adjusted for the calibration or beam-steering processes. The emitted beam is monitored by a virtual image plane described by the coordinate system of ${\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\rightharpoonup$}} \over r} _{img}}$ = (x, y). The image plane is parallel to the object plane with a separation distance of L. According to the Fraunhofer approximation, the observed image pattern can be expressed as:

 

Fig. 1. The imaging coordinate system.

Download Full Size | PDF

In our proposed radix-p algorithm, multi-layer processing is performed. In the initial processing layer, the N antenna elements are grouped into M sub-sets. Each sub-set comprises p distinguished antenna elements and is processed sequentially according to the index number of the sub-sets. For better illustrating the following processes, p=2 is chosen as an example. In order to distinguish the processing terms, the 1-st sub-set is separated from Eq. (2) as:

Figure 2 summarizes the whole processing of the proposed algorithm for p=2. In the beginning of the processing (1-st layer), the elements of the OPA have phase errors of ${\phi _1}$∼${\phi _N}$. Each two of the antenna elements are grouped and form M = N/2 sub-sets (each is circled by the red line in the left side of Fig. 2). First, the 1-st sub-set is selected, and their phases are optimized to get the maximum intensity for the interested point on the image plane. After the processing, both antenna elements’ phases are adjusted to ${\phi _{s,1}}$ . This ${\phi _{s,1}}$ is identical to the phase value of the overall field composed of the other N-2 antenna elements (as shown in the leftmost skin-colored block). Next, the 2-nd sub-set is selected and processed by the same processing. The final phases of the 2-nd sub-set will be ${\phi _{s,2}}$. Identical processing is sequentially and repeatedly performed on all of the M sub-sets. At the end of the 1-st layer processing, the phases of the M sub-sets become ${\phi _{s,1}}$∼${\phi _{s,M}}$. It is seen that both antenna elements of the same sub-set have the same phases. These elements are then synchronized in values for every succeeding processing. Hence, only half number of phases are needed to be adjusted for the next processing layer. The processing continues layer-by-layer until there remains a single phase ${\phi _v}$. Then all phases of the antenna elements are perfectly in-phase, and the intensity on the point of interest is guaranteed to be maximum.

 

Fig. 2. Processing flows of the proposed radix-p algorithm. Opt.: optimization.

Download Full Size | PDF

Once all antenna elements are adjusted to be in-phase, further procedures may be performed to estimate the relations between the phase states and the electrical control signals of each antenna element. The antenna elements are sequentially selected to be tuned with different electrical signal levels. At the tuning moments for a specific antenna element, all other elements else are maintained intact. Hence, the relations between the feedback signals and the electrical control signals can be detected. Then, according to the interference physics of light, the relations between the phase states and the control signals can be estimated. To steer the beams towards specific directions, the required phase states of each antenna element can be calculated according to diffraction optics principles. Then the electrical control signals are tuned to the corresponding levels for beam-steering.

It is noted that the free-space propagation delays of waves from different antenna elements are included in the overall phases of the far fields as shown in Eq. (1) and (2). Hence, the above beam-steering processing implicitly presumes the phase calibration is performed with the feedback signal generated at the central (0-degree) direction. Practically, it is difficult to accurately identify whether the feedback signals are from the central direction or not. Hence, phase calibration processing may be performed at multiple beam angles. The relative beam angles can be measured and used to estimate the additional phases generated from the propagation delays. Then, the actual on-chip phase errors can be separated from the overall phases estimated by the radix-p processing. It is also noted that the geometric positions of the antenna elements are not restricted in Eq. (1) and (2). Hence, the radix-p processing can also be applied to calibrate the phase errors of aperiodic OPAs. It can be seen that the propagation optical path lengths, d_i*sin(θ) for the i-th antenna element, become 0 at the central direction with arbitrary values of d_i. Therefore, the in-phase condition is maintained for any geometric arrangement of antennas. However, while the calibration feedback signals are detected off the central direction, the degree of freedom introduced by the geometric factors can result in further phase uncertainty. Hence, the geometric errors must be guaranteed within specific levels to ensure the final on-chip phase errors are correctly estimated and calibrated.

3. Performance

Simulation and numerical calculation are performed to estimate the performance of the proposed algorithm. The CPU equipped in the computer is i5-4590. A 2-D imaging system like Fig. 1 is numerically built in MATLAB. Wavelength of 1.55 μm is used in the simulation. The size of each antenna element is 0.75 μm and the neighbor elements are separated by 9.5 μm. The beam is detected 10 cm away from the antenna array. Initial phases between [0, 2π) are randomly allocated to the antenna elements, as an emulation of the phase errors. The performance purely relevant to the proposed algorithm is characterized. PSO is selected as the benchmark to separate the potential local extreme issues. The maximal iteration of every PSO processing is set as 500 to exclude the ambiguity between the convergence performance and the execution time. Techniques proposed in [14,15] are added to the PSO setup to enhance the convergence capability at the final iterations. The number of antenna elements and the particles are varied according to the analysis conditions. The proposed scheme is applied as the radix-p (p ∈ {2, 8, 16}) version of PSO. The performance between algorithms of the traditional PSO and our proposed radix-p scheme is analyzed.

3.1 Computation time

The PSO may suffer from the local extreme issues while the number of particles is insufficient. In order to ensure the global extreme results are correctly derived, multiple runs applying different numbers of particles are performed. For each run, the derived result is compared with the ideal solution, and the required numbers of particles for different OPA sizes are characterized as Fig. 3(a). It is observed that required particles exponentially increased with the sizes of the OPA. We believe this comes from the fact that the solution space grows exponentially with the number of antenna elements [12]. Though it is not required for the PSO to densely sweep the whole solution space, Fig. 3(a) may reveal that the particle number must be kept at specific levels to ensure sufficient probabilities for the particles to hit the optimal region. Hence, the required number of the particles increases with the solution space (antenna elements). As a comparison, a particle number as high as 200 is applied in [12]. Figure 3(b) shows that the calculation time approximately grows linearly with the number of the applied particles. We believe this mainly comes from the variable transferring latency, and could be much reduced for practical systems. The corresponding execution time of the PSO optimization processes is shown in Fig. 3(c), which also grows exponentially and meets the results of Fig. 3(a) and (b). Robust curve fitting is performed for the results shown in Fig. 3 through the toolbox of MATLAB. The fitting results of Fig. 3(b) is directly labeled on the picture. The regression results of Fig. 3(a) and (c) are discussed separately in Table 1. It is seen that two different types of curve models are applied in Table 1. Type I is the generic form of the exponential model. However, it is observed that coefficient “a” between Fig. 3(a) and (c) are similar (3.319 and 3.292). This obviously does not obey the linear relation observed in Fig. 3(b), and could influence the accuracy for prediction of the required particles/execution time for OPAs with larger sizes. Observable variations for the fitting results are also seen by changing the initial conditions. We believe the errors mainly come from the sensitive features of the exponential model. In order to reduce the errors, fitting model Type II is applied by fixing coefficient “b” of the exponential model for both curves. Different “b” is locally swept according to the results of Type 1, and the one with the highest correlation is chosen to be 0.09234. It is observed that the division between coefficient “a” of the two curves is 5.539/2.075 = 2.669, which is closer to the slope of Fig. 3(b). We believe type II model more accurately describes the results. According to the modified curve, the execution time for low OPA size is negative (−7.601), which is surely impractical. By regionally enlarging Fig. 3(c) into Fig. 3(d), it is seen the fitted curve deviates significantly with the data while the OPA has few elements. We believe at this region, the computation time is not dominated by only the PSO processing. Hence the practical curve behavior should be modified as the pink curve. The main reasons behind may mainly come from the basic computer architecture, including both the hardware and software configurations, which dominate the computation time at less data amount. According to the modified curve, the execution time for OPA size of 128 and 512 will respectively be 209 and 5.25 × 10¹⁷ hours. This is surely impractical, and will cast much challenges on the calibration processing. Though in practical systems the expended time may be much reduced (because no needs of transferring and processing a large amount of samples as the numerical emulation does), a faster algorithm is still desirable. The coordinates of Fig. 3(c) are extended and shown as Fig. 3(e), which clearly present the curve behaviors of the radix-p methods at larger size of antenna arrays. From Fig. 3(e), it is seen through the proposed radix-p (p ∈ {2, 8, 16}) algorithm, the computation time is much reduced, and the expended time grows linearly (rather than exponentially) with the size of the OPAs. The linear behaviors hold even while the number of antenna elements is as large as 4096, and meet the theoretical prediction. This shows the great performance improvement of the proposed algorithm.

 

Fig. 3. Performance of the computation results. (a) Required particles for getting the optimal solution. (b) Expended time for different number of particles. (c) Expended time for OPAs of different sizes. (d) Regionally enlarged picture of (c) for the PSO curves. (e) Extended picture of (c) for the radix-p curves. (f) Calibrated errors.

Download Full Size | PDF

Table 1. Coefficients of the Exponential Curves of Fig. 3

View Table | View all tables in this article

Table 2. Table 2 is used to show the regression curves of Fig. 3(c) and (e)

View Table | View all tables in this article

In the proposed algorithm, p-element optimization is performed for every sub-set processing. The time expended by the p-element processing is denoted as O_p-element. For the whole calibration processing, there are log_p(N) layers; and for the r-th layer, there is N/p^r sub-set processing. The total expended time is then a geometric series with log_p(N) terms, and can be derived as:

A residual (Res.) term is added in Eq. (9) to cover general cases, for which log_p(N) may be a fractional number. While log_p(N) is not an integer, the sub-set will be selected at its best possible efforts, and the remaining < p terms will be calculated independently. As the number is less than p, this will not cause much additional computation time. In order to more clearly observe the behavior of the computation time, a factoring function, K(N), is defined from Eq, (9) as:

Figure 4(a) shows the values of K for different N and p. It is observed that while N becomes larger than 4, the ascending rate of K significantly decreases. After N>16, K grows even slower and becomes almost constant (K_c) for all values of p. Hence, Eq. (9) can be approximated as:

Equation (11) shows the computation time scales linearly with the number of N. This shows great performance in computation time of the proposed algorithm for the cases of larger N, which meets the results of Fig. 3(c).

 

Fig. 4. (a) The factoring function at different p and OPA sizes. (b) the estimated time expense for different p.

Download Full Size | PDF

It is noted that O_p-element typically increases with the number of p as shown in Fig. 3(c) and (d). However, as shown in Fig. 4(a), K_c increases with the reduction of p due to the increasing layers and sub-sets at low-p cases. Hence, an optimized p for best reducing the total computation time may commonly be larger than 2. The optimized p is estimated according to the results of Fig. 3(d) and 4(a). Figure 4(b) shows the estimated time expense for N=64 at different p. It is seen that the time expense can be further reduced while p slightly increases as high as up to 11. After p>11, the time expense of the PSO dominates, and slow down the overall processing speed. Figure 3(e) shows the simulated time expense for different values of p with the array size up to 4096. Since the computation region must be extended to get finer image resolution for the narrower beam widths of larger array sizes, the computation data points of the radix-p methods are increased by 4 times. This results in about x4 more time expense, which is also observed between the results of Fig. 3(c) and Fig. 4(b). In Fig. 3(e), the lowest time expense is observed at the case of p = 8. For other cases, whether p becomes larger (p=16) or smaller (p=2), the time expense increases as predicted in Fig. 4(b). This suggests the possibilities of further time reduction with more optimized choices of p.

Figure 3(f) shows the final average calibrated phase errors for different values of p at different sizes of antenna arrays. For all values of p, higher phase errors are approximately observed while the array sizes increase. This is possibly caused by the numerical errors of each sub-set processing, which are accumulated with more processing stages at large array sizes. For the same array sizes, the final phase errors increase at higher p. This may mainly come from the performance variances of PSO at larger antenna sizes (higher p). As shown in Fig. 3(c), it is seen that the calibration results are sensitive to the number of applied particles. While the applied particles are around the values of the fitting curve, there are possibilities that the converged phase results are not fully optimized. As the array sizes increase, huge variations are observed for the optimized number of particles. Hence, while p increases, there are potentially more errors generated by PSO. This also explains for the fluctuations of phase errors for the p=8 case shown in Fig. 3(f). One insight perceived from the results of Fig. 3(f) is that: even if p>2 may result in less computation time as illustrated in the last paragraph, it will also introduce further phase errors. Hence, the optimized values of p should also be determined by the tolerable phase errors of the practical systems, and may deviate from the values around 11 shown in Fig. 4(b).

3.2 Calibration error

The calibration performance is analyzed for an OPA with 64 elements. Additional noise with normal distribution is added onto the phase term of each antenna element to emulate the possible variations of the control signals during the optimization processing. The standard deviation of the noise is swept from 0∼2π. The original and the calibrated patterns are respectively shown in Fig. 5(a) and (b), which reveal the monitored region between ±1^st order grating lobes (two FSRs) [16]. The legends for the red and green curves respectively denote for the strength of the added noise in unit of radian. It can be seen that while the control signals do not exhibit additional noise (the red curve), the phase errors are almost correctly calibrated. Hence, the intensity of the main lobe is significantly distinguished, and the SMSR is about 13.256 dB. This is a level close to the ideal level of 13.269 dB. We believe the ∼0.013 dB penalty is an inherent limitation of the proposed radix-p algorithm. The phases of the antenna elements belonging to different sub-sets are not adjusted simultaneously. Numerical errors exist after each sub-set processing. The errors will propagate to the processing of the succeeding layers, and result in additional calibrated errors compared to the algorithms not applying the radix-p scheme. The errors can also be observed at the edges of Fig. 5(b), for which the red curve is not symmetrical for the lobes ∼ 20 dB lower than the main lobe.

 

Fig. 5. Calibration performance. (a) Field pattern with wide range. (b) zoomed in plot of (a). (c) Intensity traces of the lobes and SMSR. (d) Calibrated errors. Note: values in legends of (a) and (b) are the strengths of the added noise in radian. Ori: original.

Download Full Size | PDF

Figure 5(c) shows the strengths of the main lobe, the maximal side lobe, and the SMSR at different noise levels. As the patterns with poor phase adjustment can be random, the main lobe is defined as the value observed at the point of interest. It is seen that the intensity of the main lobe and SMSR decrease at higher noise levels of the control signals. While the noise level is higher than ∼2 radian, the intensities of the main lobe and the SMSR are at their minimum levels. Hence no benefits are offered by calibration. This phenomenon can be further observed in Fig. 5(b) for the case of 2π noise level. Ripples with no significant lobes are observed as expected.

Figure 5(d) shows the average absolute phase errors after the calibration processes. At all noise levels, error floors are observed. An error around 0.043 radian is observed even while the noise is not added. This again can be explained by the numerical errors and the resultant non-ideal SMSR level of 13.256 dB shown in Fig. 5(a)-(c). The error floors increase with the strengths of the noise, and saturate at levels around 1.4255 as labeled by the green line and text. As the maximum value in Fig. 5(d) is 1.5538 (∼π/2), the saturated error levels may result from an un-correlated quadrature condition. It is seen that the saturation point is at the noise level of around 2 radian, which corresponds to the point of the SMSR floor shown in Fig. 5(c). Before the saturation point, the final average phase errors grow linearly with the noise levels of the control signals. The fitting curve is directly shown and labeled on Fig. 5(d). These values indicate the potential strategies to choose different control-signal generators. Assume the control signals are generated through ideal digital-to-analog converters (DACs), which can be adjusted exactly within the 2π range. Then, according to Fig. 5(d), specific calibration benefits can be gained with bit resolution higher than 2. And considering the error floor of 0.043 radian at 0 control noise, DACs with bit resolution around 8 bits may be sufficient to offer performance close to the best available levels.

4. Conclusions

An efficient radix-p algorithm is proposed and analyzed. The proposed radix-p algorithm can be combined with most of the existing optimization algorithms as a speedup booster. Theoretical analyses are performed and show that the computation time of the radix-p algorithm grows linearly with the size of the OPAs. This offers great computation improvement, in comparison to the conventional particle swarm optimization (PSO) algorithm exhibiting exponential computation time shown in our numerical results. The processing time and accuracy are analyzed through numerical methods. In the simulation environment, the proposed radix-p processing shows more than 12 times efficiency improvement than the traditional PSO algorithm. The radix-p algorithm may suffer from additional numerical errors, which results in about 0.013 dB SMSR degradation. Considering the existence of the numerical errors, the bit resolution with ∼8 may reach performance close to the algorithm’s limitation.