%                togrid4.m
%               for constructing ascii flat file with lon,lat,u(ind),v(ind)
%               where each of the variables is subsampled by 2

load ghfactor.cal
im=sqrt(-1);
ghval=[1 2 3 4 5 6 7 9];
umean=[];
for k=1:8
     kk=ghval(k);
	disp(kk)
	uarray=[]; varray=[];
	eval(concat('load gh',int2str(ghval(k)),'cof.tid'))
	eval(concat('load gh',int2str(ghval(k)),'nav.bin'))
	for i=1:7
		uarray(:,i)=subsamp(2,2,u0coef(:,i));
		varray(:,i)=subsamp(2,2,v0coef(:,i));
	end
	depth1=subsamp(2,2,depthbin);
	lat=subsamp(2,2,latbin(:));
	lon=subsamp(2,2,lonbin(:));
%               rotate to common phase
	[uarray(:,2),uarray(:,3)]=rot(uarray(:,2),uarray(:,3),-phase(kk));
	[uarray(:,4),uarray(:,5)]=rot(uarray(:,4),uarray(:,5),-phase(kk)*2);
	[uarray(:,6),uarray(:,7)]=rot(uarray(:,6),uarray(:,7),-phase(kk)*3);
	[varray(:,2),varray(:,3)]=rot(varray(:,2),varray(:,3),-phase(kk));
	[varray(:,4),varray(:,5)]=rot(varray(:,4),varray(:,5),-phase(kk)*2);
	[varray(:,6),varray(:,7)]=rot(varray(:,6),varray(:,7),-phase(kk)*3);
%               adjust amplitude

	 uarray=uarray/mag(kk);
      varray=varray/mag(kk);
	
	[i,j]=size(umean)
	[inew,j]=size(uarray);
     ind=[i+1:i+inew];

	latlon(ind,1)=lon;
	latlon(ind,2)=lat;
	depth(ind)=depth1;

	ucoef(ind,:)=uarray;
	vcoef(ind,:)=varray;

	umean(ind,3)=uarray(:,1);
	vmean(ind,3)=varray(:,1);

     um2(ind,3)=uarray(:,2);
     um2s(ind,3)=uarray(:,3);
     vm2(ind,3)=varray(:,2);
     vm2s(ind,3)=varray(:,3);

     um4(ind,3)=uarray(:,4);
     um4s(ind,3)=uarray(:,5);
     vm4(ind,3)=varray(:,4);
     vm4s(ind,3)=varray(:,5);

     um6(ind,3)=uarray(:,6);
     um6s(ind,3)=uarray(:,7);
     vm6(ind,3)=varray(:,6);
     vm6s(ind,3)=varray(:,7);
end

umean(:,1:2)=latlon;
vmean(:,1:2)=latlon;

um2(:,1:2)=latlon;
vm2(:,1:2)=latlon;
um4(:,1:2)=latlon;
vm4(:,1:2)=latlon;
um6(:,1:2)=latlon;
vm6(:,1:2)=latlon;

um2s(:,1:2)=latlon;
vm2s(:,1:2)=latlon;
um4s(:,1:2)=latlon;
vm4s(:,1:2)=latlon;
um6s(:,1:2)=latlon;
vm6s(:,1:2)=latlon;

%............................do transports

qumean(:,3)=umean(:,3).*depth(:);
qvmean(:,3)=vmean(:,3).*depth(:);

qum2(:,3)=um2(:,3).*depth(:);
qvm2(:,3)=vm2(:,3).*depth(:);
qum4(:,3)=um4(:,3).*depth(:);
qvm4(:,3)=vm4(:,3).*depth(:);
qum6(:,3)=um6(:,3).*depth(:);
qvm6(:,3)=vm6(:,3).*depth(:);

qum2s(:,3)=um2s(:,3).*depth(:);
qvm2s(:,3)=vm2s(:,3).*depth(:);
qum4s(:,3)=um4s(:,3).*depth(:);
qvm4s(:,3)=vm4s(:,3).*depth(:);
qum6s(:,3)=um6s(:,3).*depth(:);
qvm6s(:,3)=vm6s(:,3).*depth(:);


qumean(:,1:2)=latlon;
qvmean(:,1:2)=latlon;

qum2(:,1:2)=latlon;
qvm2(:,1:2)=latlon;
qum4(:,1:2)=latlon;
qvm4(:,1:2)=latlon;
qum6(:,1:2)=latlon;
qvm6(:,1:2)=latlon;

qum2s(:,1:2)=latlon;
qvm2s(:,1:2)=latlon;
qum4s(:,1:2)=latlon;
qvm4s(:,1:2)=latlon;
qum6s(:,1:2)=latlon;
qvm6s(:,1:2)=latlon;

%.........................................get rid of bad points
bad=find(um2(:,3)==nan);

umean(bad,:)=[];
vmean(bad,:)=[];

um2(bad,:)=[];
vm2(bad,:)=[];
um4(bad,:)=[];
vm4(bad,:)=[];
um6(bad,:)=[];
vm6(bad,:)=[];

um2s(bad,:)=[];
vm2s(bad,:)=[];
um4s(bad,:)=[];
vm4s(bad,:)=[];
um6s(bad,:)=[];
vm6s(bad,:)=[];

umean(bad,:)=[];
vmean(bad,:)=[];

qum2(bad,:)=[];
qvm2(bad,:)=[];
qum4(bad,:)=[];
qvm4(bad,:)=[];
qum6(bad,:)=[];
qvm6(bad,:)=[];

qum2s(bad,:)=[];
qvm2s(bad,:)=[];
qum4s(bad,:)=[];
qvm4s(bad,:)=[];
qum6s(bad,:)=[];
qvm6s(bad,:)=[];

save latlon.all latlon depth
save coefs.all ucoef vcoef

%...............................save velocity files

save umean.all umean
save vmean.all vmean

save um2.all um2
save vm2.all vm2
save um4.all um4
save vm4.all vm4
save um6.all um6
save vm6.all vm6
save um2s.all um2s
save vm2s.all vm2s
save um4s.all um4s
save vm4s.all vm4s
save um6s.all um6s
save vm6s.all vm6s

%...............................save all the tranport files
save qum2.all qum2
save qvm2.all qvm2
save qum4.all qum4
save qvm4.all qvm4
save qum6.all qum6
save qvm6.all qvm6
save qum2s.all qum2s
save qvm2s.all qvm2s
save qum4s.all qum4s
save qvm4s.all qvm4s
save qum6s.all qum6s
save qvm6s.all qvm6s
