Drexel dragonThe Math ForumDonate to the Math Forum



Search All of the Math Forum:

Views expressed in these public forums are not endorsed by Drexel University or The Math Forum.


Math Forum » Discussions » Software » comp.soft-sys.matlab

Topic: random time-series displacement data from PSD and back again
Replies: 5   Last Post: Jun 24, 2013 9:01 PM

Advanced Search

Back to Topic List Back to Topic List Jump to Tree View Jump to Tree View   Messages: [ Previous | Next ]
Derek Goring

Posts: 3,919
Registered: 12/7/04
Re: random time-series displacement data from PSD and back again
Posted: Jun 21, 2013 5:48 PM
  Click to see the message monospaced in plain text Plain Text   Click to reply to this topic Reply

On Saturday, June 22, 2013 8:10:26 AM UTC+12, Jim Schwendeman wrote:
> Howdy Forum,
>
>
>
> I have a project where I have a vibration PSD profile (G^2/Hz). I'd like to use this to create random time-series displacement data to feed into a dynamic model to evaluate performance in the presence of vibe.
>
>
>
> The problem occurs when I convert from G^2/Hz to displacement/time, and then back to the identical (or very near it) G^2/Hz profile. I create random time-series data, but when I convert back to the frequency domain, it doesn't match up. The re-created PSD has content up to twice the frequency domain I'm initially putting in, and the RMS values don't match up. This sanity check needs to work before I can confidently put the time-series data into the model, so I'm a bit hung up...
>
>
>
> Code below:
>
> % Create noise (centered 1), pump it through PSD, and then calc gRMS
>
> noise_range = 0.02
>
> noise = noise_range*randn(1,N)+(1-noise_range/2);
>
> PSD_noise = PSD.*noise;
>
>
>
> %Calculate gRMS for original PSD and noise-injected PSD
>
> gRMS = sqrt(max(cumtrapz(PSD0(:,1),PSD0(:,2))));
>
> gRMS_noise = sqrt(max(cumtrapz(f,PSD_noise)));
>
> %These match up, the gRMS of my noisy PSD is nearly identical to my input PSD
>
>
>
> % Create time and displacement
>
> t = linspace(0, N*dt, N);
>
> x = (1/sqrt(N))*real(ifft(sqrt(PSD_noise)));
>
> x = (x*9.81)./sqrt(std(x));
>
> disp('Parceval"s Theorem says this result should = gRMS')
>
> sum(abs(x).^2)
>
> %This answer is similar to the gRMS of my noisy/input PSD, but now off by ~10%
>
>
>
> %Now I convert time-series displacement back to frequency domain and should get the same PSD back out:
>
> N = length(x);
>
> xdft = fft(x);
>
> xdft = xdft(1:N/2+1);
>
> psdx = (1/(Fs*N)).*abs(xdft).^2;
>
> psdx(2:end-1) = 2*psdx(2:end-1);
>
> psdx = psdx*sqrt(N)*9.81;
>
> freq = 0:Fs/length(x):Fs/2;
>
>
>
> When I end up plotting freq vs psdx, I get frequency content over double of what was put in, and my gRMS is orders of magnitude below what was expected.
>
>
>
> Any thoughts?
>
> Jim


If PSD is the traditional power spectral density, then you have not prepared the Fourier transform correctly before ifft.
To see how to do this correctly, search "phase randomisation" (note spelling) on this newsgroup.



Point your RSS reader here for a feed of the latest messages in this topic.

[Privacy Policy] [Terms of Use]

© Drexel University 1994-2014. All Rights Reserved.
The Math Forum is a research and educational enterprise of the Drexel University School of Education.