integra.m

德国人开发的地震处理分析软件

function integra

global grafArange;
global nlbgraf;
global ngraf;
global A;
global fc;
global numero;
global deltat;
global archivo;
global letras;
global linutil;
global vacio;
global MaxA;
global directorio;
global directorio1;

izq=numero(nlbgraf);
der=numero(nlbgraf+1)-1;
Nn=der-izq+1;

%cd d:\matlab\toolbox\prima
cd (directorio1);

fid=fopen('Acel.dat','wt');
fprintf(fid,'%12.4e\n',A(izq:der));
fclose(fid);

dt=deltat(nlbgraf);

fid=fopen('Sdat.dat','wt');
fprintf(fid,'%8i%8.3f',Nn,dt);
fclose(fid);

!integra

load Avel.sal;

delete Acel.dat
delete Sdat.dat
delete Avel.sal

Amax = max(max(Avel),abs(min(Avel)));

lineas=0;
letras = [letras;letras(nlbgraf)];
nombre=[archivo(nlbgraf,1:letras(nlbgraf)-3) 'vel' vacio(letras(nlbgraf)+1:80)];
delta=dt;

numero   =  [numero;numero(length(numero))+Nn];
linutil  =  [linutil;lineas];
deltat   =  [deltat;delta];
archivo  =  [archivo;nombre];
MaxA     =  [MaxA;Amax];
A        =  [A;Avel];

Nfft=2^(fix(log(Nn)/log(2))+1);
df=1/deltat(nlbgraf)/Nfft;
Nceros=fix(fc/df)+1;
filtro=[zeros(Nceros+1,1);ones(Nfft-2*Nceros-1,1);zeros(Nceros,1)];

nvgraf=length(deltat);
izq=numero(nvgraf);
der=numero(nvgraf+1)-1;

x = A(izq:der);
y = detrend(x);
yfft = fft(y,Nfft);

yfft=yfft.*filtro;
z=ifft(yfft);
A(izq:der)=real(z(1:Nn));

Amax = max(max(A(numero(nvgraf):numero(nvgraf+1)-1)),abs(min(A(numero(nvgraf):numero(nvgraf+1)-1))));
MaxA(nvgraf) = Amax;

fid=fopen('Acel.dat','wt');
fprintf(fid,'%12.4e\n',A(izq:der));
fclose(fid);

fid=fopen('Sdat.dat','wt');
fprintf(fid,'%8i%8.3f',Nn,dt);
fclose(fid);

!integra

load Avel.sal;

delete Acel.dat
delete Sdat.dat
delete Avel.sal

%xi=[0:Nn-1]';
%a=sum(Avel.*xi)*6/(Nn*(Nn-1)*(2*Nn-1));
%Avel=Avel-a*xi;

Amax = max(max(Avel),abs(min(Avel)));

lineas=0;
letras = [letras;letras(nlbgraf)];
nombre=[archivo(nlbgraf,1:letras(nlbgraf)-3) 'des' vacio(letras(nlbgraf)+1:80)];
delta=dt;

numero   =  [numero;numero(length(numero))+Nn];
linutil  =  [linutil;lineas];
deltat   =  [deltat;delta];
archivo  =  [archivo;nombre];
MaxA     =  [MaxA;Amax];
A        =  [A;Avel];

nvgraf=length(deltat);
izq=numero(nvgraf);
der=numero(nvgraf+1)-1;

x = A(izq:der);
y = detrend(x);
yfft = fft(y,Nfft);

yfft=yfft.*filtro;
z=ifft(yfft);
A(izq:der)=real(z(1:Nn));

Amax = max(max(A(numero(nvgraf):numero(nvgraf+1)-1)),abs(min(A(numero(nvgraf):numero(nvgraf+1)-1))));
MaxA(nvgraf) = Amax;

ngraf=length(deltat);
grafArange=[-Amax Amax];

cd (directorio);

grafica;