function data=q_ica_filter(data,trail_t,fs,frq,kept_c_num)
% 命令格式為：data=q_ica_filter(data,trail_t,fs,frq,kept_c_num)
% 使用ica對多導數據分解后只選取frq頻率成份最高的kept_c_num個分量重建信號
% 輸入參數： data，輸入數據，每行為一導信號
%           trail_t,單次trail的持續時間
%           fs,采樣頻率
%           frq，保留頻段，格式為一維數組
%           kept_c_num，保留分量個數
% 輸出參數： data，重構信號
%
% 注：      輸入數據為n個trail順序連接而成。

l=length(data);
chan_num=size(data,1);
trail_n=floor(l/trail_t/fs);

data_new=[];
for i=1:trail_n
    i% one trail
    temp_index=[1:trail_t*fs]+(i-1)*trail_t*fs;
    data_temp=data(:,temp_index);
    % ica process
    [ww,ss]=runica(data_temp,'extended',-1);
    ica_temp=ww*ss*data_temp;
    
    % get rid of components with low 10hz psd
    for j=1:chan_num
        psd_temp(j,:)=q_psfd_cal(ica_temp(j,:),frq,fs);
    end
    psd_amp=mean(abs(psd_temp),2);  % psd mean value of each component in frq
    [a order_ind]=sort(psd_amp);    % sort from lowest to highest
    min_order_ind=order_ind(1:chan_num-kept_c_num);   % low components no
    ica_temp(min_order_ind,:)=0*ica_temp(min_order_ind,:); % kept high psd components
    
    % data reconstuction 
    data_temp=inv(ww*ss)*ica_temp;
    data_new=[data_new,data_temp];
end

data=data_new;