#include <iostream.h>
#include <fstream.h>

// D. Yvon, CE Saclay, DAPNIA/SPP, 08/99

#include <math.h>
#include "mayer_fft.h"
#include "dyffttools.h"


#include "transfelec.h"
#include "formepulse.h"

#include "ssthardware.h"

#define LTab (2048)
#define LSave (256)
#define Tpas (5.e-4)
//float PulseOptOrigin=LTab/16.*Tpas;
float PulseOptOrigin=0.;

//#define DYTestRun

/*
FormePulse::FormePulse() {
 	gainElec=0.;
	FullwidthHalfMax=0;
	ShapeArr=new double[LSave];
	for(int i=0; i<LSave; i++)  ShapeArr[i]=0.;
}
*/

FormePulse::FormePulse(TransFuncElec TF) {

	Dcomplex j(0.,1.);
	double Fpas=1./(Tpas*LTab);
	
	double *PulseInReel=new double[LTab];
	double *PulseInImg=new double[LTab];
	double *PulseOut=new double[LTab];

	Dcomplex* PulseFFT=new Dcomplex[LTab];
	

	for(long i=0; i<LTab; i++) PulseInReel[i]=pulseOptique(i*Tpas);
	for(long i=0; i<LTab; i++) PulseInImg[i]=0.;
	
/*
//	On fenetre le signal.
	for(long i=0; i<LTab; i++) PulseInReel[i]=PulseInReel[i]*DYWindow(i,LTab);
*/

	ofstream dlout; // sert a sortir les print test
	
#ifdef DYTestRun
		//On ecrit le signal entrant
		char s[32]="WaveIn.txt";
		dlout.open(s,ios::out | ios::trunc);
		if (!dlout.is_open())
		{	cerr<<"erreur a l'ouverture du fichier de resultats: "<<s<<endl;
			exit(-1);
		}
		dlout<<"temps"<<"\tPartieReelle"<<"\tPartieImag"<<endl;
		
		for (long k=0; k<LTab; k++) {
			dlout<<k*Tpas<<"\t"<<PulseInReel[k]<<"\t"<<PulseInImg[k]<<endl;
		}
		dlout.close();

#endif

// C'est la fft!
	fft(LTab, PulseInReel, PulseInImg);


// On Filtre dans l'espace des frequences
	for (int i=0; i<LTab; i++) {
		PulseFFT[i]=PulseInReel[i]+PulseInImg[i]*j;
	//On applique le filtre Iici
		PulseFFT[i]=TF.TsFunc(2*3.1416*i*Fpas)*PulseFFT[i];
		PulseInReel[i]=real(PulseFFT[i]);
		PulseInImg[i]=imag(PulseFFT[i]);
	}
	
#ifdef DYTestRun	
	// On écrit la FFT	
		char s1[32]="PulseFFT.txt";
		dlout.open(s1,ios::out | ios::trunc);
		if (!dlout.is_open())
		{	cerr<<"erreur a l'ouverture du fichier de resultats: "<<s1<<endl;
			exit(-1);
		}
		dlout<<"Frequence"<<"\tPartieReelle"<<"\tPartieImag"<<endl;
		
		for (long k=0; k<LTab; k++) {
			dlout<<k*Fpas<<"\t"<<PulseInReel[k]<<"\t"<<PulseInImg[k]<<endl;
		}
		dlout.close();
#endif

// FFT inverse
	ifft(LTab,PulseInReel,PulseInImg);
	for (int i=0; i<LTab; i++) {
		PulseInReel[i]/=LTab;
		PulseInImg[i]/=LTab;
	}
	

// On recherche le maximum d'amplitude
	double Max=0.;		//En volt
	int MaxIndex=0;
	
	for(int i=0; i<LSave; i++) 
		if(PulseInReel[i]>Max) {
			Max=PulseInReel[i];
			MaxIndex=i;
		}
		
	gainElec=Max/SSTFPulseCourant;
		//En Volt/ampere	
		
// On calcule la LMH

	int HalfMaxInf;
		for(HalfMaxInf=MaxIndex;HalfMaxInf>0;HalfMaxInf--)
			if(fabs(PulseInReel[HalfMaxInf])<(Max/2.)) break;
	int HalfMaxSup;	
		for(HalfMaxSup=MaxIndex;HalfMaxSup<LTab;HalfMaxSup++)
			if(fabs(PulseInReel[HalfMaxSup])<(Max/2.)) break;
	
	
	FullwidthHalfMax=(HalfMaxSup-HalfMaxInf)*Tpas;
		
// On normalise à un Volt. Et on ecrit en memoire le resultat
	ShapeArr=new double[LSave];
	for(int i=0; i<LSave; i++)  ShapeArr[i]=PulseInReel[i]/Max;

//#ifdef DYTestRun	
		//On ecrit le resultat dans un fichier 
		char s2[32]="Resultatfiltrage.txt";
		dlout.open(s2,ios::out | ios::trunc);
	//	if (!dlout.is_open()){	
	//		cerr<<"erreur a l'ouverture du fichier de resultats: "<<s2<<endl;
	//		exit(-1);
	//	}

		dlout<<"Impulsion courant de: "<<SSTFPulseCourant<<" Ampères"<<endl;
		dlout<<"Maximum d'amplitude: "<<Max<< " Volts"<<endl;
		dlout<<" gain chaine preampli: "<<gainElec<<endl;
		dlout<<"largueur a mi hauteur: "<<FullwidthHalfMax<<endl<<endl;
		dlout<<"Temps"<<"\tPartie Reelle"<<endl;

		for (int k=0; k<LSave; k++) {
			dlout<<k*Tpas<<"\t"<<ShapeArr[k]<<endl;
		}	
//#endif


	
//On tue les tableaux de passage....
	delete[] PulseInReel;
	delete[] PulseInImg;
	delete[] PulseOut;
	delete[] PulseFFT;

}

FormePulse& FormePulse::operator =(const FormePulse& FPlse) {
	ShapeArr=new double[LSave];
	for(int i=0; i<LSave; i++) ShapeArr[i]=FPlse.ShapeArr[i];
	return *this;
}

FormePulse::~FormePulse() {
	delete[] ShapeArr;
}

double FormePulse::TMax() {
	double returnVal=LSave*Tpas;
	return returnVal;
}

double FormePulse::PulseShape(double temps) {
	double returnValue;
	int ientier=(int)(temps/Tpas);
	if (ientier<0) returnValue=0.;
	else if (ientier>=LSave) returnValue= ShapeArr[LSave-1];
	else {
		double Min=ShapeArr[ientier];
		double Max=ShapeArr[ientier+1];
//		double test=temps/Tpas-ientier;
		returnValue=Min+(Max-Min)*(temps/Tpas-ientier);
	}
	return returnValue;
}


double FormePulse::LongPulseOptique() {
	double returnVal;
	returnVal=SSTLongPhotDiode/SSTLongFocale*PeriodRotTeles/2/3.1415/sin(SSTTelesElevation);
	//en secondes
	return returnVal;
//	return returnVal=30.e-3;			//Gros BUGG!!!!!
}
	
double FormePulse::pulseOptique(double temps) {
	double value;
	if((temps>PulseOptOrigin)&&(temps<=PulseOptOrigin+LongPulseOptique()))  
		value=SSTFPulseCourant;			 //	Amplitude: 1 nanoampère

	else value=0.;
//	for(long i=0; i<LTab; i++) value=sin(2*3.1416*temps/5.e-3); //pour test FFT
	return value;
}