Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I can increase and decrease the frequency of a signal using the combination of fft and a Fourier series expansion FOR loop in the code below but if the signal/array is to large it becomes extremely slow (an array that's 1x44100 takes about 2 mins to complete) I'm sure it has to do with the for loop but I'm not exactly sure how to vectorize it to improve performance. Please note that this will be used with audio signals that are 3 to 6 mins long. The 1x44100 array is only a second and it takes about 2 mins to complete

Any recommendations

%create signal
clear all, clc,clf,tic
x= linspace(0,2*pi,44100)';

%Used in exporting to ycalc audio file make sure in sync with above
phaseshift=(0)*pi/180 ; %can use mod to limit to 180 degrees


N=size_r; %to test make 50
phase = 0;
f0 = 1/T; % Exactly, one period

y=(y/max(abs(y))*.8)/2; %make the max amplitude here
C = fft(y)/N; % No semicolon to display output

A = real(C);
B = imag(C)*-1; %I needed to multiply by -1 to get the correct sign

% Single-Sided (f >= 0)
An = [A(1); 2*A(2:round(N/2)); A(round(N/2)+1)]; 
Bn = [B(1); 2*B(2:round(N/2)); B(round(N/2)+1)];

ycalc=zeros(N,1); %preallocating space for ycalc

for p=2:pmax
       %%1 step) re-create signal using equation
        w=w+(360/(pmax-1)); %used to create phaseshift


fprintf('\n- Completed in %4.4fsec or %4.4fmins\n',toc,toc/60);

subplot(2,1,1), plot(y),title('Orginal Signal');
subplot(2,1,2),plot(ycalc),title('FFT new signal');

Here's a pic of the plot if some one wants to see the output, which is correct the FOR loop is just really really slow

FFT output

share|cite|improve this question

I believe you can replace the main portion of the for loop with the following line

ycalc=An.*cos(freq_new*(1:(p-1))*t-phaseshift) +Bn.*sin(freq_new*(1:(p-1))*t-phaseshift)+(vertoff/pmax);

share|cite|improve this answer
unfortunately that didn't work but I'll try a couple of different variations and post back if any of them work – Rick T Apr 3 '11 at 2:09

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.