%% CCA than FFT
%input: data (time*ch)
%frequency: data fundamental frequency (in vector)
%harmonics: number of harmonics
%targetFreq: target frequency

%output: one dimension snr vector over frequency
%ex: out(8:0.2:11.8), out(i) =  (power of frequency(i)) / mean(power of freqnecy(!i))

function [power] = myCCA2_FFT_SNR(input, resolution, frequency, harmonics)
%column major
template = zeros(size(input,1),harmonics*2);
time = size(input,1);
T = 1/time * resolution;                     % Sample time
L = size(input,1);                     % Length of signal
t = (0:L-1)*T;                % Time vector
index=1;
corrVector = zeros(1,length(frequency));

for f = 1:length(frequency)
    for harm=1:(harmonics)
        template(:,index) = sin(2*harm*pi*frequency(f)*t);
        template(:,index+1) = cos(2*harm*pi*frequency(f)*t);
        index= index+2;
    end
    [A,B,r,U,V] = canoncorr(input,template);
    %corrVector(f) = r(1);
    tempPower(f,:) = abs(fft(U(:,1)));
    index=1;
    template=[];
end
%% snr power ratio of target/nonTarget
for f = 1:length(frequency)
    nonTargetIndex = int16(freq2Index_2(frequency(~ismember(frequency,frequency(f)))));
    
    power(f) = tempPower(f,int16(freq2Index(frequency(f))))/...
                           mean(tempPower(f,int16(freq2Index(frequency(nonTargetIndex)))));
 
end

end
function [index] = freq2Index_2(hz)
index = (hz-8)*5+1;
end
function [hz] = index2Freq(index)
hz = (index-1)/5;
end
function [index] = freq2Index(hz)
index = (hz)*5+1;
end