#include <shogun/base/init.h>
#include <shogun/io/AsciiFile.h>
#include <shogun/lib/common.h>
#include <shogun/lib/SGVector.h>
#include <shogun/io/SGIO.h>
#include <shogun/features/SparseFeatures.h>
#include <shogun/kernel/GaussianKernel.h>
#include <shogun/statistics/KernelMeanMatching.h>
#include <iostream>

using namespace shogun;

void print_message(FILE* target, const char* str)
{
	fprintf(target, "%s", str);
}

void print_warning(FILE* target, const char* str)
{
	fprintf(target, "%s", str);
}

void print_error(FILE* target, const char* str)
{
	fprintf(target, "%s", str);
}

int main(int argc, char **argv)
{
	init_shogun(&print_message, &print_warning,&print_error);
	sg_io->set_loglevel(MSG_DEBUG);

	// Read in source and target data, indices indicating offsets in file 
	// Create sparse matrix
	const char* filename=argv[1];
	const int sourceEnd = atoi(argv[2]);
	CAsciiFile* data = new CAsciiFile(filename,'r',"x1");
	int32_t num_feat=0, num_vec=0;
	SGSparseVector<float32_t>* dataMatrix = NULL;
	data->get_sparse_matrix(dataMatrix, num_feat, num_vec);
	std::cout << num_feat << " " << num_vec << std::endl;

	// Create Kernel
	CSparseFeatures<float32_t>* feats = new CSparseFeatures<float32_t>(dataMatrix, num_feat, num_vec, false);
	CGaussianKernel* kernel = new CGaussianKernel(feats, feats, 50, 10000);
	/*/ 
	for(int i=0; i<sourceEnd; i++)
		for(int j=sourceEnd; j<num_vec; j++) {
			float kij = kernel->kernel(i,j);
			if(kij != 0) {
				std::cout << "K(" << i << "," << j << ") = " << kij << std::endl;
			}
		} // */

	// Call KMM optimization code
	int sourceSz = sourceEnd;
	int targetSz = num_vec - sourceSz;
	SGVector<index_t> indices1(sourceSz);
	SGVector<index_t> indices2(targetSz);
	indices1.range_fill(0);		
	indices2.range_fill(sourceSz);

	double rho = 1;
	CKernelMeanMatching* kmm = new CKernelMeanMatching(kernel, indices1, indices2);
	SGVector<float64_t> betas = kmm->compute_weights(sourceSz*rho); // argument: bound on \sum_beta_i
	for(int i=0; i<betas.vlen; i++)
		std::cout << "Beta(" << i << ") = " << betas[i] << std::endl;

	data->close();
	SG_UNREF(feats);
	exit_shogun();
	return 0;
}
