/[MITgcm]/MITgcm_contrib/gmaze_pv/eg_write_bin2cdf_latlongrid_subdomain.m
ViewVC logotype

Annotation of /MITgcm_contrib/gmaze_pv/eg_write_bin2cdf_latlongrid_subdomain.m

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


Revision 1.3 - (hide annotations) (download)
Wed Sep 19 15:37:38 2007 UTC (18 years, 3 months ago) by gmaze
Branch: MAIN
CVS Tags: HEAD
Changes since 1.2: +0 -0 lines
General Update

1 gmaze 1.1 % Script to extract and write in netCDF format a subdomain
2     % from the 1.8 global run
3     %
4     clear
5     global sla
6     pv_checkpath
7    
8     % Load grid
9     GRID_125
10    
11     % Load list of all outputs
12     otab = latlon8grid_outputs_table;
13    
14     % Setup standard grid variables
15     lon_c = lon125;
16     lon_u = [lon125(1)-360+lon125(end) (lon125(2:end)+lon125(1:end-1))/2];
17     lat_c = lat125;
18     lat_v = [lat125(1)-(lat125(2)-lat125(1))/2 (lat125(1:end-1)+lat125(2:end))/2];
19     z_c = (cumsum(thk125)-thk125/2);
20     z_w = [0 cumsum(thk125(1:end-1))];
21    
22    
23     % Set subrange - Longitude given as degrees east
24     subdomain = 4;
25    
26     switch subdomain
27     case 1
28     sub_name = 'western_north_atlantic';
29     lonmin = lon125(2209)-180;
30     lonmax = lon125(2497-1)-180;
31     latmin = lat125(1225);
32     latmax = lat125(1497-1);
33     depmin = min(z_w);
34     depmax = z_c(29);
35     m_proj('mercator','long',[270 365],'lat',[0 60]);
36     %clf;hold on;m_coast;m_grid;
37     LIMITS = [lonmin+180 lonmax+180 latmin latmax depmin depmax]
38     %m_line(LIMITS([1 2 2 1 1]),LIMITS([3 3 4 4 3]),'color','r','linewidth',2);
39     %title(sub_name);
40    
41     case 3
42     sub_name = 'north_atlantic';
43     lonmin = lon125(2209)-180;
44     lonmax = lon125(2881-1)-180;
45     latmin = lat125(1157);
46     latmax = lat125(1565-1);
47     depmin = min(z_w);
48     depmax = z_c(29);
49     m_proj('mercator','long',[270 365],'lat',[0 60]);
50     clf;hold on;m_coast;m_grid;
51     LIMITS = [lonmin+180 lonmax+180 latmin latmax depmin depmax]
52     m_line(LIMITS([1 2 2 1 1]),LIMITS([3 3 4 4 3]),'color','k','linewidth',2);
53     title(sub_name);
54    
55     case 4
56     sub_name = 'global';
57     lonmin = lon125(1)-180;
58     lonmax = lon125(2881-1)-180;
59     latmin = lat125(1);
60     latmax = lat125(2177-1);
61     depmin = min(z_w);
62     depmax = z_c(29); depmax = z_w(2);
63     m_proj('mercator','long',[0 360],'lat',[-90 90]);
64     clf;hold on;m_coast;m_grid;
65     LIMITS = [lonmin+180 lonmax+180 latmin latmax depmin depmax]
66     m_line(LIMITS([1 2 2 1 1]),LIMITS([3 3 4 4 3]),'color','k','linewidth',2);
67     title(sub_name);
68    
69    
70     case 10
71     sub_name = 'KE';
72     lonmin = lon125(961)-180;
73     lonmax = lon125(1601-1)-180;
74     latmin = lat125(1140);
75     latmax = lat125(1523-1);
76     depmin = min(z_w);
77     depmax = z_c(25);
78     m_proj('mercator','long',[0 360],'lat',[-90 90]);
79     %clf;hold on;m_coast;m_grid;
80     LIMITS = [lonmin+180 lonmax+180 latmin latmax depmin depmax]
81     %m_line(LIMITS([1 2 2 1 1]),LIMITS([3 3 4 4 3]),'color','k','linewidth',2);
82     %title(sub_name);
83    
84    
85     end
86    
87     %refresh
88    
89     % Path of the directory to find input binary files:
90     pathi = 'ecco2_cycle1_bin/';
91    
92     % Path where the netcdf outputs will be stored:
93     patho = './ecco2_cycle1_netcdf/monthly/';
94     %patho = './ecco2_cycle1_netcdf/six_hourly/';
95    
96     % Variables to analyse (index into otab):
97     wvar = [];
98     % 3D fields:
99     wvar = [wvar 34]; % THETA
100     wvar = [wvar 35]; % THETASQ
101     %wvar = [wvar 33]; % SALTanom
102     %wvar = [wvar 47]; % VVEL
103     %wvar = [wvar 31]; % RHOAnoma
104    
105    
106    
107     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
108     for i = 1 : length(wvar)
109     ifield = wvar(i);
110     fil = otab{ifield,1};
111     l = dir(strcat(pathi,fil,sla));
112     if ifield == 33,
113     l = dir(strcat(pathi,'SALT',sla));
114     end
115     if ifield == 35,
116     l = dir(strcat(pathi,'THETA',sla));
117     end
118     if ifield == 31,
119     l = dir(strcat(pathi,'RHO',sla));
120     end
121     it = 0;
122     clear ll
123     for il = 1 : size(l,1)
124     if ~l(il).isdir & findstr(l(il).name,strcat(fil,'.')) % is'it the file type we want ?
125     it = it + 1;
126     ll(it).name = l(il).name;
127     end %if
128     end %for il
129    
130     if it ~= 0 % if we found any files to compute:
131    
132     % Create the timetable:
133     for il = 1 : size(ll,2)
134     filinprog = ll(il).name;
135     stepnum=str2num(filinprog(findstr(filinprog,fil)+length(fil)+1:length(filinprog)- ...
136     length('.data')));
137     TIME(il,:) = dtecco2(stepnum,0);
138     end
139    
140     % Translate files:
141     for il = 1 : size(ll,2)
142    
143     filinprog = ll(il).name;
144     stepnum=str2num(filinprog(findstr(filinprog,fil)+length(fil)+1:length(filinprog)- ...
145     length('.data')));
146     ID = datenum(1992,1,1)+stepnum*300/60/60/24;
147     dte = datestr(ID,'yyyymmddHHMM');
148     disp(strcat(fil,'->',datestr(ID,'yyyy/mm/dd/HH:MM'),'== Recorded in ==>',TIME(il,:)));
149     dirout = strcat(patho,sla,TIME(il,:));
150    
151     if 1 % Want to record ?
152     if ~exist(dirout,'dir')
153     mkdir(dirout);
154     end
155     pathname = strcat(pathi,fil,sla);
156     if ifield == 33,
157     pathname = strcat(pathi,'SALT',sla);
158     end
159     if ifield == 35,
160     pathname = strcat(pathi,'THETA',sla);
161     end
162     if ifield == 31,
163     pathname = strcat(pathi,'RHO',sla);
164     end
165     if ~exist(strcat(dirout,sla,sprintf('%s.%s.nc',otab{ifield,1},sub_name)),'file')
166     %if 1
167     latlon2ingrid_netcdf(pathname,strcat(dirout,sla),...
168     stepnum,otab{ifield,1},otab, ...
169     lon_c, lon_u, ...
170     lat_c, lat_v, ...
171     z_c, z_w, ...
172     sub_name, ...
173     lonmin,lonmax,latmin,latmax,depmin,depmax);
174     else
175     disp(strcat('##### Skip file (already done):',dirout,sla,...
176     sprintf('%s.%s.nc',otab{ifield,1},sub_name)))
177     end %if %file exist
178    
179     end %if 1/0 want to record ?
180     % if il==1,break,end;
181    
182     fclose('all');
183    
184     end %for il
185    
186     end %if it
187    
188     end %for i

  ViewVC Help
Powered by ViewVC 1.1.22