2019-04-17 18:05:07 +02:00
<?xml version="1.0" encoding="utf-8"?>
< !DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Strict//EN"
"http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd">
< html xmlns = "http://www.w3.org/1999/xhtml" lang = "en" xml:lang = "en" >
< head >
2020-11-12 10:30:56 +01:00
<!-- 2020 - 11 - 12 jeu. 10:28 -->
2019-04-17 18:05:07 +02:00
< meta http-equiv = "Content-Type" content = "text/html;charset=utf-8" / >
2019-05-02 14:12:11 +02:00
< title > Huddle Test of the L22 Geophones< / title >
2019-04-17 18:05:07 +02:00
< meta name = "generator" content = "Org mode" / >
2019-05-10 16:06:43 +02:00
< meta name = "author" content = "Dehaeze Thomas" / >
2020-11-12 10:30:56 +01:00
< link rel = "stylesheet" type = "text/css" href = "https://research.tdehaeze.xyz/css/style.css" / >
< script type = "text/javascript" src = "https://research.tdehaeze.xyz/js/script.js" > < / script >
< script > M a t h J a x = {
tex: {
tags: 'ams',
macros: {bm: ["\\boldsymbol{#1}",1],}
}
};
< / script >
< script type = "text/javascript" src = "https://cdn.jsdelivr.net/npm/mathjax@3/es5/tex-mml-chtml.js" > < / script >
2019-04-17 18:05:07 +02:00
< / head >
< body >
2019-05-10 16:06:43 +02:00
< div id = "org-div-home-and-up" >
< a accesskey = "h" href = "../index.html" > UP < / a >
|
< a accesskey = "H" href = "../index.html" > HOME < / a >
< / div > < div id = "content" >
2019-05-02 14:12:11 +02:00
< h1 class = "title" > Huddle Test of the L22 Geophones< / h1 >
2019-04-17 18:05:07 +02:00
< div id = "table-of-contents" >
< h2 > Table of Contents< / h2 >
< div id = "text-table-of-contents" >
< ul >
2020-11-12 10:30:56 +01:00
< li > < a href = "#orgb47b0b2" > 1. Experimental Setup< / a > < / li >
< li > < a href = "#org0b1e99e" > 2. Signal Processing< / a >
2019-04-18 17:02:47 +02:00
< ul >
2020-11-12 10:30:56 +01:00
< li > < a href = "#org35507da" > 2.1. Load data< / a > < / li >
< li > < a href = "#org69bec1f" > 2.2. Time Domain Data< / a > < / li >
< li > < a href = "#orgf0a72da" > 2.3. Computation of the ASD of the measured voltage< / a > < / li >
< li > < a href = "#org7acabeb" > 2.4. Scaling to take into account the sensibility of the geophone and the voltage amplifier< / a > < / li >
< li > < a href = "#org4b3b8a7" > 2.5. Computation of the ASD of the velocity< / a > < / li >
< li > < a href = "#orgaaca76b" > 2.6. Transfer function between the two geophones< / a > < / li >
< li > < a href = "#org5c32c39" > 2.7. Estimation of the sensor noise< / a > < / li >
2019-04-18 17:02:47 +02:00
< / ul >
< / li >
2020-11-12 10:30:56 +01:00
< li > < a href = "#org52fed77" > 3. Compare axis< / a >
2019-04-18 17:02:47 +02:00
< ul >
2020-11-12 10:30:56 +01:00
< li > < a href = "#org7f52c4b" > 3.1. Load data< / a > < / li >
< li > < a href = "#orgb105293" > 3.2. Compare PSD< / a > < / li >
< li > < a href = "#orgbce3b98" > 3.3. Compare TF< / a > < / li >
2019-04-18 17:02:47 +02:00
< / ul >
< / li >
2020-11-12 10:30:56 +01:00
< li > < a href = "#orgd32606c" > 4. Appendix< / a >
2019-04-17 18:20:17 +02:00
< ul >
2020-11-12 10:30:56 +01:00
< li > < a href = "#orge652a09" > 4.1. Computation of coherence from PSD and CSD< / a > < / li >
2019-04-17 18:20:17 +02:00
< / ul >
< / li >
2019-04-17 18:05:07 +02:00
< / ul >
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-orgb47b0b2" class = "outline-2" >
< h2 id = "orgb47b0b2" > < span class = "section-number-2" > 1< / span > Experimental Setup< / h2 >
2019-04-17 18:05:07 +02:00
< div class = "outline-text-2" id = "text-1" >
2019-04-17 18:20:17 +02:00
< p >
Two L22 geophones are used.
They are placed on the ID31 granite.
They are leveled.
< / p >
< p >
The signals are amplified using voltage amplifier with a gain of 60dB.
2019-05-10 16:06:43 +02:00
The voltage amplifiers includes:
2019-04-17 18:20:17 +02:00
< / p >
2019-05-10 16:06:43 +02:00
< ul class = "org-ul" >
< li > an high pass filter with a cut-off frequency at 1.5Hz (AC option)< / li >
< li > a low pass filter with a cut-off frequency at 1kHz< / li >
< / ul >
2019-04-17 18:20:17 +02:00
2020-11-12 10:30:56 +01:00
< div id = "orgfc98f42" class = "figure" >
2019-05-10 16:06:43 +02:00
< p > < img src = "./img/setup.jpg" alt = "setup.jpg" width = "500px" / >
2019-04-17 18:20:17 +02:00
< / p >
< p > < span class = "figure-number" > Figure 1: < / span > Setup< / p >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "org01be4e4" class = "figure" >
2019-05-10 16:06:43 +02:00
< p > < img src = "./img/geophones.jpg" alt = "geophones.jpg" width = "500px" / >
2019-04-17 18:20:17 +02:00
< / p >
< p > < span class = "figure-number" > Figure 2: < / span > Geophones< / p >
< / div >
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-org0b1e99e" class = "outline-2" >
< h2 id = "org0b1e99e" > < span class = "section-number-2" > 2< / span > Signal Processing< / h2 >
2019-04-17 18:20:17 +02:00
< div class = "outline-text-2" id = "text-2" >
2019-04-18 17:11:25 +02:00
< p >
2020-11-12 10:30:56 +01:00
< a id = "org28ac54a" > < / a >
2019-05-10 16:06:43 +02:00
< / p >
2020-11-12 10:30:56 +01:00
< div class = "note" id = "orgdb9b9cc" >
2019-05-10 16:06:43 +02:00
< p >
All the files (data and Matlab scripts) are accessible < a href = "data/huddle_test_signal_processing.zip" > here< / a > .
2019-04-18 17:11:25 +02:00
< / p >
2019-05-10 16:06:43 +02:00
< / div >
2019-04-17 18:20:17 +02:00
< / div >
2019-04-18 17:11:25 +02:00
2020-11-12 10:30:56 +01:00
< div id = "outline-container-org35507da" class = "outline-3" >
< h3 id = "org35507da" > < span class = "section-number-3" > 2.1< / span > Load data< / h3 >
2019-04-17 18:20:17 +02:00
< div class = "outline-text-3" id = "text-2-1" >
2019-04-18 17:02:47 +02:00
< p >
We load the data of the z axis of two geophones.
< / p >
2019-04-18 17:11:25 +02:00
2019-04-17 18:05:07 +02:00
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > load(< span class = "org-string" > 'mat/data_001.mat'< / span > , < span class = "org-string" > 't'< / span > , < span class = "org-string" > 'x1'< / span > , < span class = "org-string" > 'x2'< / span > );
dt = t(2) < span class = "org-type" > -< / span > t(1);
2019-04-17 18:05:07 +02:00
< / pre >
< / div >
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-org69bec1f" class = "outline-3" >
< h3 id = "org69bec1f" > < span class = "section-number-3" > 2.2< / span > Time Domain Data< / h3 >
2019-04-17 18:20:17 +02:00
< div class = "outline-text-3" id = "text-2-2" >
2019-04-17 18:05:07 +02:00
< div class = "org-src-container" >
2019-04-18 09:20:31 +02:00
< pre class = "src src-matlab" > < span class = "org-type" > figure< / span > ;
2019-04-17 18:05:07 +02:00
hold on;
2020-11-12 10:30:56 +01:00
plot(t, x1);
plot(t, x2);
2019-04-17 18:05:07 +02:00
hold off;
2020-11-12 10:30:56 +01:00
xlabel(< span class = "org-string" > 'Time [s]'< / span > );
ylabel(< span class = "org-string" > 'Voltage [V]'< / span > );
xlim([t(1), t(end)]);
2019-04-17 18:05:07 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "orgc871843" class = "figure" >
2019-04-17 18:05:07 +02:00
< p > < img src = "figs/data_time_domain.png" alt = "data_time_domain.png" / >
< / p >
2019-04-17 18:20:17 +02:00
< p > < span class = "figure-number" > Figure 3: < / span > Time domain Data< / p >
2019-04-17 18:05:07 +02:00
< / div >
< div class = "org-src-container" >
2019-04-18 09:20:31 +02:00
< pre class = "src src-matlab" > < span class = "org-type" > figure< / span > ;
2019-04-17 18:05:07 +02:00
hold on;
2020-11-12 10:30:56 +01:00
plot(t, x1);
plot(t, x2);
2019-04-17 18:05:07 +02:00
hold off;
2020-11-12 10:30:56 +01:00
xlabel(< span class = "org-string" > 'Time [s]'< / span > );
ylabel(< span class = "org-string" > 'Voltage [V]'< / span > );
xlim([0 1]);
2019-04-17 18:05:07 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "org30ca9f2" class = "figure" >
2019-04-17 18:05:07 +02:00
< p > < img src = "figs/data_time_domain_zoom.png" alt = "data_time_domain_zoom.png" / >
< / p >
2019-04-17 18:20:17 +02:00
< p > < span class = "figure-number" > Figure 4: < / span > Time domain Data - Zoom< / p >
2019-04-17 18:05:07 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-orgf0a72da" class = "outline-3" >
< h3 id = "orgf0a72da" > < span class = "section-number-3" > 2.3< / span > Computation of the ASD of the measured voltage< / h3 >
2019-04-17 18:20:17 +02:00
< div class = "outline-text-3" id = "text-2-3" >
2019-04-18 17:02:47 +02:00
< p >
We first define the parameters for the frequency domain analysis.
< / p >
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > Fs = 1< span class = "org-type" > /< / span > dt; < span class = "org-comment" > % [Hz]< / span >
2019-05-03 11:32:54 +02:00
2020-11-12 10:30:56 +01:00
win = hanning(ceil(10< span class = "org-type" > *< / span > Fs));
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
2019-05-03 11:32:54 +02:00
< p >
Then we compute the Power Spectral Density using < code > pwelch< / code > function.
< / p >
2019-04-17 18:05:07 +02:00
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > [pxx1, f] = pwelch(x1, win, [], [], Fs);
[pxx2, < span class = "org-type" > ~< / span > ] = pwelch(x2, win, [], [], Fs);
2019-04-17 18:05:07 +02:00
< / pre >
< / div >
2019-05-03 11:32:54 +02:00
< p >
2020-11-12 10:30:56 +01:00
And we plot the result on figure < a href = "#orgf270ed7" > 5< / a > .
2019-05-03 11:32:54 +02:00
< / p >
< div class = "org-src-container" >
< pre class = "src src-matlab" > < span class = "org-type" > figure< / span > ;
hold on;
2020-11-12 10:30:56 +01:00
plot(f, sqrt(pxx1));
plot(f, sqrt(pxx2));
2019-05-03 11:32:54 +02:00
hold off;
2020-11-12 10:30:56 +01:00
< span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'xscale'< / span > , < span class = "org-string" > 'log'< / span > );
< span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'yscale'< / span > , < span class = "org-string" > 'log'< / span > );
xlabel(< span class = "org-string" > 'Frequency [Hz]'< / span > ); ylabel(< span class = "org-string" > 'ASD of the measured Voltage $\left[\frac{V}{\sqrt{Hz}}\right]$'< / span > )
xlim([0.1, 500]);
2019-05-03 11:32:54 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "orgf270ed7" class = "figure" >
2019-05-03 11:32:54 +02:00
< p > < img src = "figs/asd_voltage.png" alt = "asd_voltage.png" / >
< / p >
< p > < span class = "figure-number" > Figure 5: < / span > Amplitude Spectral Density of the measured voltage< / p >
< / div >
2019-04-17 18:05:07 +02:00
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-org7acabeb" class = "outline-3" >
< h3 id = "org7acabeb" > < span class = "section-number-3" > 2.4< / span > Scaling to take into account the sensibility of the geophone and the voltage amplifier< / h3 >
2019-04-17 18:20:17 +02:00
< div class = "outline-text-3" id = "text-2-4" >
2019-04-17 18:05:07 +02:00
< p >
2020-11-12 10:30:56 +01:00
The Geophone used are L22. Their sensibility is shown on figure < a href = "#orgc586e06" > 6< / a > .
2019-04-17 18:05:07 +02:00
< / p >
2019-04-18 17:02:47 +02:00
2019-04-17 18:05:07 +02:00
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > S0 = 88; < span class = "org-comment" > % Sensitivity [V/(m/s)]< / span >
f0 = 2; < span class = "org-comment" > % Cut-off frequnecy [Hz]< / span >
2019-05-03 11:32:54 +02:00
2020-11-12 10:30:56 +01:00
S = S0< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > 2< span class = "org-type" > /< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > /< / span > f0)< span class = "org-type" > /< / span > (1< span class = "org-type" > +< / span > s< span class = "org-type" > /< / span > 2< span class = "org-type" > /< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > /< / span > f0);
2019-04-17 18:05:07 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "orgc586e06" class = "figure" >
2019-04-17 18:05:07 +02:00
< p > < img src = "figs/geophone_sensibility.png" alt = "geophone_sensibility.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 6: < / span > Sensibility of the Geophone< / p >
2019-04-17 18:05:07 +02:00
< / div >
< p >
2019-04-18 17:02:47 +02:00
We also take into account the gain of the electronics which is here set to be \(60dB\).
2019-04-17 18:05:07 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > G0_db = 60; < span class = "org-comment" > % [dB]< / span >
2019-04-17 18:05:07 +02:00
2020-11-12 10:30:56 +01:00
G0 = 10< span class = "org-type" > ^< / span > (G0_db< span class = "org-type" > /< / span > 20); < span class = "org-comment" > % [abs]< / span >
2019-04-17 18:05:07 +02:00
< / pre >
< / div >
2019-04-18 17:02:47 +02:00
< p >
2019-05-03 11:32:54 +02:00
We divide the ASD measured (in \(\text{V}/\sqrt{\text{Hz}}\)) by the gain of the voltage amplifier to obtain the ASD of the voltage across the geophone.
2019-04-18 17:02:47 +02:00
We further divide the result by the sensibility of the Geophone to obtain the ASD of the velocity in \(m/s/\sqrt{Hz}\).
< / p >
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > scaling = 1< span class = "org-type" > ./< / span > squeeze(abs(freqresp(G0< span class = "org-type" > *< / span > S, f, < span class = "org-string" > 'Hz'< / span > )));
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-org4b3b8a7" class = "outline-3" >
< h3 id = "org4b3b8a7" > < span class = "section-number-3" > 2.5< / span > Computation of the ASD of the velocity< / h3 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-3" id = "text-2-5" >
< p >
2020-11-12 10:30:56 +01:00
The ASD of the measured velocity is shown on figure < a href = "#orgdacf8d0" > 7< / a > .
2019-04-18 17:02:47 +02:00
< / p >
2019-04-17 18:05:07 +02:00
< div class = "org-src-container" >
2019-04-18 09:20:31 +02:00
< pre class = "src src-matlab" > < span class = "org-type" > figure< / span > ;
2019-04-17 18:05:07 +02:00
hold on;
2020-11-12 10:30:56 +01:00
plot(f, sqrt(pxx1)< span class = "org-type" > .*< / span > scaling);
plot(f, sqrt(pxx2)< span class = "org-type" > .*< / span > scaling);
2019-04-17 18:05:07 +02:00
hold off;
2020-11-12 10:30:56 +01:00
< span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'xscale'< / span > , < span class = "org-string" > 'log'< / span > );
< span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'yscale'< / span > , < span class = "org-string" > 'log'< / span > );
xlabel(< span class = "org-string" > 'Frequency [Hz]'< / span > ); ylabel(< span class = "org-string" > 'ASD of the measured Velocity $\left[\frac{m/s}{\sqrt{Hz}}\right]$'< / span > )
xlim([0.1, 500]);
2019-04-17 18:05:07 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "orgdacf8d0" class = "figure" >
2019-04-17 18:05:07 +02:00
< p > < img src = "figs/psd_velocity.png" alt = "psd_velocity.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 7: < / span > Amplitude Spectral Density of the Velocity< / p >
2019-04-17 18:20:17 +02:00
< / div >
2019-04-18 17:25:02 +02:00
< p >
2020-11-12 10:30:56 +01:00
We also plot the ASD in displacement (figure < a href = "#org47aa587" > 8< / a > );
2019-04-18 17:25:02 +02:00
< / p >
< div class = "org-src-container" >
< pre class = "src src-matlab" > < span class = "org-type" > figure< / span > ;
hold on;
2020-11-12 10:30:56 +01:00
plot(f, (sqrt(pxx1)< span class = "org-type" > .*< / span > scaling)< span class = "org-type" > ./< / span > (2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > f));
plot(f, (sqrt(pxx2)< span class = "org-type" > .*< / span > scaling)< span class = "org-type" > ./< / span > (2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > f));
2019-04-18 17:25:02 +02:00
hold off;
2020-11-12 10:30:56 +01:00
< span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'xscale'< / span > , < span class = "org-string" > 'log'< / span > ); < span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'yscale'< / span > , < span class = "org-string" > 'log'< / span > );
xlabel(< span class = "org-string" > 'Frequency [Hz]'< / span > ); ylabel(< span class = "org-string" > 'ASD of the displacement $\left[\frac{m}{\sqrt{Hz}}\right]$'< / span > )
xlim([0.1, 500]);
2019-04-18 17:25:02 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "org47aa587" class = "figure" >
2019-04-18 17:25:02 +02:00
< p > < img src = "figs/asd_displacement.png" alt = "asd_displacement.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 8: < / span > Amplitude Spectral Density of the Displacement< / p >
2019-04-18 17:25:02 +02:00
< / div >
2019-04-17 18:05:07 +02:00
< / div >
< / div >
2019-04-18 17:02:47 +02:00
2020-11-12 10:30:56 +01:00
< div id = "outline-container-orgaaca76b" class = "outline-3" >
< h3 id = "orgaaca76b" > < span class = "section-number-3" > 2.6< / span > Transfer function between the two geophones< / h3 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-3" id = "text-2-6" >
< p >
We here compute the transfer function from one geophone to the other.
2020-11-12 10:30:56 +01:00
The result is shown on figure < a href = "#org313c704" > 9< / a > .
2019-04-18 17:02:47 +02:00
< / p >
< p >
2020-11-12 10:30:56 +01:00
We also compute the coherence between the two signals (figure < a href = "#orgbed81d4" > 10< / a > ).
2019-04-18 17:02:47 +02:00
< / p >
2019-04-18 09:20:31 +02:00
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > [T12, < span class = "org-type" > ~< / span > ] = tfestimate(x1, x2, win, [], [], Fs);
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "org313c704" class = "figure" >
2019-04-18 17:02:47 +02:00
< p > < img src = "figs/tf_geophones.png" alt = "tf_geophones.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 9: < / span > Estimated transfer function between the two geophones< / p >
2019-04-18 17:02:47 +02:00
< / div >
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > [coh12, < span class = "org-type" > ~< / span > ] = mscohere(x1, x2, win, [], [], Fs);
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "orgbed81d4" class = "figure" >
2019-04-18 17:02:47 +02:00
< p > < img src = "figs/coh_geophones.png" alt = "coh_geophones.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 10: < / span > Cohererence between the signals of the two geophones< / p >
2019-04-18 17:02:47 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-org5c32c39" class = "outline-3" >
< h3 id = "org5c32c39" > < span class = "section-number-3" > 2.7< / span > Estimation of the sensor noise< / h3 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-3" id = "text-2-7" >
< p >
2020-11-12 10:30:56 +01:00
The technique to estimate the sensor noise is taken from (< a href = "#citeproc_bib_item_1" > Barzilai, VanZandt, and Kenny 1998< / a > ).
2019-04-18 17:02:47 +02:00
< / p >
< p >
The coherence between signals \(X\) and \(Y\) is defined as follow
\[ \gamma^2_{XY}(\omega) = \frac{|G_{XY}(\omega)|^2}{|G_{X}(\omega)| |G_{Y}(\omega)|} \]
where \(|G_X(\omega)|\) is the output Power Spectral Density (PSD) of signal \(X\) and \(|G_{XY}(\omega)|\) is the Cross Spectral Density (CSD) of signal \(X\) and \(Y\).
< / p >
< p >
The PSD and CSD are defined as follow:
< / p >
\begin{align}
|G_X(\omega)| & = \frac{2}{n_d T} \sum^{n_d}_{n=1} \left| X_k(\omega, T) \right|^2 \\
|G_{XY}(\omega)| & = \frac{2}{n_d T} \sum^{n_d}_{n=1} [ X_k^*(\omega, T) ] [ Y_k(\omega, T) ]
\end{align}
< p >
where:
< / p >
< ul class = "org-ul" >
< li > \(n_d\) is the number for records averaged< / li >
< li > \(T\) is the length of each record< / li >
< li > \(X_k(\omega, T)\) is the finite Fourier transform of the kth record< / li >
< li > \(X_k^*(\omega, T)\) is its complex conjugate< / li >
< / ul >
< p >
2020-11-12 10:30:56 +01:00
The < code > mscohere< / code > function is compared with this formula on Appendix (section < a href = "#org06003b1" > 4.1< / a > ), it is shown that it is identical.
2019-04-18 17:02:47 +02:00
< / p >
< p >
2020-11-12 10:30:56 +01:00
Figure < a href = "#orgae83a33" > 11< / a > illustrate a block diagram model of the system used to determine the sensor noise of the geophone.
2019-04-18 17:02:47 +02:00
< / p >
< p >
Two geophones are mounted side by side to ensure that they are exposed by the same motion input \(U\).
< / p >
< p >
Each sensor has noise \(N\) and \(M\).
< / p >
2020-11-12 10:30:56 +01:00
< div id = "orgae83a33" class = "figure" >
2019-04-18 17:02:47 +02:00
< p > < img src = "figs/huddle-test.png" alt = "huddle-test.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 11: < / span > Huddle test block diagram< / p >
2019-04-18 17:02:47 +02:00
< / div >
< p >
2019-05-15 10:08:11 +02:00
We here assume that each sensor has the same magnitude of instrumental noise: \(N = M\).
2019-05-02 14:12:11 +02:00
We also assume that \(S_1 = S_2 = 1\).
2019-04-18 17:02:47 +02:00
< / p >
< p >
We then obtain:
< / p >
\begin{equation}
2020-11-12 10:30:56 +01:00
\label{orgbacdde9}
2019-04-18 17:02:47 +02:00
\gamma_{XY}^2(\omega) = \frac{1}{1 + 2 \left( \frac{|G_N(\omega)|}{|G_U(\omega)|} \right) + \left( \frac{|G_N(\omega)|}{|G_U(\omega)|} \right)^2}
\end{equation}
< p >
Since the input signal \(U\) and the instrumental noise \(N\) are incoherent:
< / p >
\begin{equation}
2020-11-12 10:30:56 +01:00
\label{org21035bd}
2019-04-18 17:02:47 +02:00
|G_X(\omega)| = |G_N(\omega)| + |G_U(\omega)|
\end{equation}
< p >
2020-11-12 10:30:56 +01:00
From equations \eqref{orgbacdde9} and \eqref{org21035bd}, we finally obtain
2019-04-18 17:02:47 +02:00
< / p >
2020-11-12 10:30:56 +01:00
< div class = "important" id = "org594bc6f" >
2019-04-18 17:02:47 +02:00
\begin{equation}
2020-11-12 10:30:56 +01:00
\label{orga8ac5d9}
2019-04-18 17:02:47 +02:00
|G_N(\omega)| = |G_X(\omega)| \left( 1 - \sqrt{\gamma_{XY}^2(\omega)} \right)
\end{equation}
< / div >
< p >
2020-11-12 10:30:56 +01:00
The instrumental noise is computed below. The result in V^2/Hz is shown on figure < a href = "#org7943278" > 12< / a > .
2019-04-18 17:02:47 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > pxxN = pxx1< span class = "org-type" > .*< / span > (1 < span class = "org-type" > -< / span > coh12);
2019-04-18 09:20:31 +02:00
< / pre >
< / div >
< div class = "org-src-container" >
< pre class = "src src-matlab" > < span class = "org-type" > figure< / span > ;
2019-04-18 17:02:47 +02:00
hold on;
2020-11-12 10:30:56 +01:00
plot(f, pxx1, < span class = "org-string" > '-'< / span > );
plot(f, pxx2, < span class = "org-string" > '-'< / span > );
plot(f, pxxN, < span class = "org-string" > 'k--'< / span > );
2019-04-18 17:02:47 +02:00
hold off;
2020-11-12 10:30:56 +01:00
< span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'xscale'< / span > , < span class = "org-string" > 'log'< / span > ); < span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'yscale'< / span > , < span class = "org-string" > 'log'< / span > );
xlabel(< span class = "org-string" > 'Frequency [Hz]'< / span > ); ylabel(< span class = "org-string" > 'PSD of the measured Voltage $\left[\frac{V^2}{Hz}\right]$'< / span > );
xlim([0.1, 500]);
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
2019-04-18 09:20:31 +02:00
2020-11-12 10:30:56 +01:00
< div id = "org7943278" class = "figure" >
2019-04-18 17:02:47 +02:00
< p > < img src = "figs/intrumental_noise_V.png" alt = "intrumental_noise_V.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 12: < / span > Instrumental Noise and Measurement in \(V^2/Hz\)< / p >
2019-04-18 17:02:47 +02:00
< / div >
< p >
2020-11-12 10:30:56 +01:00
This is then further converted into velocity and compared with the ground velocity measurement. (figure < a href = "#org8ecfd08" > 13< / a > )
2019-04-18 17:02:47 +02:00
< / p >
< div class = "org-src-container" >
< pre class = "src src-matlab" > < span class = "org-type" > figure< / span > ;
hold on;
2020-11-12 10:30:56 +01:00
plot(f, sqrt(pxx1)< span class = "org-type" > .*< / span > scaling, < span class = "org-string" > '-'< / span > );
plot(f, sqrt(pxx2)< span class = "org-type" > .*< / span > scaling, < span class = "org-string" > '-'< / span > );
plot(f, sqrt(pxxN)< span class = "org-type" > .*< / span > scaling, < span class = "org-string" > 'k--'< / span > );
2019-04-18 17:02:47 +02:00
hold off;
2020-11-12 10:30:56 +01:00
< span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'xscale'< / span > , < span class = "org-string" > 'log'< / span > ); < span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'yscale'< / span > , < span class = "org-string" > 'log'< / span > );
xlabel(< span class = "org-string" > 'Frequency [Hz]'< / span > ); ylabel(< span class = "org-string" > 'ASD of the Velocity $\left[\frac{m/s}{\sqrt{Hz}}\right]$'< / span > );
xlim([0.1, 500]);
2019-04-18 09:20:31 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "org8ecfd08" class = "figure" >
2019-04-18 17:02:47 +02:00
< p > < img src = "figs/intrumental_noise_velocity.png" alt = "intrumental_noise_velocity.png" / >
2019-04-18 09:20:31 +02:00
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 13: < / span > Instrumental Noise and Measurement in \(m/s/\sqrt{Hz}\)< / p >
2019-04-18 17:02:47 +02:00
< / div >
< / div >
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-org52fed77" class = "outline-2" >
< h2 id = "org52fed77" > < span class = "section-number-2" > 3< / span > Compare axis< / h2 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-2" id = "text-3" >
2019-04-18 17:11:25 +02:00
< p >
2020-11-12 10:30:56 +01:00
< a id = "org4f6c263" > < / a >
2019-04-18 17:11:25 +02:00
< / p >
2019-05-10 16:06:43 +02:00
2020-11-12 10:30:56 +01:00
< div class = "note" id = "org7d50518" >
2019-05-10 16:06:43 +02:00
< p >
All the files (data and Matlab scripts) are accessible < a href = "data/huddle_test_compare_axis.zip" > here< / a > .
< / p >
< / div >
2019-04-18 17:02:47 +02:00
< / div >
2019-04-18 17:11:25 +02:00
2020-11-12 10:30:56 +01:00
< div id = "outline-container-org7f52c4b" class = "outline-3" >
< h3 id = "org7f52c4b" > < span class = "section-number-3" > 3.1< / span > Load data< / h3 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-3" id = "text-3-1" >
2019-04-18 17:11:25 +02:00
< p >
We first load the data for the three axis.
< / p >
2019-04-18 17:02:47 +02:00
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > z = load(< span class = "org-string" > 'mat/data_001.mat'< / span > , < span class = "org-string" > 't'< / span > , < span class = "org-string" > 'x1'< / span > , < span class = "org-string" > 'x2'< / span > );
east = load(< span class = "org-string" > 'mat/data_002.mat'< / span > , < span class = "org-string" > 't'< / span > , < span class = "org-string" > 'x1'< / span > , < span class = "org-string" > 'x2'< / span > );
north = load(< span class = "org-string" > 'mat/data_003.mat'< / span > , < span class = "org-string" > 't'< / span > , < span class = "org-string" > 'x1'< / span > , < span class = "org-string" > 'x2'< / span > );
2019-04-18 17:02:47 +02:00
< / pre >
2019-04-18 09:20:31 +02:00
< / div >
< / div >
< / div >
2019-04-18 17:02:47 +02:00
2020-11-12 10:30:56 +01:00
< div id = "outline-container-orgb105293" class = "outline-3" >
< h3 id = "orgb105293" > < span class = "section-number-3" > 3.2< / span > Compare PSD< / h3 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-3" id = "text-3-2" >
2019-04-18 17:11:25 +02:00
< p >
The PSD for each axis of the two geophones are computed.
< / p >
2019-04-18 17:02:47 +02:00
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > [pz1, fz] = pwelch(z.x1, hanning(ceil(length(z.x1)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (z.t(2)< span class = "org-type" > -< / span > z.t(1)));
[pz2, < span class = "org-type" > ~< / span > ] = pwelch(z.x2, hanning(ceil(length(z.x2)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (z.t(2)< span class = "org-type" > -< / span > z.t(1)));
2019-04-18 17:02:47 +02:00
2020-11-12 10:30:56 +01:00
[pe1, fe] = pwelch(east.x1, hanning(ceil(length(east.x1)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (east.t(2)< span class = "org-type" > -< / span > east.t(1)));
[pe2, < span class = "org-type" > ~< / span > ] = pwelch(east.x2, hanning(ceil(length(east.x2)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (east.t(2)< span class = "org-type" > -< / span > east.t(1)));
2019-04-18 17:02:47 +02:00
2020-11-12 10:30:56 +01:00
[pn1, fn] = pwelch(north.x1, hanning(ceil(length(north.x1)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (north.t(2)< span class = "org-type" > -< / span > north.t(1)));
[pn2, < span class = "org-type" > ~< / span > ] = pwelch(north.x2, hanning(ceil(length(north.x2)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (north.t(2)< span class = "org-type" > -< / span > north.t(1)));
2019-04-18 17:02:47 +02:00
< / pre >
2019-04-17 18:05:07 +02:00
< / div >
2019-04-18 17:02:47 +02:00
2019-04-18 17:11:25 +02:00
< p >
2020-11-12 10:30:56 +01:00
We compare them. The result is shown on figure < a href = "#orgf6d7705" > 14< / a > .
2019-04-18 17:11:25 +02:00
< / p >
2019-04-18 17:02:47 +02:00
2020-11-12 10:30:56 +01:00
< div id = "orgf6d7705" class = "figure" >
2019-04-18 17:02:47 +02:00
< p > < img src = "figs/compare_axis_psd.png" alt = "compare_axis_psd.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 14: < / span > Compare the measure PSD of the two geophones for the three axis< / p >
2019-04-18 17:02:47 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-orgbce3b98" class = "outline-3" >
< h3 id = "orgbce3b98" > < span class = "section-number-3" > 3.3< / span > Compare TF< / h3 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-3" id = "text-3-3" >
2019-04-18 17:11:25 +02:00
< p >
The transfer functions from one geophone to the other are also computed for each axis.
2020-11-12 10:30:56 +01:00
The result is shown on figure < a href = "#org1936d29" > 15< / a > .
2019-04-18 17:11:25 +02:00
< / p >
2019-04-18 17:02:47 +02:00
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > [Tz, fz] = tfestimate(z.x1, z.x2, hanning(ceil(length(z.x1)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (z.t(2)< span class = "org-type" > -< / span > z.t(1)));
[Te, fe] = tfestimate(east.x1, east.x2, hanning(ceil(length(east.x1)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (east.t(2)< span class = "org-type" > -< / span > east.t(1)));
[Tn, fn] = tfestimate(north.x1, north.x2, hanning(ceil(length(north.x1)< span class = "org-type" > /< / span > 100)), [], [], 1< span class = "org-type" > /< / span > (north.t(2)< span class = "org-type" > -< / span > north.t(1)));
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "org1936d29" class = "figure" >
2019-04-18 17:02:47 +02:00
< p > < img src = "figs/compare_tf_axis.png" alt = "compare_tf_axis.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 15: < / span > Compare the transfer function from one geophone to the other for the 3 axis< / p >
2019-04-18 17:02:47 +02:00
< / div >
< / div >
< / div >
< / div >
2019-04-18 17:11:25 +02:00
2020-11-12 10:30:56 +01:00
< div id = "outline-container-orgd32606c" class = "outline-2" >
< h2 id = "orgd32606c" > < span class = "section-number-2" > 4< / span > Appendix< / h2 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-2" id = "text-4" >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "outline-container-orge652a09" class = "outline-3" >
< h3 id = "orge652a09" > < span class = "section-number-3" > 4.1< / span > Computation of coherence from PSD and CSD< / h3 >
2019-04-18 17:02:47 +02:00
< div class = "outline-text-3" id = "text-4-1" >
< p >
2020-11-12 10:30:56 +01:00
< a id = "org06003b1" > < / a >
2019-04-18 17:02:47 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > load(< span class = "org-string" > 'mat/data_001.mat'< / span > , < span class = "org-string" > 't'< / span > , < span class = "org-string" > 'x1'< / span > , < span class = "org-string" > 'x2'< / span > );
dt = t(2) < span class = "org-type" > -< / span > t(1);
Fs = 1< span class = "org-type" > /< / span > dt;
win = hanning(ceil(length(x1)< span class = "org-type" > /< / span > 100));
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
< div class = "org-src-container" >
2020-11-12 10:30:56 +01:00
< pre class = "src src-matlab" > pxy = cpsd(x1, x2, win, [], [], Fs);
pxx = pwelch(x1, win, [], [], Fs);
pyy = pwelch(x2, win, [], [], Fs);
coh = mscohere(x1, x2, win, [], [], Fs);
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
< div class = "org-src-container" >
< pre class = "src src-matlab" > < span class = "org-type" > figure< / span > ;
hold on;
2020-11-12 10:30:56 +01:00
plot(f, abs(pxy)< span class = "org-type" > .^< / span > 2< span class = "org-type" > ./< / span > abs(pxx)< span class = "org-type" > ./< / span > abs(pyy), < span class = "org-string" > '-'< / span > );
plot(f, coh, < span class = "org-string" > '--'< / span > );
2019-04-18 17:02:47 +02:00
hold off;
2020-11-12 10:30:56 +01:00
< span class = "org-type" > set< / span > (< span class = "org-variable-name" > gca< / span > , < span class = "org-string" > 'xscale'< / span > , < span class = "org-string" > 'log'< / span > );
xlabel(< span class = "org-string" > 'Frequency'< / span > ); ylabel(< span class = "org-string" > 'Coherence'< / span > );
xlim([1, 500]);
2019-04-18 17:02:47 +02:00
< / pre >
< / div >
2020-11-12 10:30:56 +01:00
< div id = "org3ae3bcc" class = "figure" >
2019-04-18 17:02:47 +02:00
< p > < img src = "figs/comp_coherence_formula.png" alt = "comp_coherence_formula.png" / >
< / p >
2019-05-03 11:32:54 +02:00
< p > < span class = "figure-number" > Figure 16: < / span > Comparison of < code > mscohere< / code > and manual computation< / p >
2019-04-18 17:02:47 +02:00
< / div >
< / div >
< / div >
< / div >
< p >
< / p >
2020-11-12 10:30:56 +01:00
< style > . csl-entry { text-indent : -1.5 em ; margin-left : 1.5 em ; } < / style > < h2 class = 'citeproc-org-bib-h2' > Bibliography< / h2 >
< div class = "csl-bib-body" >
< div class = "csl-entry" > < a name = "citeproc_bib_item_1" > < / a > Barzilai, Aaron, Tom VanZandt, and Tom Kenny. 1998. “Technique for Measurement of the Noise of a Sensor in the Presence of Large Background Signals.” < i > Review of Scientific Instruments< / i > 69 (7):2767– < a href = "https://doi.org/10.1063/1.1149013" > https://doi.org/10.1063/1.1149013< / a > .< / div >
< / div >
2019-04-17 18:05:07 +02:00
< / div >
< div id = "postamble" class = "status" >
2019-05-10 16:06:43 +02:00
< p class = "author" > Author: Dehaeze Thomas< / p >
2020-11-12 10:30:56 +01:00
< p class = "date" > Created: 2020-11-12 jeu. 10:28< / p >
2019-04-17 18:05:07 +02:00
< / div >
< / body >
< / html >
