Dear Marc,
I’m using Geopsy 3.5.2 on Ubuntu for H/V analysis, and I’m having some trouble understanding how windowed H/V or PSD spectra are combined into an average over windows. I have been trying to replicate the H/V algorithm in Matlab using the spectrogram tool, and my problem is that when I manually average the individual Fourier/PSD spectra of a given number of windows, it is very close to, but not quite the same as when I calculate the average spectrum in Geopsy. When I export the individual windowed spectra from Geopsy and manually average them in Matlab, they deviate from the average spectra as Geopsy calculates it. Taking two single windows, I see that the difference between the Matlab- and Geopsy-average is only a couple percent. With all windows over an entire multi-hour recording, however, my Matlab-calculated average is 2-3 times larger than the one computed in Geopsy and the shape of the curve is a little different, too (see png attached). Hence my question: Could you maybe tell me how, in the background, Geopsy averages over individual windows? Is there another step (weighting, squaring?) involved in the process of calculating the average? I would greatly appreciate any help or information.
For the record, I found that Geopsy’s PSD spectrum equals sqrt(fs)*|X/L|/2 (where X is the double-sided complex spectrum for a single window calculated with Matlab’s fft() function, L is the number of elements in the window vector and fs the sampling frequency). I convert each window's spectrum according to this factor. I calculate the average in Matlab simply by taking the unweighted mean of all windows' spectra.
I have not changed countPerVolt or voltPerUnit from the default of 1, so the PSD is in counts^2/Hz. From there, I compared the PSD of a single window, calculated in Geopsy, with that in Matlab, taking this unit conversion factor into account, and it turns out that they are identical. I used linearly spaced frequencies in the fft, no window rejection, a 5% Tukey window with a length of 3000, 0% smoothing, nfft=7500 (halved in Geopsy, as it outputs only positive frequencies), and fixed window length of 30s. fs = 250 Hz. My recordings are several hours.
How are windowed segments averaged?
How are windowed segments averaged?
- Attachments
-
- P6a2006_avg_Matlab_vs_Geopsy.png
- (29.41 KiB) Not downloaded yet
Re: How are windowed segments averaged?
Hi,
The average is computed on a log scale. This is maybe the difference with your results.
Do you know that you can do the H/V computation with geopsy-hv? It is a core application (without a graphical interface) that might be useful for large datasets. It does the same job as in the graphical interface.
Best regards,
Marc
The average is computed on a log scale. This is maybe the difference with your results.
Do you know that you can do the H/V computation with geopsy-hv? It is a core application (without a graphical interface) that might be useful for large datasets. It does the same job as in the graphical interface.
Best regards,
Marc