FpOutput: flex_read_V7.m

File flex_read_V7.m, 6.4 KB (added by sureshraja1, 3 years ago)
Line 
1function [grid, drydep, wetdep, fail] = flex_read_V7(header,unit,nspec_ret,pspec_ret,age_ret,nest,time_ret,fin)
2%===========================================
3% flex_read_V7.m
4%-------------------------------------------
5% input
6%   - header object
7%   - optional: data units, (1==conc(default), 2==pptv, 3==time(bkwd) 4==footprint_total 9=column_total)
8%   - optional: which species
9%               which pointspecies
10%               which ageclasse
11%   - optional: nest (0==no(default)/1==yes)
12%               which time
13% output
14%   - data
15%   - dry deposition
16%   - wet deposition
17%   - fail indicator (-1=fail, 0=success)
18%-------------------------------------------
19% FLEXPART matlab import routines
20%-------------------------------------------
21% last changes: HSO, 15-Jun-2007
22%===========================================
23
24% preset output variables
25fail=1;
26
27% set default units
28prefix=['grid_conc_          ';'grid_pptv_          '
29       ;'grid_time_          ';'footprint_total     '; ...
30        'grid_conc_nest_     ';'grid_pptv_nest_     '; ...
31        'grid_time_nest_     ';'footprint_total_nest';
32        'column_total        '];
33
34if nargin==1
35  unit=1;
36end
37
38if nargin<=2
39  nest=0;
40end
41
42if nest==1
43  unit=unit+4;
44end
45
46if nargin<=6 || time_ret(1)<=0
47  time_ret=1:length(header.dates);
48end % if times
49
50
51
52if unit==4 | unit==9
53grid=zeros(header.numxgrid,header.numygrid,1,length(nspec_ret),length(pspec_ret),length(age_ret),1);
54else
55grid=zeros(header.numxgrid,header.numygrid,header.numzgrid,length(nspec_ret),length(pspec_ret),length(age_ret),length(time_ret));
56end
57wetdep=zeros(header.numxgrid,header.numygrid,length(nspec_ret),length(pspec_ret),length(age_ret),length(time_ret));
58drydep=zeros(header.numxgrid,header.numygrid,length(nspec_ret),length(pspec_ret),length(age_ret),length(time_ret));
59%header.numpointspec=1;
60
61dat_cnt=0;
62
63%--------------------------------------------------
64% Loop over all times, given in field header.dates
65%--------------------------------------------------
66for dat_i=time_ret,
67dat_cnt=dat_cnt+1;
68ks_cnt=0;
69
70for ks=nspec_ret,
71ks_cnt=ks_cnt+1;
72
73if nargin<8 || isempty(fin)
74    % filename to use is not given
75    if dat_i<=length(header.dates)
76       if unit~=4 & unit~=9 %not footprint or tot_col
77         fin=[header.path deblank(prefix((unit),:)) datestr(header.dates(dat_i),'YYYYmmddHHMMSS') num2str(ks,'_%03d') ];
78       else
79         fin=[header.path deblank(prefix((unit),:)) num2str(ks,'_%03d') ];
80       end
81    else
82        fprintf('Date does not exist! \n');
83    end
84end
85
86
87  fprintf('Opening %s \n',fin);
88  fid=fopen(fin);
89 
90if fid<=0
91    fprintf('******************************\n could not open file %s !\n******************************\n',fin);
92else
93     
94rl=fread(fid,1,'int32');
95itime=fread(fid,1,'int32');
96rl=fread(fid,2,'int32'); 
97
98%----------------------------------------------
99% read all fields
100%----------------------------------------------
101kp_cnt=0;
102
103for kp=1:header.numpointspec
104  %  fprintf('reading: %f \n',kp);
105    kp_read=0;
106    if (kp_cnt)<length(pspec_ret)
107      if find(pspec_ret(kp_cnt+1)==kp)
108         kp_cnt=kp_cnt+1;
109         kp_read=1;
110      end
111    end
112   nage_cnt=0;
113   for nage=1:max(1,header.nageclass)
114       nage_read=0;
115       if (nage_cnt)<length(age_ret)
116         if find(age_ret(nage_cnt+1)==nage)
117          nage_cnt=nage_cnt+1;
118          nage_read=1;
119         end
120       end
121      %----------------------------------------------
122      % read wetgrid
123      %----------------------------------------------
124       sp_count_i=fread(fid,1,'int32');
125       rl=fread(fid,2,'int32'); 
126       sparse_dump_i=fread(fid,sp_count_i,'int32');
127       rl=fread(fid,2,'int32'); 
128       sp_count_r=fread(fid,1,'int32');         
129       rl=fread(fid,2,'int32'); 
130       sparse_dump_r=fread(fid,sp_count_r,'float');   
131       rl=fread(fid,2,'int32');
132             
133       if (nage_read & kp_read)
134       fact=1;
135       ii=1;
136       for ir=1:length(sparse_dump_r)
137          if sparse_dump_r(ir)*fact>0       
138             n=sparse_dump_i(ii);
139             fact=fact*-1;
140             ii=ii+1;
141          else
142             n=n+1;
143          end
144          jy=floor(n/header.numxgrid);
145          ix=n-header.numxgrid*jy;
146          wetdep(ix+1,jy+1,ks_cnt,kp_cnt,nage_cnt,dat_cnt)=sparse_dump_r(ir)*fact*(-1);
147       end
148       end
149       
150   
151      %----------------------------------------------
152      % read drygrid
153      %----------------------------------------------
154       sp_count_i=fread(fid,1,'int32',8);
155       sparse_dump_i=fread(fid,sp_count_i,'int32');   
156       rl=fread(fid,2,'int32');
157       sp_count_r=fread(fid,1,'int32',8);         
158       sparse_dump_r=fread(fid,sp_count_r,'float');   
159       rl=fread(fid,2,'int32');
160
161       if (nage_read & kp_read)
162       fact=1;
163       ii=1;
164       for ir=1:length(sparse_dump_r)
165          if sparse_dump_r(ir)*fact>0       
166             n=sparse_dump_i(ii);
167             fact=fact*-1;
168             ii=ii+1;
169          else
170             n=n+1;
171          end
172          jy=floor(n/header.numxgrid);
173          ix=n-header.numxgrid*jy;
174          drydep(ix+1,jy+1,ks_cnt,kp_cnt,nage_cnt,dat_cnt)=sparse_dump_r(ir)*fact*(-1);
175       end
176       end
177     
178      %----------------------------------------------
179                % read data grid
180      %----------------------------------------------     
181   
182       sp_count_i=fread(fid,1,'int32');
183       rl=fread(fid,2,'int32'); 
184       sparse_dump_i=fread(fid,sp_count_i,'int32');   
185       rl=fread(fid,2,'int32');
186       sp_count_r=fread(fid,1,'int32');         
187       rl=fread(fid,2,'int32'); 
188       sparse_dump_r=fread(fid,sp_count_r,'float');   
189       rl=fread(fid,2,'int32');
190       
191       if (nage_read & kp_read)
192       fact=1;
193       ii=1;
194       for ir=1:length(sparse_dump_r)
195          if sparse_dump_r(ir)*fact>0       
196             n=sparse_dump_i(ii);
197             fact=fact*-1;
198             ii=ii+1;
199          else
200             n=n+1;
201          end
202          kz=floor(n/(header.numxgrid*header.numygrid));
203          if kz>length(header.zp)
204              fprintf('ERROR - number of vertical levels do not match!');
205          end
206          jy=floor((n-kz*header.numxgrid*header.numygrid)/header.numxgrid);
207          ix=floor(n-header.numxgrid*header.numygrid*kz-header.numxgrid*jy);     
208          grid(ix+1,jy+1,kz,ks_cnt,kp_cnt,nage_cnt,dat_cnt)=sparse_dump_r(ir)*(-1)*fact;
209         end
210       end
211
212   end % for nage
213end % for numpoint
214fclose(fid);
215end % file succesfully opened
216
217end % for species
218end % for dat_i
219
220
221
222fail=0;
223
224end % of function flex_read
225
226%fin
hosted by ZAMG