#include "cars_excitation.h"

int main(int argc, char** argv)
	{
	if(argc < 4)
		{
		cerr << "Usage: " << argv[0] << " <pump-parameters file> <stokes-parameters file> <output directory>\n";
		exit(-1);
		}
	string param_file_p = argv[1];
	string param_file_s = argv[2];
	string outputdir = argv[3];
	stringstream filename;

	focus3d_parameters pp, ps, pc, temp;
	if(!pp.load(param_file_p)) {exit (-1);} // load the parameters file

	if(!ps.load(param_file_s)) {exit (-1);} // load the parameters file

	cout << "PUMP PARAMETERS:\n";
	pp.print();

	cout << "STOKES PARAMETERS:\n";
	ps.print();

	temp = pp;
	temp.lambda = ps.lambda;

	if(temp != ps)
		{
		cerr << "Input error! The spatial parameters in the two files do not match! Exiting\n";
		exit(-1);
		}

	pc = ps;
	pc.lambda = 1.0/((2.0/pp.lambda) - (1.0/ps.lambda));
	cout << "CARS PARAMETERS:\n";
	pc.print();

	Array<complex<double>,3> E_focus_pump_x, E_focus_pump_y, E_focus_pump_z;
	Array<complex<double>,3> E_focus_stokes_x, E_focus_stokes_y, E_focus_stokes_z;
	Array<complex<double>,3> cars_pol_x, cars_pol_y, cars_pol_z;
	Array<complex<double>,2> E_inc_x, E_inc_y, E_inc_z;

	E_focus_pump_x.resize(pp.Nx,pp.Ny,pp.Nz); E_focus_pump_x = 0;
	E_focus_pump_y.resize(pp.Nx,pp.Ny,pp.Nz); E_focus_pump_y = 0;
	E_focus_pump_z.resize(pp.Nx,pp.Ny,pp.Nz); E_focus_pump_z = 0;

	E_focus_stokes_x.resize(ps.Nx,ps.Ny,ps.Nz); E_focus_stokes_x = 0;
	E_focus_stokes_y.resize(ps.Nx,ps.Ny,ps.Nz); E_focus_stokes_y = 0;
	E_focus_stokes_z.resize(ps.Nx,ps.Ny,ps.Nz); E_focus_stokes_z = 0;

	cars_pol_x.resize(ps.Nx,ps.Ny,ps.Nz); cars_pol_x = 0;
	cars_pol_y.resize(ps.Nx,ps.Ny,ps.Nz); cars_pol_y = 0;
	cars_pol_z.resize(ps.Nx,ps.Ny,ps.Nz); cars_pol_z = 0;

	E_inc_x.resize(pp.Ntheta,pp.Nphi);
	E_inc_y.resize(pp.Ntheta,pp.Nphi);
	E_inc_z.resize(pp.Ntheta,pp.Nphi);
	E_inc_x = 1.0;
	E_inc_y = 0.0;
	E_inc_z = 0.0;

	cout << "Calculating pump focal field\n~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n";
	focus3d focus_pump(pp);
	focus_pump.set_incident_field(E_inc_x,E_inc_y,E_inc_z);
	focus_pump.calculate_focus(E_focus_pump_x,E_focus_pump_y,E_focus_pump_z);
	focus_pump.save_results(outputdir+"pump");

	E_inc_x.resize(ps.Ntheta,ps.Nphi);
	E_inc_y.resize(ps.Ntheta,ps.Nphi);
	E_inc_z.resize(ps.Ntheta,ps.Nphi);
	E_inc_x = 1.0;
	E_inc_y = 0.0;
	E_inc_z = 0.0;

	cout << "Calculating stokes focal field\n~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n";
	focus3d focus_stokes(ps);
	focus_stokes.set_incident_field(E_inc_x,E_inc_y,E_inc_z);
	focus_stokes.calculate_focus(E_focus_stokes_x,E_focus_stokes_y,E_focus_stokes_z);
	focus_stokes.save_results(outputdir+"stokes");

	complex<double> chi = 2.0;

	cout << "Calculating CARS excitation\n~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n";
	cars_pol_x = chi*E_focus_pump_x*E_focus_pump_x*conj(E_focus_stokes_x);
	cars_pol_y = chi*E_focus_pump_y*E_focus_pump_y*conj(E_focus_stokes_y);
	cars_pol_z = chi*E_focus_pump_z*E_focus_pump_z*conj(E_focus_stokes_z);

	cout << "Saving CARS excitation\n~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n";
	ofstream fout;
	filename.str(string());
	filename << outputdir << "cars_excitation_x.bar";
	cout << filename.str() << endl;
	fout.open(filename.str().c_str());
	if(!fout)
		{	cerr << "File error!! Unable to write " << filename.str() << endl; exit(-1); }
	fout << cars_pol_x;
	fout.close();

	filename.str(string());
	filename << outputdir << "cars_excitation_y.bar";
	cout << filename.str() << endl;
	fout.open(filename.str().c_str());
	if(!fout)
		{	cerr << "File error!! Unable to write " << filename.str() << endl; exit(-1); }
	fout << cars_pol_y;
	fout.close();

	filename.str(string());
	filename << outputdir << "cars_excitation_z.bar";
	cout << filename.str() << endl;
	fout.open(filename.str().c_str());
	if(!fout)
		{	cerr << "File error!! Unable to write " << filename.str() << endl; exit(-1); }
	fout << cars_pol_z;
	fout.close();

	filename.str(string());
	filename << outputdir << "cars_simulation_parameters.dat";
	cout << filename.str() << endl;
	pc.save(filename.str());
	return(0);
	}




