# Overall description

I am trying to estimate a filtering system’s transfer function, given its input and output. This system takes $$x$$ as input . This signal is low pass filtered and added to a WGN by the system, which gives the output signal $$y$$.

\begin{align} x &= \text{BPSK} + b_1 \\ y &= h*x + b_2 \text, \end{align}

$$b_1$$ and $$b_2$$ being white gaussian Noise

## Hypothesis

• $$b_1\ll b_2 \ll x$$
• $$b_1, b_2$$ independent

The purpose of my manipulation is to estimate as best as possible the filter’s frequency response amplitude ($$|H|$$) and compare it with the applied filtering (which is known – red plot on next figures).

I had some trouble with this task and it would be very helpful if someone could help with an additional signal conditioning or identify an algorithmic issue.

# Algorithm

Theoretically, a simple estimator based on PSD (Power Spectral Density) of input and output signals can be performed with following equations:

$$H(\omega) = \frac{S_{xy}(\omega)}{S_x(\omega)}$$

with $$S_{xy}(\omega) = \text{FFT}\left(R_{xy}(\tau)\right) = \text{FFT}(\mathtt{xcorr}(x,y))$$ and $$S_{x}(\omega) = \text{FFT}\left(R_{x}(\tau)\right) = \text{FFT}(\mathtt{xcorr}(x))$$.

For a better signal to noise ratio, the signals ($$x$$ and $$y$$) are split into bursts to compute a mean FFT.

# Results

The output of this estimator is also given hereafter (Estimated TF = estimated transfer function). It can be seen that the estimator is noisy.

In this step I have tried to enhance the estimator by adding median filter on the transfer function estimate. However, this did not improve the estimator.

It seems that the Frequency response amplitude estimate shows spectral hollows compared to the input filtering, which are inducing an important error on the estimate.

# Questions

• What kind of filter (or another way to estimate the transfer function) can I use?
• Could a signal conditioning on $$x$$ or $$y$$ enhance the estimator accuracy?

The following Matlab function is used to compute the TF estimation. Please note that nb_means parameter is arbitrary equal to 16 in my computation, this parameter can evolve.

function [ H1, G1 ] = estimateTF_DSP( sigX, sigY, nb_means )
if nb_means > 1
sigX_I_reshaped = reshape(real(sigX),length(real(sigX))/nb_means,nb_means)';
sigX_Q_reshaped = reshape(imag(sigX),length(imag(sigX))/nb_means,nb_means)';
sigY_I_reshaped = reshape(real(sigY),length(real(sigY))/nb_means,nb_means)';
sigY_Q_reshaped = reshape(imag(sigY),length(imag(sigY))/nb_means,nb_means)';

for k=1:size(sigX_I_reshaped,1);
sig_X = sigX_I_reshaped(k,:)+1i*sigX_Q_reshaped(k,:);
sig_Y = sigY_I_reshaped(k,:)+1i*sigY_Q_reshaped(k,:);
% XY
Rxy = xcorr(sig_X,sig_Y);
Sxy(k,:) = fft(Rxy);
% X
Rx = xcorr(sig_X);
Sx(k,:) = fft(Rx);
end

Sx_m = mean(Sx);
Sxy_m = mean(Sxy);

else
% XY
Rxy = xcorr(sigX,sigY);
Sxy_m = fft(Rxy);
% X
Rx = xcorr(sigX);
Sx_m = fft(Rx);
end

H1 = (Sxy_m./Sx_m)./(2*pi);
G1 = 20*log10(abs(H1));
end

• Hi Math! Your notation is a bit strange at times; what's BPSK(1)? When you write $H(x)$, do you mean $h*x$ or do you mean $H(\omega = x)$? You might want to write something like $E(|b_1|)\ll E(|b_2|)$, because for WGN, $b_1 \ll b_2$ cannot be true in general. – Marcus Müller Jan 29 at 15:41
• And: this is a nice post! But, you forgot to actually ask a question. Make it as precise as possible without being overspecific, please. – Marcus Müller Jan 29 at 15:42
• Hi Marcus ! Thanks for your remarks ! So, here are my answers : - BPSK(1) stands for basically BPSK signal - In fact, H(x) stands for h*x - You're right about WGN, thanks ! Finally, I forgot to paste the essential part of my question... So I have edited my post with your remarks ! Thanks again @MarcusMüller! – Math Jan 29 at 16:26
• Can't you just use the tfestimate function? – fibonatic Jan 30 at 14:56
• I already tried this function, but with no satisfying results. My estimator seems accurate for signals with equal C/N0 (different noise draw are not affecting results). With different C/N0, the resulting estimation is really noisy... – Math Jan 30 at 15:59