1. Introduction

Owing to advancements in deep space exploration, photon-counting communication systems based on pulse position modulation (PPM) and single-photon detectors have been widely applied to deep space optical communication to achieve high photon efficiency [1–3]. In such systems, the arrival time of the optical pulse is detected and demodulated at the receiver to recover the information from transmitter. However, because demodulation depends on the arrival time of photons, the PPM slot clock discrepancy between communication terminals, including frequency and phase offsets, introduces symbol errors [4]. Therefore, clock synchronization is crucial in photon-counting PPM communication systems.

Various slot and frame synchronization schemes are applicable to different scenarios of optical PPM communication systems. Previously, a slot edge detection synchronization scheme based on integral circuits was proposed in [5], and a few schemes using the special PPM acquisition sequences for timing synchronization were studied [6–8]. In [8], under the data link up of 622 Mbps at the ground optical receiver of the Lunar Laser Communications Demonstration (LLCD) project, the combination of synchronization sequence detection based on autocorrelation values and early-late gate was used to verify the slot and frame synchronization performance. However, the synchronization performance of such schemes depends on the insertion period and length of the acquisition sequence, which occupies data transmission bandwidth. In [9], the two synchronization schemes based on pilot and inter-symbol guard slots were compared. In [10], a maximum a posteriori (MAP) estimator and particle filter were used to estimate random phase and frequency jitters for slot synchronization. Studies [11,12] adopted a slot synchronization method that used several gated clocks at the receiver with the same frequency as that of the transmitter and different phase delays. Thereafter, the clock with the maximum pulse counts inside the gate was selected as the receiver clock. However, the precision was determined using the number of phase divisions requiring multiple sampling clocks, which hinders system implementation. In [13], a real-time slot clock recovery based on oversampling and pulse matching filter was introduced. However, the synchronization accuracy was limited by the oversampling factor, and the sampling clock frequency restricted the communication transmission bandwidth.

Another feasible synchronization technology in PPM communication is based on the maximum likelihood estimator (MLE) [14,15]. R. Rogalin et al. utilized an inter-symbol guard time (ISGT) containing the timing information of symbol boundaries and MLE for excellent performance and low-complexity timing estimation scheme [16,17], which will be applied to the synchronization module in the optical receiver of NASA’s Deep Space Optical Communications (DSOC) project, scheduled for September 2023 [18]. Similarly, [19] applied MLE to the 2${\times} $1 multiple-in single-out (MISO) optical systems, which use the sampling frequency of the transmitter clock without requiring explicit slot synchronization sequences, thereby maximizing the data bandwidth resources.

However, the MLE model in [16–19] is based on the assumption of blocking mitigation [20] by a large scale detector array (e.g. a 64-pixel tungsten silicide SNSPD) [21], wherein the slot counts obey the ideal Poisson statistics well, whose mean varies linearly with the slot offset. Since the fabrication of large-scale array is currently facing technological challenges such as low multi-pixel yield, difficulty in controlling the uniformity of thin films and nanowires, and the complexity of parallel signal readout electronics, medium-scale arrays (8-elements, 16-elements) [22] or single-pixel detectors in deep space optical communications still have advantages in communication applications with different data rates in terms of stability and cost.

In order to achieve high data rate and slot signal-noise ratio (SNR), the slot width is reduced to the sub-dead time scale and the photon flux is increased. Furthermore, the blocking effects become significant [22,23] and the nonlinearity caused by blocking in the estimation model becomes more apparent, resulting in performance degradation of the original algorithm. Therefore, the present study extends the application of maximum likelihood synchronization algorithm to the sub-deadtime slot scenario with blocking loss and introduces the corresponding corrections to the model.

The main contributions are summarized as follows:

The remainder of the paper is organized as follows. Section 2 describes the model, wherein subsection 2.1 describes the MLE with deadtime correction, subsection 2.2 provides the CRLB, and subsection 2.3 introduces the slot frequency offset compensation. Section 3 introduces the simulation performance of the algorithm and the experiment based on a single-photon detector, including the experimental design and result analysis. Simulation based on detector arrays is also performed in this section. Finally, the conclusions are presented in Section 4.

2. Model

2.1 Estimation algorithm with deadtime correction

Essentially, the scenarios, wherein the proposed algorithm is applicable, are discussed. A modulated symbol consists of M and P number of PPM and ISGT slots, respectively, which contain no optical pulses. At the receiver, a same-frequency slot clock is used to divide the photon arrivals into slot bins. However, owing to transmission delay, crystal oscillator frequency drift, and Doppler frequency shift, a timing offset $\tau $ (slot phase offset) between the TX and RX symbol boundaries may exist. Without timing offset compensation, photons falling into adjacent slots cause erroneous demodulation and increased bit-error rate (BER).

Furthermore, for scenarios with severe blockage, the slot width ${T_s}$ is shorter than the detector deadtime ${\tau _\textrm{d}}$. For non-photon-number-resolution detectors, only one photon per slot is detected, and the duration of optical pulse is equal to the slot width when high-speed electro-optical modulation is applied. The measured ${\tau _\textrm{d}}$ was 28.5 ns, and to ensure the complete recovery of detector at the arrival of each pulsed slot, the parameters were set as M = 8, P = 4, ${T_s}$=10 ns. However, the inter-symbol interference caused by the detector’s deadtime may reduce the response probability of pulsed slots.

The impact of detector deadtime on estimation model is discussed as follows. For ${T_s} \le {\tau _\textrm{d}}$, the photon arrival distribution within the sub-deadtime slot $\textrm{f (t)}$ is attributed to the first photon detected in the slot. Thus, the first photon detection can be regarded as a renewal process [23].

For the pulsed slot with uniform intensity and duration of ${T_s}$, the immediate failure probability at time $\textrm{t}$, $\phi (t) = \lambda ,\textrm{ }t \in [0,{T_s}]$. Here, for ease of understanding, we can treat $\lambda $ as the Poisson intensity. And $\textrm{f(t)}$ is the probability density function of the survival time, given as follows:

Eq. (1) provides the theoretical distribution $\textrm{f (t)}$ of photon arrivals within the slot. Fig. 1(a) shows $\textrm{f (t)}$ with different Poisson intensities and histograms of simulated photon arrivals for comparison. Fig. 1(b) shows the time-digital converter (TDC) histograms of measured photon arrivals detected using a single photon avalanche diode (SPAD).

 

Fig. 1. (a) Normalized arrival time distribution of detected photons in the slot under different Poisson intensity $\lambda $.$f(t)$ and simulation results are compared. (b) TDC histogram of SPAD-detected photon timestamps.${K_\textrm{s}} = \lambda {T_s}$. The distribution degenerates to uniform as $\lambda \to 0$.

Download Full Size | PDF

Owing to the detector’s deadtime, the photon arrivals within a slot do not obey the uniform distribution mentioned in [16,17] when the Poisson intensity increases. Additionally, the average counts in each slot within a symbol are related to the offset between the symbol boundaries of TX/RX. Let the symbol boundary offset that RX succeeds TX be $\tau = k + \varepsilon $, where k and $\varepsilon $ are the integer and fractional parts of $\tau $, respectively. The counts of each slot within N symbols $x[n],n = 0,1, \ldots N(M + P) - 1$ are cyclically accumulated into the aggregated slot vector $\boldsymbol{y} = [{y_0},{y_1}, \cdots ,{y_{M + P - 1}}]$, wherein${y_m} = \sum\nolimits_{i = 0}^{N - 1} {x[i(M + P) + m]} ,m = 0,1, \ldots M + P - 1.$This vector provides sufficient statistic samples for estimating each slot count. As shown in Fig. 2, the percentage leakage of average photon counts from the last TX PPM slot to the first RX guard slot, $\varepsilon ^{\prime}({K_s},\varepsilon )$, is non-linear with fractional offset $\varepsilon $. Furthermore, it also depends on the average number of signal photons in the signal slots ${K_\textrm{s}}$ and background photons in each slot ${K_\textrm{b}}$.

 

Fig. 2. Mean photon numbers contained in each slot of a symbol at RX when $\tau $ and ${\tau _d}$ exist. The red slots represent the RX slots, wherein the photon counts are affected by non-integer timing offset $\varepsilon $, whereas the blue slots remain unaffected. PPM slots and ISGT are represented by semi-transparent and shaded rectangles, respectively.

Download Full Size | PDF

After correction, the average counts in each slot of a symbol, namely ${\lambda _m}^{\prime}(\tau )$, is derived as Eq. (3) with the approximation $({(1 - \varepsilon^{\prime}) + \varepsilon } )N{K_b} \approx N{K_b}$, which is reasonable for high slot SNR requirement (${K_\textrm{s}} \gg {K_b}$) essential for obtaining reliable BER. Therefore, the approximated noise term has negligible influence on estimation performance. Furthermore, only the signal photon term includes the offset. Fig. 2 depicts the relationship between ${\lambda _m}^{\prime}(\tau )$ and fractional offset $\varepsilon $.

The probability mass function of each slot count ${y_m}$ within ${y}$ can be obtained using Eq. (3), and the count process is modeled as Poisson statistics [24]

And ${\lambda _m}^{\prime}(\tau )$ could be considered as the sum of signal term and noise term, namely:

Then, ${K_{s,m}}^{\prime}(\tau )$ in the rest (m-1) slots can be written as cyclic shifts of ${K_{s,0}}^{\prime}(\tau )$:

The log-likelihood function $\ell ^{\prime}(\tau ;y)$ of ${P_{{Y_m}|\tau ({y_m}|\tau )}}$ given an offset $\tau $ is

Because the sum of the last two terms in $\ell ^{\prime}(\tau ;y)$ is independent of $\tau $, it can be discarded as a constant to compare the likelihood functions corresponding to M + P slots. When $\tau $ equals to the actual symbol boundary offset, the maximum likelihood estimation is obtained using Eq. (8). Therefore, the actual offset can be obtained by maximizing the function $\ell ^{\prime}(\tau ;y)$ with respect to $\tau $. Because this function is not differentiable at integer values, the zeros of partial derivative equations of $\ell ^{\prime}(\tau ;y)$ are obtained with respect to $\tau $ over M + P integer intervals.$\tau \in ({j,j + 1} ),j = 0,1, \ldots M + P - 1.$From Eqs. (6)–(7), for each j, the first and second order partial derivatives of ${K_{s,m}}^{\prime}(\tau )$ with respect to $\tau $ are non-zero when $m = ({j + M} )\bmod ({M + P} )$ or $m = j\bmod ({M + P} )\textrm{. }$Thus, the local optimal estimate of $\tau $ on each interval $(j,j + 1)$ is:

Furthermore, adding the constraint for ${\hat{\tau }_{ML,(j,j + 1)}}$ to each integer interval is crucial

Finally, the local optimal offset estimation with the largest likelihood function value is selected as the overall timing offset estimation

Furthermore, a set of photon arrivals were simulated with blocking loss. For comparison, the log-likelihood function $\ell (\tau ;y)$ from the original algorithm in [16] and $\ell ^{\prime}(\tau ;y)$ in Eq. (8) are plotted in Fig. 3. The estimated offsets before and after correction are 6.4858 and 6.7059 slot, respectively. Hence, the original algorithm degrades when affected by blockage.

 

Fig. 3. Log-likelihood function before and after deadtime correction (denoted by blue and red, respectively). The modified likelihood function loses symmetry on each integer subinterval. The horizontal coordinates of the open circles are the local optimal estimates of timing offset in each subinterval. The estimated ${\hat{\tau }_{ML}}$ equals to 6.71 and 6.49 slot with and without correction, respectively. Furthermore, $M = 8.\textrm{ }P = 4.\textrm{ }{K_s} = 2.\textrm{ }{K_b} = 0.0001.\textrm{ }N = 5000.\textrm{ }\tau \textrm{ = 6}\textrm{.70 slot}\textrm{. }{\tau _d} = 3{T_s}.$

Download Full Size | PDF

2.2 Cramér-Rao lower bound

In estimation theory and statistics, CRLB represents the lower limit of the variance of an unbiased estimator for a deterministic but unknown parameter. The minimum value of the variance of such an estimator is equal to the reciprocal of Fisher information that represents the upper limit of an unbiased estimator’s accuracy [25,26]. Thus, with the given parameters, the optimal estimation performance of the offset estimator is provided using its CRLB. However, its function with respect to input parameters varies owing to the deadtime correction. For the estimator in Eq. (9), the new CRLB value is derived to compare the performances of MLE algorithms in subsequent sections.

First, the Fisher information $I(\tau )$ is calculated, which equals to the second partial derivative of $\ell ^{\prime}(\tau ;y)$ with respect to $\tau $

Because $\ell ^{\prime}(\tau ;y)$ is not differentiable at integer slot offsets, the interval of $\tau $ is $\tau \in (j,j + 1)$, and the aggregated photon counts of each slot in ${y}$ obey the Poisson statistics. Hence, the mean of ${y_m}$ is ${\lambda _m}(\tau )^{\prime}$, wherein $1 - \varepsilon = j + 1 - \tau $

Therefore, the CRLB of corrected MLE is

Herein, the CRLB descents with $O({N^{ - 1}})$; however, it exhibits a non-linear relationship with the flux parameter ${K_s}$ and fractional offset $\varepsilon $.

Second, the estimation of $\tau $ using Eq. (9) is unbiased by including deadtime correction. From Eq. (13), it can be inferred that

According to the corollary of Jensen's inequality [27], for convex exponential functions, if X represents a random variable, then

Therefore, let $X = \exp [{ - ({j + 1 - {{\hat{\tau }}_{ML,(j,j + 1)}}} )({{K_s} + {K_b}} )} ]$.

Using Eqs. (15)–(16), the following is obtained:

However, the mean of estimated offsets approaches the actual offset as ${K_s}$ and N increase, which is also demonstrated in the numerical results of Section 3.1.

2.3 Non-ideal frequency offset compensation

In actual communication scenarios, the non-ideal factors affecting the MLE performance include the slot clock frequency offset $\gamma $ between the transmitter and receiver. When the accumulated symbol number N is sufficiently large, the estimated timing offset $\tau $ is significantly affected by the small $\gamma $ owing to its short slot duration of several nanoseconds. The frequency offset $\gamma $ is attributed to the frequency error of oscillator and on-chip PLL (such as the mixed-mode clock manager (MMCM) of FPGA), which is related to the device itself and other factors, such as ambient temperature and Doppler frequency shift.

From Fig. 4(a)-(b), for the same offset $\tau $, the existence of $\gamma $ has a significant influence on the aggregated slot vector. Let the actual slot duration at the transmitter and receiver be ${T_t}$ and ${T_r}$, respectively. Assuming an initial symbol boundary offset of ${\tau _0}{T_t}$, under the continuous accumulation effect of $\gamma $, the offset at the start of the $N\textrm{ - th}$ symbol at the receiver transforms to

 

Fig. 4. $N = 5 \times {10^4},\textrm{ }{K_s} = 1,\textrm{ }{K_\textrm{b}}\textrm{ = 0}\textrm{.0001, }\tau = 6.70\textrm{ slot}.$ (a) ${{y}_{m}}$ with no slot frequency offset. (b) ${{y}_{m}}$ with a slot clock frequency offset $\gamma = 5 \times {10^{ - 6}}$ between the TX and RX terminals. (c)-(e) Frequency offset estimation using linear least squares (LLS) and maximum likelihood estimator (MLE). Red and blue triangles correspond to the results of MLE with and without deadtime correction, respectively. (c) corresponds to K instances of phase offset estimates without frequency compensation. (d) demonstrates re-estimated phase offsets after compensating the frequency offset once. While the influence of frequency offset is neglectable in (e) because the compensation of frequency offset has been performed twice.

Download Full Size | PDF

The accumulation time used to estimate $\gamma $ is short compared to the required time during which, the crystal oscillator temperature drift and Doppler shifts trigger slot frequency variation (despite using tens of ${\hat{\tau }_{ML}}$ to estimate $\gamma $, the time required is several milliseconds); hence, $\gamma $ is considered to be constant in each frequency offset estimation. Thereafter, the real ${\tau _\textrm{n}}$ and estimated offsets ${\hat{\tau }_{ML}}$ increase linearly with $N.$Therefore, slot photon counts are collected in $K \times N$ consecutive symbols to generate K consecutive estimated offsets ${\hat{\tau }_{ML,1}},\textrm{ }\ldots \textrm{ ,}{\hat{\tau }_{ML,K}}$ using the MLE algorithm. Finally, $\gamma $ is estimated by leveraging the linear relationship between ${\hat{\tau }_{ML}}$ and the cumulative number of slots.

In Fig. 4(c), the variation in ${\hat{\tau }_{ML}}$ was simulated with respect to the cumulative number of slots when ${\tau _0} = 6.7,\textrm{ }\gamma ={-} 5 \times {10^{ - 6}},\textrm{ }{\tau _d} = 3{T_s}$. Here, $K = 40,\textrm{ }N = 5 \times {10^4}\textrm{, }{K_s} = 1,\textrm{ }{K_b} = 0.0001$. Because $\gamma $ causes linear variation in estimated offsets, its effect on MLE before and after deadtime correction is approximately the same. Applying the LLS method to $K$ estimated offsets for slope calculation, the estimated slot frequency offset $\hat{\gamma }$ is obtained. Figures 4(d)-(e) illustrate K re-estimated offsets ${\hat{\tau }_{ML}}$ after compensating slot frequency offset once and twice, respectively. Evidently, the influence of $\gamma $ can be eliminated after compensation, and proper values of K and N can be selected based on the photon flux and synchronization accuracy requirement for practical applications.

3. Quantitative performance

3.1 Simulation results

Before characterizing the estimation performance, the CRLB indicating the optimal estimator's performance with given parameters, is discussed.

In Figs. 5(a)-(b), the surface with blue grid represents the original form of CRLB in [14] without considering deadtime effects. The average photon counts within the slot varies linearly with the timing offset $\tau $. Herein, CRLB is symmetric about the fractional offset $\varepsilon $ and decreases with $O({K_\textrm{s}}^{ - 1})$ and $O({N^{ - 1}})$. The corrected CRLB is plotted on the surface with red grid. Owing to the exponential suppression of deadtime on the photon arrival time distribution, CRLB decreases linearly with $O({N^{ - 1}})$ after correction. Furthermore, the fractional offsets $\varepsilon $ and $1 - \varepsilon $ are dissimilar, and CRLB attains a minima as ${K_\textrm{s}}$ increases.

 

Fig. 5. Distribution of CRLB^1/2 before (denoted by blue grid) and after (denoted by red grid) deadtime correction. (a) The variation of CRLB^1/2 with ${K_s}$ and $\varepsilon $. $N = 5000,\textrm{ }{K_b} = 0.0001.$(b) The variation of CRLB^1/2 with $N$ and $\varepsilon $. ${K_s} = 1,\textrm{ }{K_b} = 0.0001$

Download Full Size | PDF

When $\varepsilon \to 0$ and ${K_s}$∼10 photons/symbol, CRLB attains its minima and continues to increase, in contrast to the situation when $\varepsilon \to 1$, in which the minima of CRLB appears at larger ${K_s}$, indicating an asymmetric CRLB with respect to $\varepsilon $. However, in DSOC, the received photon flux achieving a reliable BER is generally of the order 0.1∼10 photons/symbol. Therefore, the theoretical optimal performance of the corrected MLE can attain the level of ${10^{ - 2}}$ slot at root mean square error (RMSE), despite $\varepsilon $ approaching zero.

Furthermore, the performances of MLE algorithms before and after deadtime correction under various $\varepsilon ,\; {K_s},\; and\; N,$ have been discussed. To verify the estimation performance quantitatively, serially concatenated PPM-coded (SCPPM) frames of random binary data with a total length of 15120 (length of the interleaver) were randomly sent, and the response of single photon detectors to the pulsed slots was simulated using the Monte-Carlo method. The method for generating photon arrival timestamps is similar to [28], and the detector efficiency recovery model was changed to the step function for SPAD. In the simulation, the relative parameters were set as follows: ${\tau _\textrm{d}} = 3{T_\textrm{s}},\textrm{ }M = 8,\textrm{ }P = 4$. To verify the performance improvement of MLE after deadtime correction, the slot frequency offset was set to zero.

The MLE performance is characterized by $\textrm{RMSE}({\hat{\tau }_{ML}}) = \sqrt {\mathrm{\mathbb{E}}[{{{({{{\hat{\tau }}_{ML}} - \tau } )}^2}} ]}$, RMSE of ${\hat{\tau }_{ML}}$, and $\textrm{bias}({{{\hat{\tau }}_{ML}}} )= \mathrm{\mathbb{E}}({{{\hat{\tau }}_{ML}}} )- \tau$, representing the bias between the mean of ${\hat{\tau }_{ML}}$ and actual offset. In Fig. 6, the corrected $\textrm{CRLB}{({{K_s},\varepsilon } )^{1/2}}$ values are plotted for comparison.

 

Fig. 6. (a) $\textrm{RMSE}({\hat{\tau }_{ML}})$vs. $({{K_s},\varepsilon } )$ before (marked with blue grid) and after (marked with red grid) deadtime correction. The black grid surface is $\textrm{CRLB}{({{K_s},\varepsilon } )^{\textrm{1/2}}}$ after deadtime correction. (b)$\textrm{bias}({{{\hat{\tau }}_{ML}}} )$vs. $({{K_s},\varepsilon } )$ before (marked with blue grid) and after (marked with red grid) deadtime correction. (c) $\textrm{RMSE}({\hat{\tau }_{ML}})$ vs. ${K_s}$ at several fixed $\varepsilon .$(d) $\textrm{bias}({{{\hat{\tau }}_{ML}}} )$vs.${K_s}$ at several fixed $\varepsilon .$

Download Full Size | PDF

Figures 6(a)-(b) show the variation in estimation performance with $\varepsilon ,\textrm{ }{K_s}$ when $N = 5000,\textrm{ }{K_b} = 0.0001$. For each pair of $({{K_s},\varepsilon } )$, two-thousand ${\hat{\tau }_{ML}}$ values were estimated to calculate the corresponding RMSE and bias. Figures 6(c)-(d) show the relationship between RMSE and bias at several fixed fractional offsets $\varepsilon $.

In Eq. (2), owing to reduction of blocking effect, the percentage of photon leakage into adjacent guard slots, $\varepsilon ^{\prime}({K_s},\varepsilon )$ degrades to the fractional offset $\varepsilon $ at low photon fluxes

Therefore, the RMSE and bias values of ${\hat{\tau }_{ML}}$ are similar to those before correction at low fluxes. As the flux increases, particularly in the typical intensity range of communication ${K_\textrm{s}} \in [{0.1,\textrm{ }10} ],$ the original MLE experiences significant performance degradation, whereas the RMSE and bias values obtained by the corrected MLE are better than those before correction by 1∼2 orders of magnitude. The performance difference is further enlarged as $\varepsilon \to 1$. Further increase in the intensity of flux (e.g., ${K_\textrm{s}} > 10$ photons/symbol) causes performance degradation of the corrected MLE as $\varepsilon \to 0$, which can be attributed to the non-monotonic variation in CRLB with ${K_s}$, as shown in Fig. 5(a). Figure 6(b) and 6(d), and Eq. (18) indicate that the corrected estimate does not satisfy the unbiased estimation, and the bias approaches zero as flux intensity increases. Contrastingly, insufficient photon flux is not conducive to ensuring frame synchronization and low error rate in communication. Furthermore, excessive photon flux is difficult to realize in large loss scenarios, such as deep space link. Thus, the proposed MLE algorithm exhibits better performance within the suitable photon flux range for communication.

Figure 7 shows the variations in estimation performance with $\varepsilon ,\textrm{ }N$ when ${K_b} = 0.0001$ and ${K_s} = 1.0$. Furthermore, to calculate the RMSE and bias, 2000 ${\hat{\tau }_{ML}}$ values were estimated for each pair of $({N,\varepsilon } )$. Figures 7(a)-(d) indicate that the estimation performances before and after correction improve as N increases, and RMSE after correction is approximately equal to CLRB^1/2. However, RMSE experienced saturation as N increased, which can be attributed to the following reasons.

 

Fig. 7. (a)$RMSE({\hat{\tau }_{ML}})$ vs. $({N,\varepsilon } )$ before (marked with blue grid) and after (marked with red grid) deadtime correction. The black grid surface is $CRLB{({N,\varepsilon } )^{\textrm{1/2}}}$ after deadtime correction. (b) $bias({{{\hat{\tau }}_{ML}}} )$ vs. $({N,\varepsilon } )$ before (marked with blue grid) and after (marked with red grid) deadtime correction. (c)$RMSE({\hat{\tau }_{ML}})$ vs.$N$ at several fixed $\varepsilon .$ (d)$bias({{{\hat{\tau }}_{ML}}} )$ vs.$N$ at several fixed $\varepsilon .$

Download Full Size | PDF

The mean-square-error (MSE) of the estimated offsets ${\hat{\tau }_{ML}}$ consists of its variance $\mathrm{\mathbb{D}}[{{{\hat{\tau }}_{ML}}} ]$ and bias $bias({\hat{\tau }_{ML}})$:

When uncertain data is mapped to PPM slots, the distribution of aggregate slot counts obeys a combination of Poisson and binomial distributions [17], wherein the mean and variance are unequal. Furthermore, because ${\hat{\tau }_{ML}}$ is a biased estimate based on Eq. (18), there may be some terms that do not converge to zero as N increases in $\mathrm{\mathbb{D}}[{{{\hat{\tau }}_{ML}}} ]$.

In order to evaluate the timing estimation performance of the algorithm under the influence of other factors, we also characterize the timing estimation performance before and after correction under the scenarios of higher background noise (scenario 2-3), higher modulation order (scenario 4-5) and arrays of different scale (scenario 6-7) (8-pixels, 16-pixels). In these simulations, we mainly focus on the RMSE of offset estimation $\textrm{RMSE}({\hat{\tau }_{ML}})$ as signal flux intensity ${K_s}$ increases.

Table 1 listed other example scenarios 2-7 for which the estimation performance is simulated, while communication scenario 1 has already been simulated above in Figs. 6–7.

Table 1. The system parameters for which the timing estimation is simulated

View Table

In scenario 2-3, higher noise level was assumed (${K_b} = {10^{ - 3}},{10^{ - 2}}$ detected photons/slot) which may occur due to incident background light and photon leakage from electro-optic modulation, while other system parameters keep unchanged, the estimation performance is shown in two subfigures of Fig. 8(a). The $\textrm{RMSE}({\hat{\tau }_{ML}})$ after correction in these scenarios can still reach the order of 10⁻² slot under blocking loss, but we can see from the difference between dashed curves (${K_b} = {10^{ - 3}}$) and solid curves (${K_b} = {10^{ - 2}}$) that higher noise level can negatively affect the estimation performance, both for the corrected and original model.

 

Fig. 8. (a) Under higher noise level, simulated $\textrm{RMSE}({\hat{\tau }_{ML}})$vs.${K_s}$ at $\varepsilon = 0.3\textrm{ and 0}\textrm{.9}.$Dashed curves correspond to scenario 2 when ${K_b} = {10^{ - 3}}$ and solid curves correspond to scenario 3 when ${K_b} = {10^{ - 2}}$. (b) Use higher PPM order, simulated $\textrm{RMSE}({\hat{\tau }_{ML}})$vs.${K_s}$ at $\varepsilon = 0.3\textrm{ and 0}\textrm{.9}.$Dashed curves correspond to scenario 4 when M = 16 and solid curves correspond to scenario 5 when M = 64.

Download Full Size | PDF

In scenario 4-5, higher PPM order was assumed (M = 16, 64) and the symbol scale N is set to the number of symbols contained in an SCPPM frame. The simulation results are shown in two subfigures of Fig. 8(b). Increasing the PPM order clearly has a more negative impact on the estimation performance than the noise, which may significantly reduce the contrast between PPM slot counts and guard slot counts in the aggregate slot vector $\boldsymbol{y}$. And we can also find from Eq. (14) in this paper and Eq. (25) of [17] that CRLB will grow with $O({M^3})$, so the target RMSE of 10⁻² slot level at high PPM orders put forward higher requirements for low noise levels.

Although increasing noise and PPM order will both degrade timing estimation performance, we can mitigate the noise to 10⁻² level at slot frequency more than hundreds of MHz by various ways such as high modulation extinction ratio transmitter [29] and adaptive optical receivers [30]. Additionally, higher-order PPM modulation also introduces more loss on data bandwidth. Therefore, in high-speed deep space communication applications, the algorithm proposed here is still effective.

Another fact we aim to verify through simulation is that the corrected estimation is an extension of the original model under blocking conditions. In array applications, this algorithm can also alleviate synchronization performance losses that may occur in small to medium-scale arrays due to moderate blocking. In large-scale arrays capable of significantly mitigating blocking loss, the performance of the corrected algorithm approaches that of the original one. Scenario 6 and 7 illustrate such two cases, and the corresponding performance is shown in Fig. 9(a)-(b). It can be observed that our proposed approach is not only effective for synchronization estimation under strong blocking conditions for single detector. Moreover, when the array cannot completely eliminate detector blockage, the modified algorithm further compensates for estimation performance loss. And the algorithm also exhibits higher tolerance to noise levels in array detection scenarios (${K_b} = {10^{ - 2}}$).

 

Fig. 9. Detector array configuration with slight blocking. (a) Scenario 6, simulated $\textrm{RMSE}({\hat{\tau }_{ML}})$vs. ${K_s}$ at several fixed $\varepsilon $ with an 8-pixel array. (b) Scenario 7, simulated $\textrm{RMSE}({\hat{\tau }_{ML}})$vs. ${K_s}$ at several fixed $\varepsilon $ with a 16-pixel array.

Download Full Size | PDF

3.2 Experimental verification

In this section, an experiment is designed to verify the performance of the proposed MLE algorithm, and the communication performance after synchronization. In the experimental system, the slot phase and frequency offset, respectively denoted as $\tau \textrm{ and }\gamma $, exist between the transmitter and receiver. Furthermore, the detector experiences severe blockage owing to the insufficient slot duration. The schematic diagram of the experimental system is shown in Fig. 10.

 

Fig. 10. Schematic diagram of the experimental system for verifying the performance of MLE algorithm. The binary data frame is loaded into the FPGA (Xilinx, XC7K325T-2FFG900C) and sent in a loop by GTX transceiver. RF amplifier driver (iXBlue, DR-VE-10-MO) amplifies the signal and drives the laser (Omicron, LuxX + 450-500) to generate modulated laser pulses. The beam enters the free space through the collimator and undergoes a fixed attenuation (50 dB) and adjustable attenuation ($\le $40 dB) to simulate channel loss. Pulses attenuated to several photons are detected by an SPAD (Excelitas, SPCM-AQRH-13-FC, deadtime is 28.5 ns). Thereafter, the rising edges of the detector output pulses are converted into timestamps by TDC (qutools, quTAG-MC-32) as the input for estimation algorithms. An auxiliary frame with the same length of data frame is transmitted using another GTX channel and drives an SFP optical transceiver (Accelink, RTXM192-500). This auxiliary signal is used for obtaining actual phase offset to compare with the offset estimation and does not exist in practical communication scenarios. Furthermore, a PIN detector (Thorlabs, DET08CFC) is used for photoelectric conversion to feed the waveform into TDC for obtaining timestamps. The pulsed slot of each symbol in the auxiliary frame is fixed at the first slot, and the interval between adjacent pulses equals to 8-PPM + 4ISGT symbol duration.

Download Full Size | PDF

Herein, the random binary data are encoded with SCPPM frames of 8PPM + 4ISGT modulation and 2/3 code rate and using MATLAB. Frame synchronization is required prior to the BER calculation. Hence, an FSM32 (frame synchronization marker) marker with the length of 32 symbols is added at the beginning of each encoded frame to perform frame synchronization based on autocorrelation value peak searching. Each element of FSM32 is mapped to an 8PPM + 4ISGT symbol. The marker is obtained by connecting the FSM8 marker [31] twice in reverse to realize better autocorrelation characteristics and ensure frame synchronization robustness under various SNRs. The autocorrelation characteristics of FSM8 and FSM32 are shown in Fig. 11. The total frame consists of 60684 bits comprising 5072 symbols.

 

Fig. 11. Autocorrelation characteristics of FSM8 and FSM32.

Download Full Size | PDF

The actual slot width ${T_\textrm{s}}$ used herein is 12 ns instead of 10 ns, owing to the optical pulse stretching (as shown in the top-left subfigure of Fig. 2(b)), which is caused by the bandwidth limitation of laser modulation ($\le $100 MHz) and additional system jitter (∼1 ns FWHM, mainly from SPAD). Considering the guard slots insertion and encoding redundancy, the final pulse repetition frequency is 6.94 MHz, with an effective data rate of 13.89Mbps. The data rate could be improved to hundreds of Mbps or Gbps by applying high extinction ratio electro-optic modulation and detectors with high count rate and low temporal jitter parameters (such as SNSPD) to address the limiting factors.

$N = 5072$ symbols were used to generate each estimated offset ${\hat{\tau }_{ML}}.$ After the slot frequency offset compensation ($\gamma $ is about ${10^{ - 6}}$), 100 consecutive estimated phase offsets were selected to calculate the RMSE and bias from the actual phase offset ${\tau _{real}}$, which was obtained using the auxiliary frame signal detected by PIN and the single-photon signal detected by SPAD, as shown in Fig. 12(a)-(b). The auxiliary signal behaved as a periodic signal with period equal to the symbol duration including path delay. Using the rising edge of auxiliary signal as trigger, the histogram of SPAD output timestamps reflected the actual path delay of the SPAD signal relative to the auxiliary signal. The actual symbol boundary offset ${\tau _{real}}$ remained almost constant (equal to the path delay between the two signals) because the starting point of the receiver symbol was set at the auxiliary signal timestamp. Therefore, the estimation performance before and after correction could be characterized under different signal intensities.

 

Fig. 12. (a) The yellow signal is the auxiliary signal waveform output by the PIN detector, and whose period is equal to the symbol duration of 144 ns. The green signal is the single-photon response TTL pulse output by SPAD. (b) The histogram of the rising edge timestamp of the SPAD response pulse using that of the auxiliary signal as a trigger. The exact position of the symbol boundary at the receiver is determined using the lowest hard decision BER.

Download Full Size | PDF

A photometer was used to calibrate the signal intensity and set various attenuation values. Using the aforementioned methods, the experimental performances of the deadtime-corrected and uncorrected MLE algorithms under different intensities are shown in Fig. 13. In Fig. 13(a), the black solid lines in each subFig. indicate the corresponding ${\tau _{real}}$ at each signal intensity. Owing to the inter-symbol interference, ${\tau _{real}}$ decreases slightly with increase in attenuation. However, ${\hat{\tau }_{ML}}$ values estimated at each intensity were compared to the current ${\tau _{real}}$ values. The blue and red hollow circles represent ${\hat{\tau }_{ML}}$ before and after deadtime-correction, respectively, and the solid lines of the corresponding colors represent the mean of 100 estimated offsets. Figures 13(b)-(c) quantitatively compare the RMSE and bias of ${\hat{\tau }_{ML}}$ in Fig. 13(a) with simulation results and CRLB^1/2 under the same parameters. Figures 13(a)-(c) indicate that in the experimental and simulation results, the two algorithms exhibited similar performances at low photon flux, thereby validating Eq. (20) that the deadtime-corrected model degrades to a linear offset model at low flux. As the intensity increases, the RMSE of the deadtime-corrected offset estimate ${\tau _{ML}}_{ - deadtime}$ decreases, and the mean gradually approaches ${\tau _{real}}$, whereas the estimate ${\tau _{ML}}_{ - uncorrected}$ based on the linear offset model suffered serious performance degradation.

 

Fig. 13. Experimental results. (a) Estimate of PPM + ISGT symbol boundary offset between TX and RX (phase offset) using the two algorithms under different photon fluxes. $N = 5072$, estimate number = 100. (b) Experimental RMSE of ${\hat{\tau }_{ML}}$ using the following two algorithms: uncorrected (blue solid triangle) and corrected (red solid triangle). Black open triangles represent CRLB^1/2 in Eq. (14). Simulated RMSE under experimental parameters are plotted with open triangles of corresponding colors. (c)$\textrm{bias}({{{\hat{\tau }}_{ML}}} )$ of experimental and simulation results obtained using the two algorithms. Simulated bias under experimental parameters are plotted with open triangles of corresponding colors. (d) Experimental hard decision BERs of PPM demodulation under various photon fluxes using the MLE without correction (blue triangle) and with deadtime correction (red triangle). Black triangles are hard decision BERs at ideal synchronization. (e) SCPPM-decoded BERs using the following two synchronization algorithms: uncorrected (blue open pentagram) and corrected (red open pentagram). (f-g) Autocorrelation results of the FSM32 marker after phase and frequency offset compensation at different photon fluxes.

Download Full Size | PDF

Notably, the experimental and simulation results differ at low photon flux in Fig. 13(b)-(c) because the shape of photon arrival distribution within the slot approaches the modulated laser pulse superimposed with jitter. Owing to the limited laser modulation bandwidth, the optical pulse is similar to the Gaussian rather than rectangular pulse (see the top-left subfigure of Fig. 2(b)), which is inconsistent with the assumption of rectangular pulses in the model. Increasing the electro-optic modulation bandwidth can effectively overcome this limitation. At higher photon flux, the exponential distribution distortion caused by blocking effect gradually dominates and the impact of non-rectangular modulation pulses decreases; hence, the observed simulation and experimental performance tend to be consistent.

Because low photon intensity is not suitable for reliable communication, the synchronization performance at higher flux is investigated. Based on the experimental results, RMSE of the deadtime-corrected offset estimate can reach 0.02749 slots at ${K_s} = 2.45,{K_b} = 0.0037$ and approach the lower bound, whereas the uncorrected model has an RMSE of 0.18734 slots under the same conditions. The CRLB^1/2 in the Fig. 13(b) was calculated using Eq. (14) and the corresponding experimental parameters.

Simultaneously, Figures 13(f)-(g) show the frame synchronization results of the FSM32 marker after slot synchronization, including slot frequency and phase offset compensation. After frame synchronization, Figures 13(d)-(e) show the hard decision BER using PPM demodulation and decoded BER, respectively. Evidently, the frame synchronization performance and hard decision BER at ideal synchronization do not satisfy the requirements of reliable communication at low photon flux (see the black hollow triangles in Fig. 13(d)). They improved with increasing signal intensity; however, the performance degradation of the uncorrected MLE decreases the BER. The hard decision BER synchronized using the original MLE is 7.28%, whereas performance close to ideal synchronization (5.67%) was obtained using the corrected synchronization algorithm. The decoded BER decreased from 3.11 × 10⁻⁵ before correction to less than 10⁻⁶ after correction. These experimental results prove that the deadtime-corrected MLE synchronization scheme can effectively solve the performance degradation of the original model under higher photon flux.

4. Conclusion

This study investigated the synchronization of PPM photon-counting communication systems containing slot frequency and phase offsets between the terminals, and detector blockage owing to sub-deadtime slots. First, the deadtime-correction was introduced to MLE based on ISGT in [16,17] using the renewal process, and the variations in CRLB and expectation of the estimator were described. Second, a method to compensate slot frequency offset using LLS and MLE algorithms was proposed. Third, the synchronization performances (including the RMSE and bias of estimated offsets) of MLE under various parameters were quantitatively simulated before and after correction. The corrected MLE was observed to degenerate to the original form at low signal intensity; however, it effectively improved the original model’s performance degradation under higher photon flux. Finally, experiments were performed to validate the theory and simulation results. At a signal flux intensity of ${K_\textrm{s}} = 2.45$ photons/symbol and $N = 5072$, the synchronization performance of the corrected MLE attained 0.02749 slots @RMSE. The hard decision BER approached the performance at ideal synchronization (5.67%), and the decoded BER was less than 10⁻⁶. Contrastingly, the original MLE's estimation performance was 0.18734 slots @RMSE, corresponding to hard decision and decoded BERs of 7.28% and 3.11 × 10⁻⁵ respectively.

Our work provides a viable solution for realizing the synchronization of PPM photon counting communication systems with blocking loss. And in future work, we will verify our algorithm on higher-speed communication systems based on superconducting nanowire single-photon detectors.