%download the auditory 
%modelling toolbox (http://amtoolbox.sourceforge.net)
addpath('~/bin/amtoolbox-full-0/code') % adjust this
amt_start
%%
%code for autocorrelation-based pitch analysis in auditory models

clc;close all;clear;
%some parameters for the auditory model
fLow = 50; % the lowest characteristic frequency of the filter bank
fHigh = 20000; % and the highest
fCut = 1000; % cut-off frequency of the lowpass filter
%load a speech file, the example contains 27 ms of /a/ vowel
%[sample,fs] = audioread('kaksi.wav',[5900 7200]); % CA2015: load your audio file in whole length in here
fs=48000;


%%
% select the file you are reading in following
[sample] = audioread('02.wav',[fs*1,fs*2]); 

%sample=sample(10000:11000,1); %20 ms
sample=sample(10000:20000,1); %200 ms
% this reads one second of sound btw 1 and 2 seconds in your file
% pls adjust the position of the file if needed
%%

sampleLen = length(sample);

%create of a gammatone filter bank using a command from the auditory
%modeling toolbox (http://amtoolbox.sourceforge.net)

[b,a] = gammatone(erbspacebw(fLow,fHigh),fs,'complex');
%processing the signal through the filter bank
filterOut = real(ufilterbankz(b,a,sample));

%%

%emulation of the neural transduction with half-wave rectification and
%lowpass filtering of the filter bank output
rectified = filterOut.*(filterOut>0);
%a first-order IIR filter is used as the lowpass filter
beta = exp(-fCut*2*pi/fs);
outSig = filter(1-beta,[1 -beta],rectified);

for freqInd=1:size(outSig,2) % autocorrelation for each band
    auCorr(:,freqInd)= xcorr(outSig(:,freqInd),outSig(:,freqInd),1000,'coeff')';
end
sumAuCorr=sum(auCorr((end-1)/2:end,:)');


% You should now plot the results out
% The "auditory spectrogram", or some kind of estimation of timbre can 
% be plotted from outSig variable. For example the command "imagesc" could work. 
% Some kind of periodicity analysis for pitch and tonality studies can be
% viewed from sumAuCorr.