Man made like buildings and drives can be characterized by the endmembers EM4 and EM8. We simply add the two related abundance maps and show the result below.
%matplotlib inline
import os.path as osp
import numpy as np
import pysptools.util as util
import pysptools.eea as eea
import pysptools.abundance_maps as amp
import pysptools.classifiers as cls
def remove_bands(M):
"""
Remove the bands where there is atmospheric
scattering.
Remove:
[0..4]
[102..110]
[148..169]
[211..end]
"""
p1 = range(5,102)
p2 = range(111,148)
p3 = range(170,211)
Mp = M[:,:,p1+p2+p3]
return Mp
def display(image):
import matplotlib.pyplot as plt
img = plt.imshow(image)
img.set_cmap('gray')
plt.show()
plt.clf()
data_path = '../data1'
sample = '92AV3C.hdr'
data_file = osp.join(data_path, sample)
data, header = util.load_ENVI_file(data_file)
# Remove some bands
data_clean = remove_bands(data)
header['wavelength'] = range(1, data_clean.shape[2]+1)
# Display the Pine Creek cube
util.display(data, 75, 34, 0)
Extract twelve endmembers
ee = eea.NFINDR()
U = ee.extract(data_clean, 12, maxit=5, normalize=False, ATGP_init=True)
ee.display(header, suffix='Pine Creek')
Invert
am = amp.FCLS()
amaps = am.map(data_clean, U, normalize=False)
am.display(colorMap='jet', suffix='Pine Creek')
# apply a threshold and add the two abundance maps
man_made = (amaps[:,:,3] > 0.02)+(amaps[:,:,7] > 0.1)
print 'Roads and buildings'
display(man_made)