clc;
clear;
name='*.txt';
filename=[path,name];
[FileSTD_Name,path]=uigetfile(filename,'打開Sensor數據文件');        
filename=[path,FileSTD_Name];                                       
[fid,message]=fopen(filename,'rt');
[Sensor,countSTD]=fscanf(fid,'%f',inf);
status=fclose(fid);
%temp=zeros(50000,1);
%Sensor=[Sensor;temp];
Sc=50;
sensor=Sensor;
len=size(sensor);
len=len(1,1);
%sensor=sensor(1:5:len);
%數據的五點三次平滑處理
for a=1:40
sensor(1) = (69 * sensor(1) + 4 * sensor(2) - 6 * sensor(3) + 4 * sensor(4) - 5 * sensor(5)) / 70;
sensor(2) = (2 * sensor(1) + 27 * sensor(2) + 12 * sensor(3) - 8 * sensor(4) + 2 * sensor(5)) / 35;

for i = 3 :(len-2)
      sensor(i) = (-3 * sensor(i - 2) + 12 * sensor(i - 1) + 17 * sensor(i) + 12 * sensor (i + 1) - 3 * sensor(i + 2)) / 35;
end
sensor(len - 1) = (2 * sensor (len - 4) - 8 * sensor (len - 3) + 12 * sensor (len - 2) + 27 * sensor (len - 1) + 2 * sensor (len)) / 35;
sensor (len) = (-1 * sensor (len - 4)+ 4 * sensor (len - 3) - 6 *sensor (len - 2) + 4 * sensor (len - 1) + 69 * sensor (len )) / 70;
end
plot(sensor),grid on;
xlabel('point');ylabel('Voltage(V)');title('截取有用信號');
zoom xon;
pause;
[x1,y1]=ginput;
[x2,y2]=ginput;
x1=round(x1);
%%%%%%%%%%%%%%%%%%
%x1=1;
%%%%%%%%%%%%%%%%%%%
x2=round(x2);
len=x2-x1+1;
sensor1=sensor(x1:x2)*1;
plot(sensor1);
grid on;
xlabel('point');ylabel('Voltage(V)');title('截取基線');
%Detrend
[x1,y1]=ginput;
[x2,y2]=ginput;
x1=round(x1);
x2=round(x2);
ave=mean(sensor1(x1:x2));
sensor_Partly=sensor1-ave;
%補零
temp=zeros(50,1);
sensor_Partly=[sensor_Partly;temp];
len=size(sensor_Partly);
len=len(1,1);
%
t=(1:len)'*0.01; %Variable
sensor_Partly=sensor_Partly*1;%Variable%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
sensor_Partly=sensor_Partly;%%%%%%%%%%%%%%%%%%%%%%%%%%%
plot(t,sensor_Partly);
grid on;
xlabel('Time(us)');ylabel('Acceleration/g');
zoom xon;
pause;
%Spectrum before filtered
fre=abs(fft(sensor_Partly,len));
fre=fre/fre(1,1);
f=100000/len*(1:len); %Variable
plot(f,fre);
grid on;
xlabel('Frequency(KHz)');
ylabel('Amplitude');
pause;
plot(f(1:round(len/300)),fre(1:round(len/300))); %Variable
grid on;
xlabel('Frequency(KHz)');ylabel('Amplitude');
pause;
title('選擇濾波截止頻率');
%功率普
%Psd=fre.^2/len;
%I=0:len-1;
%plot_Psd=10*log10(Psd(I+1));
%plot(f(1:round(len/300)),plot_Psd(1:round(len/300)));
%grid on
%pause;
%Sampling
sensor_Partly=sensor_Partly(1:100:len); %Variable
%Filter design
[xp,yp]=ginput;
Fstop=xp(1,1);
fm=500; %Variable
fe=Fstop/fm;
f=[0 fe fe 1];
m=[1 1 0 0];
b=fir2(101,f,m); %Variable
sensor_Partly=filtfilt(b,1,sensor_Partly);
%The filter frequency response
[h,w]=freqz(b,1,128);
plot(f*fm,m,w/pi*fm,abs(h));
grid on;
xlabel('Frequency(KHz)');ylabel('Amplitude');title('Frequency Response of The Filter');
pause;
len=size(sensor_Partly);
len=len(1,1);
t=(1:len)'*1; %Variable
plot(t,sensor_Partly);
zoom xon;
grid on;
xlabel('Time(us)');ylabel('Acceleration/g');title('');
Fstop=num2str(Fstop);
lab='K濾波后曲線';
lab=[Fstop,lab];
title(lab);
pause;
%Detrend
title('截取基線');
[x1,y1]=ginput;
[x2,y2]=ginput;
x1=round(x1);
x2=round(x2);
ave=mean(sensor_Partly(x1:x2));
sensor_Partly=sensor_Partly-ave;
sensor_Partly=sensor_Partly;
plot(t,sensor_Partly);
grid on;
xlabel('Time(us)');ylabel('Voltage(V)');
pause;
Max_A=max(sensor_Partly);
Max_A=num2str(Max_A);
lab=['Peak Voltage: ',Max_A,'V'];
title(lab);
pause;
%Specturm after filtered
fre=abs(fft(sensor_Partly,len));
f=1000/len*(1:len); %Variable
plot(f(1:round(len/3)),fre(1:round(len/3))); %Variable
grid on;
xlabel('Frequency(KHz)');ylabel('Amplitude');title('Specturm After Filtered');
pause;
%integration Range
plot(t,sensor_Partly);
zoom xon;
grid on;
xlabel('Time(us)');ylabel('Voltage(V)');title('選擇積分區間');
pause;
zoom xon;
grid on;
[x1,y1]=ginput;
[x2,y2]=ginput;
x1=round(x1);
x2=round(x2);
sensor_Partly=sensor_Partly(x1:x2);
%Velocity
Velocity=cumsum(sensor_Partly); %Variable
len=x2-x1+1;
t=(1:len)*1;
plot(t,Velocity);
grid on;
xlabel('Time(us)');ylabel('Velocity(m/s)');
Max_V=max(Velocity);
Max_V=num2str(Max_V);
lab=['Peak Velocity: ',Max_V,'m/s'];
title(lab);
Max_V=str2num(Max_V);
Max_E=23.5789;
Se=Max_V./Max_E;
%Se=Se*9.8*0.544*2e4/1e6 %1e5代表1e5g/v電荷放大器 yadian
%Se=Se*9.8*1e6./10/1e6;
Se=Se*9.8*0.048*1e5/1e6;