[R] Recuperate Spectrum() amplitude

[Mike Marchywka]

----------------------------------------
> Date: Tue, 8 Feb 2011 17:34:15 +0000
> From: marie.guillot at bordeaux.inra.fr
> To: r-help at r-project.org
> Subject: [R] Recuperate Spectrum() amplitude
>
> Dear list,


> I put the code I use here to understand if the difference comes from the
> way that fft function are encoded int both software or from my own
> misunderstanding of the function.
> I am importing data as zoo, and verifying that they are regular spaced,
> and after, I simply use the non smoothed spectrum() function :
> SP<-spectrum(as.ts(Piezo))
> To recuperate the amplitude I calculate the square root of the density,
> divided by the number of observations at the the frequency I am
> interested in (1/86400=frequency of one cycle per day) :
> sqrt(SP$spec[which(SP$freq==1/86400)])/length(Piezo)
>
> Does anyone of you use both Matlab and R and have already find
> differences in results ? Or does my trouble comes from my own methodology ?
> Thank you very much in advance for your feedback,
>

Well, presumably your MATLAB results have some unit unless you publish
" MATLAB says 7 " and the trick here would be to run some simple controlled tests
with generated "data" to see if the R routines can produce similar results.
Probably the easiest place to start is to use Parseval and see if you can 
conserve energry between domains.  For example, is this what you would
expect? The 22.97==22.97 result and what units would you expect?  I've
never done a case with R where I cared about absolute amplitudes
but even if someone had an answer you may want to check if for yourself.

> asdf<-0:1000
> qwert=sin(asdf/200)
> fdasa<-abs(fft(qwert))
> max(fdasa)
> sqrt(sum(fdasa*fdasa))/sqrt(1001)
[1] 22.97086
> sqrt(sum(qwert*qwert))
[1] 22.97086

I just ran into a problem deconvolving a spectrum with an impulse response using
R. I divided fft's of two histograms and for the life of me couldn't figure out
why the cleaned histogram was time shifted until I finally realised that the result
is what you expect if the impulse response is time shifted ( duh). Sometimes you
just need to do simple checks...