/[MITgcm]/MITgcm_contrib/gael/comm/course-idma2016/matlab/idma_load_fields.m
ViewVC logotype

Contents of /MITgcm_contrib/gael/comm/course-idma2016/matlab/idma_load_fields.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.2 - (show annotations) (download)
Thu Jan 21 19:55:14 2016 UTC (9 years, 6 months ago) by gforget
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +0 -2 lines
- remove return statement

1
2 %%0) setup:
3
4 %add matlab path and load ECCO v4 grid:
5 p = genpath([pwd filesep 'gcmfaces/']); addpath(p);
6 p = genpath([pwd filesep 'MITprof/']); addpath(p);
7 p = genpath([pwd filesep 'm_map/']); addpath(p);
8 grid_load; gcmfaces_global;
9 m=mygrid.mskC(:,:,1);%standard land mask
10
11 %input directories:
12 dirNctiles='release1/nctiles_climatology/';
13 dirBinary='release1/binary_inputs/';
14
15 %constant parameters:
16 parms.rhoconst =1029; %sea water density [in kg/m^3]
17 parms.rhoi = 910; %sea ice density [in kg/m^3]
18 parms.rhosn = 330; %snow density [in kg/m^3]
19
20
21 %%1) sea surface temperature and related variables:
22
23 %THETA: potential temperature [in degree C]
24 %SALT: salinity [in psu]
25 %SIarea: fractional ice-covered area [0 to 1]
26 %MXLDEPTH: Mixed-Layer Depth [in m]
27 %
28 %note: MXLDEPTH is based on the Kara, Rochford, and Hurlburt 2000 formula
29
30 if isdir(dirNctiles);
31 sst_ecco_clim=read_nctiles([dirNctiles 'THETA/THETA'],'THETA',[],1);
32 sss_ecco_clim=read_nctiles([dirNctiles 'SALT/SALT'],'SALT',[],1);
33 mld_ecco_clim=read_nctiles([dirNctiles 'MXLDEPTH/MXLDEPTH'],'MXLDEPTH');
34 ice_ecco_clim=read_nctiles([dirNctiles 'SIarea/SIarea'],'SIarea');
35 end;
36
37 %sst_reynolds: Reynolds monthly mean sea surface temperature maps
38 %iicecover_nsidc: NSIDC monthly mean fractional ice-covered area [0 to 1]
39 %
40 %note: references are provided in Forget et al 2015 (GMD)
41
42 if isdir(dirBinary);
43 sst_reynolds_1992=read_bin([dirBinary 'sst_reynolds/reynolds_oiv2_r1_1992']);
44 sst_reynolds_1992(sst_reynolds_1992==0)=NaN;%0's denote missing values in these files
45 ice_nsidc_1992=read_bin([dirBinary 'icecover_nsidc/nsidc79_monthly_1992']);
46 ice_nsidc_1992(ice_nsidc_1992==-9999)=NaN;%-9999's denote missing values in these files
47 end;
48
49
50 %%2) heat and freshwater fluxes into the ocean (and ice pack)
51
52 %oceFWflx: net surface Fresh-Water flux into the ocean [in kg/m^2/s]
53 %oceQnet: net surface heat flux into the ocean [in W/m^2]
54 %SIatmFW: net freshwater flux from atmosphere & land [in kg/m^2/s]
55 %SIatmQnt: net heat flux to atmosphere & land [in W/m^2]
56 %
57 %note: the key difference between oceQnet and SIatmQnt (or oceFWflx and SIatmFW) is
58 % that the former is computed at the bottom of the ice pack when there is sea-ice
59 %note: the MITgcm sign conventions are opposite for SIatmQnt and oceQnet so below
60 % I flip the sign of SIatmQnt upon reading it from file.
61
62 if isdir(dirNctiles);
63 qnet_to_oce=read_nctiles([dirNctiles 'oceQnet/oceQnet'],'oceQnet');
64 qnet_to_surf=-read_nctiles([dirNctiles 'SIatmQnt/SIatmQnt'],'SIatmQnt');
65 qnet_to_ice=qnet_to_surf-qnet_to_oce;
66
67 fwf_to_oce=read_nctiles([dirNctiles 'oceFWflx/oceFWflx'],'oceFWflx');
68 fwf_to_surf=read_nctiles([dirNctiles 'SIatmFW/SIatmFW'],'SIatmFW');
69 fwf_to_ice=fwf_to_surf-fwf_to_oce;
70 end;
71
72
73 %%3) sea surface height and bottom pressure:
74
75 %ETAN: Surface Height Anomaly [in m]
76 %sIceLoad: sea-ice loading [in km/m^2]
77 %PHIBOT: Bottom Pressure Anomaly [in m^2/s^2]
78 %
79 %note: ETAN corresponds to the ocean-ice interface when there is seaice. The
80 % height contribution of seaice+snow (sIceLoad/parms.rhoconst) therefore
81 % needs to be accounted for to allow comparison with altimetry.
82 %note: bp is converted to equivalent sea level meters by dividing by g.
83
84 if isdir(dirNctiles);
85 ETAN=read_nctiles([dirNctiles 'ETAN/ETAN'],'ETAN');
86 sIceLoad=read_nctiles([dirNctiles 'sIceLoad/sIceLoad'],'sIceLoad');
87 PHIBOT=read_nctiles([dirNctiles 'PHIBOT/PHIBOT'],'PHIBOT');
88
89 ssh=ETAN+sIceLoad/parms.rhoconst;
90 bp=PHIBOT/9.81;
91 end;
92
93 %%4) near surface velocity and related vector variables:
94
95 %UVELMASS: first component of horizontal velocity [in m/s]
96 %VVELMASS: second component of horizontal velocity [in m/s]
97 %
98 %note: the loaded vector components are along grid line directions. To
99 % allow for comparison with e.g. current meters they need to be
100 % converted to zonal and meridional components using calc_UEVNfromUXVY.
101
102 UVELMASS=read_nctiles([dirNctiles 'UVELMASS/UVELMASS'],'UVELMASS',[],1);
103 VVELMASS=read_nctiles([dirNctiles 'VVELMASS/VVELMASS'],'VVELMASS',[],1);
104 oceTAUX=read_nctiles([dirNctiles 'oceTAUX/oceTAUX'],'oceTAUX');
105 oceTAUY=read_nctiles([dirNctiles 'oceTAUY/oceTAUY'],'oceTAUY');
106
107 [Ue,Vn]=calc_UEVNfromUXVY(UVELMASS,VVELMASS);
108 [Te,Tn]=calc_UEVNfromUXVY(oceTAUX,oceTAUY);
109

  ViewVC Help
Powered by ViewVC 1.1.22