% extract Eta/S/T/U/V/W in Kuroshio Region
% 140-160E, 30-40N

% define desired region
nx=4320;
prec='real*4';
region_name='Kuroshio3';
minlat=30;
maxlat=40;
minlon=140;
maxlon=160;
kx=1:90;

% extract indices for desired region
gdir='/nobackupp2/dmenemen/llc_4320/grid/';
fnam=[gdir 'Depth.data'];
[tmp fc ix jx] = ...
    quikread_llc(fnam,nx,1,prec,gdir,minlat,maxlat,minlon,maxlon);
m=length(ix);
n=length(jx);
quikpcolor(tmp')
thincolorbar

% get and save 2D grid information
close all
pin='~dmenemen/llc_4320/grid/';
pout=['~dmenemen/llc_4320/regions/' region_name '/grid/'];
eval(['mkdir ' pout])
eval(['cd ' pout])
for fnm={'AngleCS','AngleSN','Depth','RAC', ...
         'U2zonDir','V2zonDir','XC','XG','YC','YG'}
    fin=[pin fnm{1} '.data'];
    fout=[fnm{1} '_' int2str(m) 'x' int2str(n)];
    fld=read_llc_fkij(fin,nx,fc,1,ix,jx);
    writebin(fout,fld);
end

% Southwest center points, no direction
fnx={'DXC','RAW'};
fny={'DYC','RAS'};
for i=1:2
    finx=[pin fnx{i} '.data'];
    finy=[pin fny{i} '.data'];
    foutx=[fnx{i} '_' int2str(m) 'x' int2str(n)];
    fouty=[fny{i} '_' int2str(m) 'x' int2str(n)];
    fldx=read_llc_fkij(finy,nx,fc,1,ix,jx);
    fldy=read_llc_fkij(finx,nx,fc,1,ix,jx-1); % <<<<<<<<
    writebin(foutx,fldx);
    writebin(fouty,fldy);
end

% Southwest edges, no direction
fnx='DXG';
fny='DYG';
finx=[pin fnx '.data'];
finy=[pin fny '.data'];
foutx=[fnx '_' int2str(m) 'x' int2str(n)];
fouty=[fny '_' int2str(m) 'x' int2str(n)];
fldx=read_llc_fkij(finy,nx,fc,1,ix,jx-1); % <<<<<<<<
fldy=read_llc_fkij(finx,nx,fc,1,ix,jx);
writebin(foutx,fldx);
writebin(fouty,fldy);

% Southwest corner (vorticity) points, no direction, fill blank tiles
fnm='RAZ';
fin=[pin fnm '.data'];
finy=[pin fny '.data'];
fout=[fnm '_' int2str(m) 'x' int2str(n)];
fld=read_llc_fkij(fin,nx,fc,1,ix,jx-1); % <<<<<<<<
writebin(fout,fld);

% get and save 3D grid information
fnm={'hFacC'}
fin=[pin fnm{1} '.data'];
fout=[fnm{1} '_' int2str(m) 'x' int2str(n) 'x' int2str(length(kx))];
for k=1:length(kx); mydisp(k)
    fld=read_llc_fkij(fin,nx,fc,kx(k),ix,jx);
    writebin(fout,fld,1,'real*4',k-1);
end

% Southwest center points, no direction
fnx='hFacW';
fny='hFacS';
finx=[pin fnx '.data'];
finy=[pin fny '.data'];
foutx=[fnx '_' int2str(m) 'x' int2str(n) 'x' int2str(length(kx))];
fouty=[fny '_' int2str(m) 'x' int2str(n) 'x' int2str(length(kx))];
for k=1:length(kx); mydisp(k)
    fldx=read_llc_fkij(finy,nx,fc,kx(k),ix,jx);
    fldy=read_llc_fkij(finx,nx,fc,kx(k),ix,jx-1); % <<<<<<<<
    writebin(foutx,fldx,1,'real*4',k-1);
    writebin(fouty,fldy,1,'real*4',k-1);
end

% get and save 2D tracer field
pin='~dmenemen/llc_4320/MITgcm/run/';
pout=['~dmenemen/llc_4320/regions/' region_name '/'];
for fnm={'Eta'}
    eval(['mkdir ' pout fnm{1}])
    eval(['cd ' pout fnm{1}])
    for ts=10368:144:485568, disp(ts)
        fin=[pin fnm{1} '.' myint2str(ts,10) '.data'];
        dy=ts2dte(ts,25,2011,9,10,30);
        fout=[fnm{1} '_' int2str(m) 'x' int2str(n) '.' dy];
        fld=read_llc_fkij(fin,nx,fc,1,ix,jx);
        writebin(fout,fld);
    end
end

% get and save 3D tracer fields
for fnm={'Salt','Theta','W'}
    eval(['mkdir ' pout fnm{1}])
    eval(['cd ' pout fnm{1}])
    for ts=10368:144:485568, disp(ts)
        fin=[pin fnm{1} '.' myint2str(ts,10) '.data'];
        dy=ts2dte(ts,25,2011,9,10,30);
        fout=[fnm{1} '_' int2str(m) 'x' int2str(n) 'x' int2str(length(kx)) '.' dy];
        for k=1:length(kx); mydisp(k)
            fld=read_llc_fkij(fin,nx,fc,kx(k),ix,jx);
            writebin(fout,fld,1,'real*4',k-1);
        end
    end
end

% get and save 3D vector fields
eval(['mkdir ' pout 'U'])
eval(['mkdir ' pout 'V'])
eval(['cd ' pout])
for ts=10368:144:485568, disp(ts)
    finu=[pin 'U.' myint2str(ts,10) '.data'];
    finv=[pin 'V.' myint2str(ts,10) '.data'];
    dy=ts2dte(ts,25,2011,9,10,30);
    foutu=['U/U_' int2str(m) 'x' int2str(n) 'x' int2str(length(kx)) '.' dy];
    foutv=['V/V_' int2str(m) 'x' int2str(n) 'x' int2str(length(kx)) '.' dy];
    for k=1:length(kx);
        fldu= read_llc_fkij(finv,nx,fc,kx(k),ix,jx);
        fldv=-read_llc_fkij(finu,nx,fc,kx(k),ix,jx-1); % <<<<<<<<
        writebin(foutu,fldu,1,'real*4',k-1);
        writebin(foutv,fldv,1,'real*4',k-1);
    end
end