Warm hematite drill core

A new version of the hematite example is introduced for the 0.12.0 version. The main improvement is to use FCLS instead of NNLS. With the introduction of FCLS a higher precision is obtained and only four endmembers is needed to map the quartz and, this is new, to define the mask. You can consult the version 0.11.0 for the previous analysis.

Quartz and feldspar have absorption feature in the LWIR range. Hematite (Fe2O3) has no absorption feature in the sensor spectral range. The dark regions correspond to a quartz (SiO2) impurity in the sample (next figure). The broad spectral feature between 1050 and 1275 cm-1 is associated with the Si-O asymmetric stretching mode of quartz. The cube is acquired in the VLW range of infrared, between 867 to 1288 wavenumber (cm-1) (7.76 to 11.54 micrometer) and it have 165 bands. The instrument used to acquire the data come from the Telops compagny.

In [8]:
# Image displayed by the Reveal Viewer software at 1113.15 [cm-1].
import os
import os.path as osp
from IPython.core.display import Image
home = os.environ['HOME']
Image(filename=osp.join(home, 'dev/pysptools/doc/hem/hem_pic1.png'))
Out[8]:

First, we load the cube.

In [9]:
# Run on Python 2.7 and 3.x

from __future__ import print_function

%matplotlib inline

import pysptools.util as util
import pysptools.eea as eea
import pysptools.abundance_maps as amp
import numpy as np

# Load the cube
data_path = os.environ['PYSPTOOLS_DATA']

sample = 'hematite.hdr'

data_file = osp.join(data_path, sample)
data, header = util.load_ENVI_file(data_file)

# Telops cubes are flipped left-right
# Flipping them again restore the orientation
data = np.fliplr(data)
In [10]:
def display(image, colormap):
    import matplotlib.pyplot as plt
    img = plt.imshow(image)
    img.set_cmap(colormap)
    plt.colorbar()
    plt.show()
    plt.clf()

The analysis is made in tree steps. First, we extract four endmembers.

In [11]:
ee = eea.NFINDR()
U = ee.extract(data, 4, maxit=5, normalize=True, ATGP_init=True)
ee.display(header, suffix='Hematite example')

<matplotlib.figure.Figure at 0x7f6e9f483790>

Next, the related abundance maps are generated with FCLS.

In [12]:
fcls = amp.FCLS()
amaps = fcls.map(data, U, normalize=True)
fcls.display(colorMap='jet', columns=2, suffix='Hematite example')

<matplotlib.figure.Figure at 0x7f6e9f430850>

For the last step, we use the EM1 and EM4 spectra to create the mask and the quartz images respectively. Creation of the quartz image is straightforward. The EM1 spectra is a background spectra. Inverting it give a mask for the drill core.

In [13]:
# Print the quartz image
# EM2 == quartz on Python 2.7 and EM2 == quartz on Python 3.3
#if sys.version_info[:2] == (2,7):
#    quartz = amaps[:,:,1]
#if sys.version_info[:2] == (3,3):
quartz = amaps[:,:,3]
display(quartz, 'spectral')

<matplotlib.figure.Figure at 0x7f6e9f291f90>
In [14]:
# Print the mask
# EM1 == background, we use the backgroud to isolate the drill core
# and define the mask
mask = (amaps[:,:,0] < 0.2)
display(mask, 'gray')

<matplotlib.figure.Figure at 0x7f6e9f3fed50>

Finally, we extract some statistics.

In [15]:
# pixels sum
rock_surface = np.sum(mask)
quartz_surface = np.sum(quartz > 0.16)
print('Some statistics')
print('  Drill core surface (mask) in pixels:', rock_surface)
print('  Quartz surface in pixels:', quartz_surface)
print('  Hematite surface in pixels:', rock_surface - quartz_surface)

Some statistics
  Drill core surface (mask) in pixels: 1718
  Quartz surface in pixels: 1163
  Hematite surface in pixels: 555