hello
I tried 2 simple methods :
- sliding window averaging : not very effective
- to fit an exp model : will almost give the exact original signal
second option seems to be the best (simple) solution . i did not test the rest of the code as I miss some Tbx
t = 0:2^8;
N = length(t);
tau = 10; % time scale of relaxation
x_asym = 4; % relative increase of the observable x t -> infinity
x_r2s = x_asym*(1 - exp(-t/tau)); % uncorrupted sought-for signal
sig = 0.05; % noise strength
% Use AR(1) to corrupt the signal instead of white noise.
phi = 0.6; % AR coefficeint
xi1 = sig*randn(1,N);
xi2 = sig*randn(1,N);
pert1 = filter(1,[1 -phi],xi1);
pert2 = filter(1,[1 -phi],xi2);
x_prd_1 = x_r2s+pert1;
x_prd_2 = x_r2s+pert2;
% denoising using sliding avg method
N = 10;
x_prd_1savg = myslidingavg(x_prd_1, N);
% denoising using exponential fit method
% code is giving good results with template equation : % y = a.*(1-exp(b.*(x-c)));
x = t;
y = x_prd_1;
f = @(a,b,c,x) a.*(1-exp(b.*(x-c)));
obj_fun = @(params) norm(f(params(1), params(2), params(3),x)-y);
sol = fminsearch(obj_fun, [y(end),0,0]);
a_sol = sol(1);
b_sol = sol(2);
c_sol = sol(3);
x_prd_1fit = f(a_sol, b_sol,c_sol, x);
figure
plot(t,x_r2s,'b',t,x_prd_1,'r',t,x_prd_1savg,'c',t,x_prd_1fit,'k');
legend('signal','signal+noise','after sliding avg','exp fit');
function out = myslidingavg(in, N)
% OUTPUT_ARRAY = MYSLIDINGAVG(INPUT_ARRAY, N)
%
% The function 'slidingavg' implements a one-dimensional filtering, applying a sliding window to a sequence. Such filtering replaces the center value in
% the window with the average value of all the points within the window. When the sliding window is exceeding the lower or upper boundaries of the input
% vector INPUT_ARRAY, the average is computed among the available points. Indicating with nx the length of the the input sequence, we note that for values
% of N larger or equal to 2*(nx - 1), each value of the output data array are identical and equal to mean(in).
%
% * The input argument INPUT_ARRAY is the numerical data array to be processed.
% * The input argument N is the number of neighboring data points to average over for each point of IN.
%
% * The output argument OUTPUT_ARRAY is the output data array.
if (isempty(in)) | (N<=0) % If the input array is empty or N is non-positive,
disp(sprintf('SlidingAvg: (Error) empty input data or N null.')); % an error is reported to the standard output and the
return; % execution of the routine is stopped.
end % if
if (N==1) % If the number of neighbouring points over which the sliding
out = in; % average will be performed is '1', then no average actually occur and
return; % OUTPUT_ARRAY will be the copy of INPUT_ARRAY and the execution of the routine
end % if % is stopped.
nx = length(in); % The length of the input data structure is acquired to later evaluate the 'mean' over the appropriate boundaries.
if (N>=(2*(nx-1))) % If the number of neighbouring points over which the sliding
out = mean(in)*ones(size(in)); % average will be performed is large enough, then the average actually covers all the points
return; % of INPUT_ARRAY, for each index of OUTPUT_ARRAY and some CPU time can be gained by such an approach.
end % if % The execution of the routine is stopped.
out = zeros(size(in)); % In all the other situations, the initialization of the output data structure is performed.
if rem(N,2)~=1 % When N is even, then we proceed in taking the half of it:
m = N/2; % m = N / 2.
else % Otherwise (N >= 3, N odd), N-1 is even ( N-1 >= 2) and we proceed taking the half of it:
m = (N-1)/2; % m = (N-1) / 2.
end % if
for i=1:nx, % For each element (i-th) contained in the input numerical array, a check must be performed:
dist2start = i-1; % index distance from current index to start index (1)
dist2end = nx-i; % index distance from current index to end index (nx)
if dist2start<m || dist2end<m % if we are close to start / end of data, reduce the mean calculation on centered data vector reduced to available samples
dd = min(dist2start,dist2end); % min of the two distance (start or end)
else
dd = m;
end % if
out(i) = mean(in(i-dd:i+dd)); % mean of centered data , reduced to available samples at both ends of the data vector
end % for i
end
