princal.m

地震资料处理的matlab源码软件包

function princal

global fileaxes;
global numero;
global deltat;
global compX;
global compY;
global compZ;
global Tini;
global Tfin;
global A;

axes(fileaxes);hold off;
fill([0 0 1 1 0],[0 1 1 0 0],'k');
axis('off');

Ni=fix(Tini/deltat(compX));
Nf=fix(Tfin/deltat(compX));

if compZ==0
   Cova=cov(A(numero(compX)+Ni:numero(compX)+Nf),A(numero(compY)+Ni:numero(compY)+Nf))
   [V,D]=eig(Cova);
   [Dmax,Ind]=max(diag(D));
   AngH=atan(V(2,Ind)/V(1,Ind))*180/pi;
   AngV=0;
else
   Cova(1,1)=sum(A(numero(compX)+Ni:numero(compX)+Nf).*A(numero(compX)+Ni:numero(compX)+Nf));
   Cova(1,2)=sum(A(numero(compX)+Ni:numero(compX)+Nf).*A(numero(compY)+Ni:numero(compY)+Nf));
   Cova(1,3)=sum(A(numero(compX)+Ni:numero(compX)+Nf).*A(numero(compZ)+Ni:numero(compZ)+Nf));
   Cova(2,2)=sum(A(numero(compY)+Ni:numero(compY)+Nf).*A(numero(compY)+Ni:numero(compY)+Nf));
   Cova(2,3)=sum(A(numero(compY)+Ni:numero(compY)+Nf).*A(numero(compZ)+Ni:numero(compZ)+Nf));
   Cova(3,3)=sum(A(numero(compZ)+Ni:numero(compZ)+Nf).*A(numero(compZ)+Ni:numero(compZ)+Nf));
   Cova(2,1)=Cova(1,2);
   Cova(3,1)=Cova(1,3);
   Cova(3,2)=Cova(2,3);

   [V,D]=eig(Cova);
   [Dmax,Ind]=max(diag(D));
   AngH=atan(V(2,Ind)/V(1,Ind))*180/pi;
   AngV=atan(V(3,Ind)/sqrt(V(1,Ind)*V(1,Ind)+V(2,Ind)*V(2,Ind)))*180/pi;
end

text('Position',[0.0 0.20],'Units','normalized','Color','y','String',num2str(AngH));
text('Position',[0.0 0.00],'Units','normalized','Color','y','String',num2str(AngV));
hold off;