Compensators#

class DCCorrector(value: float = 0.0, *, axis: int = -1, name: str = 'mean_corrector')[source]#

Remove the DC component of a signal and set it to a target value.

Signal Model#

The empirical mean \(\mu_x\) of the input is subtracted and replaced by the target value \(\alpha\):

\[y[n] = x[n] - \mu_x + \alpha, \qquad \mu_x = \frac{1}{N} \sum_{n=0}^{N-1} x[n]\]

so that the output satisfies the DC constraint \(\frac{1}{N}\sum_n y[n] = \alpha\). This is the digital counterpart of the DC-offset cancellation performed at the output of a direct-conversion (zero-IF) receiver, where the self-mixing of the local oscillator adds a constant term to the baseband signal.

Axes: declared axis – the mean \(\mu_x\) is computed along axis (default -1, the canonical serial layout (..., N)) and broadcast over the remaining axes. The estimand is per path (D49): a DC offset belongs to a converter, so the default gives each row of a batch its own offset – the old default of 0 averaged across the batch, which is the silent coupling D51 forbids.

param value:

Target mean value \(\alpha\) of the output. Default is 0.0, i.e. a zero-mean output.

type value:

float, optional

param axis:

Axis along which \(\mu_x\) is computed. Default is -1.

type axis:

int, optional, keyword-only

param name:

Name of the corrector instance. Default is "mean_corrector".

type name:

str, optional, keyword-only

References

B. Razavi, “Design considerations for direct-conversion receivers,” IEEE Transactions on Circuits and Systems II, vol. 44, no. 6, pp. 428-435, 1997 (DC offset in zero-IF receivers).

Examples

>>> x = np.array([1.0, 2.0, 3.0, 6.0])   # mean 3.0
>>> print(DCCorrector()(x))
[-2. -1.  0.  3.]
>>> print(DCCorrector(2.0)(x))
[0. 1. 2. 5.]
class Normalizer(method: Literal['amp', 'abs', 'var', 'max'] = 'amp', *, name: str = 'signal_amplifier', value: float = 1.0)[source]#

Scale a signal by a gain derived from one of its own statistics.

Signal Model#

The block is an Amplifier whose gain \(\alpha\) is measured on the input instead of being configured:

\[y[n] = \alpha \, x[n]\]

The gain depends on the selected method and on the target value \(v\):

  • ‘amp’ – fixed gain, no measurement:

    \[\alpha = v\]
  • ‘abs’ – the largest modulus becomes \(v\):

    \[\alpha = \frac{v}{\max_n |x[n]|}\]
  • ‘var’ – the output variance becomes \(v\), with \(\sigma_x^2\) the empirical variance of the input:

    \[\alpha = \sqrt{\frac{v}{\sigma_x^2}}\]
  • ‘max’ – the largest rectangular excursion becomes \(v\) (the useful one before a DAC or a clipping stage):

    \[\alpha = \frac{v}{\max\left(\max_n |\Re e(x[n])|,\; \max_n |\Im m(x[n])|\right)}\]

Axes: element-wise – the gain \(\alpha\) is a single scalar measured over the whole array in prepare(), then applied pointwise. On a batched input the rows are therefore pooled into one estimate – normalize per trial by looping when each row must reach the target on its own.

Unlike Amplifier, the gain is not a parameter here: it is measured, so it is exposed as gain_ (decision D23) and the inherited gain argument is removed from the constructor. The first positional argument is therefore method, which makes Normalizer('max') mean what it reads.

param method:

Statistic used to derive the gain \(\alpha\). Default is 'amp' (constant gain). It is the first positional argument.

type method:

{‘amp’, ‘abs’, ‘var’, ‘max’}, optional

param value:

Target value \(v\) of the selected statistic. Must be strictly positive. Default is 1.0.

type value:

float, optional, keyword-only

param name:

Name of the instance. Default is "signal_amplifier".

type name:

str, optional, keyword-only

gain_#

Gain \(\alpha\) actually applied (data-dependent, hence the trailing underscore, decision D23). Measured from the input in prepare() for every method except 'amp', and overwritten at each call.

Type:

float

raises ValueError:

If value is not strictly positive, or if method is not one of 'amp', 'abs', 'var', 'max'.

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 3 (average energy and power of a constellation, and their normalization).

Examples

>>> x = np.array([1.0, 2.0, 3.0, 4.0])
>>> print(Normalizer('max', value=2.0)(x))
[0.5 1.  1.5 2. ]
>>> y = Normalizer(method='var')(x)          # unit-variance output
>>> print(np.round(y, 4), round(float(np.var(y)), 6))
[0.8944 1.7889 2.6833 3.5777] 1.0
class BlindIQCompensator(should_fit: bool = True, *, coef: float = 1, name: str = 'iq_compensator')[source]#

Blind IQ-imbalance compensator by whitening of the real covariance matrix.

Signal Model#

A receiver with gain and quadrature imbalance delivers an improper (non-circular) signal: the in-phase and quadrature components no longer have equal power, nor are they uncorrelated. The compensator applies the widely-linear transformation

\[y[n] = c \left(\alpha \, \Re e(x[n]) + \beta \, \Im m(x[n])\right), \qquad c = \sqrt{\frac{P}{2}}\]

where the complex coefficients \(\alpha\) and \(\beta\) are read column-wise from the \(2 \times 2\) whitening matrix \(\mathbf{M}\), i.e. \(\alpha = M_{00} + j M_{10}\) and \(\beta = M_{01} + j M_{11}\). Writing \(\mathbf{x}[n] = [\Re e(x[n]),\, \Im m(x[n])]^T\), the estimator eigen-decomposes the empirical covariance matrix

\[\widehat{\mathbf{R}} = \frac{1}{N} \sum_{n=0}^{N-1} \mathbf{x}[n] \mathbf{x}^T[n] = \mathbf{U} \boldsymbol\Lambda \mathbf{U}^T, \qquad \mathbf{M} = \boldsymbol\Lambda^{-1/2} \mathbf{U}^T\]

so that \(\mathbf{M} \widehat{\mathbf{R}} \mathbf{M}^T = \mathbf{I}_2\). The compensated signal is therefore second-order circular with the requested mean power \(P\):

\[\mathbb{E}\left[\Re e^2(y[n])\right] = \mathbb{E}\left[\Im m^2(y[n])\right] = \frac{P}{2}, \qquad \mathbb{E}\left[\Re e(y[n]) \Im m(y[n])\right] = 0\]

Whitening only constrains the second-order statistics: the residual \(2 \times 2\) rotation/reflection (and hence a possible I/Q swap or sign flip) is left to a downstream phase compensator. The estimator is therefore not the Gram-Schmidt orthogonalization procedure (GSOP) of Fatadin et al., which whitens by successive projection and keeps the in-phase axis fixed.

Axes: declared axis – accepts (N,) or (..., P, N). The estimand is per path (D49): an IQ imbalance belongs to a receiver, so each path gets its own \(2 \times 2\) covariance and its own whitening. Pooling them, which is what stacking the real and imaginary parts of a multi-path array used to do, mixes the paths and returns a signal worse than the one it was given.

param should_fit:

If True (default), \(\alpha\) and \(\beta\) are re-estimated from the input at every call (per-block regime, D22). If False, the coefficients of the last fit are reused.

type should_fit:

bool, optional

param coef:

Target mean power \(P = \mathbb{E}[|y[n]|^2]\) of the compensated signal. Default is 1.

type coef:

float, optional, keyword-only

param name:

Name of the compensator instance. Default is "iq_compensator".

type name:

str, optional, keyword-only

alpha_#

Estimated coefficient \(\alpha\) applied to the in-phase component (data-dependent, hence the trailing underscore, decision D23; initialized to 1).

Type:

complex or np.ndarray

beta_#

Estimated coefficient \(\beta\) applied to the quadrature component (data-dependent, hence the trailing underscore, decision D23; initialized to 0).

Type:

complex or np.ndarray

References

  • I. Fatadin, S. J. Savory, D. Ives, “Compensation of quadrature imbalance in an optical QPSK coherent receiver,” IEEE Photonics Technology Letters, vol. 20, no. 20, pp. 1733-1735, 2008, doi: 10.1109/LPT.2008.2004630.

  • P. J. Schreier, L. L. Scharf, Statistical Signal Processing of Complex-Valued Data, Cambridge University Press, 2010, Chapter 2 (properness, circularity and widely-linear transformations).

Examples

>>> rng = np.random.default_rng(0)
>>> s = (rng.integers(0, 2, 2000) * 2 - 1) + 1j * (rng.integers(0, 2, 2000) * 2 - 1)
>>> x = 2 * s.real + 1j * (0.5 * s.real + s.imag)   # gain + quadrature imbalance
>>> y = BlindIQCompensator(coef=1.0)(x)
>>> print(round(float(np.var(y.real)), 3), round(float(np.var(y.imag)), 3))
0.499 0.5
>>> print(round(float(abs(np.mean(y.real * y.imag))), 3), round(float(np.mean(np.abs(y)**2)), 3))
0.0 1.0
class BlindCFOCompensator(w0_init: float = 0.0, *, N_iter: int = 3, should_fit: bool = True, grid_search: bool = True, save_history: bool = False, method: Literal['grad', 'newton'] = 'newton', step_size: float = 1e-08, grid_search_tuple: tuple[float, float, float] = (-0.1, 0.1, 0.0001), name: str = 'cfo_compensator')[source]#

Blind carrier frequency offset compensator (fourth-power periodogram).

Signal Model#

The input carries a residual carrier frequency offset \(\omega_0\) (in rad/sample), which the block removes:

\[x[n] = s[n] \, e^{j \omega_0 n}, \qquad y[n] = x[n] \, e^{-j \widehat{\omega}_0 n}\]

For a QPSK-like constellation, raising the signal to the fourth power wipes out the modulation (\(s^4[n]\) is constant up to a sign) and leaves a pure tone at \(4\omega_0\). The offset is therefore obtained by maximizing the periodogram of \(x^4[n]\):

\[\widehat{\omega}_0 = \frac{1}{4} \arg \max_{\omega} \frac{1}{N} \left| \sum_{n=0}^{N-1} x^4[n] \, e^{-j\omega n} \right|^2\]

The maximization is performed in two stages: an optional coarse grid search over \(\omega_0\), followed by N_iter refinements of the criterion by Newton’s method (or by plain gradient ascent). The Newton step uses the exact first and second derivatives of the periodogram with respect to \(\omega\).

Axes: declared axis – accepts (N,) or one polarization pair (2, N); the correction \(e^{-j\widehat{\omega}_0 n}\) uses the sample index \(n\) along the last axis. The estimand is shared (D49): a frequency offset is one laser beating against one local oscillator, so the periodogram is accumulated over the pair and a single \(\widehat{\omega}_0\) is applied to both rows. A wider batch is refused with ShapeError – one scalar estimate cannot serve independent trials (D51); loop over them, one compensator each. That is not only tidier than one estimate per path – it is the better estimator, since it sees all the data.

param w0_init:

Initial value of \(\omega_0\) in rad/sample, used as the starting point when grid_search is False. Default is 0.0.

type w0_init:

float, optional

param N_iter:

Number \(N_{iter}\) of local refinement steps. Default is 3.

type N_iter:

int, optional, keyword-only

param should_fit:

If True (default), \(\widehat{\omega}_0\) is re-estimated at every call (per-block regime, D22); otherwise the value of the last fit is reused.

type should_fit:

bool, optional, keyword-only

param grid_search:

If True (default), a coarse grid search initializes the local refinement.

type grid_search:

bool, optional, keyword-only

param save_history:

If True, every intermediate value of \(\omega_0\) is appended to history_. Default is False.

type save_history:

bool, optional, keyword-only

param method:

Local optimizer. Default is "newton".

type method:

{“grad”, “newton”}, optional, keyword-only

param step_size:

Gradient step used when method="grad". Default is 1e-8.

type step_size:

float, optional, keyword-only

param grid_search_tuple:

(start, stop, step) of the grid over \(\omega_0\), in rad/sample. Default is (-0.1, 0.1, 0.0001).

type grid_search_tuple:

tuple, optional, keyword-only

param name:

Name of the compensator instance. Default is "cfo_compensator".

type name:

str, optional, keyword-only

w0_#

Estimated frequency offset \(\widehat{\omega}_0\) in rad/sample (data-dependent, hence the trailing underscore, decision D23; None before the first fit).

Type:

float

history_#

Successive iterates of \(\omega_0\) when save_history=True.

Type:

list

raises NotFittedError:

If forward is called with should_fit=False before any fit (decision D23).

References

  • A. J. Viterbi, A. M. Viterbi, “Nonlinear estimation of PSK-modulated carrier phase with application to burst digital transmission,” IEEE Transactions on Information Theory, vol. 29, no. 4, pp. 543-551, 1983 (the fourth-power nonlinearity).

  • D. C. Rife, R. R. Boorstyn, “Single-tone parameter estimation from discrete-time observations,” IEEE Transactions on Information Theory, vol. 20, no. 5, pp. 591-598, 1974 (periodogram maximization and its local refinement).

  • J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 5 (carrier and symbol synchronization).

Examples

>>> alphabet = np.array([1+1j, 1-1j, -1+1j, -1-1j]) / np.sqrt(2)
>>> rng = np.random.default_rng(2)
>>> s = alphabet[rng.integers(0, 4, 500)]
>>> x = s * np.exp(1j * 0.01 * np.arange(500))   # omega_0 = 0.01 rad/sample
>>> compensator = BlindCFOCompensator()
>>> y = compensator(x)
>>> print(round(float(compensator.w0_), 5))
0.01
class BlindPhaseCompensation(alphabet: ndarray, *, theta0: float = 0.0, should_fit: bool = True, name: str = 'phase correction')[source]#

Blind phase compensator minimizing the error vector magnitude (EVM).

Signal Model#

The input carries an unknown constant phase offset, removed by

\[y[n] = x[n] \, e^{j \widehat{\theta}}\]

No reference is available (blind regime): the estimator uses the constellation itself as a soft reference and minimizes the squared distance to the nearest alphabet point, i.e. the EVM,

\[\widehat{\theta} = \arg \min_{\theta} \sum_{n=0}^{N-1} \left| x[n] e^{j\theta} - \mathcal{P}_{\mathcal{A}}\!\left(x[n] e^{j\theta}\right) \right|^2\]

where \(\mathcal{P}_{\mathcal{A}}\) is the hard projection onto the alphabet \(\mathcal{A}\). This is the decision-directed carrier phase estimator; the criterion is solved by a nonlinear least-squares search (scipy.optimize.least_squares) started at \(\theta_0\). Being decision-directed, it inherits the \(2\pi/M\) phase ambiguity of the constellation, so \(\theta_0\) must lie in the correct basin.

Axes: declared axis – accepts (N,) or (..., P, N). The estimand is per path (D49): after a butterfly equalizer each output carries its own residual rotation, so one angle is fitted per path and theta_ holds one value each. A phase common to the paths, laser noise before any equalizer, is the shared case – fit it on the flattened signal.

param alphabet:

Modulation alphabet \(\mathcal{A}\), a 1D array of constellation symbols.

type alphabet:

np.ndarray

param theta0:

Starting point \(\theta_0\) of the optimizer, in radians. Default is 0.

type theta0:

float, optional, keyword-only

param should_fit:

If True (default), \(\widehat{\theta}\) is re-estimated at every call (per-block regime, D22); otherwise the value of the last fit is reused.

type should_fit:

bool, optional, keyword-only

param name:

Name of the processor. Default is "phase correction".

type name:

str, optional, keyword-only

theta_#

Estimated phase \(\widehat{\theta}\) in radians (data-dependent, hence the trailing underscore, decision D23).

Type:

float or np.ndarray

raises NotFittedError:

If forward is called with should_fit=False before any fit (decision D23).

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 5 (decision-directed carrier phase estimation).

Examples

>>> alphabet = np.array([1+1j, 1-1j, -1+1j, -1-1j]) / np.sqrt(2)
>>> rng = np.random.default_rng(1)
>>> s = alphabet[rng.integers(0, 4, 200)]
>>> compensator = BlindPhaseCompensation(alphabet)
>>> y = compensator(s * np.exp(-1j * 0.3))    # a -0.3 rad rotation
>>> print(round(float(compensator.theta_), 4), bool(np.allclose(y, s)))
0.3 True
class BlindPhaseSearchCompensator(alphabet: Any, *, n_test: int = 32, half_window: int = 16, name: str = 'bps')[source]#

Feedforward carrier phase recovery by blind phase search (BPS).

Signal Model#

The input carries a slowly varying phase \(\phi[n]\) – laser phase noise, typically – on top of the constellation:

\[x[n] = s[n] \, e^{j \phi[n]} + b[n]\]

For each of \(B\) test phases \(\varphi_i = i \, \frac{\pi/2}{B}\), the input is rotated, decided against the alphabet, and the squared decision distance is averaged over a sliding window of \(2W + 1\) symbols. The estimate at each symbol is the test phase of smallest averaged distance, then unwrapped modulo \(\pi/2\) so the estimate can follow the walk beyond one quadrant.

The \(\pi/2\) search range is the symmetry of a square QAM (and of QPSK): the recovered constellation is aligned up to a quadrant. That residual quadrant ambiguity is not resolved here – it is the same ambiguity every blind phase estimator leaves, and a system resolves it with differential encoding or a known field. The estimated phase trajectory is exposed as phase_ (D23).

Axes: axis -1 – the window slides along the last axis; leading axes (the two rows of a polarization pair) are processed independently, each with its own trajectory.

param alphabet:

The constellation the decision is made against.

type alphabet:

np.ndarray or Constellation

param n_test:

Number of test phases \(B\) over \([0, \pi/2)\). Default is 32.

type n_test:

int, optional, keyword-only

param half_window:

Half-width \(W\) of the averaging window. Default is 16. The window trades noise rejection against tracking speed: it must be short against the phase drift and long against the additive noise.

type half_window:

int, optional, keyword-only

param name:

Block name. Default "bps".

type name:

str, optional, keyword-only

phase_#

Estimated unwrapped phase trajectory of the last call, of the same shape as the input – one trajectory per leading-axis row (data-dependent, hence the trailing underscore, D23).

Type:

np.ndarray

References

T. Pfau, S. Hoffmann, R. Noe, “Hardware-Efficient Coherent Digital Receiver Concept With Feedforward Carrier Recovery for M-QAM Constellations,” Journal of Lightwave Technology, vol. 27, no. 8, pp. 989-999, 2009.

Examples

>>> from comnumpy.core.utils import get_alphabet
>>> alphabet = get_alphabet("PSK", 4)
>>> s = alphabet[np.arange(64) % 4]
>>> rotated = s * np.exp(1j * 0.3)
>>> y = BlindPhaseSearchCompensator(alphabet, half_window=8)(rotated)
>>> bool(np.max(np.abs(y - s)) < np.pi / 2 / 32)   # within one grid step
True
class LinearEqualizer(h: ndarray, *, method: Literal['zf', 'mmse'] = 'zf', sigma2: float = 0.0, name: str = 'equalizer')[source]#

Block linear equalizer (ZF or MMSE) for a known FIR channel.

Signal Model#

The block input is the received signal of a channel with inter-symbol interference (ISI), of known impulse response \(h[l]\) of length \(L\):

\[x[n] = \sum_{l=0}^{L-1} h[l] \, s[n-l] + b[n], \qquad b[n] \sim \mathcal{CN}\left(0, \sigma^2\right)\]

Stacking the \(N\) received samples gives \(\mathbf{x} = \mathbf{H}\mathbf{s} + \mathbf{b}\), where \(\mathbf{H}\) is the \(N \times (N-L+1)\) banded Toeplitz convolution matrix built from \(h[l]\). The block returns the estimate \(\mathbf{y} = \widehat{\mathbf{s}}\) of the \(N-L+1\) transmitted symbols, by one of two closed forms:

\[\begin{split}\widehat{\mathbf{s}}_{\mathrm{ZF}} &= \mathbf{H}^{\dagger} \mathbf{x} = \left(\mathbf{H}^H \mathbf{H}\right)^{-1} \mathbf{H}^H \mathbf{x} \\ \widehat{\mathbf{s}}_{\mathrm{MMSE}} &= \left(\mathbf{H}^H \mathbf{H} + \sigma^2 \mathbf{I}\right)^{-1} \mathbf{H}^H \mathbf{x}\end{split}\]

The zero-forcing solution cancels the ISI exactly but amplifies the noise wherever \(\mathbf{H}\) is ill-conditioned; the MMSE solution trades residual ISI against noise enhancement through the loading term \(\sigma^2\) (written for unit-power symbols, \(\mathbb{E}[|s[n]|^2] = 1\)). The two coincide when \(\sigma^2 = 0\).

Axes: axis -1 – expects a serial signal (..., N) and returns (..., N-L+1); leading axes are batch (independent frames, one equalizer matrix).

param h:

Channel impulse response \(h[l]\), a 1D array of length \(L\). It is configured, not estimated: use DataAidedFIRCompensator when the channel is unknown.

type h:

np.ndarray

param method:

Closed form used. Default is "zf".

type method:

{“zf”, “mmse”}, optional, keyword-only

param sigma2:

Noise variance \(\sigma^2\), used by "mmse" only. Default is 0.0.

type sigma2:

float, optional, keyword-only

param name:

Name of the equalizer instance. Default is "equalizer".

type name:

str, optional, keyword-only

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 9 (linear equalization of band-limited channels: peak-distortion / zero-forcing and MMSE criteria).

Examples

>>> h = np.array([1.0, 0.5])
>>> s = np.array([1.0, -1.0, 1.0, 1.0])
>>> x = np.convolve(s, h)                       # noiseless ISI channel
>>> print(np.round(LinearEqualizer(h)(x).real, 6))
[ 1. -1.  1.  1.]
>>> print(np.round(LinearEqualizer(h, method="mmse", sigma2=0.1)(x).real, 3))
[ 0.862 -0.827  0.872  0.973]
class DataAidedFIRCompensator(reference: ndarray | None = None, *, should_fit: bool = True, name: str = 'data_aided_fir')[source]#

Data-aided FIR compensator: least-squares channel estimation, then deconvolution.

Signal Model#

The block observes the output of an unknown FIR channel driven by a known reference \(d[n]\) of length \(N\) (preamble, training sequence):

\[x[n] = \sum_{l} h[l] \, d[n-l]\]

The impulse response is estimated by least squares (zero-forcing criterion) from the \(N \times N\) lower-triangular Toeplitz convolution matrix \(\mathbf{D}\) built from \(d[n]\):

\[\widehat{\mathbf{h}} = \arg \min_{\mathbf{h}} \left\| \mathbf{D} \mathbf{h} - \mathbf{x} \right\|^2 = \mathbf{D}^{\dagger} \mathbf{x}\]

The correction is the deconvolution of the input by that estimate, so that \(y[n] \simeq d[n]\):

\[y[n] = \left(x * \widehat{h}^{-1}\right)[n]\]

Because \(\mathbf{D}\) is square here, the input and the reference must have the same length \(N\).

Axes: declared axis – expects a 1D serial signal (N,) of the same length as reference.

param reference:

Known reference \(d[n]\). When the reference is produced by the chain itself, leave it unset and declare the edge with Sequential(wiring={"data_aided_fir.reference": "source"}) instead of freezing an array.

type reference:

np.ndarray, optional

param should_fit:

If True (default), the channel is re-estimated at every call (per-block regime, D22); otherwise the estimate of the last fit is reused.

type should_fit:

bool, optional, keyword-only

param name:

Name of the processor. Default is "data_aided_fir".

type name:

str, optional, keyword-only

h_#

Estimated impulse response \(\widehat{\mathbf{h}}\) (data-dependent, hence the trailing underscore, decision D23). The response is never configured here: use LinearEqualizer(h) when it is already known.

Type:

np.ndarray

raises NotFittedError:

If forward is called with should_fit=False before any fit (decision D23).

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 9 (training-sequence channel estimation and zero-forcing equalization).

Examples

>>> x_ref = np.array([1.0, -1.0, 1.0, 1.0])
>>> received = np.convolve(x_ref, [1.0, 0.5])[:4]
>>> compensator = DataAidedFIRCompensator(reference=x_ref)
>>> print(np.round(compensator(received), 6))
[ 1. -1.  1.  1.]
>>> print(np.round(compensator.h_, 6) + 0.0)
[1.  0.5 0.  0. ]
class DataAidedPhaseCompensator(reference: ndarray | None = None, *, shared: bool | None = None, name: str = 'data_aided_phase')[source]#

Data-aided phase compensator (maximum-likelihood phase from a known reference).

Signal Model#

The input carries an unknown constant phase offset with respect to a known reference \(d[n]\):

\[x[n] = d[n] \, e^{-j\theta} + b[n], \qquad y[n] = x[n] \, e^{j\widehat{\theta}}\]

In additive white Gaussian noise, the maximum-likelihood estimate of \(\theta\) is the argument of the cross-correlation between the input and the reference at lag zero:

\[\widehat{\theta} = \arg \left( \sum_{n=0}^{N-1} x^*[n] \, d[n] \right)\]

Unlike the blind estimator, this one is free of the \(2\pi/M\) constellation ambiguity, since the reference fixes the absolute phase.

Axes: declared axis – estimation expects a signal aligned with reference along the last axis; the correction is element-wise. A (..., P, N) signal needs shared= to say which D49 family the phase belongs to, and is refused otherwise.

param reference:

Known reference \(d[n]\). When the reference is produced by the chain itself, leave it unset and declare the edge with Sequential(wiring={"data_aided_phase.reference": "source"}) instead of freezing an array.

type reference:

np.ndarray, optional

param shared:

Which D49 family the phase belongs to, for a multi-path signal. True fits one angle jointly over every path – laser phase noise ahead of any equalizer is common to both polarizations. False fits one per path. Default is None, which refuses a multi-path signal rather than choosing for the caller.

type shared:

bool, optional, keyword-only

param name:

Name of the processor. Default is "data_aided_phase".

type name:

str, optional, keyword-only

theta_#

Estimated phase \(\widehat{\theta}\) in radians (data-dependent, hence the trailing underscore, decision D23). Re-estimated at every call (per-block regime, D22).

Type:

float or np.ndarray

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 5 (data-aided maximum-likelihood carrier phase estimation).

Examples

>>> d = np.array([1+1j, 1-1j, -1+1j, -1-1j]) / np.sqrt(2)
>>> x = d * np.exp(-1j * np.pi / 4)          # a -pi/4 rotation
>>> compensator = DataAidedPhaseCompensator(reference=d)
>>> y = compensator(x)
>>> print(round(float(compensator.theta_), 6), round(float(np.pi / 4), 6))
0.785398 0.785398
>>> print(np.round(y, 6))
[ 0.707107+0.707107j  0.707107-0.707107j -0.707107+0.707107j
 -0.707107-0.707107j]
class DataAidedComplexGainCompensator(reference: ndarray | None = None, *, extractor: DataExtractor = <factory>, should_fit: bool = True, shared: bool | None = None, name: str = 'complex_gain_compensator')[source]#

Data-aided compensator of a flat complex gain (amplitude and phase).

Signal Model#

The channel is assumed frequency-flat over the record, so it reduces to a single complex coefficient. A preamble is extracted from the input and compared to the known reference \(d[n]\) of length \(N_d\):

\[\widehat{g} = \arg \min_{g} \sum_{n=0}^{N_d-1} \left| g \, x[n] - d[n] \right|^2 = \frac{\sum_{n} x^*[n] \, d[n]}{\sum_{n} |x[n]|^2}\]

The estimate is the compensation gain (an estimate of the inverse channel gain), applied to the whole record:

\[y[n] = \widehat{g} \, x[n]\]

It corrects amplitude and phase at once, and generalizes DataAidedPhaseCompensator, which keeps only \(\arg(\widehat{g})\).

Axes: declared axis – estimation expects a preamble aligned with reference along the last axis (after extractor); the correction is element-wise on the full record. A (..., P, N) signal needs shared= to say which D49 family the gain belongs to, and is refused otherwise.

param reference:

Known preamble \(d[n]\). When it is produced by the chain itself, leave it unset and declare the edge with Sequential(wiring={"complex_gain_compensator.reference": "source"}) instead of freezing an array.

type reference:

np.ndarray, optional

param shared:

Which D49 family the gain belongs to, for a multi-path signal. True fits one gain jointly over every path – the case of a dual-polarization link whose nonlinear phase comes from the total intensity, or of one laser feeding both. False fits one gain per path, as after a butterfly equalizer. Default is None, which refuses a multi-path signal rather than choosing for the caller.

type shared:

bool, optional, keyword-only

param extractor:

Selects the preamble samples inside the received record. Default is a pass-through (the whole record is the preamble).

type extractor:

DataExtractor, optional, keyword-only

param should_fit:

If True (default), \(\widehat{g}\) is re-estimated at every call (per-block regime, D22); if False, the gain of the last fit is reused (reused-preamble regime).

type should_fit:

bool, optional, keyword-only

param name:

Name of the processor. Default is "complex_gain_compensator".

type name:

str, optional, keyword-only

gain_#

Estimated complex gain \(\widehat{g}\) (data-dependent, hence the trailing underscore, decision D23).

Type:

complex or np.ndarray

raises NotFittedError:

If forward is called with should_fit=False before any fit (decision D23).

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 5 (data-aided estimation of the carrier amplitude and phase).

Examples

>>> d = np.array([1+1j, 1-1j, -1+1j, -1-1j])
>>> x = 2 * np.exp(1j * np.pi / 3) * d           # channel gain 2 exp(j pi/3)
>>> compensator = DataAidedComplexGainCompensator(reference=d)
>>> y = compensator(x)
>>> print(np.round(compensator.gain_, 6))        # 0.5 exp(-j pi/3)
(0.25-0.433013j)
>>> print(np.round(y, 6))
[ 1.+1.j  1.-1.j -1.+1.j -1.-1.j]
class DataAidedSimpleSynchronizer(reference: ndarray | None = None, scale_correction: bool = True, *, save_cross_correlation: bool = True, signal_len: int | None = None, name: str = 'synchronizer')[source]#

Data-aided frame synchronizer at one-sample resolution (correlation peak).

Signal Model#

The received record is a delayed and scaled copy of a known preamble \(d[n]\) of length \(N_d\), followed by the payload:

\[x[n] = a \, d[n - m_0] + b[n]\]

The delay is obtained from the peak of the cross-correlation between the input and the preamble, normalized by the preamble energy:

\[c[m] = \frac{\sum_{n} x[n+m] \, d^*[n]}{\sum_{n} |d[n]|^2}, \qquad \widehat{m}_0 = \arg \max_{m} \left| c[m] \right|^2\]

That normalization makes the peak a least-squares estimate of the channel gain itself, \(\widehat{a} = c[\widehat{m}_0] \simeq a\). The output realigns the record on the preamble and, when scale_correction is enabled, divides it by that gain, so that the block restores the transmitted signal instead of returning a doubly-distorted copy of it:

\[y[n] = \frac{x[n + \widehat{m}_0]}{\widehat{a}}\]

For a unit channel gain (\(a = 1\)), \(\widehat{a} = 1\) and the block is a pure realignment, whatever the preamble power.

Axes: declared axis – expects a 1D serial signal (N,); the output length depends on the estimated delay and on signal_len.

param reference:

Known preamble \(d[n]\). When it is produced by the chain itself, leave it unset and declare the edge with Sequential(wiring={"synchronizer.reference": "source"}) instead of freezing an array.

type reference:

np.ndarray, optional

param scale_correction:

If True (default), the output is divided by the complex peak \(\widehat{a}\); otherwise only the delay is corrected.

type scale_correction:

bool, optional

param save_cross_correlation:

If True (default), stores \(c[m]\) and its lag axis for inspection or plotting (plot).

type save_cross_correlation:

bool, optional, keyword-only

param signal_len:

Truncates the realigned signal to that length. Default is None (keep everything after the peak).

type signal_len:

int, optional, keyword-only

param name:

Name of the synchronizer instance. Default is "synchronizer".

type name:

str, optional, keyword-only

delay_#

Estimated delay \(\widehat{m}_0\) in samples (data-dependent, hence the trailing underscore, decision D23; re-estimated at every call).

Type:

int

scale_#

Compensation gain \(1/\widehat{a}\) actually applied, or 1 when scale_correction is False.

Type:

complex

cross_corr_#

Full normalized cross-correlation \(c[m]\) (when saved).

Type:

np.ndarray

n_vect_#

Lag axis \(m\) associated with cross_corr_.

Type:

np.ndarray

raises ValueError:

If the reference has zero energy, which makes the normalized correlation undefined.

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 5 (frame and symbol synchronization by correlation with a known preamble).

Examples

>>> d = np.array([1.0, -1.0, 1.0, 1.0])
>>> x = np.hstack([np.zeros(2), d, [1.0, -1.0]])   # delayed by 2 samples
>>> synchronizer = DataAidedSimpleSynchronizer(d)
>>> y = synchronizer(x)
>>> print(synchronizer.delay_, round(float(synchronizer.scale_), 6))
2 1.0
>>> print(y)
[ 1. -1.  1.  1.  1. -1.]
>>> print(np.round(DataAidedSimpleSynchronizer(d)(0.5 * x), 6))   # 6 dB attenuation
[ 1. -1.  1.  1.  1. -1.]
class DataAidedFineSynchronizer(reference: ndarray | None = None, scale_correction: bool = True, *, up_factor: int = 2, save_cross_correlation: bool = True, signal_len: int | None = None, d_max: int | None = None, name: str = 'synchronizer')[source]#

Data-aided frame synchronizer at fractional-sample resolution.

Signal Model#

Same correlation-based model as DataAidedSimpleSynchronizer, but the delay \(\tau_0\) is no longer assumed to be an integer number of samples:

\[x[n] = a \, d\!\left(nT - \tau_0\right) + b[n]\]

The input and the preamble are first interpolated by a polyphase resampler of factor \(U\) (scipy.signal.resample_poly), which brings the search grid down to \(T/U\):

\[c[m] = \frac{\sum_{n} x_{\uparrow}[n+m] \, d_{\uparrow}^*[n]}{\sum_{n} |d_{\uparrow}[n]|^2}, \qquad \widehat{m}_0 = \max\left(0,\; \arg \max_{m} \left| c[m] \right|^2 \right)\]

As in DataAidedSimpleSynchronizer, the normalization by the preamble energy makes the peak an estimate of the channel gain, \(\widehat{a} = c[\widehat{m}_0] \simeq a\). The upsampled record is realigned on that grid, decimated back by \(U\), and optionally divided by that gain:

\[\widehat{\tau}_0 = \frac{\widehat{m}_0}{U} T, \qquad y[n] = \frac{x_{\uparrow}\!\left[U n + \widehat{m}_0\right]} {\widehat{a}}\]

Negative lags are clamped to zero, so the block only corrects non-negative delays.

Axes: declared axis – expects a 1D serial signal (N,); the output length depends on the estimated delay and on signal_len.

param reference:

Known preamble \(d[n]\). When it is produced by the chain itself, leave it unset and declare the edge with Sequential(wiring={"synchronizer.reference": "source"}) instead of freezing an array.

type reference:

np.ndarray, optional

param scale_correction:

If True (default), the output is divided by the complex peak \(\widehat{a}\); otherwise only the delay is corrected.

type scale_correction:

bool, optional

param up_factor:

Interpolation factor \(U\) applied before the correlation; the timing resolution is \(T/U\). Default is 2.

type up_factor:

int, optional, keyword-only

param save_cross_correlation:

If True (default), stores \(c[m]\) and its lag axis for inspection or plotting (plot).

type save_cross_correlation:

bool, optional, keyword-only

param signal_len:

Truncates the realigned signal to that length. Default is None.

type signal_len:

int, optional, keyword-only

param d_max:

Maximum expected delay, in samples of \(x\). Restricts the correlation window and therefore the search range. Default is None (search over the whole record).

type d_max:

int, optional, keyword-only

param name:

Name of the synchronizer instance. Default is "synchronizer".

type name:

str, optional, keyword-only

delay_#

Estimated delay \(\widehat{m}_0\), expressed in upsampled samples (data-dependent, hence the trailing underscore, decision D23; divide by up_factor to get it in input sample periods).

Type:

int

scale_#

Compensation gain \(1/\widehat{a}\) actually applied, or 1 when scale_correction is False.

Type:

complex

cross_corr_#

Full normalized cross-correlation \(c[m]\) on the upsampled grid.

Type:

np.ndarray

n_vect_#

Lag axis \(m\) associated with cross_corr_.

Type:

np.ndarray

raises ValueError:

If the reference has zero energy, which makes the normalized correlation undefined.

References

  • J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Chapter 5 (frame and symbol synchronization).

  • F. M. Gardner, “Interpolation in digital modems – Part I: Fundamentals,” IEEE Transactions on Communications, vol. 41, no. 3, pp. 501-507, 1993 (interpolation-based fractional timing recovery).

Examples

>>> d = np.array([1.0, -1.0, 1.0, 1.0])
>>> x = np.hstack([np.zeros(3), d, np.zeros(3)])   # delayed by 3 samples
>>> synchronizer = DataAidedFineSynchronizer(d, up_factor=2, signal_len=4)
>>> y = synchronizer(x)
>>> print(synchronizer.delay_, synchronizer.delay_ / synchronizer.up_factor, y.shape)
6 3.0 (4,)