Digital Signal Processing (DSP)

Core digital signal processing utilities (optic.dsp.core)

sigPow(x)	Calculate the average power of x per mode.
signalPower(x)	Calculate the total power of x.
firFilter(h, x)	Perform FIR filtering and compensate for filter delay.
rrcFilterTaps(t, alpha, Ts)	Generate Root-Raised Cosine (RRC) filter coefficients.
rcFilterTaps(t, alpha, Ts)	Generate Raised Cosine (RC) filter coefficients.
pulseShape(param)	Generate a pulse shaping filter.
clockSamplingInterp(x, inFs, outFs[, jitter])	Interpolate signal to a given sampling rate.
quantizer(x[, nBits, maxV, minV])	Quantize the input signal using a uniform quantizer with the specified precision.
lowPassFIR(fc, fs, N[, typeF])	Calculate FIR coefficients of a lowpass filter.
decimate(sigIn, param)	Decimate signal.
resample(sigIn, param)	Resample signal to a desired sampling rate.
upsample(x, factor)	Upsample a signal by inserting zeros between samples.
symbolSync(rx, tx, SpS[, mode])	Symbol synchronizer.
finddelay(x, y)	Find delay between x and y.
pnorm(x)	Normalize the average power of each componennt of x.
anorm(x)	Normalize the amplitude of each componennt of x.
gaussianComplexNoise(shapeOut, σ2=None[, seed])	Generate complex circular Gaussian noise.
gaussianNoise(shapeOut, σ2=None[, seed])	Generate Gaussian noise.
phaseNoise(lw, Nsamples, Ts[, seed])	Generate realization of a random-walk phase-noise process.
movingAverage(x, N)	Calculate the sliding window moving average of a 2D NumPy array along each column.
delaySignal(sig, delay[, Fs, NFFT])	Apply a time delay to a signal sampled at fs samples per second using FFT/IFFT algorithms.
blockwiseFFTConv(x, h[, NFFT, freqDomainFilter])	Blockwise convolution in the frequency domain using the overlap-and-save FFT method.
freqShift(x, deltaF, Fs)	Frequency shift of a signal.
iqMixing(sig, param)	Add IQ mixing to a signal.
calcMZM(sigIn, Vpi, u, Vb, ER)	Fast function to calculate the Mach-Zehnder modulator (MZM) model.
calcPM(sigIn, Vpi, u)	Fast function to calculate the phase modulator (PM) model.
levinson(r, nTaps)	Levinson-Durbin algorithm
autocorr(x, nTaps)	Estimate the autocorrelation coefficients of a signal x up to lag nTaps-1.
estimateWhiteningFilter(x, nTaps)	Estimate the coefficients of a whitening filter of order nTaps using the Levinson-Durbin algorithm.

anorm(x)

Normalize the amplitude of each componennt of x.

Parameters:

x (np.array) – Signal.

Returns:

Signal x with each component normalized in amplitude.

Return type:

np.array

Notes

The signal is scaled so that its peak magnitude is equal to one,

\[y[n] = \frac{x[n]}{\max_{m} |x[m]|}, \tag{1}\]

where the maximum is taken over all the elements of \(x\).

autocorr(x, nTaps)

Estimate the autocorrelation coefficients of a signal x up to lag nTaps-1.

Parameters:

• x (array-like) – The input signal for which to estimate the autocorrelation coefficients.

• nTaps (int) – The number of autocorrelation coefficients to estimate (lags from 0 to nTaps-1).

Returns:

r – An array of length nTaps containing the estimated autocorrelation coefficients, where r[k] is the autocorrelation at lag k.

Return type:

np.array

Notes

The autocorrelation coefficients are estimated using the unbiased estimator, which normalizes the sum of products by the number of terms that contribute to each lag. This provides a more accurate estimate of the autocorrelation, especially for larger lags.

Explicitly, the estimate at lag \(k\) is

\[\hat{r}[k] = \frac{1}{N-k}\sum_{n=k}^{N-1} x[n]\, x^*[n-k], \qquad k = 0, \ldots, n_{Taps}-1, \tag{1}\]

where \(N\) is the length of the signal.

References

[1] Gallager, R. G., Introduction to Random Signals and Applied Kalman Filtering. John Wiley & Sons, 2010.

blockwiseFFTConv(x, h, NFFT=None, freqDomainFilter=False)

Blockwise convolution in the frequency domain using the overlap-and-save FFT method.

Parameters:

• x (np.array) – Input signal.

• h (np.array) – Filter impulse response.

• NFFT (int, optional) – FFT size to be used. Must be greater than the length of the filter. If None, it will be set to the next power of 2 greater than or equal to the length of the filter. Default is None.

• freqDomainFilter (bool, optional) – If True, h is assumed to be the frequency response of the filter. If False, the FFT of h will be computed. Default is False.

Returns:

y – The filtered output signal.

Return type:

np.array

Raises:

ValueError – If NFFT is not greater than the length of the filter h.

Notes

Long convolutions are computed in the frequency domain with the overlap-and-save method. The input is split into overlapping blocks of \(N_{FFT}\) samples, each one sharing its first \(K-1\) samples with the previous block, where \(K\) is the length of the filter impulse response. For each block \(b\),

\[y_b = \mathrm{IDFT}\left\{\mathrm{DFT}\{x_b\} \cdot H\right\}, \tag{1}\]

where \(H\) is the \(N_{FFT}\)-point DFT of the zero-padded impulse response. The first \(K-1\) samples of \(y_b\) are corrupted by the circular wrap-around of the DFT and are discarded, while the remaining \(N_{FFT} - K + 1\) samples are equal to the linear convolution and are concatenated to form the output. The filter delay \(\lfloor (K-1)/2 \rfloor\) is compensated. The computational cost grows as \(\mathcal{O}(N\log N_{FFT})\) instead of \(\mathcal{O}(NK)\) for the direct convolution.

calcMZM(sigIn, Vpi, u, Vb, ER)

Fast function to calculate the Mach-Zehnder modulator (MZM) model.

Parameters:

• sigIn (np.array or float) – Complex-valued optical input field.

• Vpi (float) – Half-wave voltage of the MZM.

• u (float) – RF voltage applied to the MZM.

• Vb (float) – DC bias voltage.

• ER (float) – Extinction ratio of the MZM (in dB).

Returns:

Complex-valued optical output field after modulation.

Return type:

np.array or float

Notes

A Mach-Zehnder modulator (MZM) splits the input field between two arms, applies opposite phase shifts \(\pm\theta\) to them (push-pull operation), and recombines the two fields. With the drive voltage \(u(t)\) and the bias \(V_b\), the phase shift in each arm is

\[\theta(t) = \frac{\pi}{2}\frac{u(t) + V_b}{V_\pi}. \tag{1}\]

A finite extinction ratio \(\varepsilon = 10^{ER/10}\) is modeled by an imbalance between the fields of the arms,

\[E_{out}(t) = \frac{E_{in}(t)}{2}\left[\sqrt{1+\gamma}\,e^{j\theta(t)} + \sqrt{1-\gamma}\,e^{-j\theta(t)}\right], \qquad \gamma = \frac{2\sqrt{\varepsilon}}{\varepsilon + 1}, \tag{2}\]

which can be written as

\[E_{out}(t) = E_{in}(t)\left[c_I\cos\theta(t) + jc_Q\sin\theta(t)\right], \qquad c_{I,Q} = \frac{\sqrt{1+\gamma} \pm \sqrt{1-\gamma}}{2}. \tag{3}\]

The ratio between the maximum and the minimum output powers is \(c_I^2/c_Q^2 = \varepsilon\). For an infinite extinction ratio, \(c_I = 1\) and \(c_Q = 0\), and Eq. (3) reduces to the transfer function of the ideal MZM, \(E_{out} = E_{in}\cos\theta\).

References

[1] Y. Yamaguchi, et al, “Precise Optical Modulation Using Extinction-Ratio and Chirp Tunable Single-Drive Mach–Zehnder Modulator,” Journal of Lightwave Technology, vol. 35, no. 21, pp. 4781-4788, 1 Nov.1, 2017,

[2] Seimetz, M., High-Order Modulation for Optical Fiber Transmission. Springer Series in Optical Sciences. Springer Berlin Heidelberg, 2009.

calcPM(sigIn, Vpi, u)

Fast function to calculate the phase modulator (PM) model.

Parameters:

• sigIn (np.array or float) – Complex-valued optical input field.

• Vpi (float) – Half-wave voltage of the PM.

• u (float) – Driving voltage applied to the PM.

Returns:

Complex-valued optical output field after modulation.

Return type:

np.array or float

Notes

An optical phase modulator (PM) driven by the voltage \(u(t)\) imposes a phase shift proportional to it,

\[E_{out}(t) = E_{in}(t)\exp\left[j\pi\frac{u(t)}{V_\pi}\right], \tag{1}\]

where \(V_\pi\) is the voltage that produces a phase shift of \(\pi\) rad.

References

[1] Seimetz, M., High-Order Modulation for Optical Fiber Transmission. Springer Series in Optical Sciences. Springer Berlin Heidelberg, 2009.

clockSamplingInterp(x, inFs, outFs, jitter=0)

Interpolate signal to a given sampling rate.

Parameters:

• x (np.array) – Input signal.

• inFs (float) – Sampling frequency of the input signal.

• outFs (float) – Sampling frequency of the output signal.

• jitter (float) – Standard deviation of the time jitter (jitter rms). Default is 0.

Returns:

y – Resampled signal.

Return type:

np.array

Notes

The input samples \(x[n] = x(nT_{in})\), with \(T_{in} = 1/F_{in}\), are resampled at the instants of the output sampling clock,

\[t_m = m T_{out} + \epsilon_m, \qquad T_{out} = 1/F_{out}, \tag{1}\]

where \(\epsilon_m \sim \mathcal{N}(0, \sigma_j^2)\) models a random timing jitter with standard deviation (rms jitter) \(\sigma_j\). The signal value at \(t_m\) is obtained by linear interpolation between the two nearest input samples: for \(nT_{in} \le t_m < (n+1)T_{in}\),

\[y[m] = x[n] + \frac{t_m - nT_{in}}{T_{in}}\left(x[n+1] - x[n]\right). \tag{2}\]

decimate(sigIn, param)

Decimate signal.

Parameters:

• sigIn (np.array) – Input signal.

• param (optic.utils.parameters object, optional) –

  Parameters of the decimation process.

  • param.SpSin : samples per symbol of the input signal.

  • param.SpSout : samples per symbol of the output signal.

Returns:

sigOut – Decimated signal.

Return type:

np.array

Notes

The signal has \(\mathrm{SpS}_{in}\) samples per symbol and is reduced to \(\mathrm{SpS}_{out}\) samples per symbol by keeping one out of every \(D = \mathrm{SpS}_{in}/\mathrm{SpS}_{out}\) samples. The sampling phase \(n_0\) is chosen as the position within the symbol period where the signal variance is maximum, which, for a Nyquist-shaped signal, corresponds to the center of the eye diagram:

\[n_0 = \arg\max_{k \in \{0, \ldots, \mathrm{SpS}_{in}-1\}} \mathrm{Var}\left\{x[m\,\mathrm{SpS}_{in} + k]\right\}_m, \tag{1}\]

\[y[m] = x[n_0 + mD]. \tag{2}\]

References

[1] P. S. R. Diniz, E. A. B. da Silva, e S. L. Netto, Digital Signal Processing: System Analysis and Design. Cambridge University Press, 2010.

delaySignal(sig, delay, Fs=1, NFFT=1024)

Apply a time delay to a signal sampled at fs samples per second using FFT/IFFT algorithms.

Parameters:

• sig (np.array) – The input signal.

• delay (float) – The time delay to apply to the signal (in seconds).

• Fs (float) – Sampling frequency of the signal (in samples per second). Default is 1.

• NFFT (int, optional) – FFT size to be used. Must be greater than the length of the filter. If None, it will be set to the next power of 2 greater than or equal to the length of the (zero-padded) signal. Default is 1024.

Returns:

The delayed signal, with the same length as sig. If delay is zero, a copy of sig is returned.

Return type:

np.array

Notes

A time delay \(\tau\) corresponds to a linear phase in the frequency domain,

\[y(t) = x(t - \tau) \quad \Longleftrightarrow \quad Y(f) = X(f)\, e^{-j2\pi f\tau}, \tag{1}\]

which holds for any real \(\tau\), including fractions of the sampling period. The frequency response in Eq. (1) is applied by blockwise FFT convolution (overlap-and-save, see blockwiseFFTConv), after zero padding the signal to avoid the circular wrap-around of the delayed samples.

estimateWhiteningFilter(x, nTaps)

Estimate the coefficients of a whitening filter of order nTaps using the Levinson-Durbin algorithm.

Parameters:

• x (array-like) – The input signal from which to estimate the autocorrelation.

• nTaps (int) – The order of the whitening filter (number of coefficients).

Returns:

w – The coefficients of the whitening filter of length nTaps, where w[0] is the leading coefficient (usually 1).

Return type:

np.array

Notes

The whitening filter is the prediction error filter of the signal: the autocorrelation of \(x[n]\) is estimated with autocorr, and the coefficients \(a_k\) are obtained by solving the Yule-Walker equations with the Levinson-Durbin recursion (levinson). The output of the filter,

\[e[n] = \sum_{k=0}^{n_{Taps}-1} a_k\, x[n-k], \qquad a_0 = 1, \tag{1}\]

is the prediction error, which is approximately white when the filter order is large enough to capture the correlation of \(x[n]\).

References

[1] Levinson, N., The Wiener RMS error criterion in filter design. Journal of Mathematics and Physics, 25(1-4), 261-278, 1947.

[2] Durbin, J., The fitting of time-series models. Review of the International Statistical Institute, 28(3), 233-244, 1960.

finddelay(x, y)

Find delay between x and y.

Parameters:

• x (np.array) – Signal 1.

• y (np.array) – Signal 2.

Returns:

d – Delay between x and y, in samples.

Return type:

int

Notes

The delay is estimated from the peak of the magnitude of the cross-correlation between the two sequences,

\[R_{xy}[\tau] = \sum_{n} x[n + \tau]\, y^*[n], \qquad \hat{\tau} = \arg\max_{\tau}\, |R_{xy}[\tau]|, \tag{1}\]

so that \(y[n] \approx x[n + \hat{\tau}]\). Using the magnitude of the correlation makes the estimate insensitive to a constant phase rotation between the sequences.

firFilter(h, x)

Perform FIR filtering and compensate for filter delay.

Parameters:

• h (np.array) – Coefficients of the FIR filter (impulse response, symmetric).

• x (np.array) – Input signal.

Returns:

y – Output (filtered) signal.

Return type:

np.array

Notes

The output of a finite impulse response (FIR) filter with coefficients \(h[k]\), \(k = 0, \ldots, L-1\), is given by the discrete convolution \((h * x)[n] = \sum_{k} h[k]\, x[n-k]\). Since a filter with symmetric coefficients delays the signal by \(D = \lfloor (L-1)/2 \rfloor\) samples, this delay is removed from the output, which is computed as

\[y[n] = \sum_{k=0}^{L-1} h[k]\, x[n + D - k], \qquad n = 0, \ldots, N-1, \tag{1}\]

so that \(y[n]\) has the same length as \(x[n]\) and is time-aligned with it. The convolution is evaluated with the overlap-add method, which is efficient when \(L \ll N\). Each column of \(x\) is filtered independently.

References

[1] P. S. R. Diniz, E. A. B. da Silva, e S. L. Netto, Digital Signal Processing: System Analysis and Design. Cambridge University Press, 2010.

freqShift(x, deltaF, Fs)

Frequency shift of a signal.

Parameters:

• x (np.array) – Input signal.

• deltaF (float) – Frequency shift (Hz).

• Fs (float) – Sampling frequency (Hz).

Returns:

y – Frequency shifted signal.

Return type:

np.array

Notes

A frequency shift \(\Delta f\) is the multiplication by a complex exponential (modulation theorem),

\[y[n] = x[n]\, e^{j2\pi\Delta f\, n/F_s} \quad \Longleftrightarrow \quad Y(f) = X(f - \Delta f), \tag{1}\]

where \(F_s\) is the sampling frequency.

gaussianComplexNoise(shapeOut, σ2=1.0, seed=None)

Generate complex circular Gaussian noise.

Parameters:

• shapeOut (tuple of int) – Shape of np.array to be generated.

• σ2 (float, optional) – Variance of the noise (default is 1).

• seed (int, optional) – Seed for the random number generator.

Returns:

noise – Generated complex circular Gaussian noise.

Return type:

np.array

Notes

The samples are drawn from a zero-mean circularly-symmetric complex Gaussian distribution, \(n = n_I + j n_Q\), where the in-phase and quadrature components are independent and identically distributed, \(n_I, n_Q \sim \mathcal{N}(0, \sigma^2/2)\). Hence \(\mathbb{E}\left[|n|^2\right] = \sigma^2\) and the probability density function is

\[p(n) = \frac{1}{\pi\sigma^2}\exp\left(-\frac{|n|^2}{\sigma^2}\right). \tag{1}\]

This is the standard model for additive white Gaussian noise (AWGN) in complex baseband.

gaussianNoise(shapeOut, σ2=1.0, seed=None)

Generate Gaussian noise.

Parameters:

• shapeOut (tuple of int) – Shape of np.array to be generated.

• σ2 (float, optional) – Variance of the noise (default is 1).

• seed (int, optional) – Seed for the random number generator.

Returns:

noise – Generated Gaussian noise.

Return type:

np.array

Notes

The samples are drawn from a zero-mean real Gaussian distribution with variance \(\sigma^2\),

\[p(n) = \frac{1}{\sqrt{2\pi\sigma^2}}\exp\left(-\frac{n^2}{2\sigma^2}\right). \tag{1}\]

iqMixing(sig, param)

Add IQ mixing to a signal.

Parameters:

• sig (np.array) – Input signal.

• param (optic.utils.parameters object) –

  Parameters of IQ mixing.

  param.ampImb : Amplitude imbalance parameter in dB.[default: 0 dB] param.phaseImb : Phase imbalance parameter (in radians).[default: 0 rad] param.timeSkew : delay between I and Q components. [default: 0 s] param.Fs : simulation sampling frequency. [default: None]

Returns:

IQ-mixed signal.

Return type:

np.array

Notes

Imperfections of the in-phase (I) and quadrature (Q) branches of a receiver front-end are modeled as follows. For an input \(s = I + jQ\), an amplitude imbalance \(\epsilon\) and a phase imbalance \(\phi\) between the branches result in

\[y_I = (1-\epsilon)\left[I\cos\frac{\phi}{2} - Q\sin\frac{\phi}{2}\right], \qquad y_Q = (1+\epsilon)\left[Q\cos\frac{\phi}{2} - I\sin\frac{\phi}{2}\right], \tag{1}\]

where \(\epsilon = 10^{A_{dB}/20} - 1\) is obtained from the amplitude imbalance in dB. Equivalently, in complex notation, the IQ imbalance adds an image of the complex conjugate of the signal,

\[y = k_1 s + k_2 s^*, \tag{2}\]

with \(k_1 = \left[(1-\epsilon)e^{j\phi/2} + (1+\epsilon)e^{-j\phi/2}\right]/2\) and \(k_2 = \left[(1-\epsilon)e^{-j\phi/2} - (1+\epsilon)e^{j\phi/2}\right]/2\). Finally, a time skew \(\tau\) between the branches is applied by advancing the I component and delaying the Q component by \(\tau/2\).

levinson(r, nTaps)

Levinson-Durbin algorithm

Parameters:

• r (array-like) – Autocorrelation coefficients of the signal, where r[0] is the zero-lag autocorrelation and r[k] is the autocorrelation at lag k.

• nTaps (int) – The order of the whitening filter (number of coefficients to estimate).

Returns:

a – The coefficients of the whitening filter of length nTaps, where a[0] is the leading coefficient (usually 1) and a[1], a[2], …, a[nTaps-1] are the estimated filter coefficients.

Return type:

np.np.array

Notes

The Levinson-Durbin algorithm is an efficient method for solving the Toeplitz system of equations that arises in linear prediction and filter design. The resulting coefficients can be used to design a whitening filter that decorrelates the input signal.

In linear prediction of order \(p = n_{Taps} - 1\), the prediction error filter \(A(z) = \sum_{k=0}^{p} a_k z^{-k}\), with \(a_0 = 1\), minimizes the power of \(e[n] = \sum_{k=0}^{p} a_k x[n-k]\). Its coefficients satisfy the Yule-Walker (normal) equations, whose matrix is Toeplitz and Hermitian,

\[\sum_{k=1}^{p} a_k\, r[i-k] = -r[i], \qquad i = 1, \ldots, p, \tag{1}\]

where \(r[k]\) is the autocorrelation of \(x[n]\). The Levinson-Durbin recursion solves Eq. (1) in \(\mathcal{O}(p^2)\) operations, increasing the order one step at a time. Starting from \(E_0 = r[0]\), for \(i = 1, \ldots, p\):

\[k_i = -\frac{r[i] + \sum_{j=1}^{i-1} a_j^{(i-1)}\, r[i-j]}{E_{i-1}}, \tag{2}\]

\[a_j^{(i)} = a_j^{(i-1)} + k_i\, a_{i-j}^{(i-1)*}, \quad j = 1, \ldots, i-1, \qquad a_i^{(i)} = k_i, \tag{3}\]

\[E_i = \left(1 - |k_i|^2\right) E_{i-1}, \tag{4}\]

where \(k_i\) are the reflection coefficients and \(E_i\) is the prediction error power of order \(i\).

References

[1] Levinson, N., The Wiener RMS error criterion in filter design. Journal of Mathematics and Physics, 25(1-4), 261-278, 1947.

[2] Durbin, J., The fitting of time-series models. Review of the International Statistical Institute, 28(3), 233-244, 1960.

lowPassFIR(fc, fs, N, typeF='rect')

Calculate FIR coefficients of a lowpass filter.

Parameters:

• fc (float) – Cutoff frequency.

• fs (float) – Sampling frequency.

• N (int) – Number of filter coefficients.

• typeF (string, optional) – Type of response (‘rect’, ‘gauss’). The default is “rect”.

Returns:

h – Filter coefficients.

Return type:

np.array

Notes

With the normalized cutoff frequency \(f_u = f_c/f_s\) and the filter delay \(d = (N-1)/2\), two responses are available:

• 'rect': truncated impulse response of the ideal lowpass filter, whose frequency response is rectangular with cutoff \(f_c\),

  \[h[n] = 2f_u\,\mathrm{sinc}\left(2f_u(n-d)\right), \tag{1}\]

  where \(\mathrm{sinc}(x) = \sin(\pi x)/(\pi x)\).

• 'gauss': Gaussian filter, with frequency response \(H(f) = \exp\left[-\frac{\ln 2}{2}\left(f/f_c\right)^2\right]\), so that \(|H(f_c)|^2 = 1/2\) (3-dB cutoff at \(f_c\)), and impulse response

  \[h[n] = \sqrt{\frac{2\pi}{\ln 2}}\, f_u \exp\left[-\frac{2}{\ln 2}\left(\pi f_u (n-d)\right)^2\right]. \tag{2}\]

In both cases the coefficients are then normalized to unit sum (unit gain at DC).

References

[1] P. S. R. Diniz, E. A. B. da Silva, e S. L. Netto, Digital Signal Processing: System Analysis and Design. Cambridge University Press, 2010.

movingAverage(x, N)

Calculate the sliding window moving average of a 2D NumPy array along each column.

Parameters:

• x (np.array) – Input 2D array with shape (M, N), where M is the number of samples and N is the number of columns.

• N (int) – Size of the sliding window.

Returns:

2D array containing the sliding window moving averages along each column.

Return type:

np.array

Notes

The function pads the signal with zeros at both ends to compensate for the lag between the output of the moving average and the original signal.

The moving average over a window of \(N\) samples centered at each output sample is the FIR filter with \(N\) equal coefficients \(1/N\),

\[y[n] = \frac{1}{N}\sum_{k=-\lfloor N/2 \rfloor}^{N - 1 - \lfloor N/2 \rfloor} x[n+k], \tag{1}\]

with \(x[n] = 0\) outside the signal (zero padding at the edges).

phaseNoise(lw, Nsamples, Ts, seed=None)

Generate realization of a random-walk phase-noise process.

Parameters:

• lw (scalar) – laser linewidth.

• Nsamples (scalar) – number of samples to be draw.

• Ts (scalar) – sampling period.

• seed (int, optional) – Seed for the random number generator.

Returns:

phi – realization of the phase noise process.

Return type:

np.array

Notes

The phase noise of a laser with linewidth \(\Delta\nu\) is modeled as a Wiener process (random walk), whose increments over a sampling period \(T_s\) are independent zero-mean Gaussian random variables:

\[\phi[k+1] = \phi[k] + \Delta_k, \qquad \Delta_k \sim \mathcal{N}\left(0,\, \sigma^2\right), \qquad \sigma^2 = 2\pi\Delta\nu T_s, \tag{1}\]

with \(\phi[0] = 0\). The variance of the phase difference accumulated over an interval \(\tau\) grows linearly with it, \(2\pi\Delta\nu|\tau|\), and the corresponding optical field \(e^{j\phi(t)}\) has a Lorentzian power spectral density with full width at half maximum equal to \(\Delta\nu\).

References

[1] M. Seimetz, High-Order Modulation for Optical Fiber Transmission. em Springer Series in Optical Sciences. Springer Berlin Heidelberg, 2009.

pnorm(x)

Normalize the average power of each componennt of x.

Parameters:

x (np.array) – Signal.

Returns:

Signal x with each component normalized in power.

Return type:

np.array

Notes

The signal is scaled to unit average power,

\[y[n] = \frac{x[n]}{\sqrt{P_x}}, \qquad P_x = \frac{1}{N}\sum_{n} |x[n]|^2, \tag{1}\]

where the average power \(P_x\) is computed over all the elements of \(x\).

pulseShape(param)

Generate a pulse shaping filter.

Parameters:

param (optic.utils.parameters object, optional) –

Parameters of the pulse shaping filter:

• param.pulseType : Type of pulse shaping filter (‘rect’,’nrz’,’rrc’,’rc’, ‘doubinary’). [default: ‘rrc’]

• param.SpS : Number of samples per symbol of input signal.[default: 2]

• param.nFilterTaps : Number of filter coefficients. [default: 256]

• param.rollOff : Rolloff of RRC/RC filters. [default: 0.1]

Returns:

filterCoeffs – Array of filter coefficients (normalized).

Return type:

np.array

Notes

The pulse \(p[n]\) is sampled with \(\mathrm{SpS}\) samples per symbol period (\(T_s = 1\)), and the available shapes are:

• 'rect': rectangular pulse of duration \(T_s\).

• 'nrz': non-return-to-zero pulse, obtained by smoothing the rectangular pulse with a Gaussian window, which emulates the finite rise and fall times of a real NRZ signal.

• 'rrc' and 'rc': root-raised cosine and raised cosine pulses (see rrcFilterTaps and rcFilterTaps), sampled at \(t = n/\mathrm{SpS}\) over nFilterTaps samples.

• 'duobinary': sum of two sinc pulses one symbol period apart, \(p(t) = \mathrm{sinc}(t/T_s) + \mathrm{sinc}\left((t - T_s)/T_s\right)\).

In all cases the coefficients are normalized to unit sum,

\[\sum_{n} p[n] = 1, \tag{1}\]

that is, the filter has unit gain at DC.

quantizer(x, nBits=16, maxV=1, minV=-1)

Quantize the input signal using a uniform quantizer with the specified precision.

Parameters:

• x (np.array) – The input signal to be quantized.

• nBits (int) – Number of bits used for quantization. The quantizer will have 2^nBits levels.

• maxV (float, optional) – Maximum value for the quantizer’s full-scale range (default is 1).

• minV (float, optional) – Minimum value for the quantizer’s full-scale range (default is -1).

Returns:

The quantized output signal with the same shape as ‘x’, quantized using ‘nBits’ levels.

Return type:

np.array

Notes

A uniform quantizer with \(b\) bits maps its input to one of \(2^b\) levels equally spaced over the full-scale range \([V_{min}, V_{max}]\),

\[d_k = V_{min} + k\Delta, \qquad \Delta = \frac{V_{max} - V_{min}}{2^b - 1}, \qquad k = 0, \ldots, 2^b - 1. \tag{1}\]

Each sample is replaced by the closest level,

\[y[n] = d_{\hat{k}}, \qquad \hat{k} = \arg\min_k\, |x[n] - d_k|. \tag{2}\]

For inputs within the full-scale range, the quantization error \(e[n] = y[n] - x[n]\) is bounded by \(|e[n]| \le \Delta/2\) and, for signals that are busy enough, it is well approximated by a uniformly distributed noise with power \(\Delta^2/12\). Inputs outside the range are mapped to the closest end level.

rcFilterTaps(t, alpha, Ts)

Generate Raised Cosine (RC) filter coefficients.

Parameters:

• t (np.array) – Time values.

• alpha (float) – RC roll-off factor.

• Ts (float) – Symbol period.

Returns:

coeffs – RC filter coefficients.

Return type:

np.array

Notes

The raised cosine (RC) pulse with symbol period \(T_s\) and roll-off factor \(0 \le \alpha \le 1\) is

\[p(t) = \frac{1}{T_s}\,\mathrm{sinc}\!\left(\frac{t}{T_s}\right) \frac{\cos\left(\pi\alpha \frac{t}{T_s}\right)} {1-\left(2\alpha\frac{t}{T_s}\right)^2}, \tag{1}\]

where \(\mathrm{sinc}(x) = \sin(\pi x)/(\pi x)\). At the points \(t = \pm T_s/(2\alpha)\), where the denominator of Eq. (1) vanishes, the limit

\[p\left(\pm\frac{T_s}{2\alpha}\right) = \frac{\pi}{4T_s}\, \mathrm{sinc}\!\left(\frac{1}{2\alpha}\right) \tag{2}\]

is used. The RC pulse satisfies the Nyquist criterion for zero intersymbol interference, \(p(kT_s) = 0\) for every integer \(k \neq 0\). Its spectrum is flat up to \((1-\alpha)/(2T_s)\) and decays with a cosine-shaped transition down to zero at \((1+\alpha)/(2T_s)\).

References

[1] Proakis, J. G., & Salehi, M. (2008). Digital Communications (5th Edition). McGraw-Hill Education.

resample(sigIn, param)

Resample signal to a desired sampling rate.

Parameters:

• sigIn (np.array) – Input signal.

• param (optic.utils.parameters object, optional) –

  Parameters of the resampling process.

  • param.inFs : sampling rate of the input signal [default: 2].

  • param.outFs : sampling rate of the output signal [default: 2].

  • param.N : order of anti-aliasing filter [default: 501].

Returns:

sigOut – Resampled signal.

Return type:

np.array

Notes

The signal is converted from the sampling rate \(F_{in}\) to \(F_{out}\) by linear interpolation (see clockSamplingInterp). To avoid aliasing when \(F_{out} < F_{in}\), the signal is first filtered by a lowpass filter with cutoff \(F_{out}/2\). When \(F_{out} > F_{in}\), the interpolated signal is lowpass filtered with cutoff \(F_{in}/2\), which removes the spectral images created by the interpolation. Both filters are windowed-sinc filters with N coefficients (see lowPassFIR).

References

[1] P. S. R. Diniz, E. A. B. da Silva, e S. L. Netto, Digital Signal Processing: System Analysis and Design. Cambridge University Press, 2010.

rrcFilterTaps(t, alpha, Ts)

Generate Root-Raised Cosine (RRC) filter coefficients.

Parameters:

• t (np.array) – Time values.

• alpha (float) – RRC roll-off factor.

• Ts (float) – Symbol period.

Returns:

coeffs – RRC filter coefficients.

Return type:

np.array

Notes

The root-raised cosine (RRC) pulse with symbol period \(T_s\) and roll-off factor \(0 \le \alpha \le 1\) is the pulse whose squared magnitude spectrum is the raised cosine spectrum, so that a matched pair of RRC filters (transmitter and receiver) produces a Nyquist pulse, free of intersymbol interference. Its impulse response is

\[p(t) = \frac{1}{T_s}\, \frac{\sin\!\left[\pi \frac{t}{T_s}(1-\alpha)\right] + 4\alpha\frac{t}{T_s}\cos\!\left[\pi \frac{t}{T_s}(1+\alpha)\right]} {\pi \frac{t}{T_s}\left[1-\left(4\alpha\frac{t}{T_s}\right)^2\right]}, \tag{1}\]

with the limiting values

\[p(0) = \frac{1}{T_s}\left[1 + \alpha\left(\frac{4}{\pi}-1\right)\right], \tag{2}\]

\[p\left(\pm\frac{T_s}{4\alpha}\right) = \frac{\alpha}{T_s\sqrt{2}} \left[\left(1+\frac{2}{\pi}\right)\sin\frac{\pi}{4\alpha} + \left(1-\frac{2}{\pi}\right)\cos\frac{\pi}{4\alpha}\right]. \tag{3}\]

The bandwidth occupied by the pulse is \((1+\alpha)/(2T_s)\).

References

[1] Proakis, J. G., & Salehi, M. (2008). Digital Communications (5th Edition). McGraw-Hill Education.

sigPow(x)

Calculate the average power of x per mode.

Parameters:

x (np.array) – Signal.

Returns:

Average power of x: P = mean(abs(x)**2).

Return type:

scalar

Notes

For a discrete-time signal \(x[n]\), \(n = 0, \ldots, N-1\), the average power is the mean squared magnitude of its samples,

\[P_x = \frac{1}{N}\sum_{n=0}^{N-1} |x[n]|^2. \tag{1}\]

If \(x\) has several columns, the average in Eq. (1) is taken over all of its elements.

signalPower(x)

Calculate the total power of x.

Parameters:

x (np.array) – Signal.

Returns:

Total power of x: P = sum(abs(x)**2).

Return type:

scalar

Notes

Let \(x_k[n]\), \(k = 1, \ldots, K\), denote the \(K\) components (columns) of the signal, e.g. the polarization or spatial modes of an optical field. The total power is the sum of the average powers of the components,

\[P = \sum_{k=1}^{K} \frac{1}{N}\sum_{n=0}^{N-1} |x_k[n]|^2. \tag{1}\]

symbolSync(rx, tx, SpS, mode='amp')

Symbol synchronizer.

Parameters:

• rx (np.array) – Received symbol sequence.

• tx (np.array) – Transmitted symbol sequence.

• SpS (int) – Samples per symbol of the received signal.

• mode (string, optional) – Synchronization mode: “amp” (amplitude) or “real” (real part). The default is “amp”.

Returns:

tx – Transmitted sequence synchronized to rx.

Return type:

np.array

upsample(x, factor)

Upsample a signal by inserting zeros between samples.

Parameters:

• x (np.array) – Input signal to upsample.

• factor (int) – Upsampling factor. The signal will be upsampled by inserting factor - 1 zeros between each original sample.

Returns:

xUp – The upsampled signal with zeros inserted between samples.

Return type:

np.array

Notes

This function inserts zeros between the samples of the input signal to increase its sampling rate. The upsampling factor determines how many zeros are inserted between each original sample.

If the input signal is a 2D array, the upsampling is performed column-wise.

Formally, upsampling by an integer factor \(L\) is defined as

\begin{equation} x_{\uparrow}[n] = \begin{cases} x[n/L], & n = 0, \pm L, \pm 2L, \ldots \\ 0, & \text{otherwise.} \end{cases} \tag{1} \end{equation}

In the frequency domain, \(X_{\uparrow}(e^{j\omega}) = X(e^{j\omega L})\): the spectrum is compressed by \(L\) and \(L-1\) spectral images appear in \([-\pi, \pi)\). These images are removed by a subsequent interpolation (e.g. pulse shaping) filter.

References

[1] P. S. R. Diniz, E. A. B. da Silva, e S. L. Netto, Digital Signal Processing: System Analysis and Design. Cambridge University Press, 2010.

DSP algorithms for adaptive filtering (optic.dsp.adaptiveFiltering)

coreAdaptEqBlockTD(sigIn, symbRef, SpS, H, ...)	Adaptive equalizer core processing function (block-wise, time-domain).
coreAdaptEqBlockFD(sigIn, symbRef, SpS, H, ...)	Adaptive equalizer core processing function (block-wise, frequency-domain)
realValuedDFECore(sigIn, symbRef[, nTapsFF, ...])	Decision feedback equalizer (DFE) core implementation.
complexValuedDFECore(sigIn, symbRef[, ...])	Decision feedback equalizer (DFE) core implementation for complex-valued signals.
realValuedFFECore(sigIn, symbRef[, nTaps, ...])	Decision-directed feedforward equalizer (FFE) core implementation.
complexValuedFFECore(sigIn, symbRef[, ...])	Decision-directed feedforward equalizer (FFE) core implementation for complex-valued signals.
volterraCore(sigIn, symbRef[, order, SpS, ...])	Decision-directed Volterra equalizer core implementation

complexValuedDFECore(sigIn, symbRef, nTapsFF=5, nTapsFB=5, SpS=1, mu=0.0001, nTrain=1000, prec=<class 'numpy.complex64'>, constSymb=None, f=None, b=None, trainingMode='data-aided', preconvIters=1)

Decision feedback equalizer (DFE) core implementation for complex-valued signals.

Parameters:

• sigIn (np.array) – Input signal to be equalized.

• symbRef (np.array) – Desired (reference) signal.

• nTapsFF (int) – Number of feedforward taps

• nTapsFB (int) – Number of feedback taps

• SpS (int) – Samples per symbol

• mu (float) – Step size

• nTrain (int) – Number of training symbols

• prec (data type) – Precision

• constSymb (np.array) – Constellation symbols

• f (np.array) – Initial feedforward filter coefficients.

• b (np.array) – Initial feedback filter coefficients.

• trainingMode (str) – Operation mode (‘data-aided’, ‘fulltime’)

• preconvIters (int) – Number of pre-convergence iterations

Returns:

• sigOut (np.array) – Equalized output signal.

• f (np.array) – Final feedforward filter coefficients.

• b (np.array) – Final feedback filter coefficients.

References

[1] Proakis, J. G., & Salehi, M. (2008). Digital Communications (5th Edition). McGraw-Hill Education.

complexValuedFFECore(sigIn, symbRef, nTaps=5, SpS=1, mu=0.0001, nTrain=1000, prec=<class 'numpy.complex64'>, constSymb=None, f=None, trainingMode='data-aided', preconvIters=1)

Decision-directed feedforward equalizer (FFE) core implementation for complex-valued signals.

Parameters:

• sigIn (np.array) – Input signal to be equalized.

• symbRef (np.array) – Desired (reference) signal.

• nTaps (int) – Number of feedforward taps

• SpS (int) – Samples per symbol

• mu (float) – Step size

• nTrain (int) – Number of training symbols

• prec (data type) – Precision

• constSymb (np.array) – Constellation symbols

• f (np.array) – Initial feedforward filter coefficients

• trainingMode (str) – Operation mode (‘data-aided’, ‘fulltime’)

• preconvIters (int) – Number of pre-convergence iterations

Returns:

• sigOut (np.array) – Equalized output signal.

• f (np.array) – Final feedforward filter coefficients.

References

[1] Proakis, J. G., & Salehi, M. (2008). Digital Communications (5th Edition). McGraw-Hill Education.

coreAdaptEqBlockFD(sigIn, symbRef, SpS, H, H_, L, mu, lambdaRLS, nTaps, storeCoeff, runWL, alg, constSymb, prec, Nfft)

Adaptive equalizer core processing function (block-wise, frequency-domain)

Parameters:

• sigIn (np.array) – Input signal array.

• symbRef (np.array) – Reference symbol sequence.

• SpS (int) – Samples per symbol.

• H (np.array) – Coefficient matrix.

• H – Augmented coefficient matrix.

• L (int) – Length of the output.

• mu (float) – Step size parameter.

• lambdaRLS (float) – RLS forgetting factor.

• nTaps (int) – Number of taps.

• storeCoeff (bool) – Flag indicating whether to store coefficient matrices.

• runWL (bool) – Run widely-linear mode

• alg (str) – Equalizer algorithm.

• constSymb (np.array) – Constellation symbols.

• prec (data type) – Precision of the computations.

• Nfft (int) – FFT size used for the overlap-and-save frequency-domain filtering. Must satisfy Nfft >= nTaps. A power of two is recommended for FFT efficiency. Larger Nfft amortizes the FFT cost over more valid output samples per block, at the cost of a longer coefficient “freeze” interval (see Notes).

Returns:

• sigOut (np.array) – Equalized output array.

• H (np.array) – Coefficient matrix.

• H_ (np.array) – Augmented coefficient matrix.

• errSq (np.array) – Squared absolute error array.

• Hiter (np.array) – History of coefficient matrices.

Notes

Filtering — overlap-and-save. For a block starting at symbol index bStart, the algorithm reads a window of Nfft consecutive input samples starting at sample bStart*SpS, computes its Nfft-point circular convolution (via FFT) with each zero-padded, time-reversed row of the frozen coefficient matrix, and discards the first nTaps - 1 samples of the result (corrupted by circular wrap-around, as usual in overlap-and-save). The rows are time-reversed (rather than used as-is) because, unlike a conventional causal FIR filter, this equalizer has no group delay: output symbol ind is formed from input samples ind*SpS up to ind*SpS + nTaps - 1 (coreAdaptEqBlockTD’s indIn = indTaps + ind*SpS), i.e. it looks forward rather than backward in time.

The remaining Nvalid = Nfft - nTaps + 1 samples of the discarded-first result are exact, uncorrupted linear-convolution outputs at the full (un-decimated) sample rate; every SpS-th one of them is a valid equalizer output symbol. This yields Lb = Nvalid // SpS output symbols per block. Any leftover valid samples (when Nvalid is not a multiple of SpS) are simply discarded and transparently recomputed as part of the next block’s window; consecutive windows overlap whenever nTaps > SpS, as expected in overlap-and-save.

With nModes propagation modes, the MIMO filter is a set of nModes^2 FIR filters, one per (input mode, output mode) pair: H[m + N*nModes, :] is the response from input mode N to output mode m. In the frequency domain this is nothing more than the classic MIMO channel equation evaluated bin by bin, Y(f) = H(f) X(f): for every FFT bin f, a small nModes x nModes matrix H(f) multiplies the nModes-vector X(f). This is computed for every bin at once with a single np.einsum, instead of a nModes x nModes nested Python loop.

Adaptation — FFT correlation theorem (LMS family only). The tap update of any LMS-type algorithm (nlms, cma, dd-lms, rde, da-rde) is, at its core, a cross-correlation between an “error-like” signal and the equalizer input, evaluated at lags t = 0..nTaps-1. The FFT correlation theorem lets that correlation be computed for all nTaps lags at once from two spectra instead of a Lb x nTaps time-domain sum: corr(a, b)[t] = sum_n a[n] conj(b[n+t]) = conj(IFFT(conj(FFT(a)) * FFT(b)))[t]. Since the error only exists at the (sparsely spaced) symbol instants, it is first “upsampled” — zero-inserted between symbols — before taking its FFT; the input spectrum needed here is the very same Xf already computed for the filtering step above, so it is reused rather than recomputed. This reproduces nlmsUpBlock, cmaUpBlock, ddlmsUpBlock, rdeUpBlock and dardeUpBlock’s gradient exactly, without calling those functions (their body is the time-domain sum) or building sigInBlock at all for these algorithms.

rls and dd-rls have no such shortcut: their state Sd (the inverse input-correlation matrix) is updated through a recursive matrix-inversion-lemma step that is inherently sequential and has no simple frequency-domain form, so these two algorithms still build sigInBlock and call the shared rlsUpBlock/ddrlsUpBlock, exactly as coreAdaptEqBlockTD does. static performs no adaptation at all.

coreAdaptEqBlockTD(sigIn, symbRef, SpS, H, H_, L, mu, lambdaRLS, nTaps, storeCoeff, runWL, alg, constSymb, prec, blockLen=1)

Adaptive equalizer core processing function (block-wise, time-domain).

Parameters:

• sigIn (np.array) – Input signal array.

• symbRef (np.array) – Reference symbol sequence.

• SpS (int) – Samples per symbol.

• H (np.array) – Coefficient matrix.

• H – Augmented coefficient matrix.

• L (int) – Length of the output.

• mu (float) – Step size parameter.

• lambdaRLS (float) – RLS forgetting factor.

• nTaps (int) – Number of taps.

• storeCoeff (bool) – Flag indicating whether to store coefficient matrices.

• runWL (bool) – Run widely-linear mode

• alg (str) – Equalizer algorithm.

• constSymb (np.array) – Constellation symbols.

• prec (data type) – Precision of the computations.

• blockLen (int) – Length of the processing block.

Returns:

• sigOut (np.array) – Equalized output array.

• H (np.array) – Coefficient matrix.

• H_ (np.array) – Augmented coefficient matrix.

• errSq (np.array) – Squared absolute error array.

• Hiter (np.array) – History of coefficient matrices.

realValuedDFECore(sigIn, symbRef, nTapsFF=5, nTapsFB=5, SpS=1, mu=0.0001, nTrain=1000, prec=<class 'numpy.float32'>, constSymb=None, f=None, b=None, trainingMode='data-aided', preconvIters=1)

Decision feedback equalizer (DFE) core implementation.

Parameters:

• sigIn (np.array) – Input signal to be equalized.

• symbRef (np.array) – Desired (reference) signal.

• nTapsFF (int) – Number of feedforward taps

• nTapsFB (int) – Number of feedback taps

• SpS (int) – Samples per symbol

• mu (float) – Step size

• nTrain (int) – Number of training symbols

• prec (data type) – Precision

• constSymb (np.array) – Array of constellation symbols used for symbol decisions.

• f (np.array) – Initial feedforward coeffs

• b (np.array) – Initial feedback coeffs

• trainingMode (str) – Operation mode (‘data-aided’, ‘fulltime’)

• preconvIters (int) – Number of pre-convergence iterations

Returns:

• sigOut (np.array) – Equalized output signal.

• f (np.array) – Final feedforward filter coefficients.

• b (np.array) – Final feedback filter coefficients.

References

[1] Proakis, J. G., & Salehi, M. (2008). Digital Communications (5th Edition). McGraw-Hill Education.

realValuedFFECore(sigIn, symbRef, nTaps=5, SpS=1, mu=0.0001, nTrain=1000, prec=<class 'numpy.float32'>, constSymb=None, f=None, trainingMode='data-aided', preconvIters=1)

Decision-directed feedforward equalizer (FFE) core implementation.

Parameters:

• sigIn (np.array) – Input signal to be equalized.

• symbRef (np.array) – Desired (reference) signal.

• nTaps (int) – Number of feedforward taps

• SpS (int) – Samples per symbol

• mu (float) – Step size

• nTrain (int) – Number of training symbols

• prec (data type) – Precision

• constSymb (np.array) – Constellation symbols

• f (np.array) – Initial feedforward filter coefficients

• trainingMode (str) – Operation mode (‘data-aided’, ‘fulltime’)

• preconvIters (int) – Number of pre-convergence iterations

Returns:

• sigOut (np.array) – Equalized output signal.

• f (np.array) – Final feedforward filter coefficients.

References

[1] Proakis, J. G., & Salehi, M. (2008). Digital Communications (5th Edition). McGraw-Hill Education.

volterraCore(sigIn, symbRef, order=2, SpS=1, mu=0.0001, nTrain=1000, h1=None, h2=None, h3=None, prec=<class 'numpy.float32'>, constSymb=None, trainingMode='data-aided', preconvIters=1)

Decision-directed Volterra equalizer core implementation

Parameters:

• sigIn (np.array) – Input signal to be equalized.

• symbRef (np.array) – Desired (reference) signal.

• order (int) – Volterra series order (2 for quadratic, 3 for cubic)

• SpS (int) – Samples per symbol

• mu (float) – Step size

• nTrain (int) – Number of training symbols

• h1 (np.array) – Initial linear filter coefficients.

• h2 (np.array) – Initial quadratic filter coefficients.

• h3 (np.array) – Initial cubic filter coefficients.

• prec (data type) – Precision

• constSymb (np.array) – Constellation symbols

• trainingMode (str) – Operation mode (‘data-aided’, ‘fulltime’)

• preconvIters (int) – Number of pre-convergence iterations

Returns:

• sigOut (np.array) – Equalized output signal.

• h1 (np.array) – Final linear filter coefficients.

• h2 (np.array) – Final quadratic filter coefficients.

• h3 (np.array) – Final cubic filter coefficients.

References

[1] Diniz, P. R., da Silva, E. A. B., & Netto, S. L. Adaptive Filtering: Algorithms and Practical Implementation. Springer Science & Business Media, 2010.

DSP algorithms for equalization (optic.dsp.equalization)

edc(sigIn, param)	Electronic chromatic dispersion compensation (EDC).
mimoAdaptEqualizer(sigIn[, param, symbRef])	General \(N \times N\) MIMO adaptive equalizer with several adaptive filtering algorithms available.
manakovDBP(Ei, param)	Run the Manakov SSF digital backpropagation (symmetric, dual-pol.).
dfe(sigIn, symbRef, param)	Decision feedback adaptive equalizer (DFE) for SISO receivers.
ffe(sigIn, symbRef, param)	Decision-directed feedforward adaptive equalizer (FFE) for SISO receivers.
volterra(sigIn, symbRef, param)	Decision-directed Volterra equalizer implementation up to 3rd order for SISO receivers

dfe(sigIn, symbRef, param)

Decision feedback adaptive equalizer (DFE) for SISO receivers.

Parameters:

• sigIn (np.array) – Input signal to be equalized.

• symbRef (np.array) – Desired (reference) signal.

• param (optic.utils.parameters object) –

  DFE parameters:

  • param.nTapsFF : number of feedforward taps [default: 5]

  • param.nTapsFB : number of feedback taps [default: 5]

  • param.SpS : samples per symbol [default: 1]

  • param.mu : step size [default: 0.0001]

  • param.nTrain : number of training symbols [default: 1000]

  • param.prec : precision [default: np.float32]

  • param.M : modulation order [default: 4]

  • param.constType : constellation type (‘pam’, ‘qam’, etc.) [default: ‘pam’]

  • param.f : initial feedforward coeffs [default: None]

  • param.b : initial feedback coeffs [default: None]

  • param.trainingMode : operation mode (‘data-aided’, ‘fulltime’) [default: ‘data-aided’]

  • param.preconvIters : number of pre-convergence iterations [default: 1]

Returns:

• sigOut (np.array) – Equalized output signal.

• f (np.array) – Final feedforward filter coefficients.

• b (np.array) – Final feedback filter coefficients.

Notes

• Training mode ‘data-aided’ uses the known training symbols for adaptation, while ‘fulltime’ continues to adapt using decision-directed mode even after the training phase.

• Pre-convergence iterations can help the algorithm to converge better by restarting the adaptation process after the initial training phase.

References

[1] Proakis, J. G., & Salehi, M. (2008). Digital Communications (5th Edition). McGraw-Hill Education.

edc(sigIn, param)

Electronic chromatic dispersion compensation (EDC).

Parameters:

• sigIn (np.array) – Dispersed input signal.

• param (optic.utils.parameters object) –

  Parameters of the optical channel.

  • param.L : total fiber length [km][default: 50 km]

  • param.D : chromatic dispersion parameter [ps/nm/km][default: 16 ps/nm/km]

  • param.Fc : carrier frequency [Hz] [default: 193.1e12 Hz]

  • param.Fs : sampling frequency [Hz] [default: []]

  • param.Rs : symbol rate [baud] [default: 32e9]

  • param.NfilterCoeffs : number of filter coefficients [default: []]

  • param.Nfft : FFT size [default: []]

Returns:

sigOut – Dispersion compensated output signal.

Return type:

np.array

Notes

Chromatic dispersion is a linear effect, described in the frequency domain by the all-pass transfer function \(H_{CD}(\omega) = \exp\left(j\frac{\beta_2}{2}\omega^2 L\right)\) (see optic.models.channels.linearFiberChannel), where \(\beta_2 = -D\lambda^2/(2\pi c)\) and \(L\) is the fiber length. It is compensated by the inverse filter,

\[H_{EDC}(\omega) = \exp\left(-j\frac{\beta_2}{2}\omega^2 L\right), \tag{1}\]

which is applied by blockwise FFT convolution (see optic.dsp.core.blockwiseFFTConv). The dispersion broadens the impulse response of the channel proportionally to \(|\beta_2|L\) and to the signal bandwidth; the default number of filter coefficients follows the rule

\[N_{taps} = 2\left\lceil 6.67\,|\beta_2|\,L\,R_s^2\,\frac{F_s}{R_s}\right\rceil, \tag{2}\]

where \(R_s\) is the symbol rate and \(F_s\) the sampling rate.

References

[1] S. J. Savory, “Digital coherent optical receivers: Algorithms and subsystems”, IEEE Journal on Selected Topics in Quantum Electronics, vol. 16, nº 5, p. 1164–1179, set. 2010, doi: 10.1109/JSTQE.2010.2044751.

[2] K. Kikuchi, “Fundamentals of Coherent Optical Fiber Communications”, J. Lightwave Technol., JLT, vol. 34, nº 1, p. 157–179, jan. 2016.

ffe(sigIn, symbRef, param)

Decision-directed feedforward adaptive equalizer (FFE) for SISO receivers.

Parameters:

• sigIn (np.array) – Input signal to be equalized.

• symbRef (np.array) – Desired (reference) signal.

• param (optic.utils.parameters object) –

  FFE parameters:

  • param.nTaps : number of feedforward taps [default: 5]

  • param.mu : step size [default: 0.0001]

  • param.SpS : samples per symbol [default: 1]

  • param.nTrain : number of training symbols [default: 1000]

  • param.prec : precision [default: np.float32]

  • param.M : modulation order [default: 4]

  • param.constType : constellation type (‘pam’, ‘qam’, etc.) [default: ‘pam’]

  • param.f : initial feedforward coeffs [default: None]

  • param.trainingMode : operation mode (‘data-aided’, ‘fulltime’) [default: ‘data-aided’]

  • param.preconvIters : number of pre-convergence iterations [default: 1]

Returns:

• sigOut (np.array) – Equalized output signal.

• f (np.array) – Final feedforward filter coefficients.

Notes

• Training mode ‘data-aided’ uses the known training symbols for adaptation, while ‘fulltime’ continues to adapt using decision-directed mode even after the training phase.

• Pre-convergence iterations can help the algorithm to converge better by restarting the adaptation process after the initial training phase.

References

[1] S. Haykin, “Adaptive Filter Theory,” 5th ed., Pearson, 2013.

manakovDBP(Ei, param)

Run the Manakov SSF digital backpropagation (symmetric, dual-pol.).

Parameters:

• Ei (np.array) – Input optical signal field.

• param (optic.utils.parameters object) –

  Physical/simulation parameters of the optical channel.

  • param.Ltotal : total fiber length [km][default: 400 km]

  • param.Lspan : span length [km][default: 80 km]

  • param.hz : step-size for the split-step Fourier method [km][default: 0.5 km]

  • param.alpha : fiber attenuation parameter [dB/km][default: 0.2 dB/km]

  • param.D : chromatic dispersion parameter [ps/nm/km][default: 16 ps/nm/km]

  • param.gamma : fiber nonlinear parameter [1/W/km][default: 1.3 1/W/km]

  • param.Fc : carrier frequency [Hz] [default: 193.1e12 Hz]

  • param.Fs : simulation sampling frequency [samples/second][default: None]

  • param.prec : numerical precision [default: np.complex128]

  • param.amp : ‘edfa’, ‘ideal’, or ‘None. [default:’edfa’]

  • param.maxIter : max number of iter. in the trap. integration [default: 10]

  • param.tol : convergence tol. of the trap. integration.[default: 1e-5]

  • param.nlprMethod : adap step-size based on nonl. phase rot. [default: True]

  • param.maxNlinPhaseRot : max nonl. phase rot. tolerance [rad][default: 2e-2]

  • param.prgsBar : display progress bar? bolean variable [default:True]

  • param.saveSpanN : specify the span indexes to be outputted [default:[]]

  • param.returnParameters : bool, return channel parameters [default: False]

Returns:

• Ech (np.array) – Optical signal after nonlinear backward propagation.

• param (parameter object (struct)) – Object with physical/simulation parameters used in the split-step alg.

Notes

Digital backpropagation (DBP) compensates the deterministic linear and nonlinear impairments of the fiber by numerically solving the Manakov equations (see optic.models.channels.manakovSSF) in the reverse direction, i.e. with the signs of the attenuation, dispersion and nonlinear parameters inverted,

\[\frac{\partial A_{x,y}}{\partial z} = +\frac{\alpha}{2}A_{x,y} + j\frac{\beta_2}{2}\frac{\partial^2 A_{x,y}}{\partial t^2} - j\frac{8}{9}\gamma\left(|A_x|^2 + |A_y|^2\right)A_{x,y}, \tag{1}\]

starting from the received field and propagating it back to the transmitter over the same spans, with the split-step Fourier method. Since the ASE noise added along the link is not deterministic, it cannot be removed, which limits the performance gain of DBP.

References

[1] E. Ip e J. M. Kahn, “Compensation of dispersion and nonlinear impairments using digital backpropagation”, Journal of Lightwave Technology, vol. 26, nº 20, p. 3416–3425, 2008, doi: 10.1109/JLT.2008.927791.

[2] E. Ip, “Nonlinear compensation using backpropagation for polarization-multiplexed transmission”, Journal of Lightwave Technology, vol. 28, nº 6, p. 939–951, mar. 2010, doi: 10.1109/JLT.2010.2040135.

mimoAdaptEqualizer(sigIn, param=None, symbRef=None)

General \(N \times N\) MIMO adaptive equalizer with several adaptive filtering algorithms available.

Parameters:

• sigIn (np.array) – Input signal array.

• symbRef (np.array, optional) – Reference symbol sequence synchronized to sigIn.

• param (optic.utils.parameters object, optional) –

  Parameter object containing the following attributes:

  • param.numIter : int, number of pre-convergence iterations [default: 1]

  • param.nTaps : int, number of filter taps [default: 15]

  • param.mu : float or list of floats, step size parameter(s) [default: [1e-3]]

  • param.lambdaRLS : float, RLS forgetting factor [default: 0.99]

  • param.SpS : int, samples per symbol [default: 2]

  • param.H : np.array, coefficient matrix [default: []]

  • param.L : int or list of ints, length of the output of the training section [default: []]

  • param.Hiter : list, history of coefficient matrices [default: []]

  • param.storeCoeff : bool, flag indicating whether to store coefficient matrices [default: False]

  • param.runWL: bool, flag indicating whether to run the equalizer in the widely-linear mode [default: False]

  • param.alg : str or list of strs, specifying the equalizer algorithm(s) [default: [‘nlms’]]

  • param.constType : str, constellation type [default: ‘qam’]

  • param.M : int, modulation order [default: 4]

  • param.prgsBar : bool, flag indicating whether to display progress bar [default: True]

  • param.returnResults : bool, flag indicating whether to return all results [default: False]

  • param.prec : data type, precision of the computations [default: np.complex64]

  • param.domain : str, domain of the equalizer (‘time’ or ‘freq’) [default: ‘freq’]

  • param.Nfft : int, FFT size for frequency domain equalization [default: 128]

  • param.blockSize : int, block size for time domain equalization [default: 1]

Returns:

• sigOut (np.array) – Equalized output array.

• H (np.array) – Coefficient matrix.

• errSq (np.array) – Squared absolute error array.

• Hiter (list) – History of coefficient matrices.

Notes

Algorithms available: ‘cma’, ‘rde’, ‘nlms’, ‘dd-lms’, ‘da-rde’, ‘rls’, ‘dd-rls’, ‘static’.

References

[1] P. S. R. Diniz, Adaptive Filtering: Algorithms and Practical Implementation. Springer US, 2012.

[2] S. J. Savory, “Digital coherent optical receivers: Algorithms and subsystems”, IEEE Journal on Selected Topics in Quantum Electronics, vol. 16, nº 5, p. 1164–1179, set. 2010, doi: 10.1109/JSTQE.2010.2044751.

[3] K. Kikuchi, “Fundamentals of Coherent Optical Fiber Communications”, J. Lightwave Technol., JLT, vol. 34, nº 1, p. 157–179, jan. 2016.

[4] E. P. Da Silva e D. Zibar, “Widely Linear Equalization for IQ Imbalance and Skew Compensation in Optical Coherent Receivers”, Journal of Lightwave Technology, vol. 34, nº 15, p. 3577–3586, ago. 2016, doi: 10.1109/JLT.2016.2577716.

volterra(sigIn, symbRef, param)

Decision-directed Volterra equalizer implementation up to 3rd order for SISO receivers

Parameters:

• sigIn (np.array) – Input signal to be equalized.

• symbRef (np.array) – Desired (reference) signal.

• param (optic.utils.parameters object) –

  Volterra equalizer parameters:

  • param.n1Taps : number of taps of linear part [default: 5]

  • param.n2Taps : number of taps of quadratic part [default: 3]

  • param.n3Taps : number of taps of cubic part [default: 2]

  • param.h : list of initial filter coefficients [default: None]

  • param.SpS : samples per symbol [default: 1]

  • param.mu : step size [default: 0.001]

  • param.nTrain : number of training symbols [default: 1000]

  • param.order : Volterra series order (2 for quadratic) [default: 2]

  • param.prec : precision [default: np.float32]

  • param.M : modulation order [default: 4]

  • param.constType : constellation type (‘pam’, ‘qam’, etc.) [default: ‘pam’]

  • param.trainingMode : operation mode (‘data-aided’, ‘fulltime’) [default: ‘data-aided’]

  • param.preconvIters : number of pre-convergence iterations [default: 1]

Returns:

• sigOut (np.array) – Equalized output signal.

• h (list of np.array) – Final Volterra filter coefficients [h1, h2, h3].

Notes

• Training mode ‘data-aided’ uses the known training symbols for adaptation, while ‘fulltime’ continues to adapt using decision-directed mode even after the training phase.

• Pre-convergence iterations can help the algorithm to converge better by restarting the adaptation process after the initial training phase.

References

[1] Diniz, P. R., da Silva, E. A. B., & Netto, S. L. (2010). Adaptive Filtering: Algorithms and Practical Implementation. Springer Science & Business Media.

DSP algorithms for carrier phase and frequency recovery (optic.dsp.carrierRecovery)

bps	Blind phase search (BPS) algorithm
ddpll	Decision-directed Phase-locked Loop (DDPLL) algorithm
viterbi	Viterbi & Viterbi carrier phase recovery algorithm.
fourthPowerFOE	Estimate the frequency offset (FO) with the 4th-power method.
cpr	Carrier phase recovery function (CPR)

bps(sigIn, N, constSymb, B)

Blind phase search (BPS) algorithm

Parameters:

• sigIn (complex-valued np.array) – Received constellation symbols.

• N (int) – Half of the 2*N+1 average window.

• constSymb (complex-valued np.array) – Complex-valued constellation.

• B (int) – number of test phases.

Returns:

phaseEst – Time-varying estimated phase-shifts.

Return type:

real-valued np.array

Notes

The blind phase search (BPS) algorithm is a feedforward carrier phase estimator that works with any constellation. Due to the rotational symmetry of square QAM constellations, the phase is only identifiable modulo \(\pi/2\), so \(B\) test phases are distributed within this interval,

\[\varphi_b = \frac{b}{B}\,\frac{\pi}{2}, \qquad b = 0, 1, \ldots, B-1. \tag{1}\]

For each received symbol \(y[k]\) and test phase, the squared distance between the rotated symbol and the closest constellation point is computed,

\[d_b[k] = \min_{x \in \mathcal{X}}\left|y[k]e^{j\varphi_b} - x\right|^2. \tag{2}\]

To reduce the influence of the noise, the distances are summed over a window of \(2N+1\) consecutive symbols, and the test phase with the smallest sum is selected,

\[\hat{\varphi}[k] = \arg\min_{\varphi_b} \sum_{n=-N}^{N} d_b[k+n]. \tag{3}\]

The phase-corrected symbols are \(y[k]e^{j\hat{\varphi}[k]}\). The \(\pi/2\) ambiguity of the estimates is removed later by phase unwrapping (see cpr).

References

[1] T. Pfau, S. Hoffmann, e R. Noé, “Hardware-efficient coherent digital receiver concept with feedforward carrier recovery for M-QAM constellations”, Journal of Lightwave Technology, vol. 27, nº 8, p. 989–999, 2009, doi: 10.1109/JLT.2008.2010511.

cpr(sigIn, param=None, symbTx=None)

Carrier phase recovery function (CPR)

Parameters:

• sigIn (complex-valued np.array) – received constellation symbols.

• param (optic.utils.parameter object, optional) –

  Configuration parameters [default: None].

  • param.alg : CPR algorithm to be used [‘bps’, ‘bpsGPU’, ‘ddpll’, or ‘viterbi’] [default: ‘bps’].

  • param.shapingFactor : shaping factor, for probabilistic shaped QAM with MB dististribution.[default: 0]

  • param.constType : constellation type [‘qam’ or ‘psk’]. [default: ‘qam’]

  • param.M : constellation order. [default: 4]

  • param.returnPhases : whether to return the estimated phase shifts along with the output signal. [default: False]

  • param.runFOE : whether to run the Mth-power frequency offset estimation and compensation before CPR. [default: True]

  BPS params:

  • param.N : length of BPS the moving average window. [default: 35]

  • param.B : number of BPS test phases. [default: 64]

  DDPLL params:

  • param.tau1 : DDPLL loop filter param. 1. [default: 1/2*pi*10e6]

  • param.tau2 : DDPLL loop filter param. 2. [default: 1/2*pi*10e6]

  • param.Kv : DDPLL loop filter gain. [default: 0.1]

  • param.Ts : symbol period. [default: 1/32e9]

  • param.pilotInd : indexes of pilot-symbol locations.

  Viterbi params:

  • param.N : length of the moving average window. [default: 35]

• symbTx (complex-valued np.array, optional) – Transmitted symbol sequence. [default: None]

Returns:

• sigOut (complex-valued np.array) – Phase-compensated signal.

• phaseEst (real-valued np.array) – Time-varying estimated phase-shifts.

References

[1] T. Pfau, S. Hoffmann, e R. Noé, “Hardware-efficient coherent digital receiver concept with feedforward carrier recovery for M-QAM constellations”, Journal of Lightwave Technology, vol. 27, nº 8, p. 989–999, 2009, doi: 10.1109/JLT.2008.2010511.

[2] S. J. Savory, “Digital coherent optical receivers: Algorithms and subsystems”, IEEE Journal on Selected Topics in Quantum Electronics, vol. 16, nº 5, p. 1164–1179, set. 2010, doi: 10.1109/JSTQE.2010.2044751.

[3] H. Meyer, Digital Communication Receivers: Synchronization, Channel estimation, and Signal Processing, Wiley 1998. Section 5.8 and 5.9.

ddpll(sigIn, Ts, Kv, tau1, tau2, constSymb, symbTx, pilotInd)

Decision-directed Phase-locked Loop (DDPLL) algorithm

Parameters:

• sigIn (complex-valued np.array) – Received constellation symbols.

• Ts (float scalar) – Symbol period.

• Kv (float scalar) – Loop filter gain.

• tau1 (float scalar) – Loop filter parameter 1.

• tau2 (float scalar) – Loop filter parameter 2.

• constSymb (complex-valued np.array) – Complex-valued ideal constellation symbols.

• symbTx (complex-valued np.array) – Transmitted symbol sequence.

• pilotInd (int np.array) – Indexes of pilot-symbol locations.

Returns:

phaseEst – Time-varying estimated phase-shifts.

Return type:

real-valued np.array

Notes

The decision-directed phase-locked loop (DD-PLL) tracks the carrier phase symbol by symbol. The received symbol \(y[k]\) is first corrected by the current phase estimate \(\hat{\theta}[k]\), and a phase error signal is generated by comparing it with the decided symbol \(\hat{x}[k]\) (or with the known pilot symbol, at pilot positions),

\[u_d[k] = \mathrm{Im}\left\{y[k]e^{j\hat{\theta}[k]}\,\hat{x}^*[k]\right\} \approx |\hat{x}[k]|^2 \sin\left(\theta[k] + \hat{\theta}[k]\right), \tag{1}\]

which, for small errors, is proportional to the residual phase error. The error signal is smoothed by a second-order (proportional-integral) loop filter,

\[u_f[k] = u_f[k-1] + b_1 u_d[k-1] + b_2 u_d[k], \tag{2}\]

with \(b_{1,2} = \frac{T_s}{2\tau_1}\left[1 \mp \cot\left(\frac{T_s}{2\tau_2}\right)\right]\), given by the loop filter time constants \(\tau_1\) and \(\tau_2\), and the phase estimate for the next symbol is updated as

\[\hat{\theta}[k+1] = \hat{\theta}[k] - K_v u_f[k], \tag{3}\]

where \(K_v\) is the loop gain.

References

[1] H. Meyer, Digital Communication Receivers: Synchronization, Channel estimation, and Signal Processing, Wiley 1998. Section 5.8 and 5.9.

fourthPowerFOE(sigIn, Fs, M=4)

Estimate the frequency offset (FO) with the 4th-power method.

Parameters:

• sigIn (np.array) – Input signal.

• Fs (float) – Sampling frequency.

• M (int, optional) – M-th power order. Default is 4.

Returns:

• The output signal after applying frequency offset correction.

• The estimated frequency offset.

Return type:

np.array, float

Notes

A frequency offset \(\Delta f\) between the transmitter laser and the local oscillator rotates the received symbols as \(y[k] = x[k]e^{j2\pi\Delta f kT_s}\). Raising the signal to the \(M\)-th power (\(M = 4\) for QPSK and square QAM) largely removes the modulation, producing a strong spectral line at \(M\Delta f\). The frequency offset is estimated from the location of the peak of the spectrum of \(y^M[k]\),

\[\Delta\hat{f} = \frac{1}{M}\arg\max_{f} \left|\mathrm{DFT}\left\{y^M[k]\right\}(f)\right|, \tag{1}\]

and compensated as

\[y_c[k] = y[k]\,e^{-j2\pi\Delta\hat{f}\,kT_s}. \tag{2}\]

The estimation range is \(|\Delta f| < F_s/(2M)\), and the resolution is \(F_s/(MN)\), where \(F_s\) is the sampling rate and \(N\) the number of samples.

References

[1] S. J. Savory, “Digital coherent optical receivers: Algorithms and subsystems”, IEEE Journal on Selected Topics in Quantum Electronics, vol. 16, nº 5, p. 1164–1179, set. 2010, doi: 10.1109/JSTQE.2010.2044751.

viterbi(sigIn, N=35, M=4)

Viterbi & Viterbi carrier phase recovery algorithm.

Parameters:

• sigIn (np.array) – Input signal.

• N (int, optional) – Size of the moving average window.

• M (int, optional) – M-th power order.

Returns:

Estimated phase error.

Return type:

np.array, float

Notes

The Viterbi & Viterbi algorithm is a feedforward carrier phase estimator for constellations with \(M\)-fold rotational symmetry. Raising the received symbols \(y[k] = x[k]e^{j\theta[k]} + n[k]\) to the \(M\)-th power removes the phase modulation, since \(x^M[k]\) has a constant argument (for QPSK with points at \(\pm\pi/4\) and \(\pm 3\pi/4\), \(x^4 = -1\)), leaving \(M\theta[k]\). The noise is reduced by averaging over a window of \(N\) symbols, and the phase correction returned is

\[\hat{\varphi}[k] = -\frac{1}{M}\arg\left\{\sum_{n} y^M[k+n]\right\} - \frac{\pi}{4}, \tag{1}\]

which is unwrapped with period \(2\pi/M\) along the sequence. The corrected symbols \(y[k]e^{j\hat{\varphi}[k]}\) are aligned with the constellation up to the \(2\pi/M\) ambiguity inherent to its rotational symmetry. For square QAM constellations, the fourth power is used as well, since the average of \(x^4\) is also a negative real number (e.g. \(-0.68E_s^2\) for 16-QAM).

References

[1] S. J. Savory, “Digital coherent optical receivers: Algorithms and subsystems”, IEEE Journal on Selected Topics in Quantum Electronics, vol. 16, nº 5, p. 1164–1179, set. 2010, doi: 10.1109/JSTQE.2010.2044751.

Functions adapted to run with GPU (CuPy) processing (optic.dsp.carrierRecoveryGPU)

bpsGPU

Blind phase search (BPS) algorithm

bpsGPU(sigIn, N, constSymb, B)

Blind phase search (BPS) algorithm

Parameters:

• sigIn (complex-valued np.array) – Received constellation symbols.

• N (int) – Half of the 2*N+1 average window.

• constSymb (complex-valued np.array) – Complex-valued constellation.

• B (int) – number of test phases.

Returns:

θ – Time-varying estimated phase-shifts.

Return type:

real-valued np.array

Notes

This is the GPU (CuPy) implementation of the blind phase search algorithm: for \(B\) test phases \(\varphi_b = \frac{b}{B}\frac{\pi}{2}\), the squared distances between the rotated symbols \(y[k]e^{j\varphi_b}\) and the closest constellation points are summed over a window of \(2N+1\) symbols, and the test phase with the smallest sum is selected. See optic.dsp.carrierRecovery.bps for the description of the algorithm.

References

[1] T. Pfau, S. Hoffmann, e R. Noé, “Hardware-efficient coherent digital receiver concept with feedforward carrier recovery for M-QAM constellations”, Journal of Lightwave Technology, vol. 27, nº 8, p. 989–999, 2009, doi: 10.1109/JLT.2008.2010511.

DSP algorithms for clock and timming recovery (optic.dsp.clockRecovery)

gardnerTED	Calculate the timing error using the Gardner timing error detector.
gardnerTEDnyquist	Modified Gardner timing error detector for Nyquist pulses.
interpolator	Perform cubic interpolation using the Farrow structure.
gardnerClockRecovery	Perform clock recovery using Gardner's algorithm with a loop PI filter.
calcClockDrift	Calculate the clock drift in parts per million (ppm) from t_nco values.

calcClockDrift(t_nco_values)

Calculate the clock drift in parts per million (ppm) from t_nco values.

Parameters:

t_nco_values (np.array) – An array containing the relative time delay values provided to the NCO.

Returns:

The clock deviation in parts per million (ppm).

Return type:

float

gardnerClockRecovery(sigIn, param=None)

Perform clock recovery using Gardner’s algorithm with a loop PI filter.

Parameters:

• sigIn (numpy.np.array) – Input array representing the received signal.

• param (core.parameter) –

  Clock recovery parameters:

  • param.kp : Proportional gain for the loop filter. [default: 1e-3]

  • param.ki : Integral gain for the loop filter. [default: 1e-6]

  • param.isNyquist : is the pulse shape a Nyquist pulse? [default: True]

  • param.returnTiming : return estimated timing values. [default: False]

  • param.lpad : length of zero padding at the end of the input vector. [default: 1]

  • param.maxPPM : maximum clock rate expected deviation in PPM. [default: 500]

Returns:

Tuple containing the recovered signal (sigOut) and the timing values.

Return type:

tuple

Notes

The clock recovery is a feedback loop composed of an interpolator, a timing error detector (TED), a loop filter and a numerically controlled oscillator (NCO). At each output sample, the interpolator computes the signal at the fractional delay \(\tau\) given by the NCO (see interpolator). Once per symbol, the TED computes the timing error \(e_k\) (see gardnerTED and gardnerTEDnyquist), which is filtered by a proportional-integral (PI) loop filter,

\[v_k = k_p e_k + k_i\sum_{l \le k} e_l, \tag{1}\]

and the NCO updates the fractional delay as

\[\tau \leftarrow \tau - v_k. \tag{2}\]

Whenever \(\tau\) crosses \(\pm 1\), it is wrapped back into \([-1, 1]\) and one input sample is skipped or repeated, which allows the loop to track a clock frequency offset (drift) between the transmitter and the receiver. The integral term removes the steady-state timing error caused by such a drift.

gardnerTED(x)

Calculate the timing error using the Gardner timing error detector.

Parameters:

x (numpy.np.array) – Input array of size 3 representing a segment of the received signal.

Returns:

Gardner timing error detector (TED) value.

Return type:

float

Notes

The Gardner timing error detector (TED) works with two samples per symbol. Given the samples at the current and the previous symbol instants, \(y_k\) and \(y_{k-1}\), and the sample halfway between them, \(y_{k-1/2}\), the timing error is

\[e_k = \mathrm{Re}\left\{y_{k-1/2}^*\left(y_k - y_{k-1}\right)\right\}. \tag{1}\]

When a symbol transition occurs, the midpoint sample is close to zero if the sampling instants are correct; otherwise, its sign relative to the slope of the transition indicates whether the sampling is early or late. The TED does not depend on symbol decisions and is insensitive to the carrier phase.

gardnerTEDnyquist(x)

Modified Gardner timing error detector for Nyquist pulses.

Parameters:

x (numpy.np.array) – Input array of size 3 representing a segment of the received signal.

Returns:

Gardner timing error detector (TED) value.

Return type:

float

Notes

For Nyquist pulses with small roll-off factors, the S-curve of the classical Gardner TED (see gardnerTED) vanishes. This modified detector uses the power of the samples instead,

\[e_k = |y_{k-1/2}|^2\left(|y_{k-1}|^2 - |y_k|^2\right), \tag{1}\]

where \(y_{k-1}\) and \(y_k\) are consecutive symbol-spaced samples and \(y_{k-1/2}\) is the sample halfway between them.

interpolator(x, t)

Perform cubic interpolation using the Farrow structure.

Parameters:

• x (numpy.np.array) – Input array of size 4 representing the values for cubic interpolation.

• t (float) – Interpolation parameter.

Returns:

y – Interpolated signal value.

Return type:

float

Notes

The interpolator computes the value of the signal at a fractional position \(t \in [-1, 1]\) with respect to the sample \(x[2]\), by cubic Lagrange interpolation of the four samples \(x[0], \ldots, x[3]\), located at the positions \(-2, -1, 0\) and \(1\),

\[y(t) = \sum_{i=0}^{3} x[i]\,\ell_i(t), \qquad \ell_i(t) = \prod_{j \neq i}\frac{t - t_j}{t_i - t_j}, \tag{1}\]

which is implemented with the Farrow structure, i.e. with the polynomial coefficients

\begin{equation} \begin{aligned} \ell_0(t) &= -\tfrac{1}{6}t^3 + \tfrac{1}{6}t, & \ell_1(t) &= \tfrac{1}{2}t^3 + \tfrac{1}{2}t^2 - t, \\ \ell_2(t) &= -\tfrac{1}{2}t^3 - t^2 + \tfrac{1}{2}t + 1, & \ell_3(t) &= \tfrac{1}{6}t^3 + \tfrac{1}{2}t^2 + \tfrac{1}{3}t. \end{aligned} \tag{2} \end{equation}

In this form, the fractional delay \(t\) can be changed at every output sample without recomputing filter coefficients.

Signal synchronization functions (optic.dsp.synchronization)

syncDataSequences(rx, tx, param)

Synchronize data sequences with a given reference sequence.

syncDataSequences(rx, tx, param)

Synchronize data sequences with a given reference sequence.

Parameters:

• rx (np.array) – Received signal.

• tx (np.array) – Reference signal (or symbols).

• param (optic.utils.parameters object, optional) –

  Parameters of the synchronization process.

  • param.SpS : samples per symbol of the received signal. [default: 1]

  • param.reference : type of reference used for synchronization (‘signal’,’symbols’) [default: ‘signal’]

  • param.syncMode : use either the real part or the amplitude of the signal to syncronize [detault: ‘amp’]

  • param.pulseType : type of pulse shaping filter. [default: ‘rrc’]

  • param.rollOff : rolloff of RRC filter. [default: 0.01]

  • param.nFilterTaps : number of filter coefficients. [default: 1024]

  • param.constType : type of constellation [Default: ‘pam’]

  • param.M : modulation order. [default: 4]

Returns:

• tx_ (np.array) – Synchronized transmitted signal.

• symb (np.array) – Detected transmitted symbols.

Notes

• Signals rx and tx must have the same number of columns (modes).

• If param.reference is set to ‘signal’, rx and tx should be sampled at the same rate.