function [data_filtered] = MAP_filter_LOUIS(data0,J)
% This a function from Louis Marie's source code of the Maximum A
% Posteriori filtering
% It will be faster and easier to work with it than insert its function in
% our pipeline.

%  Input:
%data0: data for filtering
%J: coefficient for the MAP filtering

%  Output:
% data_filtered: the data filtered using the MAP filter



data0_new=data0(~isnan(data0));


NRA=length(data0_new);
SMM=sparse(NRA,NRA);
SMM=spdiags([ones(NRA,1),-4*ones(NRA,1),(6+J)*ones(NRA,1),-4*ones(NRA,1),ones(NRA,1)],[-2,-1,0,1,2],SMM);

SMM(1,1)=SMM(1,1)+SMM(1,2)+SMM(1,3);
SMM(2,1)=SMM(2,1)+SMM(3,1);

SMM(end  ,end)=SMM(end,end)+SMM(end-1,end)+SMM(end-2,end);
SMM(end-1,end)=SMM(end-1,end)+SMM(end-2,end);

data_filtered=SMM\data0_new*J;

end

