[S,fs] = audioread('/home/kuortti/work/matlab/PS-4-scramble.wav');
N = length(S);
fS = fft(S);
% Nyquist frequency is half the sample rate: each frequency bin is of width
% fs/N, where N is the number of bins
df = [0,(1:N/2-1)*(fs/N)];
% take care of the complex conjugate
dff = [-fliplr(df),df];
plot(dff,abs(fS))
ffS = fS(1:N/2);
bandA = df<1000;
bandB = df>=1000&df<2000;
bandC = df>=2000&df<=3000;
bandD = df>=3000&df<4000;
fbandA = ffS(bandA);
fbandB = ffS(bandB);
fbandC = ffS(bandC);
fbandD = ffS(bandD);
ffS(bandA) = fbandC;
ffS(bandB) = fbandB;
ffS(bandC) = fbandD;
ffS(bandD) = fbandA;