/[MITgcm]/MITgcm_contrib/gael/profilesMatlabProcessing/profiles_misc/model_interp.m
ViewVC logotype

Annotation of /MITgcm_contrib/gael/profilesMatlabProcessing/profiles_misc/model_interp.m

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


Revision 1.1 - (hide annotations) (download)
Wed Apr 13 21:16:01 2011 UTC (14 years, 3 months ago) by roquet
Branch: MAIN
new capability: interpolation and storage of collocated model profiles
in a MITprof netcdf file. Successfully tested with an annual T/S climatology of SOSE59.

1 roquet 1.1 function [MITprofCur]=model_interp(MITprofCur,model,varargin)
2     % function [MITprofCur]=model_interp(MITprofCur,model,varargin)
3     % interpolate model hydrographic fields at the location of given profiles
4     %
5     % MITprofCur: structure containing profile information
6     % model is a string used to select the model to be loaded
7     % 'OCCA' : ECCOv4 grid + OCCA atlas
8     % 'SOSE59' : SOSE59 grid + atlas
9    
10     % load model
11     [grid,atlas]=model_load(model);
12    
13     % protect global grid variables
14     global mygrid mytri MYBASININDEX
15     mygrid1=mygrid; mytri1=mytri; MYBASININDEX1=MYBASININDEX;
16     mygrid=grid.mygrid; mytri=grid.mytri; MYBASININDEX=grid.MYBASININDEX;
17    
18     % compute index of closest profile
19     MITprofCur.prof_point=gcmfaces_bindata(MITprofCur.prof_lon,MITprofCur.prof_lat);
20     [MITprofCur.ii,MITprofCur.jj]=gcmfaces_bindata(MITprofCur.prof_lon,MITprofCur.prof_lat);
21    
22     % temporal index
23     switch model
24     case 'SOSE59'
25     MITprofCur.imonth=MITprofCur.ii*0+1;
26     case 'OCCA'
27     MITprofCur.imonth=floor(MITprofCur.prof_YYYYMMDD/1e2)-1e2*floor(MITprofCur.prof_YYYYMMDD/1e4);
28     otherwise
29     error('not a valid model string');
30     end
31    
32     % extract profiles
33     nk=length(mygrid.RC);kk=ones(1,nk);
34     np=length(MITprofCur.ii);pp=ones(np,1);
35     ind2prof=sub2ind(size(atlas.T{1}),MITprofCur.ii*kk,MITprofCur.jj*kk,pp*[1:nk],MITprofCur.imonth*kk);
36     Tatlas=atlas.T{1}(ind2prof); Tatlas=interp1(-mygrid.RC',Tatlas',MITprofCur.prof_depth)';
37     Satlas=atlas.S{1}(ind2prof); Satlas=interp1(-mygrid.RC',Satlas',MITprofCur.prof_depth)';
38    
39     % record values
40     MITprofCur=setfield(MITprofCur,['prof_T_' model],Tatlas');
41     MITprofCur=setfield(MITprofCur,['prof_S_' model],Satlas');
42    
43     % put back global grid values
44     mygrid=mygrid1; mytri=mytri1; MYBASININDEX=MYBASININDEX1;
45    

  ViewVC Help
Powered by ViewVC 1.1.22