Perhaps this question is inappropriate or posted to the wrong list? Any guidance would be appreciated. Thanks.
--- On Fri, 9/23/11, Joseph Park <jpark...@att.net> wrote: > From: Joseph Park <jpark...@att.net> > Subject: Cross Spectrum : Conversion of 2-D spectrum into a single complex > array > To: r-help@r-project.org > Date: Friday, September 23, 2011, 5:02 PM > Hi, I'm wondering why the spectrum() > phase of quadrature > couple isn't purely +/-pi. > > But mostly, I'm looking for a recommended way to take a > 2-D > spectrum and convert it into a single complex array. > > Kindly consider: > > # 10 Hz sine wave 10 seconds long sampled at 50 Hz > deltaT = 1/50 > t = seq(0, 10, deltaT) > w = 2 * pi * 10 > > x = ts( sin( w * t ), deltat = deltaT ) > y = ts( sin( w * t ), deltat = deltaT ) > > # Want the cross spectrum between x and y > Sxy = spectrum( cbind( x, y ), plot = F ) > > # The phase is perfectly zero > plot( Sxy $ freq, Sxy $ phase[ ,1], type = 'l' ) > > # Make y a quadrature component of x > y = ts( cos( w * t ), deltat = deltaT ) > Sxy = spectrum( cbind( x, y ), plot = F ) > > # The phase should be either +pi or -pi > # since exp(i pi) = exp(-i pi) = -1 + 0i > # Why isn't it? Sampling issues? Nyquist has > been satisfied. > plot( Sxy $ freq, Sxy $ phase[ ,1], type = 'l' ) > > # The real question (limited to a 2-D input): > # How to combine the amplitude/phase into a single > # complex valued cross spectrum? > mySxy = complex( real = Sxy $ > spec[,1] * Sxy$ spec[,2], > > imaginary = Sxy $ phase[ ,1] ) > > # This does give something affine to a plot of Sxy > # Compare > plot( Sxy, log="dB" ) > # to: > plot( Sxy$freq, 10*log10( Re(mySxy) ), type='l') > > ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
