utils/cs_grid/rotate_csCg_EN.m

function [uE,vN] = rotate_csCg_EN(u,v)
% [uE,vN] = rotate_csCg_EN(u,v)
%
% Rotate U,V vector components of cs-grid (C-grid) to East,North direction
%  u,v   is a 2-D or 3-D horizontal components of model flow fields.
%  uE,vN are the component of the rotated flow field
% assume that the 1rst 2 dimensions of (u,v) are [6*nc nc] for a CS(nc x nc) grid
%
% >> u=rdmds('uVeltave.0000513360');
% >> v=rdmds('vVeltave.0000513360');
% >> [uE,vN]=rotate_csCg_EN(u,v);
%
% $Header: /u/gcmpack/MITgcm/utils/cs_grid/rotate_csCg_EN.m,v 1.2 2005/08/18 18:26:15 jmc Exp $

%Rac='/home/jmc/grid_cs32/' ;
Rac='grid_files/';

NN=size(u);
if length(NN) < 3, nz=1;
else nz=prod(NN(3:end));
end 
nnx=NN(1); nc=NN(2); nPg=nnx*nc;

if nnx ~= 6*nc,
 fprintf('Error in CS-dim: %i %i %i \n',NN);
 uE=zeros(NN); vN=zeros(NN);
 return
end

u=reshape(u,[nnx nc nz]);
v=reshape(v,[nnx nc nz]);

%- do simple average to put u,v at the center (A-grid):
[uu,vv] = split_UV_cub(u,v);
uu=reshape(uu,[nc+1 nc nz 6]);
vv=reshape(vv,[nc nc+1 nz 6]);

u=(uu(1:nc,:,:,:)+uu(2:nc+1,:,:,:))/2;
v=(vv(:,1:nc,:,:)+vv(:,2:nc+1,:,:))/2;

u=reshape(permute(u,[1 4 2 3]),[nPg nz]);
v=reshape(permute(v,[1 4 2 3]),[nPg nz]);

%- rotate toward E,N (lon,lat) directions :
% load COS & SIN of rotation angle:
namfil=['proj_cs',int2str(nc),'_2uEvN.bin'];
fid=fopen([Rac,namfil],'r','b'); uvEN=fread(fid,nPg*2,'real*8'); fclose(fid);
uvEN=reshape(uvEN,[nPg 2]);

uE=zeros(nPg,nz); vN=zeros(nPg,nz);
for k=1:nz;
 uE(:,k)=uvEN(:,1).*u(:,k)-uvEN(:,2).*v(:,k);
 vN(:,k)=uvEN(:,2).*u(:,k)+uvEN(:,1).*v(:,k);
end

uE=reshape(uE,NN);
vN=reshape(vN,NN);

return
1	function [uE,vN] = rotate_csCg_EN(u,v)
2	% [uE,vN] = rotate_csCg_EN(u,v)
3	%
4	% Rotate U,V vector components of cs-grid (C-grid) to East,North direction
5	% u,v is a 2-D or 3-D horizontal components of model flow fields.
6	% uE,vN are the component of the rotated flow field
7	% assume that the 1rst 2 dimensions of (u,v) are [6*nc nc] for a CS(nc x nc) grid
8	%
9	% >> u=rdmds('uVeltave.0000513360');
10	% >> v=rdmds('vVeltave.0000513360');
11	% >> [uE,vN]=rotate_csCg_EN(u,v);
12	%
13	% $Header: /u/gcmpack/MITgcm/utils/cs_grid/rotate_csCg_EN.m,v 1.2 2005/08/18 18:26:15 jmc Exp $
14
15	%Rac='/home/jmc/grid_cs32/' ;
16	Rac='grid_files/';
17
18	NN=size(u);
19	if length(NN) < 3, nz=1;
20	else nz=prod(NN(3:end));
21	end
22	nnx=NN(1); nc=NN(2); nPg=nnx*nc;
23
24	if nnx ~= 6*nc,
25	fprintf('Error in CS-dim: %i %i %i \n',NN);
26	uE=zeros(NN); vN=zeros(NN);
27	return
28	end
29
30	u=reshape(u,[nnx nc nz]);
31	v=reshape(v,[nnx nc nz]);
32
33	%- do simple average to put u,v at the center (A-grid):
34	[uu,vv] = split_UV_cub(u,v);
35	uu=reshape(uu,[nc+1 nc nz 6]);
36	vv=reshape(vv,[nc nc+1 nz 6]);
37
38	u=(uu(1:nc,:,:,:)+uu(2:nc+1,:,:,:))/2;
39	v=(vv(:,1:nc,:,:)+vv(:,2:nc+1,:,:))/2;
40
41	u=reshape(permute(u,[1 4 2 3]),[nPg nz]);
42	v=reshape(permute(v,[1 4 2 3]),[nPg nz]);
43
44	%- rotate toward E,N (lon,lat) directions :
45	% load COS & SIN of rotation angle:
46	namfil=['proj_cs',int2str(nc),'_2uEvN.bin'];
47	fid=fopen([Rac,namfil],'r','b'); uvEN=fread(fid,nPg2,'real8'); fclose(fid);
48	uvEN=reshape(uvEN,[nPg 2]);
49
50	uE=zeros(nPg,nz); vN=zeros(nPg,nz);
51	for k=1:nz;
52	uE(:,k)=uvEN(:,1).u(:,k)-uvEN(:,2).v(:,k);
53	vN(:,k)=uvEN(:,2).u(:,k)+uvEN(:,1).v(:,k);
54	end
55
56	uE=reshape(uE,NN);
57	vN=reshape(vN,NN);
58
59	return