Reply by louis June 29, 20052005-06-29
Okay tahnk you eveyrone for your responses;

I've been trying to figure out if the  signal I'm trying to synthesize
would be symmetric in any way, even, or odd. I know I'm trying to create a
'real' signal.. (I think.. the whole real vs imag. concept is still
confusing for me in practice, thought I've done plenty of reading on it)

So I have the code that performs an FFT on an audio file I have which
plays the type of sound I am trying to simulate. I output all 1024 real
and imag values after doing the FFT at a couple of time points and look
through these values to see if I can find a pattern (i.e. real , even odd
etc)

So the FFT produces a complex output, both real and imaginaries are
nonzero, but theere is no symmetrical pattern that I can see? 

The reason I want to know this is I thought that maybe I can input my
spectral values into my inverse FFT in such a way to trick it into
thinking it was say real and even, and then the inverse FFT would produce
the correct output (I thought this might be a way to get around teh fact
that i have no phase info?)

Do you think this could make sense?
Thanks!
LD
		
This message was sent using the Comp.DSP web interface on
www.DSPRelated.com
Reply by Ron June 24, 20052005-06-24
>... i thought nobody was coding their own FFT no more. Do like >anybody, use FFTW, you'll never code anything as good as that
good for what purpose? If you code your own FFT, you understand much better how it works (or doesn't work), and you own the copyright. If you design your own hardware (a custorm VLIW inside an FPGA for instance) you might be able to optimize for your specific constraints better than by just letting a compiler do its job on canned code. IMHO. YMMV. -- rhn Ron Nicholson
Reply by Michel Rouzic June 24, 20052005-06-24
damn, i thought nobody was coding their own FFT no more. Do like
anybody, use FFTW, you'll never code anything as good as that

louis wrote:
> Hello, > > I am using C++ on a linux little endian machine, and am interfacing to > OpenSound system to produce my audio sound. So my goal is to be able to > reprduce the correct raw sound data corresponding to spectral data I > already have (for example, I know the frequency versus power information > for at least 128 bins, and over about 30 time points per second) > > I am using the four1 function found on the Numerical Recipes in C pages > online (their FFT code is posted online). > > So before I do this, I am first trying to output a single frequency tone > to make sure the IFFT is working - but low and behold it doesnt work. I am > sure the numerical recipes code is reliable as many people have experience > with it, so maybe I am setting it up wrong? Maybe I don't have the > imaginary and real components set right? or maybe its even a type casting > thing? Who knows! Below is the main code: > (The last few lines outputs the raw sound data, and I later play this back > using the correct channel and sampling rate info, however I don't hear a > single frequency tone like I should) > > static unsigned char buf[30*1024]; > static unsigned char large_buffer_audio[30*1024]; > static unsigned char buf_audio[30*1024]; > > void four1(float [], unsigned long, int ); > //the third parameter in four1 corresponds to forward or > // backward fft, if set to -1 then do inverse fft > > int main() > { > const int fft_size = 1024; > float f1[2048]; > int i, j, k, t; > FILE *fd_audio; > > for (t = 0; t<30; t++ ) > { > for (k=0; k<fft_size; k++) > { > f1[ 100 ] = 1000; > f1[fft_size-100] = 1000.0; > f1[k] = 0; > } > > four1(f1-1, fft_size, -1); > > for (i=0; i< 2*fft_size; i++) > { > buf_audio[ i ] = (unsigned char)(f1[(2*i)]); > } > > > for ( i = (t*1024), j=0; j< 2*fft_size;i++,j++ ) > { > large_buffer_audio[ i ] = buf_audio[j]; > large_buffer_audio[ i ] = buf_audio[j]; > } > } > > > FILE *audiofile = fopen( "myfile_num_rec_fft", "wb" ); > fwrite(large_buffer_audio,sizeof(unsigned > char),sizeof(large_buffer_audio), audiofile ); > > fclose(fd_audio); > return 0; > } > > > Thank you for your assistance!!! Its very appreciated! > LD > > > > > This message was sent using the Comp.DSP web interface on > www.DSPRelated.com
Reply by rhnl...@yahoo.com June 23, 20052005-06-23
louis allegedly wrote:
> so just to clarify - when i assign a positive amplitude to a bin that > corresponds to say frequency 1000hz, then do the ifft, i hear a clear > sounding 1khz tone. Thne when i assign two frequency bins, say bins 3 and > 4 (corresponding to frequencies 129 and 172 hz, since each bin is 43 hz) > then i get some waveform. I used matlab to plot the summation of two sine > waves 129 and 172 hz, and the waveform looks different.
Depending on whether you assign to the real or imaginary bins (in a complex IFFT, there are 2 bins for each frequency), you could get the sum of 2 sine waves, 2 cosine waves, or a sin+cos waveform (etc.), which will all look different. IMHO. YMMV. -- rhn Ron Nicholson
Reply by Andor June 22, 20052005-06-22
louis wrote:
> okay thank you again for your post, but can you explain why theoretically > it can't work??
Because you need magnitude and phase, or equivalently real and imaginary part, to fully specify the DFT X of a time domain signal x. For example sin(t), -sin(t) and cos(t) all have the same magnitude but different phase. Regards, Andor
Reply by louis June 21, 20052005-06-21
okay thank you again for your post, but can you explain why theoretically
it can't work??

(i have taken a dsp course, but i feel my signal processing background is
still weak)

so just to clarify - when i assign a positive amplitude to a bin that
corresponds to say frequency 1000hz, then do the ifft, i hear a clear
sounding 1khz tone. Thne when i assign two frequency bins, say bins 3 and
4 (corresponding to frequencies 129 and 172 hz, since each bin is 43 hz)
then i get some waveform. I used matlab to plot the summation of two sine
waves 129 and 172 hz, and the waveform looks different. My plot doesnt
look noisy or anything, but its not right when i compare it to matlab..

So i dont know what to make of all this

Is it me? is it Numerical recipes?? My boss is insistent that it shoudl be
possible :s
(i am ready to have a nervous breakdown!!)

okay thank you kindly again for your assistance!!


		
This message was sent using the Comp.DSP web interface on
www.DSPRelated.com
Reply by Andor June 21, 20052005-06-21
louis wrote:

> theoretically shouldnt this work?
>From your original post it seems that you only have the power
(magnitude) information in spectral domain. In this case, no, theoretically it shouldn't work. Regards, Andor
Reply by louis June 20, 20052005-06-20
OKAY 
so it turns out it actually works, it was a mistake in the way i was
outputting the sound data

anyways, now ive applied this to spectral data that i have, but it doesnt
produce the right sound

theoretically shouldnt this work?
		
This message was sent using the Comp.DSP web interface on
www.DSPRelated.com
Reply by louis June 20, 20052005-06-20
Hi thanks very much for your response,

Okay below I have corrected my previously posted code - 
Here I outputting a 1khz tone, using the Numerical Recipes four1 code.
Note that I am using a 1024 point FFT, where I specify that my frequencies
range from +22050 to -22050, so each frequency bin is 43 Hz.  To output a 1
Khz tone, then I need to assign a positive amplitude to bin 1000/43 = 23

I noticed that the code below only works for when I have even numbered
bins -  so for example if I change the '23' to a '24' , then this will
output a beautiful sounding 1khz tone.  However if i leave it at bin 23,
then there is some noise.  I tried this for numerous bins, and the pattern
holds - even bins work, odd ones don't. 

I feel like i'm really close - and I am probably not assigning something 
correctly? But I've looked at this a million times and don't know for the
life of me what the problem could be! four1 really should work..!!

Is a c++ version of an ifft that might work better?? 
Thank you in advance!!


int main()
{
	const int fft_size = 1024;
	float f1[2*fft_size];
	float f2[2*fft_size];
	float *ff1, *ff2;
	
	int i, j, k, t;
	

	ff1 = f1-1;
	ff2 = f2-1;

	
	for (  t=0; t<79; t++)
	{
		for(  k=0; k<2*fft_size; k++)
		{
			
			f1[ k ] = 0.0;
			f2[ k ] = 0.0;
			
			f1[ (23*2) ]= 1000.0; //1000Hz/43Hz = bin 23
			
		}
	
		four1( ff1, (unsigned long)1024, -1);
		four1( ff2, (unsigned long)1024, -1);			
					
	
		for (i=0; i< fft_size; i++)
		{
			buf_audio[ 2*i ] 	=   (short int)(f1[(2*i)]);	// left channel
			buf_audio[ (2*i) + 1 ]  =   (short int)(f2[(2*i)]);	// right channel
		}
		
		
		for ( i = (t*1024), j=0; j < fft_size ; i++, j++ ) 
		{
			large_buffer_audio[ i ] = buf_audio[j];
			large_buffer_audio[ i ] = buf_audio[j];
		}
	}

	
	FILE *audiofile = fopen( "debug1", "wb" );
	fwrite( large_buffer_audio, sizeof(short int),
sizeof(large_buffer_audio), audiofile );

	return 0;
}
		
This message was sent using the Comp.DSP web interface on
www.DSPRelated.com
Reply by Mark Borgerding June 14, 20052005-06-14
louis wrote:
> I am > sure the numerical recipes code is reliable
I wouldn't be so quick to assume that. I'd recommend a little googling. I, for one, will never trust NR again. There are other choices. My kissfft library might be a good match for your needs. http://www.sourceforge.net/projects/kissfft It will probably be fast enough, but if you need faster and can tolerate GPL in your code -- FFTW is a great library (fftw.org). There are also C-calllable routines in GSL for FFT (also GPL). If I were doing your project, I might use Python. The Numeric module has FFT support built in. Maybe I don't have the
> imaginary and real components set right? or maybe its even a type casting > thing? Who knows! Below is the main code: > (The last few lines outputs the raw sound data, and I later play this back > using the correct channel and sampling rate info, however I don't hear a > single frequency tone like I should) > > static unsigned char buf[30*1024]; > static unsigned char large_buffer_audio[30*1024]; > static unsigned char buf_audio[30*1024]; > > void four1(float [], unsigned long, int ); > //the third parameter in four1 corresponds to forward or > // backward fft, if set to -1 then do inverse fft > > int main() > { > const int fft_size = 1024; > float f1[2048]; > int i, j, k, t; > FILE *fd_audio; > > for (t = 0; t<30; t++ ) > { > for (k=0; k<fft_size; k++) > { > f1[ 100 ] = 1000; > f1[fft_size-100] = 1000.0; > f1[k] = 0;
> } ^^^^^^^^^^ I don't think this is what you meant to do. Actually it looks like you never even initialize the latter half of your array.
> > four1(f1-1, fft_size, -1);
^^^^ why f1-1 ?? Oh I think I know -- is this b/c the NR people couldn't be bothered to actually convert their Fortran code from one-based indexing to zero-based? ( Are you starting to get the idea that I don't think much of Numerical Recipes? )