clc;clear all;close all;

pathi='/data3/yiwen/cesm2/JGR_ocean/Model_output/CTRL/';
file=[pathi,'TSRyr_000101-006112.nc']
LAT=ncread(file,'TLAT');lat=LAT(1,:);        %(-79~72)
LON=ncread(file,'TLONG');lon=LON(:,1);       %(320~319)
depth=ncread(file,'z_t');depth_m=depth./100; %unit=cm

salt=ncread(file,'SALT');temp=ncread(file,'TEMP');
%========== Model output ===========%
filename="/data3/yiwen/cesm2/JGR_ocean/Model_output/G_JRA55_test1.pop.h.0031-04.nc"
depth=ncread(filename,'z_t');depth_km=depth./1e5; 
%------- 1*1 -------%
lat_real=-80:1:90;
lon_real=0:1:360;
[xq,yq] = meshgrid(lon_real,lat_real);

n=52:61; %last-10yrs data
salt_M=squeeze(nanmean(salt(:,:,:,n),4));salt_M(isnan(salt_M))=0;
temp_M=squeeze(nanmean(temp(:,:,:,n),4));

V={salt_M temp_M}
for ii=1:size(V,2)
 for k=1:60
  V_N{ii}(:,:,k)=griddata(LON,LAT,V{ii}(:,:,k),xq,yq);
 end
  V_N2{ii}=circshift(V_N{ii},-180,2);
  V_N2{ii}(:,181:183,:)=[];
end
%====(1) W to E ====%
lon2=linspace(-180,180,size(V_N2{1},2));
lonIN=find(lon2>=-60 & lon2<=20);
% mask=squeeze(nanmean(V_N2{2}(:,:,:),2));   %Global zonal mean 
mask=squeeze(nanmean(V_N2{2}(:,lonIN,:),2)); %GIN region
mask(~isnan(mask))=1;


patho='/data3/yiwen/cesm2/JGR_ocean/Model_output/';
fn2 = [patho 'TSD_lat_mask_1X1'];
save(fn2,'mask')

