# FFT and changing frequency and vectorizing for loop

5 messages
Open this post in threaded view
|

## FFT and changing frequency and vectorizing for loop

 FFT and changing frequency and vectorizing for loopGreetings AllI 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 recommendationstia 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 ropmax=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'); -- - _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 --- On Sat, 4/2/11, Rick T <[hidden email]> wrote: From: Rick T <[hidden email]> Subject: FFT and changing frequency and vectorizing for loop To: [hidden email] 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 [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octaveI 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. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 It works it's just slow due to the for loop here's an example of the output which is correcthttp://dl.dropbox.com/u/6576402/questions/q1.png On Sat, Apr 2, 2011 at 2:41 PM, Sergei Steshenko wrote: --- On Sat, 4/2/11, Rick T <[hidden email]> wrote: From: Rick T <[hidden email]> Subject: FFT and changing frequency and vectorizing for loop To: [hidden email] 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 [hidden email] 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. -- _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave