% ------------------------------------------------------------------------
% Script: Test script for machine learning on coral using augmented PCA-SVM
% Author:   Dr. Alan Li
% Date:     Feb 28, 2017
% Project:  Coral Classification with Machine Learning
% NASA Ames Laboratory for Advanced Sensing (LAS)
% NASA Ames Research Center, February 2017
%
% ALL Dependencies:
%   MATLAB Image processing toolbox
%   VLFeat 0.9.20           (can be downloaded at http://www.vlfeat.org/index.html)
%   sift_match.m                        (Included in .zip file)
%   DBSCAN.m                (can be downloaded at http://yarpiz.com/255/ypml110-dbscan-clustering)
%   SIFT_register_down.m                (Included in .zip file)
%   coral_learn_remap.m                 (Included in .zip file)
%   coral_predict_remap.m               (Included in .zip file)
%
% ALL required images and .mat files:
%   NASA Coral Transect 1 - Fluid Lensed - Ofu American Samoa.tiff (Download from https://www.nasa.gov/ames/las/coral-classification)
%   Im_air.tif                          (Included in .zip file)
%   MAP Morphology ID Transect1.mat     (Included in .zip file)
%   TruthImage.mat                      (Included in .zip file)
%
% Description: Script to test coral functions. Currently applied to aerial data.
% Note that areas where no boundaries exist tends to produce mixed results.
% ------------------------------------------------------------------------

%% STARTUP
clear
clc
close all;

customMap = [0 0 0; 1 0.5 0 ; 0 0.5 0 ; 0 0 0.8];       % Colormap for rock, branching, mounding, and sand

% Cells for storing multiple images (only 1 image in our version)
c_ID = cell(1,1);
c_test = cell(1,1);
c_truth = cell(1,1);
c_ori = cell(1,1);      % Original image

% UAV Image
Transect1_UAV_filestr = 'NASA Coral Transect 1 - Fluid Lensed - Ofu American Samoa.tiff'; % Download from https://www.nasa.gov/ames/las/coral-classification
Im_ori = imread(Transect1_UAV_filestr);
c_ori{1} = Im_ori(:,:,1:3);

% Aerial Image
Transect1_Air_filestr = 'Im_air.tif';   % Aerial imagery (Download online)
Im_air = imread(Transect1_Air_filestr);

% Coral Morphology by MAP Estimate
load('MAP Morphology ID Transect1.mat')        
% Change values to: 1 = rock,  2 = branching, 3 = mounding, 4 = sand
morphID1_med(morphID1_med == 0) = 4; morphID1_med(morphID1_med == 3) = 0; morphID1_med(morphID1_med == 2) = 3;
morphID1_med(morphID1_med == 1) = 2; morphID1_med(morphID1_med == 0) = 1;
c_ID_all{1} = uint8(morphID1_med); 
c_ID_cover{1} = uint8(morphID1_med);
c_ID_cover{1}(c_ID_cover{1} == 4) = 1; c_ID_cover{1}(c_ID_cover{1} == 3) = 2;    % 1 = no coral cover, 2 = coral cover

% Truth
load('TruthImage.mat')         
TruthImage = fliplr(TruthImage);
% Change values to: 1 = rock,  2 = branching, 3 = mounding, 4 = sand
TruthImage(TruthImage == 0) = 4; TruthImage(TruthImage == 3) = 0; TruthImage(TruthImage == 2) = 3;
TruthImage(TruthImage == 1) = 2; TruthImage(TruthImage == 0) = 1;
TruthCover = TruthImage;
TruthCover(TruthCover == 4) = 1; TruthCover(TruthCover == 3) = 2;    % 1 = no coral cover, 2 = coral cover
%% Image Registration
[c_Air_down H_down] = SIFT_register_down(Transect1_UAV_filestr,Transect1_Air_filestr);     % Image Registration with low resolution output

%% Learn on MAP Estimate using 3rd order kernel SVM
cl_multiclass = coral_learn_remap(c_ID_all,c_ori,c_ori,4000,0);
cl_2class = coral_learn_remap(c_ID_cover,c_ori,c_ori,4000,0);

%% Predict on aerial image
[c_airpredict,c_airmulti_predict] = coral_predict_remap(cl_2class, c_Air_down, 30, 0, cl_multiclass);

%% Figures

% Truth/Reference Data
figure
ax1 = subplot(1,2,1);
colormap(ax1,gray(2))
imagesc(TruthCover)
axis image; axis off;
title('Coral Cover')
ax2 = subplot(1,2,2);
colormap(ax2,customMap)
imagesc(TruthImage)
axis image; axis off;
title('Morphology')

% MAP Estimate
figure
ax1 = subplot(1,2,1);
colormap(ax1,gray(2))
imagesc(c_ID_cover{1})
axis image; axis off;
title('Coral Cover')
ax2 = subplot(1,2,2);
colormap(ax2,customMap)
imagesc(c_ID_all{1})
axis image; axis off;
title('Morphology')

% Predictions by Augmented PCA-SVM
figure
ax1 = subplot(1,2,1);
colormap(ax1,gray(2))
imagesc(c_airpredict{1})
axis image; axis off;
title('Coral Cover')
ax2 = subplot(1,2,2);
colormap(ax2,customMap)
imagesc(ax2,c_airmulti_predict{1})
axis image; axis off;
title('Morphology')
