Index: /trunk/src/calc_field_tps.m =================================================================== --- /trunk/src/calc_field_tps.m (revision 866) +++ /trunk/src/calc_field_tps.m (revision 867) @@ -105,6 +105,6 @@ VarAttribute=[VarAttribute VarAttributeNew]; end -Attr_FF.Role='errorflag'; -VarAttribute=[VarAttribute {Attr_FF}]; +%Attr_FF.Role='errorflag'; +%VarAttribute=[VarAttribute {Attr_FF}]; FieldName(check_remove)=[]; @@ -160,6 +160,6 @@ end end -DataOut.FF=nbval==0; %put errorflag to 1 for points outside the interpolation rang -nbval(nbval==0)=1;% to avoid division by zero for averaging +%DataOut.FF=nbval==0; %put errorflag to 1 for points outside the interpolation rang +nbval(nbval==0)=NaN;% to avoid division by zero for averaging ListFieldOut=fieldnames(DataOut); for ifield=1:numel(ListFieldOut) Index: /trunk/src/find_field_bounds.m =================================================================== --- /trunk/src/find_field_bounds.m (revision 866) +++ /trunk/src/find_field_bounds.m (revision 867) @@ -41,7 +41,7 @@ NbDim=max(NbDimArray);% spatial dimension of the input field imax=find(NbDimArray==NbDim);% indices of field cells to consider -if isfield(Field,'NbDim') - NbDim=double(Field.NbDim);% deal with plane fields containing z coordinates -end +% if isfield(Field,'NbDim') +% NbDim=double(Field.NbDim);% deal with plane fields containing z coordinates +% end FieldOut.NbDim=NbDim; if NbDim<=1; return; end% stop here for 1D fields Index: /trunk/src/find_file_series.m =================================================================== --- /trunk/src/find_file_series.m (revision 866) +++ /trunk/src/find_file_series.m (revision 867) @@ -68,4 +68,8 @@ i1_series=[];i2_series=[];j1_series=[];j2_series=[]; i1_input=1;i2_input=[];j1_input=[];j2_input=[]; + if ~exist(fullfileinput,'file') + RootFile=''; + return + end end if ~exist(FilePath,'dir') @@ -74,5 +78,5 @@ if checkfileindexing NomTypePref=''; -if isempty(NomType) +if isempty(NomType)||strcmp(NomType,'*') if exist(fullfileinput,'file') [tild,RootFile]=fileparts(fileinput);% case of constant name (no indexing), get the filename without its extension Index: /trunk/src/proj_field.m =================================================================== --- /trunk/src/proj_field.m (revision 866) +++ /trunk/src/proj_field.m (revision 867) @@ -666,5 +666,5 @@ %% identify vector components - testU=isfield(CellInfo{icell},'VarIndex_vector_x') &&isfield(CellInfo{icell},'VarIndex_vector_y') ;% test for vectors + %testU=isfield(CellInfo{icell},'VarIndex_vector_x') &&isfield(CellInfo{icell},'VarIndex_vector_y') ;% test for vectors % if testU % UName=FieldData.ListVarName{CellInfo{icell}.VarIndex_vector_x}; @@ -700,6 +700,6 @@ return end - ProjData.ListVarName=[ProjData.ListVarName {FieldData.ListVarName{CellInfo{icell}.CoordIndex(end)}}]; - ProjData.VarDimName=[ProjData.VarDimName {FieldData.ListVarName{CellInfo{icell}.CoordIndex(end-1)}}]; + ProjData.ListVarName=[ProjData.ListVarName FieldData.ListVarName(CellInfo{icell}.CoordIndex(end))]; + ProjData.VarDimName=[ProjData.VarDimName FieldData.ListVarName(CellInfo{icell}.CoordIndex(end))]; nbvar=numel(ProjData.ListVarName); ProjData.VarAttribute{nbvar}.long_name='abscissa along line'; @@ -717,14 +717,14 @@ ProjData.(FieldData.ListVarName{CellInfo{icell}.CoordIndex(end)})=Xproj; for ivar=1:numel(VarIndex) - ProjData.(VarName{ivar})=FieldData.(VarName{ivar})(flagsel);% restrict vrtibles to the projection band + ProjData.(VarName{ivar})=FieldData.(VarName{ivar})(flagsel);% restrict variables to the projection band ProjData.(VarName{ivar})=ProjData.(VarName{ivar})(indsort);% sort by absissa ProjData.ListVarName=[ProjData.ListVarName VarName{ivar}]; - ProjData.VarDimName=[ProjData.VarDimName {FieldData.ListVarName{CellInfo{icell}.CoordIndex(end)}}]; + ProjData.VarDimName=[ProjData.VarDimName FieldData.ListVarName(CellInfo{icell}.CoordIndex(end))]; ProjData.VarAttribute{nbvar+ivar}=FieldData.VarAttribute{VarIndex(ivar)};%reproduce var attribute if isfield(ProjData.VarAttribute{nbvar+ivar},'Role') if strcmp(ProjData.VarAttribute{nbvar+ivar}.Role,'vector_x'); - ivar_U=ivar; + ivar_U=nbvar+ivar; elseif strcmp(ProjData.VarAttribute{nbvar+ivar}.Role,'vector_y'); - ivar_V=ivar; + ivar_V=nbvar+ivar; end end @@ -777,4 +777,6 @@ end [DataOut,VarAttribute,errormsg]=calc_field_tps(Coord,NbCentres,SubRange,FieldVar,CellInfo{icell}.FieldName,cat(3,XI,YI)); + ProjData.ListVarName=[ProjData.ListVarName {'X'}]; + ProjData.VarDimName=[ProjData.VarDimName {'X'}]; ProjData.X=Xproj; nbvar=numel(ProjData.ListVarName); @@ -784,5 +786,5 @@ ProjData.VarAttribute=[ProjData.VarAttribute VarAttribute]; for ivar=1:numel(VarAttribute) - ProjData.VarDimName=[ProjData.VarDimName {XName}]; + ProjData.VarDimName=[ProjData.VarDimName {'X'}]; if isfield(VarAttribute{ivar},'Role') if strcmp(VarAttribute{ivar}.Role,'vector_x'); Index: /trunk/src/read_field.m =================================================================== --- /trunk/src/read_field.m (revision 866) +++ /trunk/src/read_field.m (revision 867) @@ -136,6 +136,6 @@ Role=[Role {'vector_y'}]; ProjModeRequest=[ProjModeRequest {ProjModeRequestVar}]; - ListInputField=[ListInputField {''}]; - %ListOperator=[ListOperator {[r.Operator '_V']}]; + ListInputField=[ListInputField InputField(ilist)]; + %ListInputField=[ListInputField {''}]; else checkV=1; Index: /trunk/src/series.m =================================================================== --- /trunk/src/series.m (revision 866) +++ /trunk/src/series.m (revision 867) @@ -135,5 +135,5 @@ [path_series,name,ext]=fileparts(which('series'));% path to the GUI series path_series_fct=fullfile(path_series,'series');%path of the functions in subdirectroy 'series' -ActionExtList={'.m';'.sh'};% default choice of extensions (Matlab fct .m or compiled version .sh +ActionExtList={'.m';'.sh';'.py (in dev.)'};% default choice of extensions (Matlab fct .m or compiled version .sh ActionPathList=cell(NbBuiltinAction,1);%initiate the cell matrix of Action fct paths ActionPathList(:)={path_series_fct}; %set the default path to series fcts to all list members @@ -620,5 +620,6 @@ if ~exist(RootPath,'dir') i1_series=[]; - RootPath=fileparts(RootPath); %will try the upper folder + %RootPath=fileparts(RootPath); %will try the upper folder + RootFile=''; else %scan the input folder [RootPath,SubDir,RootFile,i1_series,i2_series,j1_series,j2_series,tild,FileInfo,MovieObject]=... @@ -626,5 +627,5 @@ end % if no file is found, open a browser - if isempty(i1_series) + if isempty(RootFile)&& isempty(i1_series) fileinput=uigetfile_uvmat(['wrong input at line ' num2str(iview) ':pick a new input file'],RootPath); if isempty(fileinput) @@ -722,9 +723,9 @@ if strcmp(iview,'append') % display the input data as a new line in the table iview=size(InputTable,1)+1;% the next line in InputTable becomes the current line - InputTable(iview,:)=[{RootPath},{SubDir},{RootFile},{NomType},{FileExt}]; +% InputTable(iview,:)=[{RootPath},{SubDir},{RootFile},{NomType},{FileExt}]; elseif strcmp(iview,'one') % refresh the list of input file series iview=1; %the first line in InputTable becomes the current line InputTable={'','','','',''}; - InputTable(iview,:)=[{RootPath},{SubDir},{RootFile},{NomType},{FileExt}]; +% InputTable(iview,:)=[{RootPath},{SubDir},{RootFile},{NomType},{FileExt}]; set(handles.TimeTable,'Data',[{''},{[]},{[]},{[]},{[]}]) set(handles.MinIndex_i,'Data',[]) @@ -742,4 +743,12 @@ SeriesData.Time={}; end + SeriesData.i1_series(iview+1:end)=[]; + SeriesData.i2_series(iview+1:end)=[]; + SeriesData.j1_series(iview+1:end)=[]; + SeriesData.j2_series(iview+1:end)=[]; + SeriesData.FileType(iview+1:end)=[]; + SeriesData.FileInfo(iview+1:end)=[]; + SeriesData.Time(iview+1:end)=[]; + InputTable(iview,:)=[{RootPath},{SubDir},{RootFile},{NomType},{FileExt}]; if iview >1 set(handles.InputLine,'String',num2str(iview)) @@ -1396,5 +1405,5 @@ ActionExt='.m';%default if isfield(Param.Action,'ActionExt') - ActionExt=Param.Action.ActionExt;% '.m' or '.sh' (compiled) + ActionExt=Param.Action.ActionExt;% '.m', '.sh' (compiled) or '.py' (Python) Param.Action=rmfield(Param.Action,'ActionExt');%remove from the recorded xml file to avoid interference during ImportConfig end @@ -1647,4 +1656,22 @@ set(handles.status,'UserData',StatusData) + +if strcmp(ActionExt, '.py (in dev.)') + fprintf([ + '\n' ... + '!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!\n' ... + 'The option .py is used. It is still in development.\n' ... + 'Do not use it unless you really know what you do!\n' ... + 'To try it, first install Pyp and the most recent version of FluidDyn.\n' ... + '!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!\n']) + RunMode = 'python'; + t = struct2xml(Param); + t = set(t, 1, 'name', 'Series'); + filexml = fullfile_uvmat(DirXml, '', Param.InputTable{1,3}, '.xml', OutputNomType,... + Param.IndexRange.first_i, Param.IndexRange.last_i, first_j, last_j); + save(t, filexml); +end + + %% direct processing on the current Matlab session if strcmp (RunMode,'local') @@ -1682,5 +1709,5 @@ end end -else +elseif ~strcmp(RunMode,'python') %% processing on a different session of the same computer (background) or cluster, create executable files batch_file_list=cell(NbProcess,1);% initiate the list of executable files @@ -1837,4 +1864,10 @@ system(oar_command); msgbox_uvmat('CONFIRMATION',[ActionName ' launched in cluster: press STATUS to see results']) + case 'python' + command = [ + 'LD_LIBRARY_PATH=$(echo $LD_LIBRARY_PATH | pyp "p.split('':'') | [s for s in p if ''matlab'' not in s] | '':''.join(p)") ' ... + 'python -m fluiddyn.postproc.uvmat ' filexml]; + % fprintf(['command:\n' command '\n\n']) + system(command, '-echo'); end @@ -1859,4 +1892,10 @@ if isfield(Param,'Pairs') Param=rmfield(Param,'Pairs'); %info Pairs not needed for output +end +if isfield(Param,'InputLine') + Param=rmfield(Param,'InputLine'); +end +if isfield(Param,'EditObject') + Param=rmfield(Param,'EditObject'); end Param.IndexRange.TimeSource=Param.IndexRange.TimeTable{end,1}; @@ -2328,5 +2367,5 @@ SeriesData=get(handles.series,'UserData'); % input line for which the field choice is relevant - iview=find(strcmp('netcdf',SeriesData.FileType)|strcmp('civx',SeriesData.FileType)|strcmp('civdata',SeriesData.FileType));% all nc files, icluding civ + iview=find(ismember(SeriesData.FileType,{'netcdf','civx','civdata'}));% all nc files, icluding civ hget_field=findobj(allchild(0),'name','get_field'); if ~isempty(hget_field) @@ -2337,5 +2376,5 @@ % check the existence of the first file in the series first_j=[];last_j=[];MinIndex_j=1;MaxIndex_j=1;%default setting for index j - if isfield(Param.IndexRange,'first_j');% if index j is used + if isfield(Param.IndexRange,'first_j');% if index j is used first_j=Param.IndexRange.first_j; last_j=Param.IndexRange.last_j; @@ -2346,6 +2385,15 @@ if isfield(Param.IndexRange,'PairString'); PairString=Param.IndexRange.PairString{iview}; end [i1,i2,j1,j2] = get_file_index(Param.IndexRange.first_i,first_j,PairString); - FirstFileName=fullfile_uvmat(InputTable{1},InputTable{2},InputTable{3},... - InputTable{5},InputTable{4},i1,i2,j1,j2); + LineIndex=iview(1); + if numel(iview)>1 + answer=msgbox_uvmat('INPUT_TXT',['select the line of the input table:' num2str(iview)] ,num2str(iview(1))); + LineIndex=str2num(answer); +% InputLine=str2num(get(handles.InputLine,'String')); +% if ismember(InputLine,iview) +% LineIndex=InputLine; +% end + end + FirstFileName=fullfile_uvmat(InputTable{LineIndex,1},InputTable{LineIndex,2},InputTable{LineIndex,3},... + InputTable{LineIndex,5},InputTable{LineIndex,4},i1,i2,j1,j2); if exist(FirstFileName,'file') ParamIn.Title='get_field: pick input variables and coordinates for series processing'; @@ -2372,13 +2420,13 @@ XName='X'; YName='y'; - set(handles.VelType,'visible','on') - set(handles.VelType_title,'visible','on') end set(handles.FieldName,'Value',1) set(handles.FieldName,'String',[FieldList; {'get_field...'}]); if ~strcmp(GetFieldData.FieldOption,'civdata...') - set(handles.FieldName,'Value',1:numel(FieldList))%select all input fields by default - set(handles.VelType,'visible','off') - set(handles.VelType_title,'visible','off') + if ~isempty(regexp(FieldList{1},'^vec')) + set(handles.FieldName,'Value',1) + else + set(handles.FieldName,'Value',1:numel(FieldList))%select all input fields by default + end XName=GetFieldData.Coordinates.Coord_x; YName=GetFieldData.Coordinates.Coord_y; @@ -2392,8 +2440,8 @@ TimeName=['att:' GetFieldData.Time.TimeName]; % update the time table - TimeTable{iview,2}=get_time(Param.IndexRange.MinIndex_i(iview),MinIndex_j,PairString,InputTable,SeriesData.FileInfo{iview},GetFieldData.Time.TimeName); % Min time - TimeTable{iview,3}=get_time(Param.IndexRange.first_i,first_j,PairString,InputTable,SeriesData.FileInfo{iview},GetFieldData.Time.TimeName); % first time - TimeTable{iview,4}=get_time(Param.IndexRange.last_i,last_j,PairString,InputTable,SeriesData.FileInfo{iview},GetFieldData.Time.TimeName); % last time - TimeTable{iview,5}=get_time(Param.IndexRange.MaxIndex_i(iview),MaxIndex_j,PairString,InputTable,SeriesData.FileInfo{iview},GetFieldData.Time.TimeName); % Max time + TimeTable{LineIndex,2}=get_time(Param.IndexRange.MinIndex_i(LineIndex),MinIndex_j,PairString,InputTable,SeriesData.FileInfo{LineIndex},GetFieldData.Time.TimeName); % Min time + TimeTable{LineIndex,3}=get_time(Param.IndexRange.first_i,first_j,PairString,InputTable,SeriesData.FileInfo{LineIndex},GetFieldData.Time.TimeName); % first time + TimeTable{LineIndex,4}=get_time(Param.IndexRange.last_i,last_j,PairString,InputTable,SeriesData.FileInfo{LineIndex},GetFieldData.Time.TimeName); % last time + TimeTable{LineIndex,5}=get_time(Param.IndexRange.MaxIndex_i(LineIndex),MaxIndex_j,PairString,InputTable,SeriesData.FileInfo{LineIndex},GetFieldData.Time.TimeName); % Max time case 'variable' set(handles.TimeName,'String',['var:' GetFieldData.Time.TimeName]) @@ -2409,9 +2457,11 @@ ParamIn.TimeDimName=GetFieldData.Time.TimeName; end - TimeTable{iview,1}=TimeName; + TimeTable{LineIndex,1}=TimeName; set(handles.TimeTable,'Data',TimeTable); end set(handles.Coord_x,'String',XName) set(handles.Coord_y,'String',YName) + set(handles.Coord_x,'Visible','on') + set(handles.Coord_y,'Visible','on') else msgbox_uvmat('ERROR',[FirstFileName ' does not exist']) @@ -3106,6 +3156,9 @@ ActionExtList=get(handles.ActionExt,'String'); ActionExt=ActionExtList{get(handles.ActionExt,'Value')}; -ActionList=get(handles.ActionName,'String'); -ActionName=ActionList{get(handles.ActionName,'Value')}; +if strcmp(ActionExt,'.py (in dev.)') + set(handles.RunMode,'Value',2) +end +% ActionList=get(handles.ActionName,'String'); +% ActionName=ActionList{get(handles.ActionName,'Value')}; function num_NbProcess_Callback(hObject, eventdata, handles) Index: /trunk/src/series/aver_stat.m =================================================================== --- /trunk/src/series/aver_stat.m (revision 866) +++ /trunk/src/series/aver_stat.m (revision 867) @@ -136,5 +136,5 @@ for iview=1:NbView if ~exist(filecell{iview,1}','file') - msgbox_uvmat('ERROR',['the first input file ' filecell{iview,1} ' does not exist']) + disp_uvmat('ERROR',['the first input file ' filecell{iview,1} ' does not exist'],checkrun) return end @@ -155,5 +155,5 @@ diff_time=max(max(diff(time))); if diff_time>0 - msgbox_uvmat('WARNING',['times of series differ by (max) ' num2str(diff_time)]) + disp_uvmat('WARNING',['times of series differ by (max) ' num2str(diff_time)],checkrun) end end @@ -176,9 +176,9 @@ FileExtOut='.nc';% write result as .nc files for netcdf inputs else - msgbox_uvmat('ERROR',['invalid file type input ' FileType{1}]) + disp_uvmat('ERROR',['invalid file type input ' FileType{1}],checkrun) return end if NbView==2 && ~isequal(CheckImage{1},CheckImage{2}) - msgbox_uvmat('ERROR','input must be two image series or two netcdf file series') + disp_uvmat('ERROR','input must be two image series or two netcdf file series',checkrun) return end @@ -216,5 +216,5 @@ % index_slice=i_slice:NbSlice:nbfield;% select file indices of the slice nbfiles=0; -nbmissing=0; +% nbmissing=0; %%%%%%%%%%%%%%%% loop on field indices %%%%%%%%%%%%%%%% @@ -275,5 +275,5 @@ [Field,errormsg]=proj_field(Field,Param.ProjObject); if ~isempty(errormsg) - msgbox_uvmat('ERROR',['error in aver_stat/proj_field:' errormsg]) + disp_uvmat('ERROR',['error in aver_stat/proj_field:' errormsg],checkrun) return end @@ -289,19 +289,36 @@ DataOut=Field;%default DataOut.Conventions='uvmat'; %suppress Conventions='uvmat/civdata' for civ input files - for ivar=1:length(Field.ListVarName) + errorvar=zeros(numel(Field.ListVarName));%index of errorflag associated to each variable + for ivar=1:numel(Field.ListVarName) VarName=Field.ListVarName{ivar}; - DataOut.(VarName)=double(DataOut.(VarName)); - end - else %current field - for ivar=1:length(Field.ListVarName) - VarName=Field.ListVarName{ivar}; - sizmean=size(DataOut.(VarName)); - siz=size(Field.(VarName)); - if ~isequal(DataOut.(VarName),0)&& ~isequal(siz,sizmean) - disp_uvmat('ERROR',['unequal size of input field ' VarName ', need to project on a grid'],checkrun) - return + DataOut.(VarName)=zeros(size(DataOut.(VarName)));% initiate each field to zero + NbData.(VarName)=zeros(size(DataOut.(VarName)));% initiate the nbre of good data to zero + for iivar=1:length(Field.VarAttribute) + if isequal(Field.VarDimName{iivar},Field.VarDimName{ivar})&& isfield(Field.VarAttribute{iivar},'Role')... + && strcmp(Field.VarAttribute{iivar}.Role,'errorflag') + errorvar(ivar)=iivar; % index of the errorflag variable corresponding to ivar + end + end + end + DataOut.ListVarName(errorvar(errorvar~=0))=[]; %remove errorflag from result + DataOut.VarDimName(errorvar(errorvar~=0))=[]; %remove errorflag from result + DataOut.VarAttribute(errorvar(errorvar~=0))=[]; %remove errorflag from result + end %current field + for ivar=1:length(DataOut.ListVarName) + VarName=Field.ListVarName{ivar}; + sizmean=size(DataOut.(VarName)); + siz=size(Field.(VarName)); + if ~isequal(DataOut.(VarName),0)&& ~isequal(siz,sizmean) + disp_uvmat('ERROR',['unequal size of input field ' VarName ', need to project on a grid'],checkrun) + return + else + if errorvar(ivar)==0 + check_bad=isnan(Field.(VarName));%=0 for NaN data values, 1 else else - DataOut.(VarName)=DataOut.(VarName)+ double(Field.(VarName)); % update the sum + check_bad=isnan(Field.(VarName)) | Field.(Field.ListVarName{errorvar(ivar)})~=0;%=0 for NaN or error flagged data values, 1 else end + Field.(VarName)(check_bad)=0; %set to zero NaN or data marked by error flag + DataOut.(VarName)=DataOut.(VarName)+ double(Field.(VarName)); % update the sum + NbData.(VarName)=NbData.(VarName)+ ~check_bad;% records the number of data for each point end end @@ -315,8 +332,9 @@ for ivar=1:length(Field.ListVarName) VarName=Field.ListVarName{ivar}; - DataOut.(VarName)=DataOut.(VarName)/nbfiles; % normalize the mean -end + DataOut.(VarName)=DataOut.(VarName)./NbData.(VarName); % normalize the mean +end +nbmissing=NbField-nbfiles; if nbmissing~=0 - msgbox_uvmat('WARNING',[num2str(nbmissing) ' input files are missing or skipted']) + disp_uvmat('WARNING',[num2str(nbmissing) ' input files are missing or skipted'],checkrun) end if isempty(time) % time is read from files Index: /trunk/src/series/histo_diff.m =================================================================== --- /trunk/src/series/histo_diff.m (revision 867) +++ /trunk/src/series/histo_diff.m (revision 867) @@ -0,0 +1,361 @@ +%'histo_diff': calculate the histogram of the difference between two fields +%------------------------------------------------------------------------ +% function ParamOut=histo_diff(Param) +% +%%%%%%%%%%% GENERAL TO ALL SERIES ACTION FCTS %%%%%%%%%%%%%%%%%%%%%%%%%%% +% +%OUTPUT +% ParamOut: sets options in the GUI series.fig needed for the function +% +%INPUT: +% In run mode, the input parameters are given as a Matlab structure Param copied from the GUI series. +% In batch mode, Param is the name of the corresponding xml file containing the same information +% when Param.Action.RUN=0 (as activated when the current Action is selected +% in series), the function ouput paramOut set the activation of the needed GUI elements +% +% Param contains the elements:(use the menu bar command 'export/GUI config' in series to +% see the current structure Param) +% .InputTable: cell of input file names, (several lines for multiple input) +% each line decomposed as {RootPath,SubDir,Rootfile,NomType,Extension} +% .OutputSubDir: name of the subdirectory for data outputs +% .OutputDirExt: directory extension for data outputs +% .Action: .ActionName: name of the current activated function +% .ActionPath: path of the current activated function +% .ActionExt: fct extension ('.m', Matlab fct, '.sh', compiled Matlab fct +% .RUN =0 for GUI input, =1 for function activation +% .RunMode='local','background', 'cluster': type of function use +% +% .IndexRange: set the file or frame indices on which the action must be performed +% .FieldTransform: .TransformName: name of the selected transform function +% .TransformPath: path of the selected transform function +% .InputFields: sub structure describing the input fields withfields +% .FieldName: name(s) of the field +% .VelType: velocity type +% .FieldName_1: name of the second field in case of two input series +% .VelType_1: velocity type of the second field in case of two input series +% .Coord_y: name of y coordinate variable +% .Coord_x: name of x coordinate variable +% .ProjObject: %sub structure describing a projection object (read from ancillary GUI set_object) +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%======================================================================= +% Copyright 2008-2014, LEGI UMR 5519 / CNRS UJF G-INP, Grenoble, France +% http://www.legi.grenoble-inp.fr +% Joel.Sommeria - Joel.Sommeria (A) legi.cnrs.fr +% +% This file is part of the toolbox UVMAT. +% +% UVMAT is free software; you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published +% by the Free Software Foundation; either version 2 of the license, +% or (at your option) any later version. +% +% UVMAT is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License (see LICENSE.txt) for more details. +%======================================================================= + +function ParamOut=histo_diff(Param) + +%% set the input elements needed on the GUI series when the action is selected in the menu ActionName +if isstruct(Param) && isequal(Param.Action.RUN,0)% function activated from the GUI series but not RUN + ParamOut.AllowInputSort='off';% allow alphabetic sorting of the list of input file SubDir (options 'off'/'on', 'off' by default) + ParamOut.WholeIndexRange='off';% prescribes the file index ranges from min to max (options 'off'/'on', 'off' by default) + ParamOut.NbSlice='on'; %nbre of slices ('off' by default) + ParamOut.VelType='two';% menu for selecting the velocity type (options 'off'/'one'/'two', 'off' by default) + ParamOut.FieldName='two';% menu for selecting the field (s) in the input file(options 'off'/'one'/'two', 'off' by default) + ParamOut.FieldTransform = 'on';%can use a transform function + ParamOut.ProjObject='on';%can use projection object(option 'off'/'on', + ParamOut.Mask='off';%can use mask option (option 'off'/'on', 'off' by default) + ParamOut.OutputDirExt='.histodif';%set the output dir extension + ParamOut.OutputFileMode='NbSlice';% '=NbInput': 1 output file per input file index, '=NbInput_i': 1 file per input file index i, '=NbSlice': 1 file per slice + % check the existence of the first file in the series + first_j=[]; + if isfield(Param.IndexRange,'first_j'); first_j=Param.IndexRange.first_j; end + last_j=[]; + if isfield(Param.IndexRange,'last_j'); last_j=Param.IndexRange.last_j; end + PairString=''; + if isfield(Param.IndexRange,'PairString'); PairString=Param.IndexRange.PairString; end + [i1,i2,j1,j2] = get_file_index(Param.IndexRange.first_i,first_j,PairString); + FirstFileName=fullfile_uvmat(Param.InputTable{1,1},Param.InputTable{1,2},Param.InputTable{1,3},... + Param.InputTable{1,5},Param.InputTable{1,4},i1,i2,j1,j2); + if ~exist(FirstFileName,'file') + msgbox_uvmat('WARNING',['the first input file ' FirstFileName ' does not exist']) + else + [i1,i2,j1,j2] = get_file_index(Param.IndexRange.last_i,last_j,PairString); + LastFileName=fullfile_uvmat(Param.InputTable{1,1},Param.InputTable{1,2},Param.InputTable{1,3},... + Param.InputTable{1,5},Param.InputTable{1,4},i1,i2,j1,j2); + if ~exist(FirstFileName,'file') + msgbox_uvmat('WARNING',['the last input file ' LastFileName ' does not exist']) + end + end + return +end + +%%%%%%%%%%%% STANDARD PART %%%%%%%%%%%% +ParamOut=[];%default output +%% read input parameters from an xml file if input is a file name (batch mode) +checkrun=1; +if ischar(Param) + Param=xml2struct(Param);% read Param as input file (batch case) + checkrun=0; +end +hseries=findobj(allchild(0),'Tag','series'); +RUNHandle=findobj(hseries,'Tag','RUN');%handle of RUN button in GUI series +WaitbarHandle=findobj(hseries,'Tag','Waitbar');%handle of waitbar in GUI series + +%% define the directory for result file (with path=RootPath{1}) +OutputDir=[Param.OutputSubDir Param.OutputDirExt]; + +%% root input file(s) name, type and index series +RootPath=Param.InputTable(:,1); +RootFile=Param.InputTable(:,3); +SubDir=Param.InputTable(:,2); +NomType=Param.InputTable(:,4); +FileExt=Param.InputTable(:,5); +hdisp=disp_uvmat('WAITING...','checking the file series',checkrun); +[filecell,i1_series,i2_series,j1_series,j2_series]=get_file_series(Param); +if ~isempty(hdisp),delete(hdisp),end; +%%%%%%%%%%%% +% The cell array filecell is the list of input file names, while +% filecell{iview,fileindex}: +% iview: line in the table corresponding to a given file series +% fileindex: file index within the file series, +% 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 +% i1_series(iview,fileindex) expresses the same indices as a 1D array in file indices +%%%%%%%%%%%% +NbView=numel(i1_series);%number of input file series (lines in InputTable) +NbField_j=size(i1_series{1},1); %nb of fields for the j index (bursts or volume slices) +NbField_i=size(i1_series{1},2); %nb of fields for the i index +NbField=NbField_j*NbField_i; %total number of fields + +%% determine the file type on each line from the first input file +ImageTypeOptions={'image','multimage','mmreader','video'}; +NcTypeOptions={'netcdf','civx','civdata'}; +for iview=1:NbView + if ~exist(filecell{iview,1}','file') + msgbox_uvmat('ERROR',['the first input file ' filecell{iview,1} ' does not exist']) + return + end + [FileInfo{iview},MovieObject{iview}]=get_file_info(filecell{iview,1}); + FileType{iview}=FileInfo{iview}.FileType; + CheckImage{iview}=~isempty(find(strcmp(FileType{iview},ImageTypeOptions)));% =1 for images + CheckNc{iview}=~isempty(find(strcmp(FileType{iview},NcTypeOptions)));% =1 for netcdf files + if ~isempty(j1_series{iview}) + frame_index{iview}=j1_series{iview}; + else + frame_index{iview}=i1_series{iview}; + end +end + +%% calibration data and timing: read the ImaDoc files +[XmlData,NbSlice_calib,time,errormsg]=read_multimadoc(RootPath,SubDir,RootFile,FileExt,i1_series,i2_series,j1_series,j2_series); +if size(time,1)>1 + diff_time=max(max(diff(time))); + if diff_time>0 + msgbox_uvmat('WARNING',['times of series differ by (max) ' num2str(diff_time)]) + end +end + +%% coordinate transform or other user defined transform +transform_fct='';%default +if isfield(Param,'FieldTransform')&&~isempty(Param.FieldTransform.TransformName) + addpath(Param.FieldTransform.TransformPath) + transform_fct=str2func(Param.FieldTransform.TransformName); + rmpath(Param.FieldTransform.TransformPath) +end +if ~(isfield(Param,'ProjObject')&& ismember(Param.ProjObject.ProjMode,{'polygon','rectangle','ellipse'})) + disp_uvmat('ERROR','projection in a closed domain: polygon, rectangle or ellipse is needed for histograms',checkrun) +end +%%%%%%%%%%%% END STANDARD PART %%%%%%%%%%%% + % EDIT FROM HERE + +%% check the validity of input file types +if CheckImage{1} + FileExtOut='.png'; % write result as .png images for image inputs +elseif CheckNc{1} + FileExtOut='.nc';% write result as .nc files for netcdf inputs +else + disp_uvmat('ERROR',['invalid file type input ' FileType{1}],checkrun) + return +end +if NbView==2 && ~isequal(CheckImage{1},CheckImage{2}) + msgbox_uvmat('ERROR','input must be two image series or two netcdf file series') + return +end + +%% settings for the output file +NomTypeOut=nomtype2pair(NomType{1});% determine the index nomenclature type for the output file +first_i=i1_series{1}(1); +last_i=i1_series{1}(end); +if isempty(j1_series{1})% if there is no second index j + first_j=1;last_j=1; +else + first_j=j1_series{1}(1); + last_j=j1_series{1}(end); +end + +%% Set field names and velocity types +InputFields{1}=[];%default (case of images) +if isfield(Param,'InputFields') + InputFields{1}=Param.InputFields; +end +if NbView==2 + InputFields{2}=[];%default (case of images) + if isfield(Param,'InputFields') + InputFields{2}=Param.InputFields{1};%default + if isfield(Param.InputFields,'FieldName_1') + InputFields{2}.FieldName=Param.InputFields.FieldName_1; + if isfield(Param.InputFields,'VelType_1') + InputFields{2}.VelType=Param.InputFields.VelType_1; + end + end + end +end + +% for i_slice=1:NbSlice +% index_slice=i_slice:NbSlice:nbfield;% select file indices of the slice +% nbfiles=0; +% nbmissing=0; + +%%%%%%%%%%%%%%%% loop on field indices %%%%%%%%%%%%%%%% +for index=1:NbField + update_waitbar(WaitbarHandle,index/NbField) + if ~isempty(RUNHandle)&& ~strcmp(get(RUNHandle,'BusyAction'),'queue') + disp('program stopped by user') + break + end + + %%%%%%%%%%%%%%%% loop on views (input lines) %%%%%%%%%%%%%%%% + for iview=1:NbView + % reading input file(s) + [Data{iview},tild,errormsg] = read_field(filecell{iview,index},FileType{iview},InputFields{iview},frame_index{iview}(index)); + if ~isempty(errormsg) + errormsg=['error of input reading: ' errormsg]; + break + end + if ~isempty(NbSlice_calib) + Data{iview}.ZIndex=mod(i1_series{iview}(index)-1,NbSlice_calib{iview})+1;%Zindex for phys transform + end + end + %%%%%%%%%%%%%%%% end loop on views (input lines) %%%%%%%%%%%%%%%% + %%%%%%%%%%%% END STANDARD PART %%%%%%%%%%%% + % EDIT FROM HERE + + if isempty(errormsg) + Field=Data{1}; % default input field structure + %% coordinate transform (or other user defined transform) + if ~isempty(transform_fct) + switch nargin(transform_fct) + case 4 + if length(Data)==2 + Field=transform_fct(Data{1},XmlData{1},Data{2},XmlData{2}); + else + Field=transform_fct(Data{1},XmlData{1}); + end + case 3 + if length(Data)==2 + Field=transform_fct(Data{1},XmlData{1},Data{2}); + else + Field=transform_fct(Data{1},XmlData{1}); + end + case 2 + Field=transform_fct(Data{1},XmlData{1}); + case 1 + Field=transform_fct(Data{1}); + end + end + + %% calculate tps coefficients if needed +% if isfield(Param,'ProjObject')&&isfield(Param.ProjObject,'ProjMode')&& strcmp(Param.ProjObject.ProjMode,'interp_tps') +% Field=tps_coeff_field(Field,check_proj_tps); +% end + + %field projection on an object + [Field,errormsg]=proj_field(Field,Param.ProjObject); + if ~isempty(errormsg) + disp_uvmat('ERROR',['error in histo_diff/proj_field:' errormsg],checkrun) + return + end +% nbfiles=nbfiles+1; + + %%%%%%%%%%%% MAIN RUNNING OPERATIONS %%%%%%%%%%%% + if index==1 %first field + time_1=[]; + if isfield(Field,'Time') + time_1=Field.Time(1); + end + DataOut=Field;%default + DataOut.Conventions='uvmat'; %suppress Conventions='uvmat/civdata' for civ input files + errorvar=zeros(numel(Field.ListVarName));%index of errorflag associated to each variable + for ivar=1:numel(Field.ListVarName) + VarName=Field.ListVarName{ivar}; + DataOut.(VarName)=zeros(size(DataOut.(VarName)));% initiate each field to zero + NbData.(VarName)=zeros(size(DataOut.(VarName)));% initiate the nbre of good data to zero + for iivar=1:length(Field.VarAttribute) + if isequal(Field.VarDimName{iivar},Field.VarDimName{ivar})&& isfield(Field.VarAttribute{iivar},'Role')... + && strcmp(Field.VarAttribute{iivar}.Role,'errorflag') + errorvar(ivar)=iivar; % index of the errorflag variable corresponding to ivar + end + end + end + DataOut.ListVarName(errorvar(errorvar~=0))=[]; %remove errorflag from result + DataOut.VarDimName(errorvar(errorvar~=0))=[]; %remove errorflag from result + DataOut.VarAttribute(errorvar(errorvar~=0))=[]; %remove errorflag from result + end %current field + for ivar=1:length(DataOut.ListVarName) + VarName=Field.ListVarName{ivar}; + sizmean=size(DataOut.(VarName)); + siz=size(Field.(VarName)); + if ~isequal(DataOut.(VarName),0)&& ~isequal(siz,sizmean) + disp_uvmat('ERROR',['unequal size of input field ' VarName ', need to project on a grid'],checkrun) + return + else + if errorvar(ivar)==0 + check_bad=isnan(Field.(VarName));%=0 for NaN data values, 1 else + else + check_bad=isnan(Field.(VarName)) | Field.(Field.ListVarName{errorvar(ivar)})~=0;%=0 for NaN or error flagged data values, 1 else + end + Field.(VarName)(check_bad)=0; %set to zero NaN or data marked by error flag + DataOut.(VarName)=DataOut.(VarName)+ double(Field.(VarName)); % update the sum + NbData.(VarName)=NbData.(VarName)+ ~check_bad;% records the number of data for each point + end + end + %%%%%%%%%%%% END MAIN RUNNING OPERATIONS %%%%%%%%%%%% + else + disp(errormsg) + end +end +%%%%%%%%%%%%%%%% end loop on field indices %%%%%%%%%%%%%%%% + +for ivar=1:length(Field.ListVarName) + VarName=Field.ListVarName{ivar}; + DataOut.(VarName)=DataOut.(VarName)./NbData.(VarName); % normalize the mean +end +if isempty(time) % time is read from files + if isfield(Field,'Time') + time_end=Field.Time(1);%last time read + if ~isempty(time_1) + DataOut.Time=time_1; + DataOut.Time_end=time_end; + end + end +else % time from ImaDoc prevails if it exists + DataOut.Time=time(1); + DataOut.Time_end=time(end); +end + +%% writing the result file +OutputFile=fullfile_uvmat(RootPath{1},OutputDir,RootFile{1},FileExtOut,NomTypeOut,first_i,last_i,first_j,last_j); +errormsg=struct2nc(OutputFile,DataOut); %save result file +if isempty(errormsg) + disp([OutputFile ' written']); +else + disp(['error in writting result file: ' errormsg]) +end + +%% open the result file with uvmat (in RUN mode) +if checkrun + uvmat(OutputFile)% open the last result file with uvmat +end Index: /trunk/src/series/merge_proj.m =================================================================== --- /trunk/src/series/merge_proj.m (revision 866) +++ /trunk/src/series/merge_proj.m (revision 867) @@ -196,5 +196,5 @@ FileExtOut='.nc'; %netcdf output end -NomTypeOut=NomType;% output file index will indicate the first and last ref index in the series +%NomTypeOut=NomType;% output file index will indicate the first and last ref index in the series RootFileOut=RootFile{1}; for iview=2:NbView @@ -403,8 +403,8 @@ end errormsg=''; -MergeData=Data{1};% merged field= first field by default, reproduces the glabal attributes of the first field +MergeData=Data{1};% merged field= first field by default, reproduces the global attributes of the first field NbView=length(Data); -if NbView==1 - return +if NbView==1% if there is only one field, just reproduce it in MergeData + return end @@ -419,8 +419,7 @@ if NbDim(icell)~=1 % skip field cells which are of dim 1 switch CellInfo{icell}.CoordType - case 'scattered' %case of input fields with unstructured coordinates: just concacene data + case 'scattered' %case of input fields with unstructured coordinates: just concatene data for ivar=CellInfo{icell}.VarIndex % indices of the selected variables in the list FieldData.ListVarName VarName=Data{1}.ListVarName{ivar}; - %MergeData=Data{1};% merged field= first field by default, reproduces the glabal attributes of the first field for iview=2:NbView MergeData.(VarName)=[MergeData.(VarName); Data{iview}.(VarName)]; @@ -431,4 +430,7 @@ if isfield(CellInfo{icell},'VarIndex_errorflag') && ~isempty(CellInfo{icell}.VarIndex_errorflag) FFName=Data{1}.ListVarName{CellInfo{icell}.VarIndex_errorflag};% name of errorflag variable + MergeData.ListVarName(CellInfo{icell}.VarIndex_errorflag)=[];%remove error flag variable in MergeData (will use NaN instead) + MergeData.VarDimName(CellInfo{icell}.VarIndex_errorflag)=[]; + MergeData.VarAttribute(CellInfo{icell}.VarIndex_errorflag)=[]; end % select good data on each view @@ -441,25 +443,25 @@ check_bad=isnan(Data{iview}.(VarName)) | Data{iview}.(FFName)~=0;%=0 for NaN or error flagged data values, 1 else end - Data{iview}.(VarName)(check_bad)=0; %set to zero NaN or masked data + Data{iview}.(VarName)(check_bad)=0; %set to zero NaN or data marked by error flag if iview==1 - MergeData.(VarName)=Data{1}.(VarName);% correct the field of MergeData + %MergeData.(VarName)=Data{1}.(VarName);% initiate MergeData with the first field NbAver=~check_bad;% initiate NbAver: the nbre of good data for each point elseif size(Data{iview}.(VarName))~=size(MergeData.(VarName)) errormsg='sizes of the input matrices do not agree, need to interpolate on a common grid using a projection object'; return - else - - MergeData.(VarName)=MergeData.(VarName) + Data{iview}.(VarName);%add data + else + MergeData.(VarName)=MergeData.(VarName) +double(Data{iview}.(VarName));%add data NbAver=NbAver + ~check_bad;% add 1 for good data, 0 else end end MergeData.(VarName)(NbAver~=0)=MergeData.(VarName)(NbAver~=0)./NbAver(NbAver~=0);% take average of defined data at each point + MergeData.(VarName)(NbAver==0)=NaN;% set to NaN the points with no good data end end - if isempty(FFName) - FFName='FF'; - end - MergeData.(FFName)(NbAver~=0)=0;% flag to 1 undefined summed data - MergeData.(FFName)(NbAver==0)=1;% flag to 1 undefined summed data +% if isempty(FFName) +% FFName='FF'; +% end +% MergeData.(FFName)(NbAver~=0)=0;% flag to 1 undefined summed data +% MergeData.(FFName)(NbAver==0)=1;% flag to 1 undefined summed data end end Index: /trunk/src/series/time_series.m =================================================================== --- /trunk/src/series/time_series.m (revision 866) +++ /trunk/src/series/time_series.m (revision 867) @@ -207,11 +207,11 @@ %% Set field names and velocity types InputFields{1}=[];%default (case of images) +if nbview==2 + InputFields{2}=[];%default (case of images) +end if isfield(Param,'InputFields') InputFields{1}=Param.InputFields; -end -if nbview==2 - InputFields{2}=[];%default (case of images) - if isfield(Param,'InputFields') - InputFields{2}=Param.InputFields{1};%default + if nbview==2 + InputFields{2}=Param.InputFields;%default if isfield(Param.InputFields,'FieldName_1') InputFields{2}.FieldName=Param.InputFields.FieldName_1; Index: /trunk/src/uvmat.m =================================================================== --- /trunk/src/uvmat.m (revision 866) +++ /trunk/src/uvmat.m (revision 867) @@ -2818,7 +2818,7 @@ elseif ref_j_1<0 errormsg='minimum j index reached'; - elseif ref_i_1+1>size(UvData.i1_series{2},3) + elseif ref_i_1+1>size(UvData.i1_series{2},3)&&~isempty(InputFile.NomType_1) errormsg='maximum i index reached for the second series (reload the input file to update the index bound)'; - elseif ref_j_1+1>size(UvData.i1_series{2},2) + elseif ref_j_1+1>size(UvData.i1_series{2},2)&&~isempty(InputFile.NomType_1) errormsg='maximum j index reached for the second series(reload the input file to update the index bound)'; end @@ -2826,20 +2826,22 @@ siz=size(UvData.i1_series{2}); ref_indices=ref_i_1*siz(1)*siz(2)+ref_j_1*siz(1)+1:ref_i_1*siz(1)*siz(2)+(ref_j_1+1)*siz(1); - i1_subseries=UvData.i1_series{2}(ref_indices); - ref_indices=ref_indices(i1_subseries>0); - if isempty(ref_indices)% case of pairs (free index i) - ref_indices=ref_i_1*siz(1)*siz(2)+1:(ref_i_1+1)*siz(1)*siz(2); + if ~isempty(InputFile.NomType_1) i1_subseries=UvData.i1_series{2}(ref_indices); ref_indices=ref_indices(i1_subseries>0); - end - i1_1=UvData.i1_series{2}(ref_indices(end)); - if ~isempty(UvData.i2_series{2}) - i2_1=UvData.i2_series{2}(ref_indices(end)); - end - if ~isempty(UvData.j1_series{2}) - j1_1=UvData.j1_series{2}(ref_indices(end)); - end - if ~isempty(UvData.j2_series{2}) - j2_1=UvData.j2_series{1}(ref_indices(end)); + if isempty(ref_indices)% case of pairs (free index i) + ref_indices=ref_i_1*siz(1)*siz(2)+1:(ref_i_1+1)*siz(1)*siz(2); + i1_subseries=UvData.i1_series{2}(ref_indices); + ref_indices=ref_indices(i1_subseries>0); + end + i1_1=UvData.i1_series{2}(ref_indices(end)); + if ~isempty(UvData.i2_series{2}) + i2_1=UvData.i2_series{2}(ref_indices(end)); + end + if ~isempty(UvData.j1_series{2}) + j1_1=UvData.j1_series{2}(ref_indices(end)); + end + if ~isempty(UvData.j2_series{2}) + j2_1=UvData.j2_series{1}(ref_indices(end)); + end end else% the second series (if needed) is the same file as the first @@ -3600,5 +3602,4 @@ %% get bounds and dimensions of the input field UvData.Field=find_field_bounds(UvData.Field); - testnewseries=UvData.NewSeries; UvData.NewSeries=0;% put to 0 the test for a new field series (set by RootPath_callback)