Scyllarus: C++ Hyperspectral Processing Library
Hyperspectral Image Processing Pipeline
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Pages
Classes | Enumerations | Functions | Variables
scyl Namespace Reference

Classes

class  component_input
 component_input is an object that facilitates reading in of files saved by component_output (or assembled in the correct way manually). More...
 
class  data_input
 data_input is an object that allows for manual data input into the pipeline (more 'managed' than using the pipelines 'set_I()' function. More...
 
class  hdr_input
 HDR_Input is an object that facilitates reading in of HDR files. More...
 
class  hsz_input
 HSZ_Input is an object that facilitates reading in of HSZ files. More...
 
class  input
 Input is a base class from which input classes for the pipeline are derived.. More...
 
class  tif_input
 TIF_Input is an object that facilitates reading in of GTIFF/TIFF files. More...
 
class  component_output
 component_output is a class that is used to output processed image components from a pipeline object as separate HDR files in a folder. More...
 
class  hdr_output
 HDR Output is a class that is used to output HDR (and associated .fla, .raw) files. More...
 
class  hsz_output
 HSZ Output is a class that is used to output HSZ files. More...
 
class  output
 Output is a pure virtual base class from which output classes for the pipeline are derived.. More...
 
class  pipeline
 Pipeline is an object that holds configuration properties and data items, and facilitates the processing of this data in a pipeline like fashion. More...
 
class  spectral_library
 Spectral Library is a class that loads and saves SLZ files and holds information from them for use with the Pipeline Object (or independently) More...
 
class  illuminant_options
 illuminant_options is used to give options for the Illuminant recovery option chosen. More...
 
class  material_options
 material_options is used to give options for the material recovery option chosen. More...
 
class  op_node
 op_node is used by band_maths() to perform mathematical operations. More...
 
class  logger
 Logger is used to record messages and other processing details for the scyl library functions and access them during program execution. More...
 

Enumerations

enum  INPUT_TYPE {
  HDR_INPUT, HSZ_INPUT, DATA_INPUT, COMPONENT_INPUT,
  TIF_INPUT, HUH
}
 
enum  OUTPUT_TYPE { HDR_OUTPUT, HSZ_OUTPUT, COMPONENT_OUTPUT, WHA }
 
enum  PIPELINE_STAGE {
  RECOVER_ILLUMINANT, RECOVER_DICHROMATIC, RECOVER_MATERIAL, PCA,
  SAM, SVM, UNMIX, FILTER
}
 
enum  ILLUMINANT_METHOD {
  ILLUMINANT_HRK, ILLUMINANT_FS, ILLUMINANT_GE, ILLUMINANT_GW,
  ILLUMINANT_SG, ILLUMINANT_WP
}
 
enum  DICHROMATIC_METHOD { DICHROMATIC_LS }
 
enum  MATERIAL_METHOD { MATERIAL_DA, MATERIAL_DAQ, MATERIAL_DASAM }
 
enum  CLASSIFICATION_METHOD { CLASSIFICATION_SAM, CLASSIFICATION_SVM, CLASSIFICATION_LSU }
 
enum  FILTER_METHOD { FILTER_SAVITZKY, FILTER_WIENER, FILTER_MVAVG }
 
enum  COMPONENTS : unsigned int {
  I = 1, ILLUMINANT = 2, K = 4, G = 8,
  S = 16, ELEMENTS = 32, ABUNDANCES = 64, INDEXES = 128,
  ENDABUNDANCES = 256, ENDINDEXES = 512, SPECTRAL_LIBRARY = 1024, ALL = 0xFFFF
}
 

Functions

void recover_dichromatic_parameters (const arma::fcube &I, const arma::fvec &illuminant, arma::fmat &k, arma::fmat &g, float &k_factor, arma::fcube &s, scyl::DICHROMATIC_METHOD method=scyl::DICHROMATIC_LS, int size=-1, int gray_threshold=2, int debug=0)
 recover_dichromatic_parameters() takes an input image I and illuminant and determines the factors g, k, K and S. More...
 
void recover_dichromatic_parameters_LS (const arma::fcube &I, const arma::fcube &illuminant, arma::fmat &k, arma::fmat &g, arma::fcube &st, int size=5, int gray_threshold=2, int debug=0)
 recover_dichromatic_parameters_LS() takes an input image I and illuminant and determines the factors g, k, and S. More...
 
template<class T >
void clip_lower (T &a, float f_lim, float f_set)
 Clip arma matrix values below a lower limit and set clipped values to an arbitrary value. More...
 
template<class T >
void clip_lower (T &a, float f_lim)
 Clip arma matrix values below a lower limit and set clipped values to that limit. More...
 
arma::fmat conv (const arma::fmat &a, arma::fvec kernel)
 Convolve a matrix with a 1D filter kernel (based on matlab implementation) More...
 
arma::fmat conv2 (const arma::fmat &a, arma::fmat kernel)
 Convolve a matrix with a 2D filter kernel (based on matlab implementation) More...
 
arma::fmat wiener2 (const arma::fmat &a, arma::fvec nhood)
 Apply a de-blurring wiener filter to an arma matrix (based on matlab implementation) More...
 
arma::fmat median (const arma::fmat &a, arma::u32 iSize)
 Apply a median filter to an arma matrix (based on matlab implementation) More...
 
void savitzky_golay (arma::fcube &in, int window=5)
 
void moving_average (arma::fcube &in, int window=3)
 
arma::umat select_smooth_patches (const arma::fcube &I, arma::umat &patches, arma::fmat &light_estimation_mask, arma::fmat &contrast_mask, arma::fmat &highlight_mask, int patch_size=20, int num_patches=50, int debug=0)
 select_smooth_patches() returns a list of the first num_patches smooth patches it finds while inspecting input image I. More...
 
int patch_dichromatic_decompose (const arma::fcube &patch, const arma::fmat &illuminant, arma::fmat &g_out, arma::fmat &k_out, arma::fvec &s_out, float &dichromatic_error, float &smooth_error, float alpha=50, int debug=0)
 patch_dichromatic_decompose() looks at a image patch and determines if the patch is suitable for use in calculating the greater image illuminant. More...
 
arma::fvec recover_global_illuminant (const arma::fcube &I, scyl::ILLUMINANT_METHOD method=scyl::ILLUMINANT_HRK, scyl::illuminant_options op=scyl::illuminant_options())
 recover_global_illuminant() takes an input image I and determines the illuminant. More...
 
arma::fvec recover_illuminant_HRK (const arma::fcube &I, arma::umat &patches, float alpha, int patch_size, bool filtered_already=false, bool fast_50=true, int debug=0)
 recover_illuminant_HRK() takes an input image I and determines the illuminant. More...
 
arma::fvec recover_illuminant_FS (const arma::fcube &I, arma::umat &patches, int patch_size=20, bool filtered_already=false, int debug=0)
 recover_illuminant_FS() takes an input image I and determines the illuminant. More...
 
arma::fvec recover_illuminant_GE (const arma::fcube &I, int order=1, int debug=0)
 recover_illuminant_GE() takes an input image I and determines the illuminant. More...
 
arma::fvec recover_illuminant_GW (const arma::fcube &I, int debug=0)
 recover_illuminant_GW() takes an input image I and determines the illuminant. More...
 
arma::fvec recover_illuminant_SG (const arma::fcube &I, int order=2, int debug=0)
 recover_illuminant_SG() takes an input image I and determines the illuminant. More...
 
arma::fvec recover_illuminant_WP (const arma::fcube &I, int debug=0)
 recover_illuminant_WP() takes an input image I and determines the illuminant. More...
 
arma::fmat recover_materials (const arma::fcube &s, const arma::fvec &wavelengths, arma::fcube &abundances, arma::ucube &indices, scyl::MATERIAL_METHOD method=scyl::MATERIAL_DA, scyl::material_options op=scyl::material_options())
 recover_materials() takes an input reflectance cube and returns a list of determined material clusters. More...
 
arma::fmat recover_materials_DA (const arma::fcube &s, arma::fcube &material_abundancy, int max_clusters=20, float temperature_max=0.02, float temperature_min=0.00025, float cooling_rate=0.8, float split_threshold=-1, int debug=0)
 recover_materials_DA() takes an input reflectance cube and returns a list of determined material clusters. More...
 
void generate_material_elements (const arma::fcube &s, const arma::fmat &elements, const arma::fcube &material_abundancy, arma::fcube &element_abundances, arma::ucube &element_abundance_indexes, int num_materials=5, int debug=0)
 generate_material_elements() uses the reflectance cube and abundancy cube from recover_materials_DA to calculate the element abundances and abundance indices. More...
 
void generate_nurbs (const arma::fcube &I, const arma::fvec &wavelengths, arma::fvec &knot_vector, arma::fcube &control_points, arma::fvec &wavelength_cp, int degree=2, int knot_threshold=-1, float alpha=0.1f, int iterations=10, int debug=0)
 generate_nurbs() encodes an array of discrete data into a NURBS curve representation. More...
 
void nurbs_global_interpolation (const arma::fcube &I, const arma::fvec &wavelengths, arma::fvec &knot_vector, arma::fvec &parametric_points, arma::fcube &control_points, arma::fvec &wavelength_cp, int degree=2, int debug=0)
 nurbs_global_interpolation() gives an equivalent NURBS encoding for data based on the given image and wavelengths. More...
 
int nurbs_find_span (const arma::fvec &knot_vector, int degree, float point)
 nurbs_find_span finds the span of 'point' in the vector 'knot_vector'. More...
 
arma::fvec nurbs_get_basis (const arma::fvec &knot_vector, int span, int degree, float point)
 nurbs_get_basis() finds the basis of the sub-vector of 'knot_vector' for the given point and span. More...
 
arma::fvec nurbs_minimal_knots (const arma::fcube &I, const arma::fvec &wavelengths, const arma::fvec &knot_vector, const arma::fcube &control_points, const arma::fvec &wavelength_cp, arma::fvec &min_knot_vector, arma::fcube &min_control_points, arma::fvec &min_wavelength_cp, float alpha, int degree, int target_knot_num, int debug)
 nurbs_minimal_knots() takes a set of NURBS curves, and returns a minimised NURBS representation of the original data, More...
 
int nurbs_remove_curve_knot (const arma::fvec &knot_vector, const arma::fcube &control_points, const arma::fvec &wavelength_cp, arma::fvec &new_knot_vector, arma::fcube &new_control_points, arma::fvec &new_wavelength_cp, float tolerance, int degree, int index, int multiplicity, int num_removes)
 nurbs_remove_curve_knot() takes a NURBS encoded data set and a specified knot index, and attempts to remove the knot from the set. More...
 
arma::fvec nurbs_find_parapoints (const arma::fvec &wavelengths, const arma::fvec &knot_vector, const arma::fvec &control_points, float error_tolerance, int degree)
 nurbs_find_parapoints() calculates the parametric points for the given NURBS set. More...
 
arma::fvec nurbs_evaluate_univariate (const arma::fvec &knot_vector, const arma::fvec &control_points, const arma::fvec &parameters, int degree)
 nurbs_evaluate_univariate() evaluates a B-Spline curve according to the given parameters. More...
 
arma::fcube nurbs_reconstruct_curve (const arma::fvec &knot_vector, const arma::fvec &t, const arma::fcube &control_points, const arma::fvec &wavelength_cp, arma::fvec &wavelengths, int degree)
 nurbs_reconstruct_curve() takes a NURBS representation of a data set and a set of points to evaluate it at and reconstructs the original data. More...
 
arma::fcube nurbs_reconstruct_image (const arma::fvec &knot_vector, const arma::fcube &control_points, const arma::fvec &wavelength_cp, const arma::fvec &wavelengths, arma::fvec &wavelengths_new, int degree)
 nurbs_reconstruct_image() takes a NURBS encoded image and returns the reconstructed image cube. More...
 
arma::fmat mix_spectra_vec (const arma::fvec &abundances, const arma::fvec &end_members, bool normalise=true)
 mix_spectra_vec() seems to merely replicate abundances across endmembers many bands. More...
 
arma::fvec wavelength_subset (const arma::fvec &set, const arma::fvec &bounds)
 wavelength_subset() removes elements from set that are outside of the range of elements in bounds. More...
 
void pca (const arma::fcube &I, arma::fcube &principal_components, arma::fvec &eigenvalues, int num_components=10, float threshold=0.001)
 
void unmix_by_spectra (const arma::fcube &I, const arma::fvec &wavelengths, const arma::fmat &spectra, const arma::fvec &spectra_waves, arma::fcube &abundances, arma::ucube &indices, int num_materials=5)
 unmix_by_spectra() takes an input data structure and a list of spectra and calculates abundances for each pixel in the data according to the given spectra. More...
 
void sam (const arma::fcube &I, const arma::fvec &wavelengths, const arma::fmat &spectra, const arma::fvec &spectra_waves, arma::fcube &angles, arma::ucube &indices, int num_materials=5)
 sam() takes an input data structure and a list of spectra and calculates angles for each pixel to each spectra given, returning the num_materials closest for each pixel. More...
 
void svm (const arma::fcube &I, const arma::fvec &wavelengths, const arma::fmat &spectra, const arma::fvec &spectra_waves, const arma::uvec &spectra_labels, arma::umat &indices)
 svm() takes an input data structure and a list of spectra and classifies each pixel as one of the given spectra, returning the index map. More...
 
arma::fcube spectral_derivative (const arma::fcube &I, const arma::fvec &wavelengths)
 
arma::fmat band_maths (std::vector< const arma::fcube * > images, std::string expression)
 
arma::fmat compute_greyscale (const arma::fcube &source, const arma::fvec &wave)
 compute_grayscale returns a single channel grey scale image from a given hyperspectral image cube, approximating the scene's visual appearance to a human eye. More...
 
float linear_interp (const float wavelength, const int band)
 linear_interp looks up a wavelength for a band (R, G or B) in the scyl::CAMERA_SENSITIVITY table, returning the corresponding float value from the table, or otherwise interpolating one based on the wavelengths position between elements. More...
 
arma::fcube compute_rgb (const arma::fcube &I, const arma::fvec &illuminant, const arma::fvec &wave)
 compute_rgb creates a pseudo-colour representation of a given hyperspectral image cube. More...
 
arma::fcube reconstruct_I (const arma::fvec &illuminant, const arma::fmat &k, const arma::fmat &g, const arma::fcube &s, float k_factor)
 reconstruct_I() creates a synthetic version of I based on it's constituents illuminant, k, g, s and big_k. More...
 
arma::fcube reconstruct_s (const arma::fvec &illuminant, const arma::fmat &k, const arma::fmat &g, const arma::fcube &I, float k_factor)
 reconstruct_s() creates a synthetic version of s based on it's constituents illuminant, k, g, I and k_factor. More...
 
arma::fcube reconstruct_s_elements (const arma::fmat &elements, const arma::fcube &abundances, const arma::ucube &indexes)
 
arma::fvec reconstruct_illuminant (const arma::fmat &k, const arma::fmat &g, const arma::fcube &I, const arma::fcube &s, float k_factor)
 reconstruct_illuminant() creates a synthetic version of the illuminant based on given image components I, s, k g and k_factor. More...
 
arma::fcube colour_jet (const arma::fmat &in)
 
arma::fcube resample_image (const arma::fcube &in_image, const arma::fvec &in_waves, const arma::fvec &in_waves_new)
 resample_image() is used to change the number of bands an image matrix has using NURBS. More...
 
cv::Mat arma_to_ocv (const arma::fmat &in)
 arma_to_ocv() converts a Armadillo mat of floats to a OpenCV Mat of floats. More...
 
arma::fmat ocv_to_arma (const cv::Mat &in)
 ocv_to_arma() converts a OpenCV Mat of floats to a Armadillo mat of floats. More...
 
template<class T >
void version (T &in)
 Get the version number of Scyllarus (Casted from Float (Major.MinorSub) More...
 
template<class T >
void version (arma::Mat< T > &in)
 Get the version number of Scyllarus (Casted from Float (Major.MinorSub) More...
 
std::string version_string ()
 Get the version number formatted as a string [Major.Minor.Sub(release)]. More...
 
void print_version ()
 Print Scyllarus information/copyright banner with version information. More...
 
template<class T >
void nan_to_zero (T &in)
 nan_to_zero removes NaN values from a Armadillo data structure, replacing them with 0. More...
 
template<class T >
void zero_to_one (T &in)
 zero_to_one replaces 0's in an Armadillo data structure with 1's. More...
 
template<class T >
bool arma_cmp (const T &a, const T &b, const bool verbose=false, float eps=arma::Datum< float >::eps)
 arma_cmp compares two Armadillo objects for equality +- epsilon. More...
 
template<class T >
int arma_cmp_num (const T &a, const T &b, float eps)
 
template<class T >
float arma_cmp_sig_figs (const T &a, const T &b, int num_sig_figs=6)
 arma_cmp compares two Armadillo objects for equality. Instead of using a tolerance, test number of significant figures. More...
 
arma::fvec non_negative_least_squares_2 (arma::fmat a, arma::fvec b)
 non_negative_least_squares_2() gives a non negative solution for x to the problem of Ax = b. More...
 
void remove_bands (arma::fcube &in, arma::fvec &waves, arma::uvec remove)
 remove_bands() removes image bands from the given image cube. More...
 
void scale_cube (arma::fcube &in, int height, int width, int cv_interpolate_flag=cv::INTER_LINEAR)
 scale_cube() resamples a cube spatially using OpenCV resize methods. More...
 
float angle (const arma::fvec &a, const arma::fvec &b)
 angle returns the angle between two vectors in degrees. More...
 
void image_gradient (const arma::fmat &in, arma::fmat &gx, arma::fmat &gy)
 Calculates the gradients in the x and y directions for a given image. More...
 
arma::fmat arma_slice_mean (const arma::fcube &in)
 arma_slice_mean() takes an image cube and returns a mat which is the mean image of all the slices. More...
 
void register_image_bands (const arma::fcube &in, arma::fvec &coloff, arma::fvec &rowoff, int debug=0)
 register_image_bands() takes a image cube and attempts to automatically find the alignment of the bands in the image. More...
 
void align_image_bands (arma::fcube &in, arma::fvec &coloff, arma::fvec &rowoff, int debug=0)
 align_image_bands() aligns bands in an image according to the input row and column offsets. More...
 
void sort_bands (arma::fcube &I, arma::fvec &wavelengths)
 sort_bands() manipulates a Image cube and wavelengths vector. More...
 
arma::fmat resample_endmembers (const arma::fmat &in_endmembers, const arma::fvec &in_waves, const arma::fvec &in_waves_new)
 resample_endmembers() is used to change the number of bands an endmember matrix has using NURBS. More...
 
template<class T >
int arma_cmp_num (const T &a, const T &b, float eps=arma::Datum< float >::eps, bool verbose=false)
 arma_cmp_num compares two Armadillo objects for equality +- epsilon and returns the number of elements exceeding the condition. More...
 

Variables

const int PHOTOPIC_START = 390
 
const int PHOTOPIC_LENGTH = 441
 Starting wavelength of PHOTOPIC_VALUES. More...
 
const int PHOTOPIC_END = 830
 Length of PHOTOPIC_VALUES. More...
 
const float PHOTOPIC_VALUES [441]
 Ending wavelength of PHOTOPIC_VALUES. More...
 
const int CAMERA_SENSITIVITY_LENGTH = 35
 
const float CAMERA_SENSITIVITY_WAVELENGTHS [35]
 
const float CAMERA_SENSITIVITY [35][3]
 
const float COLOUR_JET [100][3]
 
const int VERSION_MAJOR = 1
 
const int VERSION_MINOR = 0
 
const int VERSION_SUB = 3
 
const float VERSION_FLOAT = 1.03f
 

Enumeration Type Documentation

Method used by process_classification()

Enumerator
CLASSIFICATION_SAM 

SAM - Spectral Angle Mapping.

CLASSIFICATION_SVM 

SVM - Support Vector Machines.

CLASSIFICATION_LSU 

LSU - Linear Spectral Unmixing.

enum scyl::COMPONENTS : unsigned int

Flags for loading and storing components with component_input and component_output classes.

Enumerator
I 

I.

ILLUMINANT 

ILLUMINANT.

K 

K.

G 

G.

S 

S.

ELEMENTS 

ELEMENTS.

ABUNDANCES 

ABUNDANCES.

INDEXES 

INDEXES.

ENDABUNDANCES 

ENDABUNDANCES.

ENDINDEXES 

ENDINDEXES.

SPECTRAL_LIBRARY 

SPECTRAL_LIBRARY.

ALL 

Method used by recover_dichromatic_parameters()

Enumerator
DICHROMATIC_LS 

LS - Use least squares method of dichromatic recovery.

Method used by filter()

Enumerator
FILTER_SAVITZKY 

SAVITZKY - Savitzky Golay filter (spectral)

FILTER_WIENER 

WIENER - Wiener filter (spatial)

FILTER_MVAVG 

MVAVG - Moving Average filter (spectral)

Method used by recover_global_illuminant()

Enumerator
ILLUMINANT_HRK 

HRK - Use Huynh Robles-Kelly method for illuminant recovery.

ILLUMINANT_FS 

FS - Use Finlayson method for illuminant recovery.

ILLUMINANT_GE 

GE - Use Grey Edge method for illuminant recovery.

ILLUMINANT_GW 

GW - Use Grey World method for illuminant recovery.

ILLUMINANT_SG 

SG - Use Shade of Grey method for illuminant recovery.

ILLUMINANT_WP 

WP - Use White Patch method for illuminant recovery.

Returned by input type() method

Enumerator
HDR_INPUT 

HDR - HDR/FLA input type.

HSZ_INPUT 

HSZ - HSZ input type.

DATA_INPUT 

DATA_INPUT - Data Input Type (For cameras/custom)

COMPONENT_INPUT 

COMPONENT_INPUT - Component Input type (Files)

TIF_INPUT 

TIF - Tiff / GTiff input type.

HUH 

HUH - Error type.

Method used by recover_materials()

Enumerator
MATERIAL_DA 

DA - Deterministic Annealing Material Recovery.

MATERIAL_DAQ 

DAQ - Quick Deterministic Annealing Material Recovery (DA method on a subsampled image, then spectral unmixing to classify each pixel based on the recovered material list)

MATERIAL_DASAM 

DASAM - SAM + Deterministic Annealing Material Recovery (DA method on a subsampled image, then Spectral Angle Mapping to classify each pixel based on the recovered material list)

Returned by output type() method (Not implemented yet) TODO: This.

Enumerator
HDR_OUTPUT 

HDR - HDR/FLA output type.

HSZ_OUTPUT 

HSZ - HSZ output type.

COMPONENT_OUTPUT 

COMPONENT_OUTPUT - Component output type (Files)

WHA 

WHA - Error type.

Stage the pipeline is currently processing reported by logger

Enumerator
RECOVER_ILLUMINANT 

RECOVER_ILLUMINANT stage.

RECOVER_DICHROMATIC 

RECOVER_DICHROMATIC stage.

RECOVER_MATERIAL 

RECOVER_MATERIAL stage.

PCA 

PCA stage.

SAM 

SAM stage.

SVM 

SVM stage.

UNMIX 

UNMIX stage.

FILTER 

FILTER stage.

Function Documentation

void scyl::align_image_bands ( arma::fcube &  in,
arma::fvec &  rowoff,
arma::fvec &  coloff,
int  debug = 0 
)

align_image_bands() aligns bands in an image according to the input row and column offsets.

align_image_bands() will align each band in the image according to the given offsets. the function register_image_bands() can be used to generate the rowoff and coloff vectors. The resultant image will be smaller than the original; areas where there is not data present for all bands will be cropped off.

Parameters
in- The input image cube to be aligned (image_height x image_width x image_bands)
rowoff- Row offsets for each band (image_bands)
coloff- Column offsets for each band (image_bands)
debug- Control the amount of output that is generated (0 - 5, default 0)
float scyl::angle ( const arma::fvec &  a,
const arma::fvec &  b 
)

angle returns the angle between two vectors in degrees.

Parameters
athe first vector
bthe other vector
Returns
angle in degrees.
template<class T >
bool scyl::arma_cmp ( const T &  a,
const T &  b,
const bool  verbose = false,
float  eps = arma::Datum< float >::eps 
)

arma_cmp compares two Armadillo objects for equality +- epsilon.

Parameters
aobject a
bobject b
epsThe tolerance difference between the values.
Returns
boolean comparison result
template<class T >
int scyl::arma_cmp_num ( const T &  a,
const T &  b,
float  eps 
)
template<class T >
int scyl::arma_cmp_num ( const T &  a,
const T &  b,
float  eps = arma::Datum<float>::eps,
bool  verbose = false 
)

arma_cmp_num compares two Armadillo objects for equality +- epsilon and returns the number of elements exceeding the condition.

Parameters
aobject a
bobject b
epsThe tolerance difference between the values.
Returns
boolean comparison result
template<class T >
float scyl::arma_cmp_sig_figs ( const T &  a,
const T &  b,
int  num_sig_figs = 6 
)

arma_cmp compares two Armadillo objects for equality. Instead of using a tolerance, test number of significant figures.

Parameters
aobject a
bobject b
num_sig_figsThe number of significant figures that must be the same.
Returns
float percentage of pixels that failed.
arma::fmat scyl::arma_slice_mean ( const arma::fcube &  in)

arma_slice_mean() takes an image cube and returns a mat which is the mean image of all the slices.

Parameters
ininput cube to be mean to.
Returns
a mat which is the mean of all the slices of the input cube.
cv::Mat scyl::arma_to_ocv ( const arma::fmat &  in)

arma_to_ocv() converts a Armadillo mat of floats to a OpenCV Mat of floats.

Parameters
inArmadillo mat to be converted
Returns
OpenCV Mat (Not transposed)
Examples:
cli_pipeline.cpp.
arma::fmat scyl::band_maths ( std::vector< const arma::fcube * >  images,
std::string  expression 
)

band_maths() allows users to perform mathematical operations on images by supplying string formula specifying the operations they would like to perform

Images are specified in the following format: Ix[y] where x is the image number (starting from 1) and y is the band number (starting from 0).

Example Operations:

(-I1[3] + 7)*4 + 40
I2[3]^2 + I3[4]^2 > 0.5

Available operators:

 +
 - (Situationally used for negation)
 /
 *
 ^
 < (results in binary map)
 > (results in binary map)

Unknown expressions and other 'junk' characters will likely result in errors.

Parameters
images- A list of pointers of images that will be used in the operation, with the first in the list being I1.
expression- The string expression to evaluate using the image list given
Returns
the result as an image
template<class T >
void scyl::clip_lower ( T &  a,
float  f_lim,
float  f_set 
)
inline

Clip arma matrix values below a lower limit and set clipped values to an arbitrary value.

Parameters
aarma vector/matrix/cube
f_limlower limit
f_setvalue to set clipped values to
template<class T >
void scyl::clip_lower ( T &  a,
float  f_lim 
)
inline

Clip arma matrix values below a lower limit and set clipped values to that limit.

Parameters
aarma vector/matrix/cube
f_limlower limit
arma::fcube scyl::colour_jet ( const arma::fmat &  in)

colour_jet() applies a cool to warm colourisation function to a given single channel image.

The image returned by colour_jet() is three channels representing (R, G, B). The image given can be any scale, the returned image will be [0, 1]. There is some rudimentary outlier removal (To try and prevent 'blank' results being produced. if some value is much much larger than the average then rather than scaling to the max the image will be scaled to two times the mean instead.

Parameters
in- Image to be colourised with a jet
Returns
- The 'jet' RGB image.
arma::fmat scyl::compute_greyscale ( const arma::fcube &  I,
const arma::fvec &  wavelengths 
)

compute_grayscale returns a single channel grey scale image from a given hyperspectral image cube, approximating the scene's visual appearance to a human eye.

Parameters
I- The input image
wavelengths- The corresponding wavelengths
Returns
The resultant grey image representation.
arma::fcube scyl::compute_rgb ( const arma::fcube &  I,
const arma::fvec &  illuminant,
const arma::fvec &  wavelengths 
)

compute_rgb creates a pseudo-colour representation of a given hyperspectral image cube.

If the bands correspending to 'R', 'G', and 'B' can be found within the image, they will be used, or otherwise interpolated. For images that do not cover the entire visible spectrum, an alternate approximation is used (choosing a band from near each end and the middle one to represent the three channels).

Parameters
I- The image to be represented
illuminant- Illumination power function for the image
wavelengths- The corresponding wavelength vector
Returns
arma::fmat scyl::conv ( const arma::fmat &  a,
arma::fvec  kernel 
)

Convolve a matrix with a 1D filter kernel (based on matlab implementation)

Parameters
aarma matrix
kernelarma vec with data values (dimension must be an odd number greater than 1)
arma::fmat scyl::conv2 ( const arma::fmat &  a,
arma::fmat  kernel 
)

Convolve a matrix with a 2D filter kernel (based on matlab implementation)

Parameters
aarma matrix
kernelsquare matrix with data values (dimension must be an odd number greater than 1)
void scyl::generate_material_elements ( const arma::fcube &  s,
const arma::fmat &  elements,
const arma::fcube &  material_abundancy,
arma::fcube &  element_abundances,
arma::ucube &  element_abundance_indexes,
int  num_materials = 5,
int  debug = 0 
)

generate_material_elements() uses the reflectance cube and abundancy cube from recover_materials_DA to calculate the element abundances and abundance indices.

generate_material_elements() is required to be called to process the output of recover_materials_DA before saving a file as HSZ format (as it generates the ratios of materials contained in each pixel).

Parameters
s- Reflectance cube (height x width x bands)
elements- Mat containing elements (Materials) (num_materials x bands)
material_abundancy- Cube containing material abundancy (from recover materials) (height x width x bands)
element_abundances- Element Abundance Cube (pass an empty, uninitialised mat) (height x width x 5)
element_abundance_indexes- Element Abundance Indexes Cube (pass an empty, uninitialised mat) (height x width x 5)
num_materials- Number of materials to use (per sample) in abundance mat. (Compression level) ( > 1, default 5)
debug- Control the amount of output that is generated (0 - 5, default 0)
Exceptions
ExceptionError
void scyl::generate_nurbs ( const arma::fcube &  I,
const arma::fvec &  wavelengths,
arma::fvec &  knot_vector,
arma::fcube &  control_points,
arma::fvec &  wavelength_cp,
int  degree = 2,
int  knot_threshold = -1,
float  alpha = 0.1f,
int  iterations = 10,
int  debug = 0 
)

generate_nurbs() encodes an array of discrete data into a NURBS curve representation.

generate_nurbs() takes an array 'I' and a set of indices 'wavelengths' and encodes the data using NURBS. For the purposes of encoding, each 'pixel' in I is considered to be one curve, that is, NURBS curves are calculated along the 3rd dimension of the data (Thus resulting in image_width x image_height NURBS curves). All of the curves will share the same set of control points, which are returned in the two arrays 'wavelength_cp' and 'knot_vector'. The NURBS control points will be given in 'control_points'. This function makes calls to nurbs_global_interpolation(), nurbs_minimal_knots() and nurbs_reconstruct_curve().

To encode a single curve, pass I as 1x1xN. Each row+col in the image is encoded as a NURB curve.

Parameters
I- Input 'Image' to be encoded. (image_width x image_height x image_bands)
wavelengths- Wavelength indices of each band in the image (image bands)
knot_vector- Resultant knot vector for the NURBS encoding (num_cp + degree + 1) (pass an empty, uninitialised vec)
control_points- Resultant control points for the NURBS encoding (image_width x image_height x num_cp) (pass an empty, uninitialised cube)
wavelength_cp- Resultant knot vector for the NURBS encoding (num_cp) (pass an empty, uninitialised vec)
degree- Polynomial degree for NURBS encoding (Default 2)
knot_threshold- Target number of knots for NURBS encoding. (Default image_bands - 2)
alpha- Alpha value used during knot removal (Default 0.1)
iterations- Max number of iterations when calculating NURBS encoding (Default 10)
debug- Control the amount of output that is generated (0 - 5, default 0)
Exceptions
ExceptionError
void scyl::image_gradient ( const arma::fmat &  in,
arma::fmat &  gx,
arma::fmat &  gy 
)

Calculates the gradients in the x and y directions for a given image.

Parameters
ininput image
gxresultant gradient in x direction.
gyresultant gradient in y direction.
float scyl::linear_interp ( const float  wavelength,
const int  band 
)

linear_interp looks up a wavelength for a band (R, G or B) in the scyl::CAMERA_SENSITIVITY table, returning the corresponding float value from the table, or otherwise interpolating one based on the wavelengths position between elements.

Parameters
wavelength- Wavelength to lookup
band- Colour band
Returns
Interpolated result or otherwise a max/min.
arma::fmat scyl::median ( const arma::fmat &  a,
arma::u32  iSize 
)

Apply a median filter to an arma matrix (based on matlab implementation)

Parameters
aarma matrix
iSizesize of the kernel (must be an odd number greater than 1)
arma::fmat scyl::mix_spectra_vec ( const arma::fvec &  abundances,
const arma::fvec &  end_members,
bool  normalise = true 
)

mix_spectra_vec() seems to merely replicate abundances across endmembers many bands.

Parameters
abundances- Input vectors of abundances
end_members- Input vectors of end_members
normalise- Normalise the result?
Returns
The resulting 'mixed' spectra.
void scyl::moving_average ( arma::fcube &  in,
int  window = 3 
)

moving_average applies the moving average filter in the spectral domain. Available window sizes are 3+.

The function works in-place and uses at most (window / 2) of additional memory to perform the filter.

Parameters
in- The input image (rows x cols x bands)
window- The window size (3 or greater and odd)
template<class T >
void scyl::nan_to_zero ( T &  in)

nan_to_zero removes NaN values from a Armadillo data structure, replacing them with 0.

Parameters
inthe data structure to be manipulated.
Examples:
cli_pipeline.cpp.
arma::fvec scyl::non_negative_least_squares_2 ( arma::fmat  a,
arma::fvec  b 
)

non_negative_least_squares_2() gives a non negative solution for x to the problem of Ax = b.

non_negative_least_squares_2() follows the implementation seen in Matlab, based on Lawson and Hanson, "Solving Least Squares Problems", Prentice-Hall, 1974. The algorithm makes use of the Armadillo 'solve()' function to iteratively find a non negative solution.

Parameters
a- A, the input matrix.
b- b, the input vector
Returns
x, the non negative least squares solution. [-1 -1] indicates failure.
arma::fvec scyl::nurbs_evaluate_univariate ( const arma::fvec &  knot_vector,
const arma::fvec &  control_points,
const arma::fvec &  parameters,
int  degree 
)

nurbs_evaluate_univariate() evaluates a B-Spline curve according to the given parameters.

nurbs_evaluate_univariate() uses a set of parameters (give by nurbs_find_parapoints()) and returns the points along the curve evaluated at those points. This function makes calls to nurbs_find_span() and nurbs_get_basis().

Parameters
knot_vector- NURBS knot vector (cp_num + degree + 1)
control_points- NURBS control points (image_width x image_height x cp_num)
parameters- Points to evaluate the curve at
degree- The degree of the knot vector
Returns
The list of points as evaluated along the given curve.
arma::fvec scyl::nurbs_find_parapoints ( const arma::fvec &  wavelengths,
const arma::fvec &  knot_vector,
const arma::fvec &  control_points,
float  error_tolerance,
int  degree 
)

nurbs_find_parapoints() calculates the parametric points for the given NURBS set.

nurbs_find_parapoints() uses a binary search method to calculate an independent parametric evaluation of the given data, to within the given tolerance. If the function does not converge after a number of iterations (50) the error tolerance will be increased and the function will continue, repeating if necessary until convergence. This function makes calls to nurbs_evaluate_univariate().

Parameters
wavelengths- Original wavelengths (image_bands)
control_points- NURBS control points (image_width x image_height x cp_num)
knot_vector- NURBS knot vector (cp_num)
error_tolerance- Maximum sum error (as a distance from the given wavelengths)
degree- The degree of the knot vector
Returns
The vector of calculated points.
int scyl::nurbs_find_span ( const arma::fvec &  knot_vector,
int  degree,
float  point 
)

nurbs_find_span finds the span of 'point' in the vector 'knot_vector'.

The 'span' is considered to be the position (index) the point 'point' would reside at if it was placed sequentially into the vector 'knot_vector'.

Parameters
knot_vector- The vector to search
degree- The degree of the knot vector
point- The point to find the span for
Returns
The span of 'point' in 'knot_vector'
arma::fvec scyl::nurbs_get_basis ( const arma::fvec &  knot_vector,
int  span,
int  degree,
float  point 
)

nurbs_get_basis() finds the basis of the sub-vector of 'knot_vector' for the given point and span.

Parameters
knot_vector- The vector to search
span- The knot span of point (given by function nurbs_find_span())
degree- The degree of the knot vector
point- The point to find the basis around
Returns
The basis vector (size degree + 1)
void scyl::nurbs_global_interpolation ( const arma::fcube &  I,
const arma::fvec &  wavelengths,
arma::fvec &  knot_vector,
arma::fvec &  parametric_points,
arma::fcube &  control_points,
arma::fvec &  wavelength_cp,
int  degree = 2,
int  debug = 0 
)

nurbs_global_interpolation() gives an equivalent NURBS encoding for data based on the given image and wavelengths.

The outputs from this function give a NURBS encoded version of the input image. The encoding is not minimised at all by this function. This function makes calls to nurbs_find_span() and nurbs_get_basis().

Parameters
I- Input 'Image' to be encoded. (image_width x image_height x image_bands)
wavelengths- Wavelength indices of each band in the image (image_bands)
knot_vector- Resultant knot vector for the NURBS encoding (image_bands + degree + 1) (pass an empty, uninitialised vec)
parametric_points- Resultant parametric points for the NURBS encoding (image_bands) (pass an empty, uninitialised vec)
control_points- Resultant control points for the NURBS encoding (image_width x image_height x image_bands) (pass an empty, uninitialised cube)
wavelength_cp- Resultant knot vector for the NURBS encoding (image_bands) (pass an empty, uninitialised vec)
degree- Polynomial degree for NURBS encoding (Default 2)
debug- Control the amount of output that is generated (0 - 5, default 0)
arma::fvec scyl::nurbs_minimal_knots ( const arma::fcube &  I,
const arma::fvec &  wavelengths,
const arma::fvec &  knot_vector,
const arma::fcube &  control_points,
const arma::fvec &  wavelength_cp,
arma::fvec &  min_knot_vector,
arma::fcube &  min_control_points,
arma::fvec &  min_wavelength_cp,
float  alpha,
int  degree,
int  target_knot_num,
int  debug 
)

nurbs_minimal_knots() takes a set of NURBS curves, and returns a minimised NURBS representation of the original data,

nurbs_minimal_knots() attempts to minimise the number of control points in the output data by performing iterative checks on candidate knots in the data for removal. If it finds a satisfactory action to perform, it will remove a knot and continue until no more knots can be removed. The function makes calls to nurbs_remove_curve_knot(), nurbs_find_parapoints() and nurbs_reconstruct_curve().

Parameters
I- Input 'Image' to be encoded. (image_width x image_height x image_bands)
wavelengths- Wavelength indices of each band in the image (image_bands)
knot_vector- NURBS knot vector (cp_num + degree + 1)
control_points- NURBS control points (image_width x image_height x cp_num)
wavelength_cp- NURBS wavelength control points (cp_num)
min_knot_vector- Resultant knot vector after removal (pass an empty, uninitialised vec)
min_control_points- Resultant control points after removal (pass an empty, uninitialised cube)
min_wavelength_cp- Resultant wavelength control points after removal (pass an empty, uninitialised vec)
alpha- Alpha value used during knot removal
degree- Polynomial degree for NURBS encoding
target_knot_num- Target number of knots for NURBS encoding.
debug- Control the amount of output that is generated
Returns
A vector of some significance.
arma::fcube scyl::nurbs_reconstruct_curve ( const arma::fvec &  knot_vector,
const arma::fvec &  t,
const arma::fcube &  control_points,
const arma::fvec &  wavelength_cp,
arma::fvec &  wavelengths,
int  degree 
)

nurbs_reconstruct_curve() takes a NURBS representation of a data set and a set of points to evaluate it at and reconstructs the original data.

nurbs_reconstruct_curve() evaluates the NURBS curves given in 'knot_vector' and 'control_points' at each point given in 't' and returns a reconstructed cube accordingly. This function makes calls to nurbs_find_span() and nurbs_get_basis().

Parameters
knot_vector- NURBS knot vector (cp_num + degree + 1)
t- Points to evaluate the curve at
control_points- NURBS control points (image_width x image_height x cp_num)
wavelength_cp- NURBS wavelength control points (cp_num)
wavelengths- The resultant wavelengths of the output cube (ideally identical to t)
degree- The degree of the knot vector
Returns
The reconstructed data cube
arma::fcube scyl::nurbs_reconstruct_image ( const arma::fvec &  knot_vector,
const arma::fcube &  control_points,
const arma::fvec &  wavelength_cp,
const arma::fvec &  wavelengths,
arma::fvec &  wavelengths_new,
int  degree 
)

nurbs_reconstruct_image() takes a NURBS encoded image and returns the reconstructed image cube.

nurbs_reconstruct_image() evaluates the given NURBS encoded image (knot_vector, control_points, wavelength_cp) at each wavelength given in the 'wavelengths' vector. It will return an Image cube with as many bands as there are elements in 'wavelengths', along with 'wavelengths_new', which is ideally identical to 'wavelengths' but may differ a small amount. This function makes calls to nurbs_find_parapoints() and nurbs_reconstruct_curve().

Parameters
knot_vector- NURBS knot vector (cp_num + degree + 1)
control_points- NURBS control points (image_width x image_height x cp_num)
wavelength_cp- NURBS wavelength control points (cp_num)
wavelengths- Wavelengths to evaluate the curve at
wavelengths_new- The resultant wavelengths of the output cube (ideally identical to t)
degree- The degree of the knot vector
Returns
The reconstructed image cube
int scyl::nurbs_remove_curve_knot ( const arma::fvec &  knot_vector,
const arma::fcube &  control_points,
const arma::fvec &  wavelength_cp,
arma::fvec &  new_knot_vector,
arma::fcube &  new_control_points,
arma::fvec &  new_wavelength_cp,
float  tolerance,
int  degree,
int  index,
int  multiplicity,
int  num_removes 
)

nurbs_remove_curve_knot() takes a NURBS encoded data set and a specified knot index, and attempts to remove the knot from the set.

nurbs_remove_curve_knot() first calculates if the knot is removable from the set before removing it if possible. If it cannot make a change, it will return the original data.

Parameters
knot_vector- NURBS knot vector (cp_num + degree + 1)
control_points- NURBS control points (image_width x image_height x cp_num)
wavelength_cp- NURBS wavelength control points (cp_num)
new_knot_vector- New knot vector (with knot(s) removed) (pass an empty, uninitialised vec)
new_control_points- New control points (with point(s) removed) (pass an empty, uninitialised cube)
new_wavelength_cp- New wavelength control points (with point(s) removed) (pass an empty, uninitialised vec)
tolerance- Tolerance for knot removal (Distance of new solution from old solution, use inf to always remove knots)
degree- The degree of the knot vector
index- Index of the knot to be removed (in knot_vector)
multiplicity- Multiplicity of the knot to be removed
num_removes- Number of times to attempt to remove a knot (will try subsequent knots after failure).
Returns
The number of knots removed by nurbs_remove_curve_knot()
arma::fmat scyl::ocv_to_arma ( const cv::Mat &  in)

ocv_to_arma() converts a OpenCV Mat of floats to a Armadillo mat of floats.

Parameters
inOpenCV Mat to be converted
Returns
Armadillo mat (Not transposed)
int scyl::patch_dichromatic_decompose ( const arma::fcube &  patch,
const arma::fmat &  illuminant,
arma::fmat &  g_out,
arma::fmat &  k_out,
arma::fvec &  s_out,
float &  dichromatic_error,
float &  smooth_error,
float  alpha = 50,
int  debug = 0 
)

patch_dichromatic_decompose() looks at a image patch and determines if the patch is suitable for use in calculating the greater image illuminant.

Parameters
patch- Input patch cube (patch_height x patch_width x image_bands)
illuminant- The estimate for the greater image illuminant
g_out- Specular mat output (pass an empty, uninitialised mat)
k_out- Diffuse mat output (pass an empty, uninitialised mat)
s_out- Reflectance vec output (pass an empty, uninitialised vec)
dichromatic_error- Dichromatic error output (pass an empty, uninitialised double)
smooth_error- Smooth error output (pass an empty, uninitialised double)
alpha- The alpha value used by the illuminant recovery process
debug- Control the amount of output that is generated
Returns
0 or 1, patch is or isn't dichromatic.
Exceptions
ExceptionError
void scyl::pca ( const arma::fcube &  I,
arma::fcube &  principal_components,
arma::fvec &  eigenvalues,
int  num_components = 10,
float  threshold = 0.001 
)

pca() performs Principal Component Analysis upon a Hyperspectral Image.

the pca() function treats each band of the given data as a sample, with each pixel with in the image being a variable. The PCA method used is based on 'Turk and Pentland''s 1991 paper 'Eigenfaces For Recognition' which uses a fast approach to calculate the principal components.

Parameters
I- The input image. (either reflectance, or irradiance)
principal_components- The resulting principal component images. Pass an empty uninitialised cube.
eigenvalues- Corresponding Eigenvalues for the principal components
num_components- Desired number of components to output. The actual number may be less depending on the threshold value and number of bands in the image.
threshold
void scyl::print_version ( )
inline

Print Scyllarus information/copyright banner with version information.

arma::fcube scyl::reconstruct_I ( const arma::fvec &  illuminant,
const arma::fmat &  k,
const arma::fmat &  g,
const arma::fcube &  s,
float  k_factor 
)

reconstruct_I() creates a synthetic version of I based on it's constituents illuminant, k, g, s and big_k.

reconstruct_I() rebuilds I according to the following equation:

I = illuminant * ((g * s) + (k * k_factor))
Parameters
illuminant- Illuminant vector (image_bands)
k- Diffuse mat (image_height x image_width)
g- Specular mat (image_height x image_width)
s- Reflectance cube (image_height x image_width x image_bands)
k_factor- Factor of K
Returns
The reconstructed image cube. (image_height x image_width x image_bands)
Exceptions
ExceptionError
Examples:
cli_pipeline.cpp.
arma::fvec scyl::reconstruct_illuminant ( const arma::fmat &  k,
const arma::fmat &  g,
const arma::fcube &  I,
const arma::fcube &  s,
float  k_factor 
)

reconstruct_illuminant() creates a synthetic version of the illuminant based on given image components I, s, k g and k_factor.

reconstruct_illuminant() rebuilds the illuminant according to the following equation:

illuminant = I / ((g * s) + (k * k_factor))

Where in the above equation returns a cube, this function takes the median value of each band and returns the illuminant as a vector.

Parameters
k- Diffuse mat (image_height x image_width)
g- Specular mat (image_height x image_width)
I- Image cube (image_height x image_width x image_bands)
s- Reflectance cube (image_height x image_width x image_bands)
k_factor- Factor of K
Returns
The reconstructed illuminant (image_bands)
arma::fcube scyl::reconstruct_s ( const arma::fvec &  illuminant,
const arma::fmat &  k,
const arma::fmat &  g,
const arma::fcube &  I,
float  k_factor 
)

reconstruct_s() creates a synthetic version of s based on it's constituents illuminant, k, g, I and k_factor.

reconstruct_s() rebuilds s according to the following equation:

s = ((I / illuminant) - (k * k_factor))) / g
Parameters
illuminant- Illuminant vector (image_bands)
k- Diffuse mat (image_height x image_width)
g- Specular mat (image_height x image_width)
I- Image cube (image_height x image_width x image_bands)
k_factor- Factor of K
Returns
The reconstructed reflectance cube (image_height x image_width x image_bands)
arma::fcube scyl::reconstruct_s_elements ( const arma::fmat &  elements,
const arma::fcube &  abundances,
const arma::ucube &  indexes 
)
Examples:
cli_pipeline.cpp.
void scyl::recover_dichromatic_parameters ( const arma::fcube &  I,
const arma::fvec &  illuminant,
arma::fmat &  k,
arma::fmat &  g,
float &  k_factor,
arma::fcube &  s,
scyl::DICHROMATIC_METHOD  method = scyl::DICHROMATIC_LS,
int  size = -1,
int  gray_threshold = 2,
int  debug = 0 
)

recover_dichromatic_parameters() takes an input image I and illuminant and determines the factors g, k, K and S.

The recover_dichromatic_parameters() function is a wrapper function for subsequent methods of dichromatic parameter recovery. Based on the method selected, the function will check the input arguments and then call the corresponding sub-function. Returned values are normalised.

Parameters
I- Input image cube (image_width x image_height x image_bands)
illuminant- Illuminant vector (image_bands)
k- Specular mat output (pass an empty, uninitialised mat)
g- Shading mat output (pass an empty, uninitialised mat)
k_factor- K factor output (pass an empty, uninitialised float)
s- Reflectance cube output (pass an empty, uninitialised cube)
method- scyl::DICHROMATIC_METHOD method to use for dichromatic parameter recovery. (Default LS).
size- Neighbourhood size of dichromatic parameter recovery. ( > 0, if -1, set according to method)
gray_threshold- Grey threshold for dichromatic parameter recovery. ( > 0, default 2)
debug- Control the amount of output that is generated (0 - 5, default 0)
void scyl::recover_dichromatic_parameters_LS ( const arma::fcube &  I,
const arma::fcube &  illuminant,
arma::fmat &  k,
arma::fmat &  g,
arma::fcube &  st,
int  size = 5,
int  gray_threshold = 2,
int  debug = 0 
)

recover_dichromatic_parameters_LS() takes an input image I and illuminant and determines the factors g, k, and S.

The recover_dichromatic_parameters_LS() uses a least squares based method to recover the Diffuse, Specular and Reflectance components from the input image I.

Parameters
I- Input image cube (image_width x image_height x image_bands)
illuminant- Illuminant vector (image_bands)
k- Specular mat output (pass an empty, uninitialised mat)
g- Shading mat output (pass an empty, uninitialised mat)
st- Reflectance cube output (pass an empty, uninitialised cube)
size- Neighbourhood size of dichromatic parameter recovery. ( > 0, if -1, set according to method)
gray_threshold- Grey threshold for dichromatic parameter recovery. ( > 0, default 2)
debug- Control the amount of output that is generated (0 - 5, default 0)
Exceptions
ExceptionError
arma::fvec scyl::recover_global_illuminant ( const arma::fcube &  I,
scyl::ILLUMINANT_METHOD  method = scyl::ILLUMINANT_HRK,
scyl::illuminant_options  op = scyl::illuminant_options() 
)

recover_global_illuminant() takes an input image I and determines the illuminant.

The recover_global_illuminant() function is a wrapper function for subsequent methods of illuminant recovery. Based on the method selected, the function will check the input arguments and then call the corresponding sub-function. Options for the method chosen are given via the illuminant_options object, if no object (with set options) is given, the default object will be used.

There are several Illumination recovery methods available: ILLUMINANT_HRK - Huynh Robles-Kelly method for illuminant recovery. ILLUMINANT_FS - Finlayson method for illuminant recovery. ILLUMINANT_GE - Grey Edge method for illuminant recovery. ILLUMINANT_GW - Grey World method for illuminant recovery. ILLUMINANT_SG - Shade of Grey method for illuminant recovery. ILLUMINANT_WP - White Patch method for illuminant recovery.

Parameters
I- Input image cube (image_width x image_height x image_bands)
method- The method used to determine the illuminant.
op- Illuminant Options object (see class definition for details) with options set for the method used.
Returns
Illuminant vector (length: image_bands)
arma::fvec scyl::recover_illuminant_FS ( const arma::fcube &  I,
arma::umat &  patches,
int  patch_size = 20,
bool  filtered_already = false,
int  debug = 0 
)

recover_illuminant_FS() takes an input image I and determines the illuminant.

recover_illuminant_FS() uses the Finlayson and Schaefer method for estimating the illuminant from an image. It estimates the dichromatic plane and optimises to find the illuminant.

Parameters
I- Input image cube (image_width x image_height x image_bands)
patches- A mat of pre-selected patches, or a zero mat of size nx1, where n is the number of patches you wish to select, or (default) an empty mat (will be populated with the selected patches)
patch_size- The size (patch_size x patch_size) for the patches used in illuminant recovery
filtered_already- Determines if the image has been prefiltered (and thus do not filter it in this function)
debug- Control the amount of output that is generated (0 - 5, default 0)
Returns
Illuminant vector (length: image_bands)
arma::fvec scyl::recover_illuminant_GE ( const arma::fcube &  I,
int  order = 1,
int  debug = 0 
)

recover_illuminant_GE() takes an input image I and determines the illuminant.

recover_illuminant_GE() uses the Grey Edge method for estimating the illuminant from an image. It uses the Minkowski mean of the image's gradient functions to estimate the illuminant.

Parameters
I- Input image cube (image_width x image_height x image_bands)
order- The order of the Minkowski mean used in the calculation
debug- Control the amount of output that is generated (0 - 5, default 0)
Returns
arma::fvec scyl::recover_illuminant_GW ( const arma::fcube &  I,
int  debug = 0 
)

recover_illuminant_GW() takes an input image I and determines the illuminant.

recover_illuminant_GW() uses the Grey World method for estimating the illuminant from an image. It uses the mean brightness of each band as the illumination power.

Parameters
I- Input image cube (image_width x image_height x image_bands)
debug- Control the amount of output that is generated (0 - 5, default 0)
Returns
arma::fvec scyl::recover_illuminant_HRK ( const arma::fcube &  I,
arma::umat &  patches,
float  alpha,
int  patch_size,
bool  filtered_already = false,
bool  fast_50 = true,
int  debug = 0 
)

recover_illuminant_HRK() takes an input image I and determines the illuminant.

The recover_illuminant_HRK() implements the Huynh Robles-Kelly method for illuminant recovery. The number of patches used in the illuminant recovery is determined based on the input size of patches, if it is empty then the default of 50 will be used, if it is a single column mat, then the number of rows determines the number of patches chosen, otherwise if the number of columns is 4 it is assumed patch locations are provided (and these will be used)

Parameters
I- Input image cube (image_width x image_height x image_bands)
patches- A mat of pre-selected patches, or a zero mat of size nx1, where n is the number of patches you wish to select, or (default) an empty mat (will be populated with the selected patches)
alpha- The alpha value used by the illuminant recovery process ( > 0, default 50)
patch_size- The size (patch_size x patch_size) for the patches used in illuminant recovery
filtered_already- Determines if the image has been prefiltered (and thus do not filter it in this function)
fast_50- Perform a faster Illuminant recovery for images with more than 50 bands
debug- Control the amount of output that is generated (0 - 5, default 0)
Returns
Illuminant vector (length: image_bands)
Exceptions
ExceptionError
arma::fvec scyl::recover_illuminant_SG ( const arma::fcube &  I,
int  order = 2,
int  debug = 0 
)

recover_illuminant_SG() takes an input image I and determines the illuminant.

recover_illuminant_SG() uses the Shade of Grey method for estimating the illuminant from an image. It calculates the p mean of the image brightness for each band, where p is the 'order'.

Parameters
I- Input image cube (image_width x image_height x image_bands)
order- The order of the mean used in the calculation
debug- Control the amount of output that is generated (0 - 5, default 0)
Returns
arma::fvec scyl::recover_illuminant_WP ( const arma::fcube &  I,
int  debug = 0 
)

recover_illuminant_WP() takes an input image I and determines the illuminant.

recover_illuminant_WP() uses the White Patch method for estimating the illuminant from an image. It uses the max brightness of each band as the illumination power.

Parameters
I- Input image cube (image_width x image_height x image_bands)
debug- Control the amount of output that is generated (0 - 5, default 0)
Returns
arma::fmat scyl::recover_materials ( const arma::fcube &  s,
const arma::fvec &  wavelengths,
arma::fcube &  abundances,
arma::ucube &  indices,
scyl::MATERIAL_METHOD  method = scyl::MATERIAL_DA,
scyl::material_options  op = scyl::material_options() 
)

recover_materials() takes an input reflectance cube and returns a list of determined material clusters.

The recover_materials() function is a wrapper function for subsequent methods of material recovery. Based on the method selected, the function will check the input arguments and then call the corresponding sub-function. Returned values are normalised.

Using the DA method will cluster the whole image using deterministic annealing. Using the DAQ/DASAM method will first use DA to cluster a subsampled version of the image (8 times is default) to recover a representative materials list, and then assign each pixel material clusters using the 'unmix_by_spectra'/'sam' method. The DA method is significantly slower but ultimately more accurate (no possibility to miss materials) while the DAQ method is much faster at the expense of potentially lower quality results. In practice, minimal differences are observed.

Parameters
s- Input reflectance cube (image_width x image_height x image_bands)
wavelengths- The wavelength vector corresponding to the bands in s (image_bands)
abundances- Material abundances per pixel output (Pass an empty, uninitialised cube)
indices- Material indices per pixel output (Pass an empty, uninitialised cube)
method- scyl::MATERIAL_METHOD method to use for material recovery. (Default DA, DAQ).
op- Material Options object (see class definition for details) with options set for the method used.
Returns
A mat containing the resultant material signatures. (num_materials x bands)
Examples:
cli_pipeline.cpp.
arma::fmat scyl::recover_materials_DA ( const arma::fcube &  s,
arma::fcube &  material_abundancy,
int  max_clusters = 20,
float  temperature_max = 0.02,
float  temperature_min = 0.00025,
float  cooling_rate = 0.8,
float  split_threshold = -1,
int  debug = 0 
)

recover_materials_DA() takes an input reflectance cube and returns a list of determined material clusters.

The recover_materials_DA() uses an iterative deterministic annealing process to cluster the data points in the input reflectance cube into groups of distinct materials.

Parameters
s- Input reflectance cube (image_width x image_height x image_bands)
material_abundancy- Material abundancy cube output (pass an empty, uninitialised cube)
max_clusters- Max number of materials to find during clustering
temperature_max- Maximum temperature of the DA method.
temperature_min- Minimum temperature of the DA method.
cooling_rate- Cooling rate for DA process
split_threshold- Threshold for splitting clusters in the DS process
debug- Control the amount of output that is generated (0 - 5, default 0)
Exceptions
ExceptionError
Returns
A mat containing the resultant material signatures.
void scyl::register_image_bands ( const arma::fcube &  in,
arma::fvec &  rowoff,
arma::fvec &  coloff,
int  debug = 0 
)

register_image_bands() takes a image cube and attempts to automatically find the alignment of the bands in the image.

register_image_bands() aligns image bands by choosing a 'reference' band which the other bands are then aligned to. The alignment is achieved by transforming each layer of the image into the Fourier domain and searching for the best correlation between the two signals.

Parameters
in- The input image cube to be aligned (inage_height x image_width x image_bands)
rowoff- Resultant row offsets for each band (image_bands)
coloff- Resultant column offsets for each band (image_bands)
debug- Control the amount of output that is generated (0 - 5, default 0)
void scyl::remove_bands ( arma::fcube &  in,
arma::fvec &  waves,
arma::uvec  remove 
)

remove_bands() removes image bands from the given image cube.

remove bands() removed images bands from a image cube based on the input vector band indices. The function assumes bands are in ascending order (i.e 1,3,4,6). It will also remove the corresponding items from the supplied wavelengths vector.

Parameters
in- The image cube to be manipulated
waves- Corresponding wavelength vector to be manipulated
remove- Vector of band indices to remove.
arma::fmat scyl::resample_endmembers ( const arma::fmat &  in_endmembers,
const arma::fvec &  in_waves,
const arma::fvec &  in_waves_new 
)

resample_endmembers() is used to change the number of bands an endmember matrix has using NURBS.

Parameters
in_endmembers- The input matrix of endmembers (num_endmembers x bands)
in_waves- The corresponding wavelengths for the input endmembers (bands)
in_waves_new- The wavelengths to resample to (new_bands)
Returns
The resampled matrix of endmembers (num_endmembers x new_bands)
arma::fcube scyl::resample_image ( const arma::fcube &  in_image,
const arma::fvec &  in_waves,
const arma::fvec &  in_waves_new 
)

resample_image() is used to change the number of bands an image matrix has using NURBS.

Parameters
in_image- The input image (rows x cols x bands)
in_waves- The corresponding wavelengths for the input image (bands)
in_waves_new- The wavelengths to resample to (new_bands)
Returns
The resampled image (rows x cols x new_bands)
void scyl::sam ( const arma::fcube &  I,
const arma::fvec &  wavelengths,
const arma::fmat &  spectra,
const arma::fvec &  spectra_waves,
arma::fcube &  angles,
arma::ucube &  indices,
int  num_materials = 5 
)

sam() takes an input data structure and a list of spectra and calculates angles for each pixel to each spectra given, returning the num_materials closest for each pixel.

Given an image 'I' and a list of spectra, first the spectra and Image will be resized to be compatible sizes (same number of bands) using NURBS and band removal. Then, based on the angle of each spectra in the given list to each pixel in the given image, a list of num_materials angles and the corresponding indexes (of the spectra) will be returned. This is calculated on a pixelwise basis.

Parameters
I- Input image to be unmixed (height x width x bands)
wavelengths- Corresponding wavelengths for I (bands)
spectra- Spectra to unmix I by (num_spectra x spectra_bands) Note this will be resized if spectra_bands is different to bands
spectra_waves- set to arma::fvec() if same as wavelengths to skip resize processing. Corresponding wavelengths for spectra.
angles- Calculated angles, pass an empty, uninitialised arma::fcube
indices- Calculated indexes, pass an empty, uninitialised arma::ucube
num_materials- The number of spectra to classify per pixel (determines the number of bands in abundances and indexes)
void scyl::savitzky_golay ( arma::fcube &  in,
int  window = 5 
)

savitzky_golay applies the Savitzky Golay filter in the spectral domain. Available window sizes are 5,7 and 9.

The function works in-place and uses at most 5 bands of additional memory to perform the filter.

Parameters
in- The input image (rows x cols x bands)
window- The window size (5, 7 and 9 only)
void scyl::scale_cube ( arma::fcube &  in,
int  height,
int  width,
int  cv_interpolate_flag = cv::INTER_LINEAR 
)

scale_cube() resamples a cube spatially using OpenCV resize methods.

scale_cube() will perform the image scaling in place using at most 1 band of temporary memory overhead. Should only be used to scale proportionally.

Parameters
inThe image cube to be resampled
heightReturned image cube height
widthReturned image cube width
cv_interpolate_flagOptional OpenCV interpolation method. Default is INTER_LINEAR, see (http://docs.opencv.org/modules/imgproc/doc/geometric_transformations.html#resize)
arma::umat scyl::select_smooth_patches ( const arma::fcube &  I,
arma::umat &  patches,
arma::fmat &  light_estimation_mask,
arma::fmat &  contrast_mask,
arma::fmat &  highlight_mask,
int  patch_size = 20,
int  num_patches = 50,
int  debug = 0 
)

select_smooth_patches() returns a list of the first num_patches smooth patches it finds while inspecting input image I.

Parameters
I- Input image cube (image_width x image_height x image_bands)
patches- Selected patches mat output (pass an empty, uninitialised mat)
light_estimation_mask- Light estimation mask mat output (pass an empty, uninitialised mat)
contrast_mask- Contrast mask mat output (pass an empty, uninitialised mat)
highlight_mask- Highlight mask mat output (pass an empty, uninitialised mat)
patch_size- Size of each patch (patch_size x patch_size), pixels
num_patches- Number of patches to select
debug- Control the amount of output that is generated (0 - 5, default 0)
Returns
a patch map showing the locations of selected patches.
Exceptions
ExceptionError
void scyl::sort_bands ( arma::fcube &  I,
arma::fvec &  wavelengths 
)

sort_bands() manipulates a Image cube and wavelengths vector.

Based on the wavelengths vector, bands in the image are sorted into order of ascending wavelengths. The operation is done in place (with one band of temporary memory) to reduce the memory overhead of sorting.

Parameters
I- The image cube to be sorted
wavelengths- The unsorted wavelengths (returned sorted)
arma::fcube scyl::spectral_derivative ( const arma::fcube &  I,
const arma::fvec &  wavelengths 
)

spectral_derivative calculates the derivative of the image in the spectral domain.

Each band in the result is given as the change from the previous band to the next.

Parameters
I- The input image (rows x cols x bands)
wavelengths- The corresponding wavelengths for the input image (bands)
Returns
The spectral derivative cube.
void scyl::svm ( const arma::fcube &  I,
const arma::fvec &  wavelengths,
const arma::fmat &  spectra,
const arma::fvec &  spectra_waves,
const arma::uvec &  spectra_labels,
arma::umat &  indices 
)

svm() takes an input data structure and a list of spectra and classifies each pixel as one of the given spectra, returning the index map.

Given an image 'I' and a list of spectra, first the spectra and Image will be resized to be compatible sizes (same number of bands) using NURBS and band removal. Then, using a Support Vector Machine, each pixel will be classified as one of the given spectra. This function uses OpenCV for SVM.

If you are using data from an endmember library, this code can be used to get the uvec of labels required for the function:

spectra_labels = arma::zeros<arma::uvec>(spectral_library->num_endmembers());
label_set = spectral_library->labels_mean();
label_list = spectral_library->labels();
for (int i = 0; i < spectra_labels.n_elem; i++)
{
  for (int j = 0; j < label_set.size(); j++)
  {
    if (label_list.at(i).compare(label_set.at(j)) == 0)
    {
      spectra_labels(i) = j + 1;
      break;
    }
  }
}
Parameters
I- Input image to be unmixed (height x width x bands)
wavelengths- Corresponding wavelengths for I (bands)
spectra- Spectra to unmix I by (num_spectra x spectra_bands) Note this will be resized if spectra_bands is different to bands
spectra_waves- set to arma::fvec() if same as wavelengths to skip resize processing. Corresponding wavelengths for spectra
spectra_labels- The list of labels corresponding to each spectra given (see above for how to make this).
indices- Calculated indexes, pass an empty, uninitialised arma::umat
void scyl::unmix_by_spectra ( const arma::fcube &  I,
const arma::fvec &  wavelengths,
const arma::fmat &  spectra,
const arma::fvec &  spectra_waves,
arma::fcube &  abundances,
arma::ucube &  indices,
int  num_materials = 5 
)

unmix_by_spectra() takes an input data structure and a list of spectra and calculates abundances for each pixel in the data according to the given spectra.

Given an image 'I' and a list of spectra, first the spectra and Image will be resized to be compatible sizes (same number of bands) using NURBS and band removal. Then, based on the distance of each spectra in the given list to each pixel in the given image, a list of num_materials abundances and the corresponding indexes (of the spectra) will be returned. This is calculated on a pixelwise basis.

To use this function for a problem such as spectral unmixing (for example, if you have a spectral signatures and a list of base element spectra), use the routine as follows: Rotate the signatures so that each one represents a 'pixel' in an image cube, and pass the data to the function. the resulting data structures will indicate how much of each 'base element' comprises each pixel (abundances being proportions corresponding to the index at each pixel location)

Parameters
I- Input image to be unmixed (height x width x bands) Note function is designed to work with reflectances, so it is advised that I and spectra are both reflectance.
wavelengths- Corresponding wavelengths for I (bands)
spectra- Spectra to unmix I by (num_spectra x spectra_bands) Note this will be resized if spectra_bands is different to bands
spectra_waves- set to arma::fvec() if same as wavelengths to skip resize processing. Corresponding wavelengths for spectra.
abundances- Calculated abundances, pass an empty, uninitialised arma::fcube
indices- Calculated indexes, pass an empty, uninitialised arma::ucube
num_materials- The number of spectra to classify per pixel (determines the number of bands in abundances and indexes)
template<class T >
void scyl::version ( T &  in)

Get the version number of Scyllarus (Casted from Float (Major.MinorSub)

Parameters
inData to be populated with Version Number
template<class T >
void scyl::version ( arma::Mat< T > &  in)

Get the version number of Scyllarus (Casted from Float (Major.MinorSub)

Parameters
inData to be populated with Version Numbers
std::string scyl::version_string ( )
inline

Get the version number formatted as a string [Major.Minor.Sub(release)].

Returns
Version number as formatted string
arma::fvec scyl::wavelength_subset ( const arma::fvec &  set,
const arma::fvec &  bounds 
)

wavelength_subset() removes elements from set that are outside of the range of elements in bounds.

wavelength_subset() will determine the minimum and maximum values in the given bounds array and remove any elements from the given set that lie outside of these values.

Parameters
set- Input vector containing set to be 'clipped'
bounds- Input vector describing range to 'clip' set to
Returns
The 'clipped' result.
arma::fmat scyl::wiener2 ( const arma::fmat &  a,
arma::fvec  nhood 
)

Apply a de-blurring wiener filter to an arma matrix (based on matlab implementation)

Parameters
aarma matrix
nhoodsize of the kernel (vector with height and width)
template<class T >
void scyl::zero_to_one ( T &  in)

zero_to_one replaces 0's in an Armadillo data structure with 1's.

Parameters
inthe data structure to be manipulated.

Variable Documentation

const float scyl::CAMERA_SENSITIVITY[35][3]
const int scyl::CAMERA_SENSITIVITY_LENGTH = 35
const float scyl::CAMERA_SENSITIVITY_WAVELENGTHS[35]
Initial value:
=
{
380, 390, 400, 410, 420, 430,
440, 450, 460, 470, 480, 490,
500, 510, 520, 530, 540, 550,
560, 570, 580, 590, 600, 610,
620, 630, 640, 650, 660, 670,
680, 690, 700, 710, 720
}
const float scyl::COLOUR_JET[100][3]
const int scyl::PHOTOPIC_END = 830

Length of PHOTOPIC_VALUES.

const int scyl::PHOTOPIC_LENGTH = 441

Starting wavelength of PHOTOPIC_VALUES.

const int scyl::PHOTOPIC_START = 390
const float scyl::PHOTOPIC_VALUES[441]

Ending wavelength of PHOTOPIC_VALUES.

const float scyl::VERSION_FLOAT = 1.03f
const int scyl::VERSION_MAJOR = 1
const int scyl::VERSION_MINOR = 0
const int scyl::VERSION_SUB = 3