help-octave
[Top][All Lists]

## Re: FFT and changing frequency and vectorizing for loop

 From: Sergei Steshenko Subject: Re: FFT and changing frequency and vectorizing for loop Date: Sat, 2 Apr 2011 17:41:42 -0700 (PDT)

--- On Sat, 4/2/11, Rick T <address@hidden> wrote:

Subject: FFT and changing frequency and vectorizing for loop
Date: Saturday, April 2, 2011, 3:02 PM

FFT and changing frequency and vectorizing for loop
Greetings All
I can increase and decrease the frequency of thesignal using the combination of
fft and a series expansion for loop in the code below

but if the signal/array is to large it becomes extremelyslow (an array that's
1x44100 takes about 2 mins to complete) I'm sure it has to do with the for loop
butI'm not exactly sure how to vectorize them to improve performance

Any recommendations
tia sal22
%create signalclear all, clc,clf,ticx= linspace(0,2*pi,44100)';
%Used in exporting to ycalc audio file make sure in sync with above

freq_orig=1;freq_new=4vertoff=0;vertoffConj=0;vertoffInv=0;vertoffInvConj=0;phaseshift=(0)*pi/180
; %can use mod to limit to 180 degrees

y=sin(freq_orig*(x)); [size_r,size_c]=size(y);
N=size_r; %to test make 50T=2*pi;dt=T/N;t=linspace(0,T-dt,N)';phase = 0;

f0 = 1/T; % Exactly, one period
y=(y/max(abs(y))*.8)/2; %make the max amplitude hereC = 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)];
%needed to put the ' to get vaules in rows

Bn = [B(1); 2*B(2:round(N/2)); B(round(N/2)+1)]; %needed to put the ' to get
vaules in ro
pmax=N/2;ycalc=zeros(N,1); %preallocating space for ycalcw=0;for p=2:pmax

%        %%1 step) re-create signal using equation
ycalc=ycalc+An(p)*cos(freq_new*(p-1).*t-phaseshift)+Bn(p)*sin(freq_new*(p-1).*t-phaseshift)+(vertoff/pmax);

w=w+(360/(pmax-1)); %used to create phaseshift         phaseshift=w;
end;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');

--
-

-----Inline Attachment Follows-----

_______________________________________________
Help-octave mailing list
https://mailman.cae.wisc.edu/listinfo/help-octave

I don't think you are doing a right thing - not WRT loop/vectorization,
but WRT final result.

If I understand correctly your code, you multiply by 4 each individual
frequency. If so, your signal will appear in your output 4 times in a row.

Try a simple example with, say, just 16 samples.

There is a package for 'octave' doing polyphase resampling:

http://octave.sourceforge.net/signal/function/resample.html
.

Regards,
Sergei.