gael/bulkMatlab/ncep_load_fields.m

function [mystruct_out]=ncep_load_fields(ycur,hcur);
%loads the fields and does the interpolation to the 1x1 ECCO grid
domaine_global_def;

mask=squeeze(tmask3D(:,:,1));

rep_in='/net/altix3700/raid4/king/data_1x1_92-03/bulk/';

files_ncep=strvcat('NCEP_RG_spfh2m_','NCEP_RG_tmp2m_degC_','NCEP_RG_u10m_','NCEP_RG_v10m_');
files_ncep=strvcat(files_ncep,'NCEP_RG_dsw_','NCEP_RG_dlw_','NCEP_RG_rain_','NCEP_RG_pres_');

var_ecco=strvcat('aqh','atemp','uwind','vwind','swdn','lwdn','rain','slp');

ncep_grid; jpig=length(x); jpjg=length(y);
lon2D_g=x'*ones(1,jpjg); lat2D_g=ones(jpig,1)*y; mask_g=ones(jpig,jpjg);
recl=jpig*jpjg*4;

[x1,y1] = meshgrid([lon2D_g(end,1)-360  lon2D_g(:,1)' lon2D_g(1:2,1)'+360],lat2D_g(1,:));
[x2,y2] = meshgrid(lon2D_t(:,1),lat2D_t(1,:));

for fcur=1:size(files_ncep,1)

%[rep_in deblank(files_ncep(fcur,:)) num2str(ycur)]
fid=fopen([rep_in deblank(files_ncep(fcur,:)) num2str(ycur)],'r','b');
recl=jpig*jpjg*4; position0=recl*(hcur-1); status=fseek(fid,position0,'bof');
field_mod=fread(fid,[jpig jpjg],'float32'); fclose(fid);
%field_mod=fread(fid,jpig*jpjg*4,'float32'); fclose(fid); field_mod=squeeze(mean(reshape(field_mod,[jpig jpjg 4]),3));

%pressure fix needed: last point missing (Charmaine interpolation 2.5x2.5->ncep grid)
if ~isempty(findstr(var_ecco(fcur,:),'slp'));
field_mod(end,:)=(field_mod(end-1,:)+field_mod(1,:))/2;
end;

tmp1=[field_mod(end,:);field_mod;field_mod(1,:);field_mod(2,:)]';
field_int=interp2(x1,y1,tmp1,x2,y2)';

eval([deblank(var_ecco(fcur,:)) '=field_int.*mask;']);

end


rhoConstFresh=999.8;
atemp=273.15+atemp; %conversion to Kelvin
precip=rain*rhoConstFresh; %m/s->kg/m^2/s

mystruct_out=struct('aqh',aqh,'atemp',atemp,'uwind',uwind,'vwind',vwind,'swdn',swdn,'lwdn',lwdn,'precip',precip,'slp',slp);

1	function [mystruct_out]=ncep_load_fields(ycur,hcur);
2	%loads the fields and does the interpolation to the 1x1 ECCO grid
3	domaine_global_def;
4
5	mask=squeeze(tmask3D(:,:,1));
6
7	rep_in='/net/altix3700/raid4/king/data_1x1_92-03/bulk/';
8
9	files_ncep=strvcat('NCEP_RG_spfh2m_','NCEP_RG_tmp2m_degC_','NCEP_RG_u10m_','NCEP_RG_v10m_');
10	files_ncep=strvcat(files_ncep,'NCEP_RG_dsw_','NCEP_RG_dlw_','NCEP_RG_rain_','NCEP_RG_pres_');
11
12	var_ecco=strvcat('aqh','atemp','uwind','vwind','swdn','lwdn','rain','slp');
13
14	ncep_grid; jpig=length(x); jpjg=length(y);
15	lon2D_g=x'ones(1,jpjg); lat2D_g=ones(jpig,1)y; mask_g=ones(jpig,jpjg);
16	recl=jpigjpjg4;
17
18	[x1,y1] = meshgrid([lon2D_g(end,1)-360 lon2D_g(:,1)' lon2D_g(1:2,1)'+360],lat2D_g(1,:));
19	[x2,y2] = meshgrid(lon2D_t(:,1),lat2D_t(1,:));
20
21	for fcur=1:size(files_ncep,1)
22
23	%[rep_in deblank(files_ncep(fcur,:)) num2str(ycur)]
24	fid=fopen([rep_in deblank(files_ncep(fcur,:)) num2str(ycur)],'r','b');
25	recl=jpigjpjg4; position0=recl*(hcur-1); status=fseek(fid,position0,'bof');
26	field_mod=fread(fid,[jpig jpjg],'float32'); fclose(fid);
27	%field_mod=fread(fid,jpigjpjg4,'float32'); fclose(fid); field_mod=squeeze(mean(reshape(field_mod,[jpig jpjg 4]),3));
28
29	%pressure fix needed: last point missing (Charmaine interpolation 2.5x2.5->ncep grid)
30	if ~isempty(findstr(var_ecco(fcur,:),'slp'));
31	field_mod(end,:)=(field_mod(end-1,:)+field_mod(1,:))/2;
32	end;
33
34	tmp1=[field_mod(end,:);field_mod;field_mod(1,:);field_mod(2,:)]';
35	field_int=interp2(x1,y1,tmp1,x2,y2)';
36
37	eval([deblank(var_ecco(fcur,:)) '=field_int.*mask;']);
38
39	end
40
41
42	rhoConstFresh=999.8;
43	atemp=273.15+atemp; %conversion to Kelvin
44	precip=rain*rhoConstFresh; %m/s->kg/m^2/s
45
46	mystruct_out=struct('aqh',aqh,'atemp',atemp,'uwind',uwind,'vwind',vwind,'swdn',swdn,'lwdn',lwdn,'precip',precip,'slp',slp);
47