%                togrid.m
%               for constructing ascii flat file with lon,lat,u(ind),v(ind)
%               where each of the variables is subsampled by 2
%               For use with multiple depth data


clear
load c:\rocky\gayhead\ghdata\ghfactor.cal
zval=input('Depth level to choose (1 to 5) or 0 for mean ');
im=sqrt(-1);
ghval=[1 2 3 4 5 6 7 9];
latlon=[];
!ncd ghdata
for k=1:8
     kk=ghval(k);
	disp(kk)
	uarray=[]; varray=[];
	if zval>0
            	eval(concat('load gh',int2str(ghval(k)),'va.tid'))
		     u0coef=r3d(uacoef,nrow,0,0,zval);
		     v0coef=r3d(vacoef,nrow,0,0,zval);
     else
   			eval(concat('load gh',int2str(ghval(k)),'cof.tid'))
     end
	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(:));
%:::::::::::::::::::::::::::::::::::   Phase and amp corrections ::::::::::::

%               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(latlon);	
	[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;

end
uout(:,1:2)=latlon;
uout(:,3:9)=ucoef;
vout(:,1:2)=latlon;
vout(:,3:9)=vcoef;

uout=denan(uout);
vout=denan(vout);

%.............................   into surf subdirectory
!ncd ghsurf

save uout uout
save vout vout

%..........................................name

% saved for use in readgrd.m

name=(['ghcoef',int2str(zval),'.grd']);
save name.grd name
disp('get out of matlab, and run gridm, and it will get you back in ')

