source: trunk/src/series/time_series.m @ 511

Last change on this file since 511 was 478, checked in by sommeria, 12 years ago

series fcts updated to fit with new waitbar fct and background run mode, and various bug repairs

File size: 21.1 KB
Line 
1%'time_series': extract a time series, used with series.fig
2% this function can be used as a template for applying a global operation on a series of input fields
3%------------------------------------------------------------------------
4% function GUI_input=time_series(Param)
5%
6%%%%%%%%%%% GENERAL TO ALL SERIES ACTION FCTS %%%%%%%%%%%%%%%%%%%%%%%%%%%
7%
8% This function is used in four modes by the GUI series:
9%           1) config GUI: with no input argument, the function determine the suitable GUI configuration
10%           2) interactive input: the function is used to interactively introduce input parameters, and then stops
11%           3) RUN: the function itself runs, when an appropriate input  structure Param has been introduced.
12%           4) BATCH: the function itself proceeds in BATCH mode, using an xml file 'Param' as input.
13%
14% This function is used in four modes by the GUI series:
15%           1) config GUI: with no input argument, the function determine the suitable GUI configuration
16%           2) interactive input: the function is used to interactively introduce input parameters, and then stops
17%           3) RUN: the function itself runs, when an appropriate input  structure Param has been introduced.
18%           4) BATCH: the function itself proceeds in BATCH mode, using an xml file 'Param' as input.
19%
20%OUTPUT
21% GUI_input=list of options in the GUI series.fig needed for the function
22%
23%INPUT:
24% In run mode, the input parameters are given as a Matlab structure Param copied from the GUI series.
25% In batch mode, Param is the name of the corresponding xml file containing the same information
26% In the absence of input (as activated when the current Action is selected
27% in series), the function ouput GUI_input set the activation of the needed GUI elements
28%
29% Param contains the elements:(use the menu bar command 'export/GUI config' in series to see the current structure Param)
30%    .InputTable: cell of input file names, (several lines for multiple input)
31%                      each line decomposed as {RootPath,SubDir,Rootfile,NomType,Extension}
32%    .OutputSubDir: name of the subdirectory for data outputs
33%    .OutputDirExt: directory extension for data outputs
34%    .Action: .ActionName: name of the current activated function
35%             .ActionPath:   path of the current activated function
36%    .IndexRange: set the file or frame indices on which the action must be performed
37%    .FieldTransform: .TransformName: name of the selected transform function
38%                     .TransformPath:   path  of the selected transform function
39%                     .TransformHandle: corresponding function handle
40%    .InputFields: sub structure describing the input fields withfields
41%              .FieldName: name of the field
42%              .VelType: velocity type
43%              .FieldName_1: name of the second field in case of two input series
44%              .VelType_1: velocity type of the second field in case of two input series
45%    .ProjObject: %sub structure describing a projection object (read from ancillary GUI set_object)
46%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
47%
48function ParamOut=time_series(Param)
49
50%% set the input elements needed on the GUI series when the action is selected in the menu ActionName
51if ~exist('Param','var') % case with no input parameter
52    ParamOut={'AllowInputSort';'off';...% allow alphabetic sorting of the list of input files (options 'off'/'on', 'off' by default)
53        'WholeIndexRange';'off';...% prescribes the file index ranges from min to max (options 'off'/'on', 'off' by default)
54        'NbSlice';'on'; ...%nbre of slices ('off' by default)
55        'VelType';'two';...% menu for selecting the velocity type (options 'off'/'one'/'two',  'off' by default)
56        'FieldName';'two';...% menu for selecting the field (s) in the input file(options 'off'/'one'/'two', 'off' by default)
57        'FieldTransform'; 'on';...%can use a transform function
58        'ProjObject';'on';...%can use projection object(option 'off'/'on',
59        'Mask';'off';...%can use mask option   (option 'off'/'on', 'off' by default)
60        'OutputDirExt';'.tseries';...%set the output dir extension
61               ''};
62        return
63end
64
65%%%%%%%%%%%% STANDARD PART (DO NOT EDIT) %%%%%%%%%%%%
66%% select different modes,  RUN, parameter input, BATCH
67% BATCH  case: read the xml file for batch case
68if ischar(Param)
69        Param=xml2struct(Param);
70        checkrun=0;
71% RUN case: parameters introduced as the input structure Param
72else
73    hseries=guidata(Param.hseries);%handles of the GUI series
74    if isfield(Param,'Specific')&& strcmp(Param.Specific,'?')
75        checkrun=1;% will only search interactive input parameters (preparation of BATCH mode)
76    else
77        checkrun=2; % indicate the RUN option is used
78    end
79end
80ParamOut=Param; %default output
81OutputDir=[Param.OutputSubDir Param.OutputDirExt];
82
83%% root input file(s) and type
84RootPath=Param.InputTable(:,1);
85RootFile=Param.InputTable(:,3);
86SubDir=Param.InputTable(:,2);
87NomType=Param.InputTable(:,4);
88FileExt=Param.InputTable(:,5);
89
90% get the set of input file names (cell array filecell), and the lists of
91% input file or frame indices i1_series,i2_series,j1_series,j2_series
92[filecell,i1_series,i2_series,j1_series,j2_series]=get_file_series(Param);
93% filecell{iview,fileindex}: cell array representing the list of file names
94%        iview: line in the table corresponding to a given file series
95%        fileindex: file index within  the file series,
96% i1_series(iview,ref_j,ref_i)... are the corresponding arrays of indices i1,i2,j1,j2, depending on the input line iview and the two reference indices ref_i,ref_j
97% i1_series(iview,fileindex) expresses the same indices as a 1D array in file indices
98% set of frame indices used for movie or multimage input
99% numbers of slices and file indices
100
101NbSlice=1;%default
102if isfield(Param.IndexRange,'NbSlice')&&~isempty(Param.IndexRange.NbSlice)
103    NbSlice=Param.IndexRange.NbSlice;
104end
105nbview=numel(i1_series);%number of input file series (lines in InputTable)
106nbfield_j=size(i1_series{1},1); %nb of fields for the j index (bursts or volume slices)
107nbfield_i=size(i1_series{1},2); %nb of fields for the i index
108nbfield=nbfield_j*nbfield_i; %total number of fields
109nbfield_i=floor(nbfield/NbSlice);%total number of  indexes in a slice (adjusted to an integer number of slices)
110nbfield=nbfield_i*NbSlice; %total number of fields after adjustement
111
112%determine the file type on each line from the first input file
113ImageTypeOptions={'image','multimage','mmreader','video'};
114NcTypeOptions={'netcdf','civx','civdata'};
115for iview=1:nbview
116    if ~exist(filecell{iview,1}','file')
117        displ_uvmat('ERROR',['the first input file ' filecell{iview,1} ' does not exist'],checkrun)
118        return
119    end
120    [FileType{iview},FileInfo{iview},MovieObject{iview}]=get_file_type(filecell{iview,1});
121    CheckImage{iview}=~isempty(find(strcmp(FileType{iview},ImageTypeOptions)));% =1 for images
122    CheckNc{iview}=~isempty(find(strcmp(FileType{iview},NcTypeOptions)));% =1 for netcdf files
123    if ~isempty(j1_series{iview})
124        frame_index{iview}=j1_series{iview};
125    else
126        frame_index{iview}=i1_series{iview};
127    end
128end
129
130%% calibration data and timing: read the ImaDoc files
131[XmlData,NbSlice_calib,time,errormsg]=read_multimadoc(RootPath,SubDir,RootFile,FileExt,i1_series,i2_series,j1_series,j2_series);
132if size(time,1)>1
133    diff_time=max(max(diff(time)));
134    if diff_time>0
135        displ_uvmat('WARNING',['times of series differ by (max) ' num2str(diff_time)],checkrun)
136    end   
137end
138
139%% coordinate transform or other user defined transform
140transform_fct='';%default
141if isfield(Param,'FieldTransform')&&~isempty(Param.FieldTransform.TransformName)
142    addpath(Param.FieldTransform.TransformPath)
143    transform_fct=str2func(Param.FieldTransform.TransformName);
144    rmpath(Param.FieldTransform.TransformPath)
145end
146
147%%%%%%%%%%%% END STANDARD PART  %%%%%%%%%%%%
148 % EDIT FROM HERE
149
150%% check the validity of  input file types
151if CheckImage{1}
152    FileExtOut='.png'; % write result as .png images for image inputs
153elseif CheckNc{1}
154    FileExtOut='.nc';% write result as .nc files for netcdf inputs
155else
156    displ_uvmat('ERROR',['invalid file type input ' FileType{1}],checkrun)
157    return
158end
159if nbview==2 && ~isequal(CheckImage{1},CheckImage{2})
160        displ_uvmat('ERROR','input must be two image series or two netcdf file series',checkrun)
161    return
162end
163NomTypeOut='_1-2_1';% output file index will indicate the first and last ref index in the series
164% if NbSlice~=nbfield_j
165%     answer=_uvmat('INPUT_Y-N',['will not average slice by slice: for so cancel and set NbSlice= ' num2str(nbfield_j)]);
166%     if ~strcmp(answer,'Yes')
167%         return
168%     end
169% end
170if checkrun==1
171    return % stop here for input checks
172end
173
174%% Set field names and velocity types
175InputFields{1}=[];%default (case of images)
176if isfield(Param,'InputFields')
177    InputFields{1}=Param.InputFields;
178end
179if nbview==2
180    InputFields{2}=[];%default (case of images)
181    if isfield(Param,'InputFields')
182        InputFields{2}=Param.InputFields{1};%default
183        if isfield(Param.InputFields,'FieldName_1')
184            InputFields{2}.FieldName=Param.InputFields.FieldName_1;
185            if isfield(Param.InputFields,'VelType_1')
186                InputFields{2}.VelType=Param.InputFields.VelType_1;
187            end
188        end
189    end
190end
191%%% TO UPDATE
192if isequal(InputFields{1},'get_field...')
193    hget_field=findobj(allchild(0),'name','get_field');%find the get_field... GUI
194    if numel(hget_field)>1
195        delete(hget_field(2:end)) % delete multiple occurerence of the GUI get_fioeld
196    elseif isempty(hget_field)
197        filename=filecell{1,1};
198      % filename=name_generator(filebase{1},i1_series{1}(1),j1_series{1}(1),FileExt{1},NomType{1},1,i2_series{1}(1),num_j2{1}(1),SubDir{1});
199       idetect(iview)=exist(filename,'file');
200       hget_field=get_field(filename);
201       return
202    end
203    SubField=read_get_field(hget_field); %read the names of the variables to plot in the get_field GUI
204    if isempty(SubField)
205        delete(hget_field)
206        filename=filecell{1,1};
207       %filename=name_generator(filebase{1},i1_series{1}(1),j1_series{1}(1),FileExt{1},NomType{1},1,i2_series{1}(1),j2_series{1}(1),SubDir{1});
208        hget_field=get_field(filename);
209        SubField=read_get_field(hget_field); %read the names of the variables to plot in the get_field GUI
210    end
211end
212%%%%%%%
213
214%% Initiate output fields
215%initiate the output structure as a copy of the first input one (reproduce fields)
216[DataOut,tild,errormsg] = read_field(filecell{1,1},FileType{1},InputFields{1},1);
217if ~isempty(errormsg)
218    displ_uvmat('ERROR',['error reading ' filecell{1,1} ': ' errormsg],checkrun)
219    return
220end
221time_1=[];
222if isfield(DataOut,'Time')
223    time_1=DataOut.Time(1);
224end
225if CheckNc{iview}
226    if isempty(strcmp('Conventions',DataOut.ListGlobalAttribute))
227        DataOut.ListGlobalAttribute=['Conventions' DataOut.ListGlobalAttribute];
228    end
229    DataOut.Conventions='uvmat';
230    DataOut.ListGlobalAttribute=[DataOut.ListGlobalAttribute {Param.Action}];
231    ActionKey='Action';
232    while isfield(DataOut,ActionKey)
233        ActionKey=[ActionKey '_1'];
234    end
235    DataOut.(ActionKey)=Param.Action;
236    DataOut.ListGlobalAttribute=[DataOut.ListGlobalAttribute {ActionKey}];
237    if isfield(DataOut,'Time')
238        DataOut.ListGlobalAttribute=[DataOut.ListGlobalAttribute {'Time','Time_end'}];
239    end
240end
241
242%% LOOP ON SLICES
243nbmissing=0; %number of undetected files
244for i_slice=1:NbSlice
245    index_slice=i_slice:NbSlice:nbfield;% select file indices of the slice
246    nbfile=0;
247    nbmissing=0;
248   
249    %%%%%%%%%%%%%%%% loop on field indices %%%%%%%%%%%%%%%%
250    for index=index_slice       
251        if checkrun
252            update_waitbar(hseries.Waitbar,index/(nbfield))
253            stopstate=get(hseries.RUN,'BusyAction');
254        else
255            stopstate='queue';
256        end
257       
258        %%%%%%%%%%%%%%%% loop on views (input lines) %%%%%%%%%%%%%%%%
259        if isequal(stopstate,'queue')% enable STOP command
260            Data=cell(1,nbview);%initiate the set Data
261            nbtime=0;
262            dt=[];
263            % loop on views (in case of multiple input series)
264            for iview=1:nbview
265                % reading input file(s)
266                [Data{iview},tild,errormsg] = read_field(filecell{iview,index},FileType{iview},Param.InputFields,frame_index{iview}(index));
267                if ~isempty(errormsg)
268                    errormsg=['time_series/read_field/' errormsg];
269                    display(errormsg)
270                    break
271                end
272                timeread(iview)=0;
273                if isfield(Data{iview},'Time')
274                    timeread(iview)=Data{iview}.Time;
275                    nbtime=nbtime+1;
276                end
277                if ~isempty(NbSlice_calib)
278                    Data{iview}.ZIndex=mod(i1_series{iview}(index)-1,NbSlice_calib{iview})+1;%Zindex for phys transform
279                end
280            end
281           
282            % coordinate transform (or other user defined transform)
283            if ~isempty(transform_fct)
284                if nbview==2
285                    [Data{1},Data{2}]=transform_fct(Data{1},XmlData{1},Data{2},XmlData{2});
286                    if isempty(Data{2})
287                        Data(2)=[];
288                    end
289                else
290                    Data{1}=transform_fct(Data{1},XmlData{1});
291                end
292            end
293           
294            % field calculation (vort, div...)
295            if strcmp(FileType{1},'civx')||strcmp(FileType{1},'civdata')
296                Data{1}=calc_field(InputFields{1}.FieldName,Data{1});%calculate field (vort..)
297            end
298           
299            % field substration (for two input file series)
300            if length(Data)==2
301                [Field,errormsg]=sub_field(Data{1},Data{2}); %substract the two fields
302            else
303                Field=Data{1};
304            end
305            if Param.CheckObject
306                [Field,errormsg]=proj_field(Field,Param.ProjObject);
307            end
308            nbfile=nbfile+1;
309           
310            % initiate the time series at the first iteration
311            if nbfile==1
312                % stop program if the first field reading is in error
313                if ~isempty(errormsg)
314                    displ_uvmat('ERROR',['error in time_series/sub_field:' errormsg],checkrun)
315                    return
316                end
317                DataOut=Field;%default
318                DataOut.NbDim=Field.NbDim+1; %add the time dimension for plots
319                nbvar=length(Field.ListVarName);
320                if nbvar==0
321                    displ_uvmat('ERROR','no input variable selected in get_field',checkrun)
322                    return
323                end
324                testsum=2*ones(1,nbvar);%initiate flag for action on each variable
325                if isfield(Field,'VarAttribute') % look for coordinate and flag variables
326                    for ivar=1:nbvar
327                        if length(Field.VarAttribute)>=ivar && isfield(Field.VarAttribute{ivar},'Role')
328                            var_role=Field.VarAttribute{ivar}.Role;%'role' of the variable
329                            if isequal(var_role,'errorflag')
330                                displ_uvmat('ERROR','do not handle error flags in time series',checkrun)
331                                return
332                            end
333                            if isequal(var_role,'warnflag')
334                                testsum(ivar)=0;  % not recorded variable
335                                eval(['DataOut=rmfield(DataOut,''' Field.ListVarName{ivar} ''');']);%remove variable
336                            end
337                            if isequal(var_role,'coord_x')| isequal(var_role,'coord_y')|...
338                                    isequal(var_role,'coord_z')|isequal(var_role,'coord')
339                                testsum(ivar)=1; %constant coordinates, record without time evolution
340                            end
341                        end
342                        % check whether the variable ivar is a dimension variable
343                        DimCell=Field.VarDimName{ivar};
344                        if ischar(DimCell)
345                            DimCell={DimCell};
346                        end
347                        if numel(DimCell)==1 && isequal(Field.ListVarName{ivar},DimCell{1})%detect dimension variables
348                            testsum(ivar)=1;
349                        end
350                    end
351                end
352                for ivar=1:nbvar
353                    if testsum(ivar)==2
354                        eval(['DataOut.' Field.ListVarName{ivar} '=[];'])
355                    end
356                end
357                DataOut.ListVarName=[{'Time'} DataOut.ListVarName];
358            end
359           
360            % add data to the current field
361            for ivar=1:length(Field.ListVarName)
362                VarName=Field.ListVarName{ivar};
363                VarVal=Field.(VarName);
364                if testsum(ivar)==2% test for recorded variable
365                    if isempty(errormsg)
366                        if isequal(Param.ProjObject.ProjMode,'inside')% take the average in the domain for 'inside' mode
367                            if isempty(VarVal)
368                                displ_uvmat('ERROR',['empty result at frame index ' num2str(i1_series{iview}(index))],checkrun)
369                                return
370                            end
371                            VarVal=mean(VarVal,1);
372                        end
373                        VarVal=shiftdim(VarVal,-1); %shift dimension
374                        DataOut.(VarName)=cat(1,DataOut.(VarName),VarVal);%concanete the current field to the time series
375                    else
376                        DataOut.(VarName)=cat(1,DataOut.(VarName),0);% put each variable to 0 in case of input reading error
377                    end
378                elseif testsum(ivar)==1% variable representing fixed coordinates
379                    eval(['VarInit=DataOut.' VarName ';']);
380                    if isempty(errormsg) && ~isequal(VarVal,VarInit)
381                        displ_uvmat('ERROR',['time series requires constant coordinates ' VarName],checkrun)
382                        return
383                    end
384                end
385            end
386           
387            % record the time:
388            if isempty(time)% time read in ncfiles
389                if isfield(Field,'Time')
390                    DataOut.Time(nbfile,1)=Field.Time;
391                else
392                    DataOut.Time(nbfile,1)=index;%default
393                end
394            else % time from ImaDoc prevails  TODO: correct
395                DataOut.Time(nbtime,1)=i1_series{1}(index);% TODO : generalise
396            end
397           
398            % record the number of missing input fields
399            if ~isempty(errormsg)
400                nbmissing=nbmissing+1;
401                display(['index=' num2str(index) ':' errormsg])
402            end
403        end
404    end
405    %%%%%%% END OF LOOP WITHIN A SLICE
406   
407    %remove time for global attributes if exists
408    Time_index=find(strcmp('Time',DataOut.ListGlobalAttribute));
409    if ~isempty(Time_index)
410        DataOut.ListGlobalAttribute(Time_index)=[];
411    end
412    DataOut.Conventions='uvmat';
413    for ivar=1:numel(DataOut.ListVarName)
414        VarName=DataOut.ListVarName{ivar};
415        eval(['DataOut.' VarName '=squeeze(DataOut.' VarName ');']) %remove singletons
416    end
417   
418    % add time dimension
419    for ivar=1:length(Field.ListVarName)
420        DimCell=Field.VarDimName(ivar);
421        if testsum(ivar)==2%variable used as time series
422            DataOut.VarDimName{ivar}=[{'Time'} DimCell];
423        elseif testsum(ivar)==1
424            DataOut.VarDimName{ivar}=DimCell;
425        end
426    end
427    indexremove=find(~testsum);
428    if ~isempty(indexremove)
429        DataOut.ListVarName(1+indexremove)=[];
430        DataOut.VarDimName(indexremove)=[];
431        if isfield(DataOut,'Role') && ~isempty(DataOut.Role{1})%generaliser aus autres attributs
432            DataOut.Role(1+indexremove)=[];
433        end
434    end
435   
436    %shift variable attributes
437    if isfield(DataOut,'VarAttribute')
438        DataOut.VarAttribute=[{[]} DataOut.VarAttribute];
439    end
440    DataOut.VarDimName=[{'Time'} DataOut.VarDimName];
441    DataOut.Action=Param.Action;%name of the processing programme
442    test_time=diff(DataOut.Time)>0;% test that the readed time is increasing (not constant)
443    if ~test_time
444        DataOut.Time=[1:filecounter];
445    end
446   
447    % display nbmissing
448    if ~isequal(nbmissing,0)
449        displ_uvmat('WARNING',[num2str(nbmissing) ' files skipped: missing files or bad input, see command window display'],checkrun)
450    end
451   
452    %name of result file
453    OutputFile=fullfile_uvmat(RootPath{1},OutputDir,RootFile{1},FileExtOut,NomTypeOut,i1_series{1}(1),i1_series{1}(end),i_slice,[]);
454    errormsg=struct2nc(OutputFile,DataOut); %save result file
455    if isempty(errormsg)
456        display([OutputFile ' written'])
457    else
458        displ_uvmat('ERROR',['error in Series/struct2nc: ' errormsg],checkrun)
459    end
460end
461
462%% plot the time series (the last one in case of multislices)
463if checkrun
464    figure
465    haxes=axes;
466    plot_field(DataOut,haxes)
467       
468    %% display the result file using the GUI get_field
469    hget_field=findobj(allchild(0),'name','get_field');
470    if ~isempty(hget_field)
471        delete(hget_field)
472    end
473    get_field(OutputFile,DataOut)
474end
475
476
Note: See TracBrowser for help on using the repository browser.