CB-Frequency-Coding-of-Motor-Kinematics / Supporting Functions / find_freqmax.m
find_freqmax.m
Raw
function [freqmax,valmax]=find_freqmax(type,data,Ch,timeframes,n_method,n_sumwindow,extract_range)
% type: 1: raw continuous .ns6 (actually loads only the 'data' matrix) or 2: data from TFplot(.mat) 
% Ch: 1 or 2 if it's from TFplot .mat, index of the channel in the data matrix if it's from raw data
% n_method: normalization method for peak amplitude
           % 0: none, ignores n_sumwindow
           % 1: normalized by its baseline defined by n_sumwindow
% n_sumwindow: freq range for self normalization 
% extract_range: range for finding peak
if type==1
    display('** should load .ns6 file (fs=30khz), downsampled to 1khz **')
    [r,c]=size(data);
    if r>c
        data=data';
    end
    timeframes=timeframes*30000;
    data=data(Ch,timeframes(1):timeframes(2));
    data=decimate(data,30);
    Y=fft(data);
    Fs = 1000;                                                
    L = length(data);
    P2 = abs(Y/L);
    P1 = P2(1:L/2+1);
    P1(2:end-1) = 2*P1(2:end-1);
    f = Fs*(0:(L/2))/L;
%     figure()
%     plot(f,P1)
%     xlim([0 60])
    if n_method==1
        [~,idx1]=min(abs(f-n_sumwindow(1)));
        [~,idx2]=min(abs(f-n_sumwindow(2)));
        base=mean(P1(idx1:idx2));
        P1=P1/base;
        [~,idx1]=min(abs(f-extract_range(1)));
        [~,idx2]=min(abs(f-extract_range(2)));
        [a,b]=max(P1(idx1:idx2));
        freqmax=f(b+idx1-1);
        valmax=P1(b+idx1-1);
    end
elseif type==2    
    display('not finished yet')
end
end