FT Funcion Gaussiana

Aqui esta la transformada de una funcion Gaussiana de la forma

$f(t)=\exp[-\frac{1}{2\sigma^2}(t-t_0)^2]$

cuya transformada es

$f(\omega)=\sigma \exp(-\frac{\omega^2 \sigma^2}{2}-i\omega t_0)$

[code language=”cpp”]

</pre>
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <iostream>
#include <fstream>
#include <string>
#include <sstream>
#include <vector>
using namespace std;

double pi = 2.0*asin(1.0);
double c_air = 331.5;
double fc = 1000.0; // escogemos una frecuencia central

/////////////////////////////////////////////////////////////////////////
///////////////////////////////////
void imprimir100(int a,double fa[])
{
ofstream myfile100;myfile100.open ("dat100.dat");
for(int i=0; i < a; i++)
{
myfile100 << i << " " << fa[i] << "\n";
}
myfile100.close();
}
///////////////////////////////////

///////////////////////////////////
void FTa(double sigma,double t0,double wc)
{
int iw,Nw=300; double w[Nw+1],wi,wf,dw;
double fwr[Nw+1],fwi[Nw+1];
wi = -1.0*wc;
wf = +1.0*wc;
dw = (wf-wi)/double(Nw);
ofstream myfile;myfile.open ("FTa.dat");
for (iw=0 ; iw<=Nw ; iw++)
{
w[iw] = wi+double(iw)*dw;
fwr[iw] = sigma*exp(-0.5*w[iw]*w[iw]*sigma*sigma)*cos(w[iw]*t0);
fwi[iw] = sigma*exp(-0.5*w[iw]*w[iw]*sigma*sigma)*sin(w[iw]*t0);
myfile << w[iw] << " " << fwr[iw] << " " << fwi[iw] << "\n";
}
myfile.close();
}
///////////////////////////////////

///////////////////////////////////
void FTn(int Nt,double t[],double ft[],double dt,double wc)
{
int iw,Nw=300; double w[Nw+1],wi,wf,dw;
double fwr[Nw+1],fwi[Nw+1];
wi = -1.0*wc;
wf = +1.0*wc;
dw = (wf-wi)/double(Nw);

ofstream myfile;myfile.open ("FTn.dat");
for (iw=0 ; iw<=Nw ; iw++)
{
w[iw] = wi+double(iw)*dw;
fwr[iw] = 0.0;
fwi[iw] = 0.0;
for (int it=0 ; it<=Nt ; it++)
{
fwr[iw] = fwr[iw]+ft[it]*cos(w[iw]*t[it]);
fwi[iw] = fwi[iw]+ft[it]*sin(w[iw]*t[it]);
}
fwr[iw] = (dt/(sqrt(2.0*pi)))*fwr[iw];
fwi[iw] = (dt/(sqrt(2.0*pi)))*fwi[iw];
myfile << w[iw] << " " << fwr[iw] << " " << fwi[iw] << "\n";
}
myfile.close();
}
///////////////////////////////////

main()
{
double tc,wc,t0,sigma;
int it,Nt=100;; double t[Nt+1],ti,tf,dt;
double ft[Nt+1];

tc = 1.0/fc;
wc = (2.0*pi)/tc;
t0 = 2.0*tc;
sigma = 0.5*tc;

ti = 0.0*tc;
tf = +4.0*tc;
dt = (tf-ti)/double(Nt);

for (it=0 ; it<=Nt ; it++)
{
t[it] = ti+double(it)*dt;
ft[it] = exp( -(0.5/(sigma*sigma))*((t[it]-t0)*(t[it]-t0)) );
}
imprimir100(Nt,ft);

FTa(sigma,t0,wc);
FTn(Nt,t,ft,dt,wc);

}
<pre>
[/code]