Problems with discrete control simulations

Dear Matlab community:

I'm not experienced with the use of Matlab. For laboratory work at the university I need to do system identification from some given data. I recieve data in four columns: time, imput with removed offset, actual input and actual output. I need do remove the offset from the actual output, and using the half of both input and output with no offset data, find the values of the model parameters using the Linear Least Squares method. Up to here I think I made things right. To check how good is the estimation, i have to check the response of the calculated system using ths filter command, and here comes my problem: I only obtain as output an impulse (from zero to infinite) which does nothing to do with the input signal. I suspect the mistake is related with the filter command, but i cannot see where is it. I really appreciate if you can help me with some advice to find the problem.

clc
clear
load Student_data.mat
vtime=Log(:,1);% takes time data from Log
per=(vtime(2)-vtime(1))/1000 %calculates the period
time=[0:per:per*length(vtime)-per]; %creates a new time vector starting from zero
not_off_input=Log(:,2); %takes input with no offset data from Log
off_input=Log(:,3);%takes input with offset data from Log
off_output=Log(:,4);%takes ooutput with offset data from Log
plot(time,off_input); hold on% 
plot(time, off_output);hold on%
avv=mean(off_output(1:168))%calculates the average of the output data till the 168th element
not_off_output=off_output-avv;% substracts the mean to obtain the output data with no offset
figure
plot(time,not_off_input);hold on
plot(time, not_off_output);hold on
col1=not_off_output(5:(length(not_off_output)/2)-1); % first column of the leas squares parameters matrix, from the 5th to the half-matrix element 
col2=not_off_output(4:(length(not_off_output)/2)-2);
col3=not_off_input(5:(length(not_off_output)/2)-1);
col4=not_off_input(4:(length(not_off_output)/2)-2);
psi=[col1 col2 col3 col4]; 
y=not_off_output(6:length(not_off_output)/2);%otputs matrix
theta=y\psi;% matrix of errors
a1=-theta(1);
a2=-theta(2);
b1=theta(3);
b2=theta(4);
figure
out_fil=filter([0 b1 b2],[1 a1 a2],not_off_input);% use of filter command to simulate the data 
plot(time,out_fil);hold on
clc
clear

Attached: code, input signal and final messy response. Thanks very much for your help.

Best regards,

Ander