Loading

Seismic Facies Identification Challenge

[Explainer] Additional Features with the Stationary Wavelet Transform (SWT)

Demonstrate the use of the Stationary Wavelet Transform (SWT) to perform a frequency based decomposition (filter banks) on the input time series data

edward_beeching

Hello everyone, I was previously a Senior Geophysicist for a Seismic Processing company and I now work in Deep Reinforcement Learning for robotics. I find this competition really interesting as it combines my previous experience in Geophysics with my current area of research, Deep Neural Networks.

Here in this Colab Notebook I demonstrate the use of the Stationary Wavelet Transform (SWT) to perform a frequency based decomposition (filter banks) on the input time series data. For my submission I have used this as input to a UNET architecture, written in PyTorch. I would be happy to share the UNET implementation if this post gets some attention.

Best of luck with the competition!

 

 

In [1]:
# **Addtional features with the Stationary Wavelet Transform**

 # By Edward Beeching, an ex Senior Geophysicist
In [7]:
!wget "https://datasets.aicrowd.com/default/aicrowd-public-datasets/seamai-facies-challenge/v0.1/public/data_train.npz"
!wget "https://datasets.aicrowd.com/default/aicrowd-public-datasets/seamai-facies-challenge/v0.1/public/labels_train.npz"
--2020-09-24 10:55:23--  https://datasets.aicrowd.com/default/aicrowd-public-datasets/seamai-facies-challenge/v0.1/public/labels_train.npz
Resolving datasets.aicrowd.com (datasets.aicrowd.com)... 35.189.208.115
Connecting to datasets.aicrowd.com (datasets.aicrowd.com)|35.189.208.115|:443... connected.
HTTP request sent, awaiting response... 302 FOUND
Location: https://s3.us-west-002.backblazeb2.com/aicrowd-public-datasets/seamai-facies-challenge/v0.1/public/labels_train.npz?X-Amz-Algorithm=AWS4-HMAC-SHA256&X-Amz-Credential=002ae2491b744be0000000002%2F20200924%2Fus-west-002%2Fs3%2Faws4_request&X-Amz-Date=20200924T105524Z&X-Amz-Expires=3600&X-Amz-SignedHeaders=host&X-Amz-Signature=2c0cee92c2b941a6ad8e58c6ad0675c083fd0f41ca0542fe5a1edab05fba6d0f [following]
--2020-09-24 10:55:24--  https://s3.us-west-002.backblazeb2.com/aicrowd-public-datasets/seamai-facies-challenge/v0.1/public/labels_train.npz?X-Amz-Algorithm=AWS4-HMAC-SHA256&X-Amz-Credential=002ae2491b744be0000000002%2F20200924%2Fus-west-002%2Fs3%2Faws4_request&X-Amz-Date=20200924T105524Z&X-Amz-Expires=3600&X-Amz-SignedHeaders=host&X-Amz-Signature=2c0cee92c2b941a6ad8e58c6ad0675c083fd0f41ca0542fe5a1edab05fba6d0f
Resolving s3.us-west-002.backblazeb2.com (s3.us-west-002.backblazeb2.com)... 206.190.215.254
Connecting to s3.us-west-002.backblazeb2.com (s3.us-west-002.backblazeb2.com)|206.190.215.254|:443... connected.
HTTP request sent, awaiting response... 200 
Length: 7160425 (6.8M) [application/octet-stream]
Saving to: ‘labels_train.npz’

labels_train.npz    100%[===================>]   6.83M  10.3MB/s    in 0.7s    

2020-09-24 10:55:25 (10.3 MB/s) - ‘labels_train.npz’ saved [7160425/7160425]

In [3]:
import numpy as np 
import matplotlib.pyplot as plt
import pywt
%matplotlib inline
In [8]:
data=np.load('data_train.npz')['data'][:,:,580]
label = np.load('labels_train.npz')['labels'][:,:,580]
plt.rcParams['figure.figsize'] = [12, 8]
plt.imshow(data,cmap='gray',interpolation='none')
plt.show()

We will now apply the stationary wavelet transform to decompose the signal into frequency bands or "scales". I use a the coiflet 4 wavelet here, but there are many, have a look on pywavelets wavelet browser http://wavelets.pybytes.com/wavelet/coif4/

In [9]:
coefs = pywt.swt(np.pad(data,((0,18), (0,0))), "coif4", level=6, trim_approx=True, axis=0)
In [16]:
plt.rcParams['figure.figsize'] = [32, 16]

plt.subplot(2,4,1)
plt.imshow(data,cmap='gray',interpolation='none')
plt.title("Full bandwidth")
plt.subplot(2,4,2)
plt.imshow(label)
plt.title("Labels")
# I assume 4ms sampling
freqs = ['62-125', ["32-62"], ["16-32"], ["8-16"], ["4-8"], ["2-4"]]


for i, (coef, freq) in enumerate(zip(coefs[::-1], freqs),3):
  plt.subplot(2,4,i)
  plt.imshow(coef,cmap='gray',interpolation='none')
  plt.title(freq)
plt.tight_layout()
plt.show()