A collection of useful (or at least interesting) Matlab functions that I've created over the years as well as some from around the place that have caught my eye.
In my earlier post (Analysis of suspected MH370 pings) I manually inspected a audio waveform to determine ping times of arrival. Whilst this is managable for the short segment I looked at, to perform this analysis over many hours or days of recordings would be a bit tedious and probably error prone.
This afternoon I put together a quick MATLAB script to automate this process and ran it over the whole Youtube video audio.
The algorithm is as follows:
Read acoustic data from file in 1s segments
For each 1s segment use a Hilbert Transform (using FFT) to calculate the envelope of the waveform
Calculate noise by smoothing the envelope with a window of length 0.1s (or 10x expected pulse width)
Envelope (red); Smoothed (black)
Calculate SNR as the envelope amplitude divided by the noise
Locate rising edges in the SNR that exceed the detection threshold (ie start of pulse)
SNR (green) with detected pulse (x)
At end of processing plot detection times of arrival module ping interval (1.106s from previous analysis)
Time of Arrival Ping Graph
In the resulting plot that there are 2 distinct sections of coherent detections. The others are false detects from other in-band transients.
Interestingly I had not noticed the 2 individual segments previously, but upon a second listen a distinct discontinuity can be heard hence the 'realignment' of the ping times of arrival. This plot shows that the first segment corresponds to a time where the TPL was going away from the source as the times of arrival are increasing. The second segment has the TPL coming closer to the source near CPA due to the decrease and stabilisation of time of arrival.
So this is a bit of a deviation from my normal MATLAB oriented posts. In-fact MATLAB was not used once in this analysis!
I've been following the MH370 disappearance since it happened. With such a lack of confirmed information there have been a lot of interesting theories discussed on multiple forums (Whirlpool and PPRUNE are a couple of interesting ones).
As I have a
knowledge in acoustics (I am an Analyst / Software Engineer for an
underwater acoustics company and I specialise in tracking things with
pingers fitted to them), I've been trying to spread facts when I know them and try and dispel myths.
Today I saw an interesting video from Sky News on a forum of the suspected MH370 pings recorded by ADV Ocean Shield using the SUBSALV Towed Pinger Locator TPL-25 (pinger audible 18 seconds in):
I thought I would see if I could glean any
information from the audio in the publicly released video. In particular I was wanting to look at the inter-ping interval to see
how steady it was to try and determine if the noise is man-made or
natural (ie marine mammal noise).
The signal is quite noisy with ambient ship noise and other transients (scuffing of shoes on deck!) but a repeated click can be heard.
Process
Extracted the audio from the Youtube video using a web-based MP3 extractor
Opened the MP3 in Adobe Audition and trimmed the file to the time segment corresponding to good audio (start time 41.787s from beginning of file, duration 22.462 seconds)
Bandpass filtered the audio between 2kHz and 10kHz to remove LF ambient noise (NOTE: the 'pinger' is audible at about 3.4kHz however it obviously has been demodulated down from the raw signal observed on the TPL-25 to bring into the audible range)
Filtered Time Domain Signal - Red ticks added at repeated 'ping' noise
A manual search of the waveform (zooming in tight where necessary) was performed to measure the start time of each ping transient to within a couple of wavelengths
Results
The detection times of the pings are listed in the table below along with some basic statistics.
Ping #
Time (s)
Interval
(s)
1
0.758
2
1.861
1.103
3
2.964
1.103
4
4.067
1.103
5
5.171
1.104
6
6.274
1.103
7
7.38
1.106
8
8.486
1.106
9
9.592
1.106
10
10.696
1.104
11
11.801
1.105
12
12.907
1.106
13
14.014
1.107
14
15.12
1.106
15
16.226
1.106
16
17.333
1.107
18
19.549
1.108
19
20.657
1.108
20
21.764
1.107
Mean
1.10544
Std Dev
0.00172
Single Ping zoom view - Selected portion is 12ms duration
As you can see, the ping interval is very stable with an average of 1.105s and standard deviation of less than 2ms.
This is incredibly stable over the 20 seconds that I've looked at and very unlikely to be a marine mammal (ie beaked or sperm whale) which DO create broadband pulses at regular intervals, but I've never seen them that regular.
This strongly suggests the recorded signal is man-made. The signal duration is consistently between 10 and 15 ms in duration. Highly probably it is from the Underwater Locator Beacon due to the very close match to the expected ping rate of 0.9Hz or 1 ping per 1.111 seconds, and the expected ping duration of 10ms.
However it COULD be from other repetitive signal sources such as a stray echo sounder nearby; although that would be a soul crushing disappointment to all onboard ADV Ocean Shield and the other involved in the search. Also.. WHAT CRAZY PERSON WOULD PUT A SOUNDER ON DURING AN ACOUSTIC SEARCH!?!
One last thing to note about the observed ping intervals is that they appear to increase slightly over time (order of 4ms). Whilst the sample set isn't long enough to show consistency, I see patterns like these all the time when tracking moving items fitted with pingers from stationary receivers (conceptually the same thing as the stationary pinger and moving receiver on ADV Ocean Shield).
It could be either:
- The decreasing battery voltage is causing the ping interval to increase. However it doesn't seem likely this phenomenon would be observed over such a short time window
- The audio recording segment corresponds to when the towed pinger locator has performed a Closest Point of Approach (CPA) and is now opening in range from the source. When this happens a effect similar to Doppler occurs in that the time between successive pings appears to increase as because between successive pings the receiver has moved further away.
For the segment analysed it corresponds to a slant range increase of 6 metres over 20 seconds, or 0.3m/s or 0.6 knots. If the towed pinger locator is transiting at 2 knots and the depth of the receiver is known, you could estimate a very rough range to the source by looking at the vector component in the horizontal direction if you assume the vessel passed over it. If you had the entire data set you could do much more interesting analysis!
On other forums I've had people ask about using Doppler to localise the pinger. Due to
the pulse length (10ms) there is an inherent 100Hz bandwidth in the
signal, which is equivalent to the frequency shift you would expect if
moving at 4m/s. In other words, Doppler measurements of a single pulse
is NOT going to be useful in attempting to localise.
However, because the pinger has a relatively stable repetition rate
(about 1.106s) you can do a cool trick to localise
it (assuming the source is stationary):
Perform a pass over the area of interest (ie where the pinger is)
Measure the ping times of arrival (either manually with waveform inspection or write an envelope detector)
Plot the ping times of arrival versus the time of arrival modulo ping interval. ie in MATLAB: plot(toa, mod(toa, 1.106), '.');
The plot indicates the relative distance from the source / receiver.
Hopefully the plot should contain a local minimum representing a
closest point of approach (ie like http://imgur.com/v1ZclaX).
This means that the source lies in a line perpendicular to the vessels
course centred at the location of the receiver (ie TPL-25) at CPA time
Perform another pass over the area of interest on a course perpendicular to the first
Measure the ping times of arrival and look for CPA as before
You should now have two intersecting lines. Guess what should be at the intersection? :)
Even if a consistent sequence of pings over CPA is not measured, a single segment can provide an estimate of a bearing to the source due to the reduction (towards) or increase (away) of the ping interval. This relies on having a stable ping repetition rate.
Edit: I have added a new post where I wrote a simple MATLAB script to automatically detect and visualise the MH370 pings: Automatic detection of MH370 ULB pings
Update 19th Apr - Added zoomed in plot of a single ping.
Completed a Physics / Computer Science double degree before joining an acoustics company as an Analyst where MATLAB was my primary algorithm development environment.
Then transitioned into a software engineering role involving developing and implementing signal processing algorithms in MATLAB, C++ and C#. The latter of which has become an increasingly viable environment for signal processing development.
