Symmetries of Second-Order Probabilistic Parameters in CSP

Do we need to consider all cycle frequencies, both positive and negative? Do we need to consider all delays and frequencies in our second-order CSP parameters?

As you progress through the various stages of learning CSP (intimidation, frustration, elucidation, puzzlement, and finally smooth operation), the symmetries of the various functions come up over and over again. Exploiting symmetries can result in lower computational costs, quicker debugging, and easier mathematical development.

What exactly do we mean by ‘symmetries of parameters?’ I’m talking primarily about the evenness or oddness of the time-domain functions in the delay $\tau$ and cycle frequency $\alpha$ variables and of the frequency-domain functions in the spectral frequency $f$ and cycle frequency $\alpha$ variables. Or a generalized version of evenness/oddness, such as $f(-x) = g(x)$ , where $f(x)$ and $g(x)$ are closely related functions. We have to consider the non-conjugate and conjugate functions separately, and we’ll also consider both the auto and cross versions of the parameters. We’ll look at higher-order cyclic moments and cumulants in a future post.

You can use this post as a resource for mathematical development because I present the symmetry equations. But also each symmetry result is illustrated using estimated parameters via the frequency smoothing method (FSM) of spectral correlation function estimation. The time-domain parameters are obtained from the inverse transforms of the FSM parameters. So you can also use this post as an extension of the second-order verification guide to ensure that your estimator works for a wide variety of input parameters.

The Non-Conjugate Cyclic Autocorrelation Function

Let’s start our symmetry odyssey with arguably the simplest CSP function: the non-conjugate cyclic autocorrelation function. Recall that this function is equal to the conventional autocorrelation function when the cycle frequency $\alpha$ is zero. The non-conjugate cyclic autocorrelation is defined by

$\displaystyle R_x^\alpha(\tau) = \left\langle x(t+\tau/2)x^*(t-\tau/2) e^{-i2\pi\alpha t} \right\rangle_t \hfill (1)$

where the notation $\left\langle \cdot \right\rangle_t$ denotes the infinite time average:

$\displaystyle \left\langle x(t) \right\rangle_t = \lim_{T\rightarrow\infty} \frac{1}{T} \int_{-T/2}^{T/2} x(t) \, dt. \hfill (2)$

The symmetry in the lag parameter $\tau$ is determined by examining $\displaystyle R_x^\alpha(\tau)$ for negative $\tau$ ,

$\displaystyle R_x^\alpha(-\tau) = \left\langle x(t-\tau/2)x^*(t+\tau/2) e^{-i2\pi\alpha t} \right\rangle_t \hfill (3)$

A little algebra yields

$\displaystyle R_x^\alpha(-\tau) = \left[ \left\langle x(t+\tau/2) x^* (t-\tau/2) e^{-i2\pi(-\alpha) t} \right\rangle_t \right]^* =R_x^{-\alpha}(\tau)^* \hfill (4a)$

Note that (4a) also reveals the symmetry in the cycle frequency parameter, since it can be rearranged as

$\displaystyle R_x^\alpha(\tau) = R_x^{-\alpha}(-\tau)^*. \hfill(4b)$

Let’s illustrate our sequence of symmetry results using our old friend, the rectangular-pulse BPSK signal. It has ten samples per bit, or $f_{bit} = 0.1$ Hz, and a carrier frequency offset of $0.05$ Hz. The carrier frequency offset renders the signal complex valued, and we also add complex-valued white Gaussian noise. The power of the signal is unity and the power of the noise is $0.1$ . The true parameters for this signal are known: frequency domain and time domain.

First, let’s take a look at the autocorrelation, which yields the symmetry relation (4) with $\alpha = 0$ ,

$\displaystyle R_x^0(-\tau) = R_x^0(\tau)^* \hfill (5)$

This implies that the real part and the magnitude are even functions of $\tau$ and the imaginary part is an odd function of $\tau$ . Figure 1 shows the estimated autocorrelation, confirming these even/odd predictions:

rxat_nc_0 — Figure 1. Estimated autocorrelation function (non-conjugate cyclic autocorrelation for a cycle frequency of 0). Note that the real part and the magnitude are even functions of the lag and the imaginary part is an odd function of the lag.

Moving on to the cycle frequency of $\alpha = f_{bit} = 0.1$ , we obtain the estimates shown in Figure 2. The result in (4) implies that the real part of the CAF for $\alpha = f_{bit}$ at $\tau$ should be equal to the real part of the CAF for $\alpha = -f_{bit}$ at $-\tau$ (verify with $\tau = 4$ ). It also implies that the imaginary part for $\alpha$ and $\tau$ should be the negative of that for $-\alpha$ and $-\tau$ (verify with $\tau = -6$ ). So the equation and the estimates are consistent.

rxat_nc — Figure 2. Estimated non-conjugate cyclic autocorrelation when the cycle frequency is equal to the bit rate or its negative.

The Conjugate Cyclic Autocorrelation Function

A similar analysis applies to the conjugate cyclic autocorrelation function. Let’s go through it. The function is defined by

$\displaystyle R_{x^*}^\alpha(\tau) = \left\langle x(t+\tau/2)x(t-\tau/2) e^{-i2\pi \alpha t} \right\rangle_t \hfill (6)$

Symmetry in $\tau$ follows easily

$\displaystyle R_{x^*}^\alpha(-\tau) = R_{x^*}^\alpha(\tau). \hfill (7)$

because we don’t have any conjugations to deal with. This simple result implies that the real and imaginary parts (and the magnitude) are even functions of $\tau$ . This symmetry is validated using estimates in Figure 3.

rxat_c — Figure 3. Estimated conjugate cyclic autocorrelation functions for three different cycle frequencies. The real part, imaginary part, and magnitude should all be even functions of the lag variable according to the mathematical result (7).

The symmetry relation (7) tells us nothing about the symmetry in cycle frequency $\alpha$ . So let’s look at $R_{x^*}^{-\alpha}(\tau)$ ,

$\displaystyle R_{x^*}^{-\alpha}(\tau) = \left\langle x(t+\tau/2) x(t-\tau/2)e^{-i2\pi(-\alpha)t} \right\rangle_t \hfill (8)$

$\displaystyle = \left\langle \left[ x^*(t+\tau/2)x^*(t-\tau/2) e^{-i2\pi\alpha t} \right]^* \right\rangle_t \hfill (9)$

$\displaystyle = \left\langle x^*(t+\tau/2)x^*(t-\tau/2) e^{-i2\pi\alpha t} \right\rangle_t^* \hfill (10)$

But (10) is just the conjugated conjugate CAF for the conjugated signal $x^*(t)$ ,

$\displaystyle R_{x^*}^{-\alpha}(\tau) = \left[ R_{(x^*)^*}^\alpha (\tau) \right]^*, \hfill (11)$

$\displaystyle R_{x^*}^\alpha (\tau) = \left[ R_{(x^*)^*}^{-\alpha} (\tau) \right]^*. \hfill (12)$

The symmetry relation (12) is validated by estimating the conjugate CAF of $x(t)$ for $\alpha = 2f_c$ and the conjugate CAF of $x^*(t)$ for $\alpha = -2f_c$ and plotting the results. Both functions should be even functions of $\tau$ , and they should be conjugates of each other. This is verified in Figure 4.

rxat_c_conj — Figure 4. Illustration of the symmetry in cycle frequency for the conjugate CAF. The conjugate CAF itself possesses no general symmetry in cycle frequency, but is the conjugate of the conjugate CAF for conjugated data and the negated cycle frequency.

The Non-Conjugate Cyclic Cross Correlation Function

Let’s stay with the time-domain parameters and introduce cross functions, starting with the non-conjugate cyclic cross correlation. The basic definition is

$\displaystyle R_{xy}^\alpha(\tau) = \left\langle x(t+\tau/2) y^*(t-\tau/2) e^{-i2\pi\alpha t} \right\rangle_t \hfill (13)$

Consider the function $R_{xy}^\alpha(-\tau)$ ,

$\displaystyle R_{xy}^\alpha(-\tau) = \left\langle x(t-\tau/2)y^*(t+\tau/2) e^{-i2\pi\alpha t} \right\rangle_t \hfill (14)$

$\displaystyle = \left\langle [y(t+\tau/2)x^*(t-\tau/2)]^* e^{-i2\pi\alpha t} \right\rangle_t \hfill (15)$

$\displaystyle = \left[ \left\langle y(t+\tau/2)x^*(t-\tau/2) e^{-i2\pi(-\alpha)t} \right\rangle_t \right]^* \hfill (16)$

$\displaystyle = R_{yx}^{-\alpha}(\tau)^* \hfill (17)$

The relation (17) also contains the symmetry relation for cycle frequency $\alpha$ and the symmetry in ‘function order,’ that is, the symmetry relating to switching the order of $x(t)$ and $y(t)$ in the integrand. The symmetry relation for the cross cyclic CAF is illustrated for $\tau$ in Figure 5, which shows the (conventional) cross correlation symmetry because $\alpha = 0$ . Here $x(t)$ is the same rectangular-pulse BPSK signal as we’ve used in the previous illustrations, and $y(t)$ is the same signal shifted by $6$ samples and subject to a different carrier phase and independent Gaussian noise. (Shift the noise-free signal making up $x(t)$ to the right by $6$ samples to get the noise-free $y(t)$ .) We therefore expect the cross correlation to peak at $\tau = -6$ in the XY cross function and at $\tau = 6$ in the YX version.

rxyat_nc_0 — Figure 5. Illustration of the symmetries in lag and function order (XY vs YX) for the non-conjugate cyclic cross correlation. See (17).

A more complicated example is shown in Figure 6, which features the cycle frequencies of $\pm f_{bit}$ .

rxyat_nc_1 — Figure 6. Illustration of the symmetries in lag, cycle frequency, and function order (XY vs YX) for the non-conjugate cyclic cross correlation. See (17).

The Conjugate Cyclic Cross Correlation Function

A similar analysis yields the symmetry relations for the conjugate cyclic cross correlation function,

$\displaystyle R_{xy^*}^\alpha(-\tau) = R_{yx^*}^\alpha(\tau) \hfill (18)$

and

$\displaystyle R_{yx^*}^\alpha(\tau) = \left[ R_{(x^*y^*)^*}^{-\alpha} (\tau) \right]^* \hfill (19)$

The symmetry in $\tau$ is validated in Figure 7, whereas the symmetry in cycle frequency is validated in Figure 8.

rxyat_c_0 — Figure 7. Illustration of the lag-variable symmetry of the conjugate cross cyclic correlation function. See (18).

rxyat_c_0_conj — Figure 8. Illustration of the cycle-frequency symmetry in the conjugate cyclic cross correlation function. See (19).

The Non-Conjugate Spectral Correlation Function

The symmetry derivations and illustrations for the cyclic correlation functions were fun, yes, but let’s move on to the important function: spectral correlation. First, our oldest and dearest friend, the non-conjugate spectral correlation function.

The SCF is the temporal correlation between the time-series obtained from two narrowband spectral components of a signal. The narrowband components can be obtained using a simple sliding Fourier transform,

$\displaystyle X_T(t, f) = \int_{t-T/2}^{t+T/2} x(u) e^{-i2\pi f u} \, du \hfill (20)$

Their idealized correlation is

$\displaystyle \lim_{T\rightarrow\infty} \left\langle \frac{1}{T} X_T(t, f+\alpha/2) X_T^*(t, f-\alpha/2) \right\rangle_t \hfill (21)$

where we first average over all time, and then let the bandwidth $1/T$ of the narrowband components approach zero. As we go along looking at the symmetry properties of the non-conjugate, conjugate, cross, and conjugate cross SCFs, we’ll need an expression for the narrowband component of the conjugate of $x(t)$ ,

$\displaystyle X_{T,*}(t,f) = \int_{t-T/2}^{t+T/2} x^*(u) e^{-i2\pi f u} \, du \hfill (22)$

$\displaystyle = \left[ \int_{t-T/2}^{t+T/2} x(u) e^{-i2\pi (-f) u} \, du \right]^* \hfill (23)$

$\displaystyle = X_T(t, -f)^*. \hfill (24)$

Let’s look at the spectral correlation for the conjugated signal,

$\displaystyle S_{(x^*)}^\alpha(f) = \lim_{T\rightarrow\infty} \left\langle \frac{1}{T} X_{T,*}(t, f+\alpha/2) X_{T,*}^*(t, f-\alpha/2) \right\rangle_t \hfill (25)$

$\displaystyle = \left\langle \frac{1}{T} X_T^*(t, -(f+\alpha/2)) [ X_T^*(t, -(f-\alpha/2)) ]^* \right\rangle \hfill (26)$

$\displaystyle = \left\langle \frac{1}{T} X_T(t, -(f)+\alpha/2) X_T^*(t, (-f)+\alpha/2) \right\rangle \hfill (27)$

$\displaystyle = S_x^\alpha(-f) \hfill (28)$

So the non-conjugate SCF doesn’t possess any symmetry in $f$ on its own, but if we include the conjugated-signal SCF, we obtain a sort of symmetry.

Symmetry in cycle frequency $\alpha$ is a bit easier,

$\displaystyle S_x^{-\alpha} (f) = \lim_{T\rightarrow\infty} \left\langle \frac{1}{T} X_T(t, f-\alpha/2) X_T^*(t, f+\alpha/2) \right\rangle_t \hfill (29)$

$\displaystyle = \lim_{T\rightarrow\infty} \left[ \left\langle \frac{1}{T} X_T(t, f+\alpha/2) X_T^*(t, f-\alpha/2) \right\rangle_t \right]^* \hfill (30)$

$\displaystyle = \left[ S_x^\alpha(f) \right]^*. \hfill (31)$

So the information in the non-conjugate SCFs for negative cycle frequencies is redundant with that for positive cycle frequencies. The two forms of symmetry for the non-conjugate SCF are illustrated in Figures 9-11.

sxaf_nc_0 — Figure 9. Illustration of symmetry in frequency f for the non-conjugate SCF for cycle frequency of zero (the PSD). The function itself is not symmetric in f, but it has a mirror-image relationship with the corresponding function for the conjugated-signal SCF with the same cycle frequency.

sxaf_nc_1 — Figure 10. Illustration of symmetry in cycle frequency for the non-conjugate spectral correlation function. The SCFs for positive cycle frequencies are complex-conjugates of those for negative cycle frequencies.

sxaf_nc_1_conj — Figure 11. Illustration of symmetry in frequency f for the non-conjugate SCF for cycle frequency of the BPSK bit rate. The function itself is not symmetric in f, but it has a mirror-image relationship with the corresponding function for the conjugated-signal SCF with the same cycle frequency.

The Conjugate Spectral Correlation Function

Using an analysis approach similar to that for the non-conjugate SCF, we obtain the symmetry relations for the conjugate SCF:

$\displaystyle S_{x^*}^\alpha(-f) = S_{x^*}^\alpha(f) \hfill (32)$

$\displaystyle S_{x^*}^{-\alpha}(f) = S_{(x^*)^*}^\alpha(f)^* \hfill (33)$

The conjugate SCF is even in $f$ , and has no cycle-frequency symmetry itself, but possesses a symmetry relation with the conjugate SCF for the conjugated signal. Figures 12-14 illustrate and validate the symmetries.

sxaf_c_-1_conj — Figure 12. Illustration of the symmetry in the frequency and cycle frequency variables for the conjugate spectral correlation function. See (32) and (33).

sxaf_c_0_conj — Figure 13. Illustration of the symmetry in the frequency and cycle frequency variables for the conjugate spectral correlation function. See (32) and (33).

The Non-Conjugate Cross Spectral Correlation Function

The non-conjugate cross SCF is defined by

$\displaystyle S_{xy}^\alpha(f) = \lim_{T\rightarrow\infty} \left\langle \frac{1}{T} X_T(t, f+\alpha/2) Y_T^*(t, f-\alpha/2) \right\rangle_t \hfill (34)$

Using the same kinds of analysis as above, we obtain a symmetry relation for frequency $f$ of

$\displaystyle S_{xy}^\alpha(f) = S_{(x^*y^*)^*}^\alpha(-f), \hfill (35)$

and for cycle frequency $\alpha$ of

$\displaystyle S_{xy}^\alpha(f) = S_{yx}^{-\alpha}(f)^*. \hfill (36)$

Equation (36) also provides the symmetry in function order. So there is no symmetry in $f$ for the function itself, but there is between the function and the non-conjugate cross SCF for conjugated $x(t)$ and $y(t)$ . Similarly, there is no symmetry in $\alpha$ for the function itself, but there is between (34) and the conjugated function-reversed version of (34) with negated cycle frequency.

The symmetries are illustrated in Figures 15 and 16.

sxyaf_nc_0_conj — Figure 15. Illustration of the symmetry relation between the non-conjugate cross SCF and the function-reversed non-conjugate cross SCF for conjugated inputs. See (35).

sxyaf_nc_1 — Figure 16. Illustration of the symmetry relation between (34) and the function-reversed version of (34) with negated cycle frequency. See (36).

The Conjugate Cross Spectral Correlation Function

If you can bear it, there is one more case. The conjugate cross spectral correlation function is defined by

$\displaystyle S_{xy^*}^\alpha(f) = \lim_{T\rightarrow\infty} \left\langle \frac{1}{T} X_T(t, f+\alpha/2) Y_T(t, \alpha/2-f) \right\rangle_t \hfill (37)$

The symmetry relations are actually relatively easy to derive here. For frequency $f$ and function order we have

$\displaystyle S_{xy^*}^\alpha(f) = S_{yx^*}^\alpha(-f) \hfill (38)$

and for cycle frequency we have

$\displaystyle S_{xy^*}^\alpha(f) = \left[ S_{(y^*x^*)^*}^{-\alpha}(f) \right]^* \hfill (39)$

Once again, no simple symmetry involving just the function itself (that is, (37)); all symmetry relations depend on reversing the function order and/or involve the conjugate cross spectral correlation function for conjugated inputs. The relations (38) and (39) are illustrated and validated in Figures 17 and 18.

sxyaf_c_0 — Figure 17. Illustration of the symmetry relation between the XY conjugate cross SCF and the YX conjugate cross SCF for the frequency variable f. See (38).

sxyaf_c_0_conj — Figure 18. Illustration of the symmetry relation between the XY conjugate cross SCF and the YX conjugated-input conjugate cross SCF for the cycle-frequency variable. See (39).

Symmetries for the Spectral Coherence Function

I can hear you saying, “Yeah, great, but what about the coherence function?” And to be honest, that makes me tired, as this has been a long post. But I’m here for you. The non-conjugate spectral coherence function, greatly useful in blind CSP, is defined by

$\displaystyle C_x^\alpha(f) = \frac{S_x^\alpha(f)}{\left[ S_x^0(f+\alpha/2) S_x^0(f-\alpha/2) \right]^{1/2}} \hfill (40)$

The denominator of (40) is symmetric in $\alpha$ : it is the same when $\alpha$ is replaced by $-\alpha$ . But it is not symmetric in $f$ . Recall that the numerator is not symmetric in $f$ , but obeys the cycle-frequency symmetry

$\displaystyle S_x^{-\alpha}(f) = S_x^\alpha(f)^*$ .

So we have for the non-conjugate coherence

$\displaystyle C_x^{-\alpha}(f) = C_x^{\alpha}(f)^*.\hfill (41)$

For the conjugate coherence,

$\displaystyle C_{x^*}^\alpha(f) = \frac{S_{x^*}^\alpha(f)}{\left[ S_x^0(f+\alpha/2) S_x^0(\alpha/2-f) \right]^{1/2}} \hfill (42)$

both numerator and denominator are symmetric in $f$ , but not in $\alpha$ . So we have the symmetry

$\displaystyle C_{x^*}^\alpha(-f) = C_{x^*}^\alpha(f) \hfill (43)$

The cross coherences don’t seem to possess much symmetry. I’ll leave those as an exercise for the reader.

Discussion

I think the main point for CSP practitioners is that the non-conjugate spectral correlation for $\alpha \ge 0$ is redundant with respect to that for $\alpha \le 0$ . If you’ve estimated the non-conjugate SCF for $\alpha \ge 0$ , you’ve got all the non-conjugate information you need. On the other hand, you have to estimate the conjugate spectral correlation (or coherence) for the entire range of cycle frequency to ensure you find all the valid cycle frequencies.

Comments, corrections, suggestions, and compliments welcome below!

Author: Chad Spooner

I'm a signal processing researcher specializing in cyclostationary signal processing (CSP) for communication signals. I hope to use this blog to help others with their cyclo-projects and to learn more about how CSP is being used and extended worldwide. View all posts by Chad Spooner

17 thoughts on “Symmetries of Second-Order Probabilistic Parameters in CSP”

Chad Spooner says:

December 12, 2019 at 3:10 pm

My favorite part of this post is the phrase “conjugated conjugate CAF for the conjugated signal.”

Reply
John Tyler says:

January 23, 2020 at 1:13 pm

Hello Dr. Spooner,

First of all I would like to thank you for this wonderful blog!

I am currently implementing the Cross-CAF (non-conjugate) and am having some trouble with the symmetry relations. I am able to reasonably replicate Fig. 5. However for Fig. 6(a) I get the correct results only when I plot the (XY) CAF for alpha = +fbit. In Fig. 6(b) I get the correct results when plotting the (YX) CAF for alpha = –fbit. Is it possible there is a typographical error in Fig. 6 and the signs for fbit were reversed? I went through a few iterations to calculate symmetry according to (17).

To make sure I am correctly using your rectangular BPSK signal, I made x[n] equal to the rectangular BPSK signal with fc = 0.05 and fbit = 0.1. If I interpreted your text correctly, I defined y[n] as a delayed version (by 6 samples) of x[n] by adding six zeros in front of it and discarding the last six samples. Is this correct?

Also, is there anyway around not having to compute Ryx (when solving for Rxy) for symmetry? It adds to the computation time.

When computing the CAF, do we average by dividing by N, the number of points in x[n], or by N-tau?

When multiplying by the exponent (before averaging) do we need to use a time vector of 1:N-tau or tau+1:N? The approach I take in my code is to only use those values of x[n] and y[n] that contribute to the average (values multiplied by zero are discarded). Therefore, I thought it would make sense to start the time vector concurrent with the non-delayed version of x[n] i.e. xn(tau+1:end).

Thank You!

-John

Reply
1. Chad Spooner says:
  
  January 23, 2020 at 5:32 pm
  
  Great questions. I believe these questions are on the minds of others, or will be someday, and so your posing them and my attempts at answers should help a lot of other people over time.
  
  Is it possible there is a typographical error in Fig. 6 and the signs for fbit were reversed? I went through a few iterations to calculate symmetry according to (17).
  
  Yes, it is possible. I checked and could not find an error though. More likely the problem is in the interpretation of the delay of $-6$ samples, considered next.
  
  If I interpreted your text correctly, I defined y[n] as a delayed version (by 6 samples) of x[n] by adding six zeros in front of it and discarding the last six samples. Is this correct?
  
  I will make the text more clear: $y(t)$ is the noise-free BPSK signal in $x(t)$ shifted to the left by $6$ samples. So you’d want to discard the first six samples of $x(t)$ and append six zeros to the end to create $y(t)$ . This is likely the reason for the difference in our plots.
  
  Also, is there anyway around not having to compute Ryx (when solving for Rxy) for symmetry?
  
  I don’t understand this question. I get the sense you are asking something about (13)-(17), but I can’t grasp it.
  
  When computing the CAF, do we average by dividing by N, the number of points in x[n], or by N-tau?
  
  When I estimate the CAF directly in the time domain, I divide by $N$ . Dividing by $N - |\tau|$ is more accurate because that represents the number of points that are actually summed over, but I don’t care much, because in all cases of interest to me, $N \gg \tau_*$ , where $\tau_*$ is the largest lag of any practical interest. In fact, typically when I create something like $x(t + \tau)$ , I circularly shift the vector of samples, so I’m always averaging $N$ samples, and I just live with the tiny edge effect.
  
  When multiplying by the exponent (before averaging) do we need to use a time vector of 1:N-tau or tau+1:N?
  
  The key thing is to maintain the proper relationship between the phase of the applied complex exponential and the samples of the lag product as indexed by $t$ . If you mess around with the lag product by truncating it or lopping off the end points or something, you just need to make sure that you are adjusting the time origin of the vector of time values inside the exponential. As noted above, I don’t typically do that because I live with the edge effects.
  
  More importantly, I don’t often directly estimate the CAF in the time domain. I inverse transform an estimate of the spectral correlation function.
  
  Reply
  1. John Tyler says:
    
    January 24, 2020 at 6:30 am
    
    Hello Dr. Spooner,
    
    Thank you for the swift response!
    
    What I mean is that, according to my understanding, in (13)–(17) it implies that we need to solve for two things in order to properly account for symmetry. We need to solve for Rxy(alpha)(tau) and Ryx(alpha)(tau). Therefore, in my code I solve for Rxy and also solve for Ryx, where the order in the integrand is reversed. Hence, the increase in computation time. Here is part of my code (Matlab) for calculating the cross-CAF (non-conjugate)
    
    % Define N, the length of x[n] and y[n]
    N = length(xn);
    
    % n is the index for x[n]
    n = 0:N-1;
    
    tauLoop = 0:tauMax;
    
    for tau = tauLoop
    
    RxyTauAlpha(tau+1,:) = (sum(xn.*[zeros(1,tau) ynConj(1:end-tau)].*exp(-1i*2*pi*Alpha’*n*Ts),2)/(N-tau));
    RyxTauAlpha(tau+1,:) = (sum(yn.*[zeros(1,tau) xnConj(1:end-tau)].*exp(-1i*2*pi*Alpha’*n*Ts),2)/(N-tau));
    
    end
    
    Now, the resulting matrices get multiplied by the correction factor exp(-1i*pi*Alpha*tau).
    
    Then, I divide the end result, which should be a matrix where each row corresponds to a specific tau and each column corresponds to a specific alpha, into 4 quadrants (with 2 already solved for) so that I can apply symmetry:
    
    % Find index where Alpha = 0
    B = Alpha == 0;
    
    % Define Ryxq3 and Ryxq4, both include Tau,Alpha = 0 row/col.
    Ryxq4 = Ryxq4q3(:,1:find(B));
    Ryxq3 = Ryxq4q3(:,find(B):end);
    
    % Pre-allocate variables for speed
    colAlpha = 1:size(Ryxq3,2);
    maxColAlf = size(Ryxq3,2);
    Rxyq1 = zeros(size(Ryxq3));
    Rxyq2 = zeros(size(Ryxq4));
    
    for col = colAlpha
    
    Rxyq1(:,col) = conj(Ryxq3(:,maxColAlf-col+1));
    
    Rxyq2(:,col) = conj(Ryxq4(:,maxColAlf-col+1));
    end
    
    % We need to flipud Rxyq1 and Rxyq2 so that the index matches -tauMax:0
    Rxyq1 = flipud(Rxyq1);
    Rxyq2 = flipud(Rxyq2);
    
    % Define xCAF
    xCAF = [Rxyq1(1:end-1,1:end-1) Rxyq2(1:end-1,:) ; Rxyq4q3];
    
    When using the above code with the correct definition of x[n] and y[n], I am not able to replicate Fig. 5 and Fig. 6. I am not sure what I am doing wrong with the implementation of the NC cross-CAF and any guidance you can provide would be greatly appreciated!
    
    Thank You!
    
    -John
    
    Reply
    1. Chad Spooner says:
      
      January 24, 2020 at 9:05 am
      
      What I mean is that, according to my understanding, in (13)–(17) it implies that we need to solve for two things in order to properly account for symmetry. We need to solve for Rxy(alpha)(tau) and Ryx(alpha)(tau).
      
      Here are the equations for all readers’ reference:
      
      I think what the equations mean is the following. Suppose we want to compute the XY non-conjugate cross cyclic autocorrelation for all lags $\tau$ and all cycle frequencies $\alpha$ . Can we get the whole function by estimating parts of it, and then using those parts to fill in all the parts where we didn’t actually perform an estimate? The answer is no, in general, you can’t. If you want the whole function, you have to estimate the whole function. If you happen to already have the other cross function (YX), then you can use that to fill in some parts of the XY cross function using the equations.
      
      If you want to verify my symmetry claims, then, yes, you need to estimate both kinds of cross functions to check that the equations indeed hold.
      
      What is $T_s$ ? Is it equal to one? I notice that your time vector $n$ is in samples, but I’m not sure about your cycle-frequency vector Alpha–is it in normalized Hz or physical Hz? Because you insert $T_s$ in
      
      RxyTauAlpha(tau+1,:) = (sum(xn.*[zeros(1,tau) ynConj(1:end-tau)].*exp(-1i*2*pi*Alpha’*n*Ts),2)/(N-tau));
      RyxTauAlpha(tau+1,:) = (sum(yn.*[zeros(1,tau) xnConj(1:end-tau)].*exp(-1i*2*pi*Alpha’*n*Ts),2)/(N-tau));
      
      but not in
      
      Now, the resulting matrices get multiplied by the correction factor exp(-1i*pi*Alpha*tau).
      
      If the first set of equations is consistent, then Alpha is in physical Hz. So the second equation is flawed because $\tau$ is in normalized units of samples, not seconds. If the second equation is consistent, then Alpha is in normalized Hz, and the first set of equations is flawed because you are multiplying a quantity in normalized Hz by one in physical Hz. All is consistent if $T_s = 1$ . But if that is the case, why include it?
      
      I appreciate what you are trying to do with the matrices: take good advantage of MATLAB’s native ability to efficiently do matrix math. But I suggest that you take a step back at this point, and write some code that focuses on a single cycle frequency at a time, just to make sure you’ve got the basic estimation correct. Then, when you do, move on to matrices and efficient estimation of the whole function.
      
      Reply
      1. John Tyler says:
        
        January 24, 2020 at 9:25 am
        
        Hello Dr. Spooner,
        
        In answer to your question(s) above, I am not using normalized frequency or time units. Ts is equal to 1/fs where fs is the sampling frequency in Hz. However, when simulating your rectangular BPSK signal the units are normalized. Either way, to the best of my knowledge, my code is consistent. Originally, I estimated the (non-cross) CAF for one Tau and Alpha at a time and used nested ‘for’ loops. Since then, my code have evolved and I was able to eliminate those nested loops.
        
        For simplicity, I took some shortcuts in my previous comment and wrote the correction factor (as shown above) but in my actual code it is defined as:
        
        for alf = 1:length(Alpha)
        exponent(:,alf) = exp(-1i*pi*Alpha(alf)*Ts*tauLoop’);
        end
        
        where tauLoop = 0:tauMax;
        
        I then do the following multiplications to get the cross (XY and YX) CAF (NC) for quadrants 3 and 4 of the final matrix:
        
        Rxyq4q3 = RxyTauAlpha.*exponent;
        Ryxq4q3 = RyxTauAlpha.*exponent;
        
        Later, I apply the symmetry relations—as explained above—to find quadrants 1 and 2 but I am not sure what I am doing wrong.
        
        Just to make sure I understand the cross CAF symmetry: In order to estimate the NC XY cross-CAF, I need to simulate the XY and YX NC cross-CAFs for Tau = 0:tauMax and Alpha = -alphaMax:dAlpha:alphaMax and then apply the symmetry relations? If so, I am having some difficulty replicating Fig. 5 and Fig. 6.
        
        Again, I would like to thank you for this blog and guidance that you provide!
        
        Thank You,
        
        -John
      2. Chad Spooner says:
        
        January 24, 2020 at 9:51 am
        
        Can you show me the values of $T_s$ and the vector Alpha?
        
        Just to make sure I understand the cross CAF symmetry: In order to estimate the NC XY cross-CAF, I need to simulate the XY and YX NC cross-CAFs for Tau = 0:tauMax and Alpha = -alphaMax:dAlpha:alphaMax and then apply the symmetry relations?
        
        I would say that the symmetry relations imply that if you want to estimate the NC XY cross-CAF, you have to estimate it for all cycle frequencies and delays (lags) of interest–there is nothing free here. If you desire, you can replace some of those calculations with appropriate calculations using the YX cross-CAF. This is unlike some of the other symmetry relations. Take (4b) for example,
        
        Suppose you estimate the CAF for a positive CF $\alpha_0$ . Then what is the CAF for $-\alpha_0$ ? You don’t have to estimate it, you can just find it from that for $\alpha_0$ by time-reversing $R_x^{\alpha_0}(\tau)$ and conjugating it.
John Tyler says:

January 24, 2020 at 10:24 am

Hello Dr. Spooner,

Here are the Taus and Alphas that I used in the simulation for the rectangular BPSK signal:

fs = 1; % Sampling frequency (Hz)
Ts = 1/fs;

tauMax = 20;

Tau = 0:tauMax;

Alpha = -fs:0.05:fs;

Then the symmetry relations in (13)–(17) are used to find the values for negative Tau.

Thank You!

-John

Reply
1. Chad Spooner says:
  
  January 24, 2020 at 3:05 pm
  
  I see you are including CFs that are not exhibited by the signal, which has (for both conjugate and non-conjugate functions) CFs of the form $k/10$ .
  
  Can you verify that your non-conjugate cross CAF for $\alpha = 0.1$ has the same magnitude as mine?
  
  Reply
  1. John Tyler says:
    
    January 27, 2020 at 7:08 am
    
    Hello Dr. Spooner,
    
    I went ahead and plotted the NC XY and YX cross CAFs for Alpha = fbit = 0.1.
    
    The XY NC Cross CAF is viewable here: https://i.postimg.cc/zvjS0DSf/XYx-CAF-alf-0p1-2.jpg
    
    The YX version is plotted here: https://i.postimg.cc/Z5yvqKBy/YXx-CAF-alf-0p1-2.jpg
    
    It looks like I am getting very similar results but the peaks are shifted in the opposite direction. These plots were made with the code posted above.
    
    x[n] is defined as your rectangular BPSK noise-free signal with a carrier of 0.05 and a bit rate of 0.1. y[n] is x[n] with its first 6 samples discarded and the last 6 samples filled in with zeros.
    
    Again, thank you for your help and assistance with CS processing!
    
    -John
    
    Reply
    1. Chad Spooner says:
      
      January 27, 2020 at 9:40 am
      
      John: Thank you for your persistence in this matter! I did find a problem with my cross computations. I think you can replicate all my results if you shift x[n] to the right by six samples to obtain y[n], instead of to the left. I’ve updated the post accordingly. Let me know if you encounter any other differences between your plots and mine.
      
      Reply
John Tyler says:

January 27, 2020 at 12:48 pm

Hello Dr. Spooner,

My pleasure! I just wanted to make sure I was understanding the material correctly.

I have another question:

The graphs from my previous post were created using the following delay product average:

for tau = 0:20

RxyTauAlpha(tau/+1,:) = (sum(xn(tau+1:end).*ynConj(1:end-tau).*exp(-1i*2*pi*Alpha’*n(1:N-tau)*Ts),2)/(N-tau));
RyxTauAlpha(tau/+1,:) = (sum(yn(tau+1:end).*xnConj(1:end-tau).*exp(-1i*2*pi*Alpha’*n(1:N-tau)*Ts),2)/(N-tau));

end

In the exponent, you will notice that n(1:N-tau) is used for the time vector, where N = length(xn) and n = 0:N-1. I am not sure if this is the correct time vector because when n = 0, my code assumes that the delayed signal, ynConj(1:end-tau), begins. However, I am only able to replicate your plots (Fig. 6) when using the code with the time vector as shown in this post. With this code I get the following plots for the NC cross CAF (XY and YX):

https://i.postimg.cc/6pNZGy7R/XYx-CAF-alf-0p1-3.jpg
https://i.postimg.cc/bvy23W9w/YXx-CAF-alf-0p1-3.jpg

When I plot the NC cross CAF (XY and YX) with what I understand to be the ‘correct’ time vector, n(tau+1:N) I get the following plots:

https://i.postimg.cc/L5Lkg7pB/XYx-CAF-alf-0p1-4.jpg
https://i.postimg.cc/8CYLKDkm/YXx-CAF-alf-0p1-4.jpg

which are different than your plots shown in Fig. 6.

Therefore, I am not sure which method is more accurate/correct, n(1:N-tau) or n(tau+1:N)?

In these plots, x[n] is the rectangular BPSK signal without noise with a carrier of 0.05 and a bit rate of 0.1. y[n] is equal to x[n] with its last 6 samples discarded and 6 zeros added to the beginning.

Thank You!

-John

Reply
1. Chad Spooner says:
  
  January 28, 2020 at 5:34 pm
  
  Two things bother me here.
  
  The first is that is appears you are computing the two cross CAFs (XY, YX) for the same cycle frequency, $\alpha = 0.1 = f_{bit}$ . But Figure 6 in the post considers XY with $\alpha = f_{bit}$ and YX with $\alpha = -f_{bit}$ . What happens when you do that? Also, how do your plots look for $\alpha=0$ , as in Figure 5?
  
  The second is the discontinuity near $\tau = 0$ in your fourth plot. I can’t find anything like that in my plots (posted or otherwise). Maybe we should figure that out first?
  
  Reply
  1. John Tyler says:
    
    January 29, 2020 at 5:46 am
    
    Hello Dr. Spooner,
    
    I apologize for flooding this topic with posts but I just want to make sure that I understand how to apply symmetry for the cross-CAFs.
    
    I am attempting to write ‘one algorithm to rule them all,’ where I can have a function that can take some signal inputs and compute the C/NC CAF/SCF and C/NC Cross-CAF/SCF.
    
    The reason for computing two NC cross-CAFs (XY and YX) in my code is so that I can apply the symmetry relations (13)–(17). Unless I am missing something, is this not the way to compute the cross CAF for Tau’s and Alpha’s of interest, especially when those values of interest are located on all 4 quadrants?
    
    In Fig. 6 you listed the XY cross-CAF with Alpha = +fbit and the YX Cross-CAF with Alpha = –fbit, is this a typo?
    
    I plotted the data you requested in your previous comment (XY Alpha = +fbit, YX Alpha = –fbit) using the time vector n(1:N-tau):
    
    https://i.postimg.cc/MKdLQFTk/XYx-CAF-alf-0p1-5.jpg
    https://i.postimg.cc/tCY3BSkG/YXx-CAF-alf-neg0p1-5.jpg
    
    The discontinuity in the 4th plot (in my previous comment) only occurs when I use a different time vector n(tau+1:N). It may be due to the concatenation of the matrices (from 4 quadrants) after applying symmetry.
    
    I have another question regarding the C Cross-CAF. I am not sure that I implemented the symmetry relations in (18) and (19) correctly because I continue to get plots similar to Figs. 7 and 8 but the real and imaginary parts are inverted.
    
    Is there a reference that you can point me to that describes the symmetrical relations for the C/NC Cross-CAFs?
    
    Again, thank you for taking the time to answer my questions and I really appreciate all that you have done with this blog along with your patience for answering questions!
    
    Thank You!
    
    -John
    
    Reply
    1. Chad Spooner says:
      
      February 2, 2020 at 12:30 pm
      
      No apology needed!
      
      I am attempting to write ‘one algorithm to rule them all,’ where I can have a function that can take some signal inputs and compute the C/NC CAF/SCF and C/NC Cross-CAF/SCF.
      
      I like that–I have a C program that does all second-order auto and cross statistics in the time and frequency domains. Very handy.
      
      The reason for computing two NC cross-CAFs (XY and YX) in my code is so that I can apply the symmetry relations (13)–(17). Unless I am missing something, is this not the way to compute the cross CAF for Tau’s and Alpha’s of interest, especially when those values of interest are located on all 4 quadrants?
      
      Yes, this is the point that you and I have been circling around for a while. (You might want to review my previous replies.) I don’t think the symmetry for the non-conjugate cross CAF is useful in an estimator–where “useful” means “saves computations.” Suppose you want to compute the XY cross CAF. Start out by estimating it for all positive cycle frequencies $\alpha$ and all $\tau$ of interest. Now we want to obtain it cheaply, if possible, for negative $\alpha$ . The formulas imply that you can obtain the XY cross CAF for negative $\alpha$ from the YX cross CAF for positive $\alpha$ . So you have to estimate the YX cross CAF. But if you are going to do that, why not just estimate the XY cross CAF for the negative $\alpha$ ? You aren’t saving anything. Do you agree?
      
      In Fig. 6 you listed the XY cross-CAF with Alpha = +fbit and the YX Cross-CAF with Alpha = –fbit, is this a typo?
      
      I haven’t been able to find any further problems. I’ll look some more in the coming week.
      
      I plotted the data you requested in your previous comment (XY Alpha = +fbit, YX Alpha = –fbit) using the time vector n(1:N-tau):
      
      I was hoping to compare your functions for Figure 5 to mine ( $\alpha=0$ ).
      
      I have another question regarding the C Cross-CAF. I am not sure that I implemented the symmetry relations in (18) and (19) correctly because I continue to get plots similar to Figs. 7 and 8 but the real and imaginary parts are inverted.
      
      Are the magnitudes the same as mine? I think maybe there is still some difference in how the delay is applied to either $x(t)$ or $y(t)$ . The important thing from my point of view is that for each of your comparisons, do the obtained estimates obey the mathematical symmetries I’ve posted, or not?
      
      Is there a reference that you can point me to that describes the symmetrical relations for the C/NC Cross-CAFs?
      
      For the NC cross CAF, I’ve posted the sequence of equations that result in the symmetry relationships (13)-(17). Do you mean to ask for the similar sequence of equations for (18) and (19)? I have them. I can add them to the post. I’m not sure what else there is to refer to, but you might want to check out Professor Napolitano’s books.
      
      Reply
      1. John Tyler says:
        
        February 3, 2020 at 8:39 am
        
        Hello Dr. Spooner,
        
        You’re right! I am not really saving computation time. Rather, I have the ability to directly estimate the CAF for all Alpha’s but only positive Tau’s (including zero). The reason I am using symmetry is calculate the –Tau portion of the CAF.
        
        With Fig. 6 I meant that your version of Fig. 6 has the XY CAF for Alpha = –fbit and the YX CAF for Alpha = +fbit. Whereas I can only reproduce Fig. 6 when solving for XY CAF for Alpha = +fbit and YX CAF for Alpha = –fbit. What I am trying to say is that it seems that the signs on your version of Fig. 6 are inverted, which leads me to think that it’s a typo.
        
        I went ahead and plotted Figs. 1–8 using my algorithm of the Cross CAF. The delay products are all of the form:
        
        (sum(xn(tau+1:N).*ynConj(1:N-tau).*exp(-1i*2*pi*Alpha’*n(1:N-tau)*Ts),2)/(N-tau));
        
        where the conjugation is changed based on the requested CAF (C vs. NC).
        
        Here are the links to my version of Figs. 1–8:
        
        https://i.postimg.cc/sXFcxqLB/NC-x-CAF-Fig1.png
        https://i.postimg.cc/3Rv3C49S/NC-x-CAF-Fig2.png
        https://i.postimg.cc/28JTg1xF/C-x-CAF-Fig3.png
        https://i.postimg.cc/FHzzkRtq/C-x-CAF-Fig4.png
        https://i.postimg.cc/fyS5bZz0/NC-x-CAF-Fig5.png
        https://i.postimg.cc/Mpp5cypS/NC-x-CAF-Fig6.png
        https://i.postimg.cc/gkSRT1mC/C-x-CAF-Fig7.png
        https://i.postimg.cc/gkQjRpY2/C-x-CAF-Fig8.png
        
        In Fig. 3 I get the conjugate for Alpha = 0 the other plots are the same.
        
        For Figs. 7 and 8, the symmetry is consistent with the mathematical symmetry relations, but the real and imaginary parts are multiplied by –1 when compared to your plots.
        
        Thank You!
        
        -John
      2. Chad Spooner says:
        
        February 3, 2020 at 2:04 pm
        
        I went ahead and plotted Figs. 1–8 using my algorithm of the Cross CAF. The delay products are all of the form:
        (sum(xn(tau+1:N).*ynConj(1:N-tau).*exp(-1i*2*pi*Alpha’*n(1:N-tau)*Ts),2)/(N-tau));
        
        The functions peak where they should, auto and cross. The shapes of the functions’ magnitudes match mine in all cases. The only differences have to do with the real and imaginary parts, which are rotated relative to mine. I think all of this boils down to the relationship between the $t$ vector in the exponential of the lag product and the symbol-clock phase in xn, and the carrier phase in xn.
        
        A straightforward implementation of the asymmetric cross CAF
        
        $\displaystyle \frac{1}{T} \int x(t) y^*(t-\tau) e^{-i2 \pi \alpha t} \, dt$
        
        would be
        
        $\displaystyle \frac{1}{N} xn(1:\mbox{\rm end}) \mbox{\rm circshift}(ynConj, \tau) e^{-i 2 \pi \alpha [1:N]}$
        
        Then you’d want to multiply by the phase factor mentioned in the Cyclic Autocorrelation post to get to the symmetric version.
        
        So at this point I think the relationship between the phase of the applied exponential in the sum and the phase of the sine wave in the quadratic functional involving $x(t)$ and $y(t)$ is simply not the same as the phase relationship the way I did it.
        
        Does that help?