gael/bulkMatlab/ecmwf_load_fields.m

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

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

if ycur<2002
rep_in='/net/altix3700/raid4/king/data_1x1_92-03/bulk_ECMWF/';
files_ecmwf=strvcat('ERA40g_spfh2m_','ERA40g_tmp2m_degC_','ERA40g_u10m_','ERA40g_v10m_');
files_ecmwf=strvcat(files_ecmwf,'ERA40g_dsw_','ERA40g_dlw_','ERA40g_rain_','../forcing_ECMWF/ERA40g_pres_r1_');
var_ecco=strvcat('aqh','atemp','uwind','vwind','swdn','lwdn','rain','slp');
else
rep_in='/net/altix3700/raid4/king/data_1x1_92-03/bulk_ECMWF/';
files_ecmwf=strvcat('EAG_spfh2m_r1_','EAG_tmp2m_degC_r1_','EAG_u10m_r1_','EAG_v10m_r1_');
files_ecmwf=strvcat(files_ecmwf,'EAG_dsw_r1_','EAG_dlw_r1_','EAG_rain_r1_','../forcing_ECMWF/EAG_pres_r1_');
var_ecco=strvcat('aqh','atemp','uwind','vwind','swdn','lwdn','rain','slp');
end;


ecmwf_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_ecmwf,1)

fid=fopen([rep_in deblank(files_ecmwf(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));

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

%if ycur<2002; slp=NaN*aqh; end; %no pressure field in data_1x1_92-03/bulk_ECMWF/
swdn=-swdn; %wrong sign compared to NCEP
lwdn=-lwdn; %wrong sign compared to NCEP

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

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


1	gforget	1.1	function [mystruct_out]=ecmwf_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			if ycur<2002
8			rep_in='/net/altix3700/raid4/king/data_1x1_92-03/bulk_ECMWF/';
9			files_ecmwf=strvcat('ERA40g_spfh2m_','ERA40g_tmp2m_degC_','ERA40g_u10m_','ERA40g_v10m_');
10			files_ecmwf=strvcat(files_ecmwf,'ERA40g_dsw_','ERA40g_dlw_','ERA40g_rain_','../forcing_ECMWF/ERA40g_pres_r1_');
11			var_ecco=strvcat('aqh','atemp','uwind','vwind','swdn','lwdn','rain','slp');
12			else
13			rep_in='/net/altix3700/raid4/king/data_1x1_92-03/bulk_ECMWF/';
14			files_ecmwf=strvcat('EAG_spfh2m_r1_','EAG_tmp2m_degC_r1_','EAG_u10m_r1_','EAG_v10m_r1_');
15			files_ecmwf=strvcat(files_ecmwf,'EAG_dsw_r1_','EAG_dlw_r1_','EAG_rain_r1_','../forcing_ECMWF/EAG_pres_r1_');
16			var_ecco=strvcat('aqh','atemp','uwind','vwind','swdn','lwdn','rain','slp');
17			end;
18
19
20			ecmwf_grid; jpig=length(x); jpjg=length(y);
21			lon2D_g=x'ones(1,jpjg); lat2D_g=ones(jpig,1)y; mask_g=ones(jpig,jpjg);
22			recl=jpigjpjg4;
23
24			[x1,y1] = meshgrid([lon2D_g(end,1)-360 lon2D_g(:,1)' lon2D_g(1:2,1)'+360],lat2D_g(1,:));
25			[x2,y2] = meshgrid(lon2D_t(:,1),lat2D_t(1,:));
26
27			for fcur=1:size(files_ecmwf,1)
28
29			fid=fopen([rep_in deblank(files_ecmwf(fcur,:)) num2str(ycur)],'r','b');
30			recl=jpigjpjg4; position0=recl*(hcur-1); status=fseek(fid,position0,'bof');
31			field_mod=fread(fid,[jpig jpjg],'float32'); fclose(fid);
32			%field_mod=fread(fid,jpigjpjg4,'float32'); fclose(fid); field_mod=squeeze(mean(reshape(field_mod,[jpig jpjg 4]),3));
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			%if ycur<2002; slp=NaN*aqh; end; %no pressure field in data_1x1_92-03/bulk_ECMWF/
42			swdn=-swdn; %wrong sign compared to NCEP
43			lwdn=-lwdn; %wrong sign compared to NCEP
44
45			rhoConstFresh=999.8;
46			atemp=273.15+atemp; %conversion to Kelvin
47			precip=rain*rhoConstFresh; %m/s->kg/m^2/s
48			slp=slp/100;
49
50			mystruct_out=struct('aqh',aqh,'atemp',atemp,'uwind',uwind,'vwind',vwind,'swdn',swdn,'lwdn',lwdn,'precip',precip,'slp',slp);
51
52