#!/usr/bin/python
#
# Copyright (c) 2013-2015, Zhouhan LIN
# All rights reserved.
#
# Redistribution and use in source and binary forms, with or without
# modification, are permitted provided that the following conditions are met:
#
# * Redistributions of source code must retain the above copyright notice, this
# list of conditions and the following disclaimer.
#
# * Redistributions in binary form must reproduce the above copyright notice,
# this list of conditions and the following disclaimer in the documentation
# and/or other materials provided with the distribution.
#
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
import pdb
import time
import sys
import scipy.io as sio
import numpy
import theano
import pylab as pl
from sklearn.decomposition import PCA
from sklearn.metrics import confusion_matrix
# the following color map is for generating wholeimage classification figures.
cmap = numpy.asarray( [[0, 0, 0],
[0, 205, 0],
[127, 255, 0],
[46, 139, 87],
[0, 139, 0],
[160, 82, 45],
[0, 255, 255],
[255, 255, 255],
[216, 191, 216],
[255, 0, 0],
[139, 0, 0],
[0, 0, 255],
[255, 255, 0],
[238, 154, 0],
[85, 26, 139],
[255, 127, 80]], dtype='int32')
def result_analysis(prediction, train_truth, valid_truth, test_truth,
verbose=False):
assert prediction.shape == test_truth.shape
print "Detailed information in each category:"
print " Number of Samples"
print "Class No. TRAIN VALID TEST RightCount RightRate"
for i in xrange(test_truth.min(), test_truth.max()+1):
right_prediction = ( (test_truth-prediction) == 0 )
right_count = numpy.sum(((test_truth==i) * right_prediction)*1)
print "%d\t\t%d\t%d\t%d\t%d\t%f" % \
(i,
numpy.sum((train_truth==i)*1),
numpy.sum((valid_truth==i)*1),
numpy.sum((test_truth==i)*1),
right_count,
right_count * 1.0 / numpy.sum((test_truth==i)*1)
)
total_right_count = numpy.sum(right_prediction*1)
print "Overall\t\t%d\t%d\t%d\t%d\t%f" % \
(train_truth.size,
valid_truth.size,
test_truth.size,
total_right_count,
total_right_count * 1.0 / test_truth.size
)
cm = confusion_matrix(test_truth, prediction)
pr_a = cm.trace()*1.0 / test_truth.size
pr_e = ((cm.sum(axis=0)*1.0/test_truth.size) * \
(cm.sum(axis=1)*1.0/test_truth.size)).sum()
k = (pr_a - pr_e) / (1 - pr_e)
print "kappa index of agreement: %f" % k
print "confusion matrix:"
print cm
# Show confusion matrix
pl.matshow(cm)
pl.title('Confusion matrix')
pl.colorbar()
if verbose:
pl.show()
else:
filename = 'conf_mtx_' + str(time.time()) + '.png'
pl.savefig(filename)
#-------------------------------------------------------------------------------
"""
The scale_to_unit_interval() and tile_raster_images() functions are from the
Deep Learning Tutorial repo:
https://github.com/lisa-lab/DeepLearningTutorials
Below are the corresponding licence.
LICENSE
=======
Copyright (c) 2010--2015, Deep Learning Tutorials Development Team
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:
* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright
notice, this list of conditions and the following disclaimer in the
documentation and/or other materials provided with the distribution.
* Neither the name of Theano nor the names of its contributors may be
used to endorse or promote products derived from this software without
specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS ''AS IS'' AND ANY
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY
DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
"""
def scale_to_unit_interval(ndar, eps=1e-8):
""" Scales all values in the ndarray ndar to be between 0 and 1 """
ndar = ndar.copy()
ndar -= ndar.min()
ndar *= 1.0 / (ndar.max() + eps)
return ndar
def tile_raster_images(X, img_shape, tile_shape, tile_spacing=(0, 0),
scale_rows_to_unit_interval=True,
output_pixel_vals=True):
"""
Transform an array with one flattened image per row, into an array in
which images are reshaped and layed out like tiles on a floor.
This function is useful for visualizing datasets whose rows are images,
and also columns of matrices for transforming those rows
(such as the first layer of a neural net).
:type X: a 2-D ndarray or a tuple of 4 channels, elements of which can
be 2-D ndarrays or None;
:param X: a 2-D array in which every row is a flattened image.
:type img_shape: tuple; (height, width)
:param img_shape: the original shape of each image
:type tile_shape: tuple; (rows, cols)
:param tile_shape: the number of images to tile (rows, cols)
:param output_pixel_vals: if output should be pixel values (i.e. int8
values) or floats
:param scale_rows_to_unit_interval: if the values need to be scaled before
being plotted to [0,1] or not
:returns: array suitable for viewing as an image.
(See:`Image.fromarray`.)
:rtype: a 2-d array with same dtype as X.
"""
assert len(img_shape) == 2
assert len(tile_shape) == 2
assert len(tile_spacing) == 2
out_shape = [
(ishp + tsp) * tshp - tsp
for ishp, tshp, tsp in zip(img_shape, tile_shape, tile_spacing)
]
if isinstance(X, tuple):
assert len(X) == 4
# Create an output numpy ndarray to store the image
if output_pixel_vals:
out_array = numpy.zeros((out_shape[0], out_shape[1], 4),
dtype='uint8')
else:
out_array = numpy.zeros((out_shape[0], out_shape[1], 4),
dtype=X.dtype)
#colors default to 0, alpha defaults to 1 (opaque)
if output_pixel_vals:
channel_defaults = [0, 0, 0, 255]
else:
channel_defaults = [0., 0., 0., 1.]
for i in xrange(4):
if X[i] is None:
# if channel is None, fill it with zeros of the correct
# dtype
dt = out_array.dtype
if output_pixel_vals:
dt = 'uint8'
out_array[:, :, i] = numpy.zeros(
out_shape,
dtype=dt
) + channel_defaults[i]
else:
# use a recurrent call to compute the channel and store it
# in the output
out_array[:, :, i] = tile_raster_images(
X[i], img_shape, tile_shape, tile_spacing,
scale_rows_to_unit_interval, output_pixel_vals)
return out_array
else:
# if we are dealing with only one channel
H, W = img_shape
Hs, Ws = tile_spacing
# generate a matrix to store the output
dt = X.dtype
if output_pixel_vals:
dt = 'uint8'
out_array = numpy.zeros(out_shape, dtype=dt)
for tile_row in xrange(tile_shape[0]):
for tile_col in xrange(tile_shape[1]):
if tile_row * tile_shape[1] + tile_col < X.shape[0]:
this_x = X[tile_row * tile_shape[1] + tile_col]
if scale_rows_to_unit_interval:
# if we should scale values to be between 0 and 1
# do this by calling the `scale_to_unit_interval`
# function
this_img = scale_to_unit_interval(
this_x.reshape(img_shape))
else:
this_img = this_x.reshape(img_shape)
# add the slice to the corresponding position in the
# output array
c = 1
if output_pixel_vals:
c = 255
out_array[
tile_row * (H + Hs): tile_row * (H + Hs) + H,
tile_col * (W + Ws): tile_col * (W + Ws) + W
] = this_img * c
return out_array
#-------------------------------------------------------------------------------
def PCA_tramsform_img(img=None, n_principle=3):
"""
This function trainsforms an HSI by 1-D PCA. PCA is fitted on the whole data
and is conducted on the spectral dimension, rendering the image from size
length * width * dim to length * width * n_principle.
Parameters:
img: initial unregularizaed HSI.
n_principle: Target number of principles we want.
Return:
reg_img: Regularized, transformed image.
WARNNING: RELATIVE ENERGY BETWEEN PRINCIPLE COMPONENTS CHANGED IN THIS
IMPLEMENTATION. YOU MAY NEED TO ADD PENALTY MULTIPLIERS IN THE HIGHER NETWORKS
TO REIMBURSE IT.
"""
length = img.shape[0]
width = img.shape[1]
dim = img.shape[2]
# reshape img, HORIZONTALLY strench the img, without changing the spectral dim.
reshaped_img = numpy.asarray(img.reshape(length*width, dim),
dtype=theano.config.floatX)
pca = PCA(n_components=n_principle)
pca_img = pca.fit_transform(reshaped_img)
# Regularization: Think about energy of each principles here.
reg_img = scale_to_unit_interval(ndar=pca_img, eps=1e-8)
reg_img = numpy.asarray(reg_img.reshape(length, width, n_principle),
dtype=theano.config.floatX)
energy_dist = pca.explained_variance_ratio_
residual = 1 - numpy.sum(energy_dist[0: n_principle])
return reg_img, energy_dist, residual
def T_pca_constructor(hsi_img=None, gnd_img=None, n_principle=3, window_size=1,
flag='supervised', merge=False):
"""
This function constructs the spectral and spatial facade for each training
pixel.
Spectral data returned are simply spectra.
spatial data returned are the former n_principle PCs of a neibor region
around each extracted pixel. Size of the neibor region is determined by
window_size. And all the values for a pixel are flattened to 1-D size. So
data_spatial is finally a 2-D numpy.array.
All the returned data are regularized to [0, 1]. Set window_size=1 to
get pure spectral returnings.
Parameters:
hsi_img=None: 3-D numpy.ndarray, dtype=float, storing initial
hyperspectral image data.
gnd_img=None: 2-D numpy.ndarray, dtype=int, containing tags for pixeles.
The size is the same to the hsi_img size, but with only
1 band.
n_principle: Target number of principles we want to consider in the
spatial infomation.
window_size: Determins the scale of spatial information incorporated.
Must be odd.
flag: For 'unsupervised', all possible pixels except marginals
are processed. For 'supervised', only pixels with
non-zero tags are processed.
Return:
data_spectral: 2-D numpy.array, sized (sample number) * (band number).
Consists of regularized spectra for all extracted pixels.
data_spatial: 2-D numpy.array, sized (sample number) * (window_size^2).
Consists the former n_principle PCs of a neibor region
around each extracted pixel. Size of the neibor region
is determined by window_size.
gndtruth: 1-D numpy.array, sized (sample number) * 1. Truth value
for each extracted pixel.
extracted_pixel_ind:2-D numpy.array, sized (length) * (width). Indicating
which pixels are selected.
"""
# PCA transformation
pca_img, _, _ = PCA_tramsform_img(img=hsi_img, n_principle=n_principle)
# Regularization
hsi_img = scale_to_unit_interval(ndar=hsi_img, eps=1e-8)
length = hsi_img.shape[0]
width = hsi_img.shape[1]
dim = hsi_img.shape[2]
# reshape img, HORIZONTALLY strench the img, without changing the spectral dim.
reshaped_img = numpy.asarray(hsi_img.reshape(length*width, dim),
dtype=theano.config.floatX)
reshaped_gnd = gnd_img.reshape(gnd_img.size)
# mask ensures marginal pixels eliminated, according to window_size
threshold = (window_size-1) / 2
if window_size >= 1 and window_size < width-1 and window_size < length-1:
mask_false = numpy.array([False, ] * width)
mask_true = numpy.hstack((numpy.array([False, ] * threshold, dtype='bool'),
numpy.array([True, ] * (width-2*threshold)),
numpy.array([False, ] * threshold, dtype='bool')))
mask = numpy.vstack((numpy.tile(mask_false, [threshold, 1]),
numpy.tile(mask_true, [length-2*threshold, 1]),
numpy.tile(mask_false, [threshold, 1])))
reshaped_mask = mask.reshape(mask.size)
else:
print >> sys.stderr, ('window_size error. choose 0 < window_size < width-1')
# construct groundtruth, and determine which pixel to process
if flag == 'supervised':
extracted_pixel_ind = (reshaped_gnd > 0) * reshaped_mask
gndtruth = reshaped_gnd[extracted_pixel_ind]
extracted_pixel_ind = numpy.arange(reshaped_gnd.size)[extracted_pixel_ind]
elif flag == 'unsupervised':
extracted_pixel_ind = numpy.arange(reshaped_gnd.size)[reshaped_mask]
gndtruth = numpy.array([], dtype='int')
else:
print >> sys.stderr, ('\"flag\" parameter error. ' +
'What type of learning you are doing?')
return
# construct data_spectral
data_spectral = reshaped_img[extracted_pixel_ind, :]
# construct data_spatial
if window_size == 1:
data_spatial = numpy.array([])
else:
data_spatial = numpy.zeros([extracted_pixel_ind.size,
window_size * window_size * n_principle],
dtype=theano.config.floatX)
i = 0
for ipixel in extracted_pixel_ind:
ipixel_h = ipixel % width
ipixel_v = ipixel / width
data_spatial[i, :] = \
pca_img[ipixel_v-threshold : ipixel_v+threshold+1,
ipixel_h-threshold : ipixel_h+threshold+1, :].reshape(
window_size*window_size*n_principle)
i += 1
# if we want to merge data, merge it
if merge:
data_spectral = numpy.hstack((data_spectral, data_spatial))
return data_spectral, data_spatial, gndtruth, extracted_pixel_ind
def train_valid_test(data, ratio=[6, 2, 2], batch_size=50, random_state=None):
"""
This function splits data into three parts, according to the "ratio" parameter
given in the lists indicating training, validating, testing data ratios.
data: a list containing:
1. A 2-D numpy.array object, with each patterns listed in
ROWs. Input data dimension MUST be larger than 1.
2. A 1-D numpy.array object, tags for each pattern.
'0' indicates that the tag for the corrresponding
pattern is unknown.
ratio: A list having 3 elements, indicating ratio of training,
validating and testing data ratios respectively.
batch_size: bathc_size helps to return an appropriate size of training
samples, which has divisibility over batch_size.
NOTE: batch_size cannot be larger than the minimal size of
all the trainin, validate and test dataset!
random_state: If we give the same random state and the same ratio on the
same data, the function will yield a same split for each
function call.
return:
[train_data_x, train_data_y]:
[valid_data_x, valid_data_y]:
[test_data_x , test_data_y ]:
Lists containing 2 numpy.array object, first for data and
second for truth. They are for training, validate and test
respectively. All the tags are integers in the range
[0, data[1].max()-1].
split_mask
"""
rand_num_generator = numpy.random.RandomState(random_state)
#---------------------------split dataset-----------------------------------
random_mask = rand_num_generator.random_integers(1, sum(ratio), data[0].shape[0])
split_mask = numpy.array(['tests', ] * data[0].shape[0])
split_mask[random_mask <= ratio[0]] = 'train'
split_mask[(random_mask <= ratio[1]+ratio[0]) * (random_mask > ratio[0])] = 'valid'
train_data_x = data[0][split_mask == 'train', :]
train_data_y = data[1][split_mask == 'train']-1
valid_data_x = data[0][split_mask == 'valid', :]
valid_data_y = data[1][split_mask == 'valid']-1
test_data_x = data[0][split_mask == 'tests', :]
test_data_y = data[1][split_mask == 'tests']-1
# tackle the batch size mismatch problem
mis_match = train_data_x.shape[0] % batch_size
if mis_match != 0:
mis_match = batch_size - mis_match
train_data_x = numpy.vstack((train_data_x, train_data_x[0:mis_match, :]))
train_data_y = numpy.hstack((train_data_y, train_data_y[0:mis_match]))
mis_match = valid_data_x.shape[0] % batch_size
if mis_match != 0:
mis_match = batch_size - mis_match
valid_data_x = numpy.vstack((valid_data_x, valid_data_x[0:mis_match, :]))
valid_data_y = numpy.hstack((valid_data_y, valid_data_y[0:mis_match]))
mis_match = test_data_x.shape[0] % batch_size
if mis_match != 0:
mis_match = batch_size - mis_match
test_data_x = numpy.vstack((test_data_x, test_data_x[0:mis_match, :]))
test_data_y = numpy.hstack((test_data_y, test_data_y[0:mis_match]))
return [train_data_x, train_data_y], \
[valid_data_x, valid_data_y], \
[test_data_x , test_data_y], split_mask
def prepare_data(hsi_img=None, gnd_img=None, window_size=7, n_principle=3,
batch_size=50, merge=False, ratio=[6, 2, 2]):
"""
Process the data from file path to splited train-valid-test sets; Binded in
dataset_spectral and dataset_spatial respectively.
Parameters
----------
hsi_img=None: 3-D numpy.ndarray, dtype=float, storing initial
hyperspectral image data.
gnd_img=None: 2-D numpy.ndarray, dtype=int, containing tags for pixeles.
The size is the same to the hsi_img size, but with only
1 band.
window_size: Size of spatial window. Pass an integer 1 if no spatial
infomation needed.
n_principle: This many principles you want to incorporate while
extracting spatial info.
merge: If merge==True, the returned dataset_spectral has
dataset_spatial stacked in the tail of it; else if
merge==False, the returned dataset_spectral and
dataset_spatial will have spectral and spatial information
only, respectively.
Return
------
dataset_spectral:
dataset_spatial:
extracted_pixel_ind:
split_mask:
"""
data_spectral, data_spatial, gndtruth, extracted_pixel_ind = \
T_pca_constructor(hsi_img=hsi_img, gnd_img=gnd_img, n_principle=n_principle,
window_size=window_size, flag='supervised')
################ separate train, valid and test spatial data ###############
[train_spatial_x, train_y], [valid_spatial_x, valid_y], [test_spatial_x, test_y], split_mask = \
train_valid_test(data=[data_spatial, gndtruth], ratio=ratio,
batch_size=batch_size, random_state=123)
# convert them to theano.shared values
train_set_x = theano.shared(value=train_spatial_x, name='train_set_x', borrow=True)
valid_set_x = theano.shared(value=valid_spatial_x, name='valid_set_x', borrow=True)
test_set_x = theano.shared(value=test_spatial_x, name='test_set_x', borrow=True)
train_set_y = theano.shared(value=train_y, name='train_set_y', borrow=True)
valid_set_y = theano.shared(value=valid_y, name='valid_set_y', borrow=True)
test_set_y = theano.shared(value=test_y, name='test_set_y', borrow=True)
dataset_spatial = [(train_set_x, train_set_y), (valid_set_x, valid_set_y),
(test_set_x, test_set_y)]
############### separate train, valid and test spectral data ###############
[train_spectral_x, train_y], [valid_spectral_x, valid_y], [test_spectral_x, test_y], split_mask = \
train_valid_test(data=[data_spectral, gndtruth], ratio=ratio,
batch_size=batch_size, random_state=123)
# if we want to merge data, merge it
if merge:
train_spectral_x = numpy.hstack((train_spectral_x, train_spatial_x))
valid_spectral_x = numpy.hstack((valid_spectral_x, valid_spatial_x))
test_spectral_x = numpy.hstack((test_spectral_x, test_spatial_x))
# convert them to theano.shared values
train_set_x = theano.shared(value=train_spectral_x, name='train_set_x', borrow=True)
valid_set_x = theano.shared(value=valid_spectral_x, name='valid_set_x', borrow=True)
test_set_x = theano.shared(value=test_spectral_x, name='test_set_x', borrow=True)
train_set_y = theano.shared(value=train_y, name='train_set_y', borrow=True)
valid_set_y = theano.shared(value=valid_y, name='valid_set_y', borrow=True)
test_set_y = theano.shared(value=test_y, name='test_set_y', borrow=True)
dataset_spectral = [(train_set_x, train_set_y), (valid_set_x, valid_set_y),
(test_set_x, test_set_y)]
return dataset_spectral, dataset_spatial, extracted_pixel_ind, split_mask
if __name__ == '__main__':
""" Sample usage. """
print '... Testing function result_analysis'
import random
from sklearn import svm, datasets
# import some data to play with
print '... loading Iris data'
iris = datasets.load_iris()
X = iris.data
y = iris.target
n_samples, n_features = X.shape
p = range(n_samples)
random.seed(0)
random.shuffle(p)
X, y = X[p], y[p]
half = int(n_samples / 2)
# Run classifier
print '... classifying'
classifier = svm.SVC(kernel='linear')
y_ = classifier.fit(X[:half], y[:half]).predict(X[half:])
result_analysis(y_, y[:half], numpy.asarray([]), y[half:])
# load .mat files
print '... loading KSC data'
hsi_file = u'/home/hantek/data/hsi_data/kennedy/Kennedy_denoise.mat'
gnd_file = u'/home/hantek/data/hsi_data/kennedy/Kennedy_groundtruth.mat'
data = sio.loadmat(hsi_file)
img = numpy.float_(data['Kennedy176'])
data = sio.loadmat(gnd_file)
gnd_img = data['Kennedy_groundtruth']
gnd_img = gnd_img.astype(numpy.int32)
print '... spliting train-valid-test sets'
dataset_spectral, dataset_spatial, extracted_pixel_ind, split_mask = \
prepare_data(hsi_img=img, gnd_img=gnd_img, window_size=7, n_principle=3, batch_size=50, merge=True)
if raw_input('Spliting finished. Do you want to check the data (Y/n)? ') == 'Y':
spectral_train_x = dataset_spectral[0][0].get_value()
spectral_train_y = dataset_spectral[0][1].get_value()
spectral_valid_x = dataset_spectral[1][0].get_value()
spectral_valid_y = dataset_spectral[1][1].get_value()
spectral_test_x = dataset_spectral[2][0].get_value()
spectral_test_y = dataset_spectral[2][1].get_value()
spatial_train_x = dataset_spatial[0][0].get_value()
spatial_train_y = dataset_spatial[0][1].get_value()
spatial_valid_x = dataset_spatial[1][0].get_value()
spatial_valid_y = dataset_spatial[1][1].get_value()
spatial_test_x = dataset_spatial[2][0].get_value()
spatial_test_y = dataset_spatial[2][1].get_value()
print 'shape of:'
print 'spectral_train_x: \t',
print spectral_train_x.shape,
print 'spectral_train_y: \t',
print spectral_train_y.shape
print 'spectral_valid_x: \t',
print spectral_valid_x.shape,
print 'spectral_valid_y: \t',
print spectral_valid_y.shape
print 'spectral_test_x: \t',
print spectral_test_x.shape,
print 'spectral_test_y: \t',
print spectral_test_y.shape
print 'spatial_train_x: \t',
print spatial_train_x.shape,
print 'spatial_train_y: \t',
print spatial_train_y.shape
print 'spatial_valid_x: \t',
print spatial_valid_x.shape,
print 'spatial_valid_y: \t',
print spatial_valid_y.shape
print 'spatial_test_x: \t',
print spatial_test_x.shape,
print 'spatial_test_y: \t',
print spatial_test_y.shape
print 'total tagged pixel number: %d' % extracted_pixel_ind.shape[0]
print 'split_mask shape: %d' % split_mask.shape
print '... checking tags in spatial and spectral data'
trainset_err = numpy.sum((spectral_train_y-spatial_train_y) ** 2)
validset_err = numpy.sum((spectral_valid_y-spatial_valid_y) ** 2)
testset_err = numpy.sum((spectral_test_y-spatial_test_y) ** 2)
if testset_err + validset_err + trainset_err == 0:
print 'Checking test PASSED.'
else:
print 'Checking test FAILED.'
if raw_input('Do you want to save results to data.mat (Y/n)? ') == 'Y':
print '... saving datasets'
sio.savemat('data.mat', {
'spectral_train_x': dataset_spectral[0][0].get_value(),
'spectral_train_y': dataset_spectral[0][1].get_value(),
'spectral_valid_x': dataset_spectral[1][0].get_value(),
'spectral_valid_y': dataset_spectral[1][1].get_value(),
'spectral_test_x': dataset_spectral[2][0].get_value(),
'spectral_test_y': dataset_spectral[2][1].get_value(),
'spatial_train_x': dataset_spatial[0][0].get_value(),
'spatial_train_y': dataset_spatial[0][1].get_value(),
'spatial_valid_x': dataset_spatial[1][0].get_value(),
'spatial_valid_y': dataset_spatial[1][1].get_value(),
'spatial_test_x': dataset_spatial[2][0].get_value(),
'spatial_test_y': dataset_spatial[2][1].get_value(),
'extracted_pixel_ind': extracted_pixel_ind,
'split_mask': split_mask})
print 'Done.'
