Deconvolution creates new peaks

Hi, I'm trying to remove some noise on an accelerometer signal, this is mounted on the back of a hammer and I would ideally position it on the front to measure the acceleration of the impact. Both the accelerometer and the load cell signals have seven impacts. To remove the influnce of the body of the hammer I tryed to estimate the transfer function using tfestimate and then performing a deconvolution of the signal. This seems to work at first but zooming on the impacts, some peaks have doubled, like the fourth one:

this is my code:

load('accCorrectedReduced.mat')
load('forceCorrectedReduced.mat')
L = length(acc);
fs=51201;

%% deconvolution

nfft = L;
[Txy, f] = tfestimate(force, acc, hann(nfft), 0.7*fs, nfft, fs);

Y = fft(acc, nfft);


TxyFull = [Txy; conj(Txy(end-1:-1:2))];

fFull = linspace(0,fs,nfft);

TxyInterp = interp1(f, Txy, fFull(1:length(f)), 'linear', 'extrap');
TxyFull = [TxyInterp, conj(TxyInterp(end-1:-1:2))];


X_estimated = Y(1:end-1) ./ TxyFull.';


x_deconvolved = ifft(X_estimated);


figure;
plot( acc);
hold on
plot(-abs(x_deconvolved));

Answer 1

[Paul]

Hi Riccardo,

It might be simpler to just use the twosided option of tfestimate

Also, the input data is single, might want to convert to double for processing.

load('accCorrectedReduced.mat')
load('forceCorrectedReduced.mat')
whos

  Name             Size              Bytes  Class     Attributes

  acc        2119721x1             8478884  single              
  ans              1x47                 94  char                
  force      2119721x1             8478884  single              

L = length(acc);
fs = 51201;
%% deconvolution
nfft = L;
[Txy, f] = tfestimate(force, acc, hann(nfft), 0.7*fs, nfft, fs,'twosided');
Y = fft(acc, nfft);
%TxyFull = [Txy; conj(Txy(end-1:-1:2))];
%fFull = linspace(0,fs,nfft);
%TxyInterp = interp1(f, Txy, fFull(1:length(f)), 'linear', 'extrap');
%TxyFull = [TxyInterp, conj(TxyInterp(end-1:-1:2))];
%X_estimated = Y(1:end-1) ./ TxyFull.';
X_estimated = Y./Txy;
x_deconvolved = ifft(X_estimated);
figure;
plot(acc);
hold on
%plot(-abs(x_deconvolved));
plot(-x_deconvolved)

Zoom in on each peak

[~,locs] = findpeaks(x_deconvolved,'MinPeakHeight',20,'MinPeakProminence',10);
numel(locs)

ans = 7

n = 1:numel(acc);
for ii = 1:numel(locs)
    figure
    index = (-1000:1000) + locs(ii);
    plot(n(index),acc(index))
    hold on
    plot(n(index),-x_deconvolved(index))
end

[Riccardo]

Thank you very much this seems to be working