#include"Pulsetrain.h"
#include"Constants.h"
#include<iostream>
#include<math.h>
#include<vector>

Pulsetrain::Pulsetrain() {}
Pulsetrain::~Pulsetrain() {}


Pulsetrain::Pulsetrain(double lambda_SI, double I0_SI, std::vector<double> t0_SI, double phi_SI, double FWHM_SI, double t0_env, double FWHM_env) {
	this->lambda = lambda_SI * 1e-9 / au_L;
	this->omega = 2 * M_PI * au_c / this->lambda;
	this->E0 = sqrt(2.0 * I0_SI * 1e4 / (const_epsilon_0 * const_c)) / au_F;
	for (auto& entry : t0_SI) {
		t0.push_back(entry / (au_T * 1e15));
	}
	this->FWHM = FWHM_SI / (au_T * 1e15);
	this->tau = this->FWHM / (2 * sqrt(log(2.0)));
	this->phi = phi_SI;
	this->t0_envelope = t0_env / (au_T * 1e15);
	this->tau_envelope = FWHM_env / (au_T * 1e15 * 2 * sqrt(log(2.0)));
}

void Pulsetrain::set(double lambda_SI, double I0_SI, std::vector<double> t0_SI, double phi_SI, double FWHM_SI, double t0_env, double FWHM_env) {
	this->lambda = lambda_SI * 1e-9 / au_L;
	this->omega = 2 * M_PI * au_c / this->lambda;
	this->E0 = sqrt(2.0 * I0_SI * 1e4 / (const_epsilon_0 * const_c)) / au_F;
	for (auto& entry : t0_SI) {
		t0.push_back(entry / (au_T * 1e15));
	}
	this->FWHM = FWHM_SI / (au_T * 1e15);
	this->tau = this->FWHM / (2 * sqrt(log(2.0)));
	this->phi = phi_SI;
	this->t0_envelope = t0_env / (au_T * 1e15);
	this->tau_envelope = FWHM_env / (au_T * 1e15 * 2 * sqrt(log(2.0)));
}

double Pulsetrain::getField(double t_au) {
	if (!t0.empty()) {
		double E = 0;
		for (auto& entry : t0) {
			E += E0 * cos(omega * (t_au - entry) + phi) * exp(-pow(t_au - entry, 2.0) / (2.0 * pow(tau, 2.0)));
		}
		return E * exp(-pow(t_au - t0_envelope, 2.0) / (2.0 * pow(tau_envelope, 2.0)));
	}
	else return 0;
}

