clear
stem='SO2';
delta=7.89474e-5;
cpnames={'X','Y','Z'};
load('/Users/eebmy/Desktop/FrcBristol/geometre.mat');
load('/Users/eebmy/Desktop/FrcBristol/geometso.mat');
geometre=geometre;
geometso=geometso;
data1=load('/Users/eebmy/Desktop/FrcBristol/S2/SO2.HAS');
[m,n]=size(data1);
st=(n-1)/3; %Number of Stations in the model stem
for rec=1:st
filename=strcat(stem,'_',num2str(rec),'.ASC');
%
%M=[data1(:,[3*rec-1 3*rec 3*rec+1])];
%Rotate the x-y-z in WAVE coordinate to Global Coordinate &
% X=x, Y=-z, Z=y
M=[data1(:,3*rec-1),-1*data1(:,3*rec+1),data1(:,3*rec) ];
% dlmwrite(filename,M)
dlmwrite(filename,M, 'delimiter', ' ','precision',6);
data=load([stem,'_',num2str(rec),'.ASC']);
for i=1:3
[sactrace] = msac_new(data(:,i),delta);
sactrace.e=sactrace.b+((sactrace.npts)-1)*delta;
sactrace.kstnm=[stem,'_',num2str(rec)];
sactrace.kevnm=stem;
sactrace.kcmpnm=cpnames{i};
fname=[stem,'_',num2str(rec),'.',cpnames{i}];
sactrace.stlo = geometre(rec,1);
sactrace.stla = geometre(rec,2);
sactrace.stdp = geometre(rec,3);
sactrace.az = geometre(rec,4);
sactrace.dist = geometre(rec,6);
sactrace.evlo = geometso(1);
sactrace.evla = geometso(2);
sactrace.evdp = geometso(3);
if strcmpi(cpnames{i},'X')
% sactrace.cmpaz = geometre(rec,4); Wrote by James
sactrace.cmpaz = 90
sactrace.cmpinc = 90;
elseif strcmpi(cpnames{i},'Y')
sactrace.cmpaz = 0;
sactrace.cmpinc = 90;
elseif strcmpi(cpnames{i},'Z')
% sactrace.cmpaz = geometre(rec,4) - 90;
sactrace.cmpaz = 0;
sactrace.cmpinc = 0;
end
sactrace.endian='b';
msac_write(fname,sactrace);
end
end