#ifndef STIM_ENVI_H
#define STIM_ENVI_H

#include "../envi/envi_header.h"
#include "../envi/bsq.h"
#include "../envi/bip.h"
#include "../envi/bil.h"
#include <iostream>
#include "../image/image.h"

namespace stim{

/** This class implements reading of ENVI hyperspectral files. These files can be stored in multiple orientations
	(including BSQ, BIP, and BIL) in order to optimize streaming speed depending on applications. Basic ENVI
	files are stored on disk as a large binary file with a corresponding header. Code for reading and processing
	ENVI header files is in the envi_header class.
*/
class envi{	

	void* file;		//void pointer to the relevant file reader (bip, bsq, or bil - with appropriate data type)

public:

	envi_header header;

	/// Allocate memory for a new ENVI file based on the current interleave format (BIP, BIL, BSQ) and data type.
	bool allocate(){

		file = NULL;	//set file to a NULL pointer

		if(header.interleave == envi_header::BSQ){
			if(header.data_type ==envi_header::float32)
				return(file = new bsq<float>());
			else if(header.data_type == envi_header::float64)
				return(file = new bsq<double>());
		}
		else if(header.interleave == envi_header::BIP){
			if(header.data_type ==envi_header::float32)
				return(file = new bip<float>());
			else if(header.data_type == envi_header::float64)
				return(file = new bip<double>());
		}
		else if(header.interleave == envi_header::BIL){
			if(header.data_type ==envi_header::float32)
				return(file = new bil<float>());
			else if(header.data_type == envi_header::float64)
				return(file = new bil<double>());
		}

		exit(1);	//if the function hasn't already returned, we don't handle this state

	}

	/// Open an existing ENVI file given the file and header names.

	/// @param filename is the name of the ENVI binary file
	/// @param headername is the name of the ENVI header file
	bool open(std::string filename, std::string headername){

		//allocate memory
		allocate();

		//load the header
		header.load(headername);

		//load the file
		if(header.interleave == envi_header::BSQ) {		//if the infile is bsq file
			if(header.data_type == envi_header::float32) {
				return ((bsq<float>*)file)->open(filename, header.samples, header.lines, header.bands, header.header_offset, header.wavelength);
			}
			else if(header.data_type == envi_header::float64) {
				return ((bsq<double>*)file)->open(filename, header.samples, header.lines, header.bands, header.header_offset, header.wavelength);
			}
			else
				return false;
		}

		else if(header.interleave == envi_header::BIL) {		//if the infile is bil file
			if(header.data_type == envi_header::float32) {
				return ((bil<float>*)file)->open(filename, header.samples, header.lines, header.bands, header.header_offset, header.wavelength);
			}
			else if(header.data_type == envi_header::float64) {
				return ((bil<double>*)file)->open(filename, header.samples, header.lines, header.bands, header.header_offset, header.wavelength);
			}
			else
				return false;
		}
				
		else if(header.interleave == envi_header::BIP) {		//if the infile is bip file
			if(header.data_type == envi_header::float32) {
				return ((bip<float>*)file)->open(filename, header.samples, header.lines, header.bands, header.header_offset, header.wavelength);
			}
			else if(header.data_type == envi_header::float64) {
				return ((bip<double>*)file)->open(filename, header.samples, header.lines, header.bands, header.header_offset, header.wavelength);
			}
			else
				return false;		
		}
		
		else{
		std::cout<<"ERROR: unidentified type file       "<<headername<<std::endl;
		exit(1);
		}

	}

	/// Normalize a hyperspectral ENVI file given a band number and threshold.

	/// @param outfile is the name of the normalized file to be output
	/// @param band is the band label to be output
	/// @param threshold is a threshold value specified such that normalization will only be done to values in the band > threshold (preventing division by small numbers)
	bool normalize(std::string outfile, double band, double threshold = 0.0){
	
		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->normalize(outfile, band, threshold);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->normalize(outfile,band, threshold);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIL){		//if the infile is bil file
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->normalize(outfile, band);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->normalize(outfile,band);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIP){		//if the infile is bip file
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->normalize(outfile, band);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->normalize(outfile,band);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else{
			std::cout<<"ERROR: unidentified file type"<<std::endl;
			exit(1);
		}
		return false;
	}

	/// Performs piecewise linear baseline correction of a hyperspectral file/

	/// @param outfile is the file name for the baseline corrected output
	/// @param w is a list of band labels to serve as baseline points (zero values)
	bool baseline(std::string outfile, std::vector<double> w){

		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->baseline(outfile, w);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->baseline(outfile,w);
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}

		else if(header.interleave == envi_header::BIL){		//if the infile is bil file
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->baseline(outfile, w);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->baseline(outfile, w);
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}

		else if(header.interleave == envi_header::BIP){		//if the infile is bip file 
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->baseline(outfile, w);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->baseline(outfile, w);
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}

		else{
			std::cout<<"ERROR: unidentified file type"<<std::endl;
			exit(1);
		}
	}

	/// Converts ENVI files between interleave types (BSQ, BIL, and BIP)

	/// @param outfile is the file name for the converted output
	/// @param interleave is the interleave format for the destination file
	bool convert(std::string outfile, stim::envi_header::interleaveType interleave){

		if(header.interleave == envi_header::BSQ){			//if the infile is bsq file

			if(header.data_type ==envi_header::float32){		//if the data type of infile is float
				if(interleave == envi_header::BSQ){
					std::cout<<"ERROR:  is already BSQ file"<<std::endl;
					exit(1);
				}
				else if(interleave == envi_header::BIL)			//if the target file is bil file
					return ((bsq<float>*)file)->bil(outfile);
				else if(interleave == envi_header::BIP)			//if the target file is bip file
					return ((bsq<float>*)file)->bip(outfile);
			}

			else if(header.data_type == envi_header::float64){		//if the data type is float
				if(interleave == envi_header::BSQ){
					std::cout<<"ERROR:  is already BSQ file"<<std::endl;
					exit(1);
				}
				else if(interleave == envi_header::BIL)
					return ((bsq<double>*)file)->bil(outfile);
				else if(interleave == envi_header::BIP)
					return ((bsq<double>*)file)->bip(outfile);
			}
			
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);				
			}
		}

		else if(header.interleave == envi_header::BIL){

			if(header.data_type ==envi_header::float32){		//if the data type of infile is float
				if(interleave == envi_header::BIL){
					std::cout<<"ERROR:  is already BIL file"<<std::endl;
					exit(1);
				}
				else if(interleave == envi_header::BSQ)			//if the target file is bsq file
					return ((bil<float>*)file)->bsq(outfile);
				else if(interleave == envi_header::BIP)			//if the target file is bip file
					return ((bil<float>*)file)->bip(outfile);
			}

			else if(header.data_type == envi_header::float64){		//if the data type is float
				if(interleave == envi_header::BIL){
					std::cout<<"ERROR:  is already BIL file"<<std::endl;
					exit(1);
				}
				else if(interleave == envi_header::BSQ)
					return ((bil<double>*)file)->bsq(outfile);
				else if(interleave == envi_header::BIP)
					return ((bil<double>*)file)->bip(outfile);
			}
			
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);				
			}
		}

		else if(header.interleave == envi_header::BIP){
		
			if(header.data_type ==envi_header::float32){		//if the data type of infile is float
				if(interleave == envi_header::BIP){
					std::cout<<"ERROR:  is already BIP file"<<std::endl;
					exit(1);
				}
				else if(interleave == envi_header::BIL)			//if the target file is bil file
					return ((bip<float>*)file)->bil(outfile);
				else if(interleave == envi_header::BSQ)			//if the target file is bsq file
					return ((bip<float>*)file)->bsq(outfile);
			}

			else if(header.data_type == envi_header::float64){		//if the data type is float
				if(interleave == envi_header::BIP){
					std::cout<<"ERROR:  is already BIP file"<<std::endl;
					exit(1);
				}
				else if(interleave == envi_header::BIL)			//if the target file is bil file
					return ((bip<double>*)file)->bil(outfile);
				else if(interleave == envi_header::BSQ)			//if the target file is bsq file
					return ((bip<double>*)file)->bsq(outfile);
			}
			
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);				
			}
		}
		
		else{
			std::cout<<"ERROR: unidentified interleave type"<<std::endl;
			exit(1);
		}
		return false;
	
	}

	/// Builds a mask from a band image and threshold value

	/// @param mask_band is the label for the band that will be used to build the mask
	/// @param threshold is a value selected such that all band values greater than threshold will have a mask value of 'true'
	/// @param p is memory of size X*Y that will store the resulting mask
	bool build_mask(double mask_band, double threshold, unsigned char* p = NULL)	{

		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->build_mask(mask_band, threshold, p);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->build_mask(mask_band, threshold, p);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIL){		//if the infile is bil file
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->build_mask(mask_band, threshold, p);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->build_mask(mask_band, threshold, p);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIP){		//if the infile is bip file
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->build_mask(mask_band, threshold, p);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->build_mask(mask_band, threshold, p);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		return false;
	}

	/// Applies a mask to the ENVI file.

	/// @param outfile is the name of the resulting masked output file
	/// @param p is memory of size X*Y containing the mask (0 = false, all other values are true)
	bool apply_mask(std::string outfile, unsigned char* p)	
	{

		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->apply_mask(outfile, p);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->apply_mask(outfile, p);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIL){		//if the infile is bil file
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->apply_mask(outfile, p);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->apply_mask(outfile, p);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIP){		//if the infile is bip file
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->apply_mask(outfile, p);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->apply_mask(outfile, p);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else{
			std::cout<<"ERROR: unidentified file type"<<std::endl;
			exit(1);
		}
		return false;
	}

	/// Compute the ratio of two baseline-corrected peaks. The result is stored in a pre-allocated array.

	/// @param lb1 is the label value for the left baseline point for the first peak (numerator)
	/// @param rb1 is the label value for the right baseline point for the first peak (numerator)
	/// @param pos1 is the label value for the first peak (numerator) position
	/// @param lb2 is the label value for the left baseline point for the second peak (denominator)
	/// @param rb2 is the label value for the right baseline point for the second peak (denominator)
	/// @param pos2 is the label value for the second peak (denominator) position
	/// @param result is a pointer to a pre-allocated array at least X * Y * sizeof(T) in size
	bool ph_to_ph(double lb1, double rb1, double pos1, double lb2, double rb2, double pos2, void * result){
		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->ph_to_ph(lb1, rb1, pos1, lb2, rb2, pos2, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->ph_to_ph(lb1, rb1, pos1, lb2, rb2, pos2, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIL){		//if the infile is bil file
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->ph_to_ph(lb1, rb1, pos1, lb2, rb2, pos2, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->ph_to_ph(lb1, rb1, pos1, lb2, rb2, pos2, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIP){		//if the infile is bip file
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->ph_to_ph(lb1, rb1, pos1, lb2, rb2, pos2, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->ph_to_ph(lb1, rb1, pos1, lb2, rb2, pos2, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else{
			std::cout<<"ERROR: unidentified file type"<<std::endl;
			exit(1);
		}
		return false;
	}

	/// Compute the ratio between a peak area and peak height.

	/// @param lb1 is the label value for the left baseline point for the first peak (numerator)
	/// @param rb1 is the label value for the right baseline point for the first peak (numerator)
	/// @param pos1 is the label value for the first peak (numerator) position
	/// @param lb2 is the label value for the left baseline point for the second peak (denominator)
	/// @param rb2 is the label value for the right baseline point for the second peak (denominator)
	/// @param pos2 is the label value for the second peak (denominator) position
	/// @param result is a pointer to a pre-allocated array at least X * Y * sizeof(T) in size
	bool pa_to_ph(double lb1, double rb1, double lab1, double rab1, double lb2, double rb2, double pos, void * result){
		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->pa_to_ph(lb1, rb1, lab1, rab1, lb2, rb2, pos, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->pa_to_ph(lb1, rb1, lab1, rab1, lb2, rb2, pos, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIL){		//if the infile is bil file
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->pa_to_ph(lb1, rb1, lab1, rab1, lb2, rb2, pos, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->pa_to_ph(lb1, rb1, lab1, rab1, lb2, rb2, pos, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIP){		//if the infile is bip file
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->pa_to_ph(lb1, rb1, lab1, rab1, lb2, rb2, pos, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->pa_to_ph(lb1, rb1, lab1, rab1, lb2, rb2, pos, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else{
			std::cout<<"ERROR: unidentified file type"<<std::endl;
			exit(1);
		}
		return false;
	}

	/// Compute the ratio between two peak areas.

	/// @param lb1 is the label value for the left baseline point for the first peak (numerator)
	/// @param rb1 is the label value for the right baseline point for the first peak (numerator)
	/// @param lab1 is the label value for the left bound (start of the integration) of the first peak (numerator)
	/// @param rab1 is the label value for the right bound (end of the integration) of the first peak (numerator)
	/// @param lb2 is the label value for the left baseline point for the second peak (denominator)
	/// @param rb2 is the label value for the right baseline point for the second peak (denominator)
	/// @param lab2 is the label value for the left bound (start of the integration) of the second peak (denominator)
	/// @param rab2 is the label value for the right bound (end of the integration) of the second peak (denominator)	
	/// @param result is a pointer to a pre-allocated array at least X * Y * sizeof(T) in size
	bool pa_to_pa(double lb1, double rb1, double lab1, double rab1,
					double lb2, double rb2, double lab2, double rab2, void* result){
		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->pa_to_pa(lb1, rb1, lab1, rab1, lb2, rb2, lab2, rab2, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->pa_to_pa(lb1, rb1, lab1, rab1, lb2, rb2, lab2, rab2, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIL){		//if the infile is bil file
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->pa_to_pa(lb1, rb1, lab1, rab1, lb2, rb2, lab2, rab2, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->pa_to_pa(lb1, rb1, lab1, rab1, lb2, rb2, lab2, rab2, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIP){		//if the infile is bip file
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->pa_to_pa(lb1, rb1, lab1, rab1, lb2, rb2, lab2, rab2, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->pa_to_pa(lb1, rb1, lab1, rab1, lb2, rb2, lab2, rab2, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else{
			std::cout<<"ERROR: unidentified file type"<<std::endl;
			exit(1);
		}
		return false;
	}

	/// Compute the centroid of a baseline corrected peak.

	/// @param lb is the label value for the left baseline point
	/// @param rb is the label value for the right baseline point
	/// @param lab is the label for the start of the peak
	/// @param rab is the label for the end of the peak
	/// @param result is a pointer to a pre-allocated array at least X * Y * sizeof(T) in size
	bool cpoint(double lb1, double rb1, double lab1, double rab1, void* result){
		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->cpoint(lb1, rb1, lab1, rab1, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->cpoint(lb1, rb1, lab1, rab1, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIL){		//if the infile is bil file
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->cpoint(lb1, rb1, lab1, rab1, (float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->cpoint(lb1, rb1, lab1, rab1, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else if(header.interleave == envi_header::BIP){		//if the infile is bip file
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->cpoint(lb1, rb1, lab1, rab1,(float*)result);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->cpoint(lb1, rb1, lab1, rab1, (double*)result);
			else
				std::cout<<"ERROR: unidentified data type"<<std::endl;
		}

		else{
			std::cout<<"ERROR: unidentified file type"<<std::endl;
			exit(1);
		}
		return false;
	}

	/// Closes the ENVI file.
	bool close(){
		if(header.interleave == envi_header::BSQ){
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->close();
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->close();
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}

		else if(header.interleave == envi_header::BIL){
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->close();
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->close();
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}

		else if(header.interleave == envi_header::BIP){
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->close();
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->close();
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}
		return false;
	}

	/// Retrieve a single pixel and stores it in pre-allocated memory.

	/// @param p is a pointer to pre-allocated memory at least sizeof(T) in size.
	/// @param n is an integer index to the pixel using linear array indexing.
	bool pixel(void * p, unsigned n){
		if(header.interleave == envi_header::BSQ){
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->pixel((float*)p, n);
			else if(header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->pixel((double*)p, n);
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}
		else if(header.interleave == envi_header::BIL){
			if(header.data_type ==envi_header::float32)
				return ((bil<float>*)file)->pixel((float*)p, n);
			else if(header.data_type == envi_header::float64)
				return ((bil<double>*)file)->pixel((double*)p, n);
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}
		else if(header.interleave == envi_header::BIP){
			if(header.data_type ==envi_header::float32)
				return ((bip<float>*)file)->pixel((float*)p, n);
			else if(header.data_type == envi_header::float64)
				return ((bip<double>*)file)->pixel((double*)p, n);
			else{
				std::cout<<"ERROR: unidentified data type"<<std::endl;
				exit(1);
			}
		}
		return false;
	}

	/// Saves a header file describing the current ENVI file parameters.
	bool save_header(std::string filename){

		//save the header file here
		header.save(filename);

		return true;
	}

	/// Retrieve a single band (by numerical label) and stores it in pre-allocated memory.

	/// @param p is a pointer to an allocated region of memory at least X * Y * sizeof(T) in size.
	/// @param wavelength is a floating point value (usually a wavelength in spectral data) used as a label for the band to be copied.
	bool band(void* ptr, double wavelength){

		if(header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if(header.data_type ==envi_header::float32)
				return ((bsq<float>*)file)->band((float*)ptr, wavelength);
			else if (header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->band((double*)ptr, wavelength);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIL){
			if (header.data_type == envi_header::float32)
				return ((bil<float>*)file)->band((float*)ptr, wavelength);
			else if (header.data_type == envi_header::float64)
				return ((bil<double>*)file)->band((double*)ptr, wavelength);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIP){
			if (header.data_type == envi_header::float32)
				return ((bip<float>*)file)->band((float*)ptr, wavelength);
			else if (header.data_type == envi_header::float64)
				return ((bip<double>*)file)->band((double*)ptr, wavelength);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		return false;
	}

	/// Retrieve a single band (based on index) and stores it in pre-allocated memory.

	/// @param p is a pointer to an allocated region of memory at least X * Y * sizeof(T) in size.
	/// @param page <= B is the integer number of the band to be copied.
	bool band_index(void* ptr, unsigned int b){
		if (header.interleave == envi_header::BSQ){		//if the infile is bsq file
			if (header.data_type == envi_header::float32)
				return ((bsq<float>*)file)->band_index((float*)ptr, b);
			else if (header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->band_index((double*)ptr, b);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIL){
			if (header.data_type == envi_header::float32)
				return ((bil<float>*)file)->band_index((float*)ptr, b);
			else if (header.data_type == envi_header::float64)
				return ((bil<double>*)file)->band_index((double*)ptr, b);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIP){
			if (header.data_type == envi_header::float32)
				return ((bip<float>*)file)->band_index((float*)ptr, b);
			else if (header.data_type == envi_header::float64)
				return ((bip<double>*)file)->band_index((double*)ptr, b);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		return false;
	}
	
	/// Helper function that loads a mask into memory given a filename.

	/// @param mask is a pointer to pre-allocated memory of size X*Y
	/// @param maskname is the file name for the image that will serve as the mask
	bool load_mask(unsigned char * mask, std::string maskname){
		//open the mask file
		stim::image<unsigned char> mask_image(maskname);
		mask_image.data_noninterleaved(mask);
		//save mask file into memory
		//memcpy(mask, mask_image.data_noninterleaved(), mask_image.size());
		//mask_image.clear();
		return true;
	}

	//p:start positon; N: number of pixels saved in X;
	bool feature_matrix(void * X,  unsigned char * mask, unsigned start, unsigned N)
	{
		//save pixels in X as floating numbers: float * p
		float * p = (float*)malloc(header.bands * sizeof(float));		//save pixel information
		unsigned pixels = header.samples * header.lines;		//calculate pixel numbers in a band
		unsigned count = 0;		//for counting use
		unsigned j = 0;		//memory the pointer location in X

		//create two indices into the mask (mask index and pixel index)
		unsigned mi = 0;	//valid pixel index
		unsigned pi = 0;	//actual pixel in the mask

		//find the actual pixel index for the mask index "start"
		while(mi < start){
			if(mask[pi])
				mi++;

			pi++;
		}

		for(unsigned i = pi; i < pixels; i++){
			if(mask[i] != 0){
				pixel(p, i);
				//copy p to X
				for(unsigned k = 0; k < header.bands; k++){
					((float *)X)[j] = p[k];
					j++;
				}
				count++;
				if(count == N)
					break;
			}
			else
				continue;
		}
		if(count < N){
			std::cout << "number of valid pixels in the mask : " << count <<"is less than N: "<< N;
			exit(1);
		}
		free(p);
		return true;
	}

	/// Calculate the mean value for all masked (or valid) pixels in a band and returns the average spectrum

	/// @param p is a pointer to pre-allocated memory of size [B * sizeof(T)] that stores the mean spectrum
	/// @param mask is a pointer to memory of size [X * Y] that stores the mask value at each pixel location
	bool avg_band(void * p, unsigned char* mask){
		if (header.interleave == envi_header::BSQ){
			if (header.data_type == envi_header::float32)
				return ((bsq<float>*)file)->avg_band((float*)p, mask);
			else if (header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->avg_band((double*)p,  mask);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIL){
			if (header.data_type == envi_header::float32)
				return ((bil<float>*)file)->avg_band((float*)p,  mask);
			else if (header.data_type == envi_header::float64)
				return ((bil<double>*)file)->avg_band((double*)p,  mask);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIP){
			if (header.data_type == envi_header::float32)
				return ((bip<float>*)file)->avg_band((float*)p, mask);
			else if (header.data_type == envi_header::float64)
				return ((bip<double>*)file)->avg_band((double*)p, mask);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		return false;
	}

	/// Calculate the covariance matrix for all masked pixels in the image.

	/// @param co is a pointer to pre-allocated memory of size [B * B] that stores the resulting covariance matrix
	/// @param avg is a pointer to memory of size B that stores the average spectrum
	/// @param mask is a pointer to memory of size [X * Y] that stores the mask value at each pixel location
	bool co_matrix(void* co, void* avg, unsigned char* mask){
		if (header.interleave == envi_header::BSQ){
			if (header.data_type == envi_header::float32)
				return ((bsq<float>*)file)->co_matrix((float*)co, (float*)avg, mask);
			else if (header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->co_matrix((double*)co, (double*)avg, mask);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIL){
			if (header.data_type == envi_header::float32)
				return ((bil<float>*)file)->co_matrix((float*)co, (float*)avg, mask);
			else if (header.data_type == envi_header::float64)
				return ((bil<double>*)file)->co_matrix((double*)co, (double*)avg, mask);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIP){
			if (header.data_type == envi_header::float32)
				return ((bip<float>*)file)->co_matrix((float*)co, (float*)avg, mask);
			else if (header.data_type == envi_header::float64)
				return ((bip<double>*)file)->co_matrix((double*)co, (double*)avg, mask);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		return false;
	}


	/// Crop a region of the image and save it to a new file.

	/// @param outfile is the file name for the new cropped image
	/// @param x0 is the lower-left x pixel coordinate to be included in the cropped image
	/// @param y0 is the lower-left y pixel coordinate to be included in the cropped image
	/// @param x1 is the upper-right x pixel coordinate to be included in the cropped image
	/// @param y1 is the upper-right y pixel coordinate to be included in the cropped image
	bool crop(std::string outfile,unsigned x0, unsigned y0, unsigned x1, unsigned y1){

		if (header.interleave == envi_header::BSQ){
			if (header.data_type == envi_header::float32)
				return ((bsq<float>*)file)->crop(outfile, x0, y0, x1, y1);
			else if (header.data_type == envi_header::float64)
				return ((bsq<double>*)file)->crop(outfile, x0, y0, x1, y1);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIL){
			if (header.data_type == envi_header::float32)
				return ((bil<float>*)file)->crop(outfile, x0, y0, x1, y1);
			else if (header.data_type == envi_header::float64)
				return ((bil<double>*)file)->crop(outfile, x0, y0, x1, y1);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		else if (header.interleave == envi_header::BIP){
			if (header.data_type == envi_header::float32)
				return ((bip<float>*)file)->crop(outfile, x0, y0, x1, y1);
			else if (header.data_type == envi_header::float64)
				return ((bip<double>*)file)->crop(outfile, x0, y0, x1, y1);
			else{
				std::cout << "ERROR: unidentified data type" << std::endl;
				exit(1);
			}
		}
		return false;
	}

};

}	//end namespace rts

#endif