Ignore:
Timestamp:
Mar 4, 2015, 12:01:38 AM (9 years ago)
Author:
sommeria
Message:

various bug fixes

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/series/time_series.m

    r875 r880  
    9090        if ismember(SeriesData.ProjObject.ProjMode,{'inside','outside'})
    9191            answer=msgbox_uvmat('INPUT_TXT','set bin size for histograms (or keep ''auto'' by default)?','auto');
    92             ParamOut.ActionInput.VarMesh=str2double(answer);
     92            ParamOut.ActionInput.VarMesh=str2num(answer);
    9393        end
    9494    end
     
    295295nbfile=0;% not used , to check
    296296nbmissing=0;
    297 VarMesh=NaN;
     297VarMesh=[];
    298298checkhisto=0;
    299299if isfield(Param,'ProjObject') && ismember(Param.ProjObject.ProjMode,{'inside','outside'})
    300300    checkhisto=1;
    301     if isfield(Param.ActionInput,'VarMesh')%case of histograms
    302     VarMesh=Param.ActionInput.VarMesh;
     301    if isfield(Param,'ActionInput') && isfield(Param.ActionInput,'VarMesh')%case of histograms
     302        VarMesh=Param.ActionInput.VarMesh;
    303303    end
    304304end
     
    320320        if ~isempty(errormsg)
    321321            errormsg=['time_series / read_field / ' errormsg];
    322             display(errormsg)
    323             break
     322            disp(errormsg)
     323            break% exit the loop on iview
    324324        end
    325325        if ~isempty(NbSlice_calib)
     
    327327        end
    328328    end
    329     if isempty(errormsg)
    330         Field=Data{1}; % default input field structure
    331         % coordinate transform (or other user defined transform)
    332         if ~isempty(transform_fct)
    333             switch nargin(transform_fct)
    334                 case 4
    335                     if length(Data)==2
    336                         Field=transform_fct(Data{1},XmlData{1},Data{2},XmlData{2});
    337                     else
    338                         Field=transform_fct(Data{1},XmlData{1});
     329    if ~isempty(errormsg)
     330        continue %in case of input error skip the current input field, go to next index
     331    end
     332    Field=Data{1}; % default input field structure
     333    % coordinate transform (or other user defined transform)
     334    if ~isempty(transform_fct)
     335        switch nargin(transform_fct)
     336            case 4
     337                if length(Data)==2
     338                    Field=transform_fct(Data{1},XmlData{1},Data{2},XmlData{2});
     339                else
     340                    Field=transform_fct(Data{1},XmlData{1});
     341                end
     342            case 3
     343                if length(Data)==2
     344                    Field=transform_fct(Data{1},XmlData{1},Data{2});
     345                else
     346                    Field=transform_fct(Data{1},XmlData{1});
     347                end
     348            case 2
     349                Field=transform_fct(Data{1},XmlData{1});
     350            case 1
     351                Field=transform_fct(Data{1});
     352        end
     353    end
     354   
     355    %field projection on an object
     356    if Param.CheckObject
     357        % calculate tps coefficients if needed
     358        if isfield(Param.ProjObject,'ProjMode')&& strcmp(Param.ProjObject.ProjMode,'interp_tps')
     359            Field=tps_coeff_field(Field,check_proj_tps);
     360        end
     361        [Field,errormsg]=proj_field(Field,Param.ProjObject,VarMesh);
     362        if ~isempty(errormsg)
     363            disp_uvmat('ERROR',['time_series / proj_field / ' errormsg],checkrun)
     364            return
     365        end
     366    end
     367    %         nbfile=nbfile+1;
     368   
     369    % initiate the time series at the first iteration
     370    if index==1
     371        % stop program if the first field reading is in error
     372        if ~isempty(errormsg)
     373            disp_uvmat('ERROR',['time_series / sub_field / ' errormsg],checkrun)
     374            return
     375        end
     376        if ~isfield(Field,'ListVarName')
     377            disp_uvmat('ERROR','no variable in the projected field',checkrun)
     378            return
     379        end
     380        DataOut=Field;%output reproduced the first projected field by default
     381        nbvar=length(Field.ListVarName);
     382        if nbvar==0
     383            disp_uvmat('ERROR','no input variable selected',checkrun)
     384            return
     385        end
     386        if checkhisto%case of histograms
     387            testsum=zeros(1,nbvar);%initiate flag for action on each variable
     388            for ivar=1:numel(Field.ListVarName)% list of variable names before projection (histogram)
     389                VarName=Field.ListVarName{ivar};
     390                if isfield(Data{1},VarName)
     391                    DataOut.ListVarName=[DataOut.ListVarName {[VarName 'Histo']}];
     392                    DataOut.VarDimName=[DataOut.VarDimName {{'Time',VarName}}];
     393%                     if isfield(DataOut.VarAttribute{ivar},'Role')
     394%                     DataOut.VarAttribute{ivar}=rmfield(DataOut.VarAttribute{ivar},'Role');
     395%                     end
     396                    StatName=pdf2stat;% get the names of statistical quantities to calcuilate at each time
     397                    for istat=1:numel(StatName)
     398                        DataOut.ListVarName=[DataOut.ListVarName {[VarName StatName{istat}]}];
     399                        DataOut.VarDimName=[DataOut.VarDimName {'Time'}];
    339400                    end
    340                 case 3
    341                     if length(Data)==2
    342                         Field=transform_fct(Data{1},XmlData{1},Data{2});
    343                     else
    344                         Field=transform_fct(Data{1},XmlData{1});
    345                     end
    346                 case 2
    347                     Field=transform_fct(Data{1},XmlData{1});
    348                 case 1
    349                     Field=transform_fct(Data{1});
     401                    testsum(ivar)=1;
     402                    DataOut.(VarName)=Field.(VarName);
     403                    DataOut.([VarName 'Histo'])=zeros([nbfield numel(DataOut.(VarName))]);
     404                    VarMesh=Field.(VarName)(2)-Field.(VarName)(1);
     405                end
    350406            end
    351         end
    352        
    353         %field projection on an object
    354         if Param.CheckObject
    355             % calculate tps coefficients if needed
    356             if isfield(Param.ProjObject,'ProjMode')&& strcmp(Param.ProjObject.ProjMode,'interp_tps')
    357                 Field=tps_coeff_field(Field,check_proj_tps);
    358             end
    359             [Field,errormsg]=proj_field(Field,Param.ProjObject,VarMesh);
    360             if ~isempty(errormsg)
    361                 msgbox_uvmat('ERROR',['time_series / proj_field / ' errormsg])
    362                 return
    363             end
    364         end
    365         %         nbfile=nbfile+1;
    366        
    367         % initiate the time series at the first iteration
    368         if index==1
    369             % stop program if the first field reading is in error
    370             if ~isempty(errormsg)
    371                 disp_uvmat('ERROR',['time_series / sub_field / ' errormsg],checkrun)
    372                 return
    373             end
    374             DataOut=Field;%default
    375             nbvar=length(Field.ListVarName);
    376             if nbvar==0
    377                 disp_uvmat('ERROR','no input variable selected',checkrun)
    378                 return
    379             end
    380             if checkhisto%case of histograms
    381                 testsum=zeros(1,nbvar);%initiate flag for action on each variable
    382                 for ivar=1:numel(Field.ListVarName)% list of variable names before projection (histogram)
    383                     VarName=Field.ListVarName{ivar};
    384                     if isfield(Data{1},VarName)
    385                         testsum(ivar)=1;
    386                         DataOut.(VarName)=Field.(VarName);
    387                         DataOut.([VarName 'Histo'])=zeros([nbfield numel(DataOut.(VarName))]);
    388                         VarMesh=Field.(VarName)(2)-Field.(VarName)(1);
    389                     end
    390                 end
    391                 disp(['mesh for histogram = ' num2str(VarMesh)])
    392             else
    393                 testsum=2*ones(1,nbvar);%initiate flag for action on each variable
    394                 if isfield(Field,'VarAttribute') % look for coordinate and flag variables
    395                     for ivar=1:nbvar
    396                         if length(Field.VarAttribute)>=ivar && isfield(Field.VarAttribute{ivar},'Role')
    397                             var_role=Field.VarAttribute{ivar}.Role;%'role' of the variable
    398                             if isequal(var_role,'errorflag')
    399                                 disp_uvmat('ERROR','do not handle error flags in time series',checkrun)
    400                                 return
    401                             end
    402                             if isequal(var_role,'warnflag')
    403                                 testsum(ivar)=0;  % not recorded variable
    404                                 eval(['DataOut=rmfield(DataOut,''' Field.ListVarName{ivar} ''');']);%remove variable
    405                             end
    406                             if strcmp(var_role,'coord_x')||strcmp(var_role,'coord_y')||strcmp(var_role,'coord_z')||strcmp(var_role,'coord')
    407                                 testsum(ivar)=1; %constant coordinates, record without time evolution
    408                             end
     407            disp(['mesh for histogram = ' num2str(VarMesh)])
     408        else
     409            testsum=2*ones(1,nbvar);%initiate flag for action on each variable
     410            if isfield(Field,'VarAttribute') % look for coordinate and flag variables
     411                for ivar=1:nbvar
     412                    if length(Field.VarAttribute)>=ivar && isfield(Field.VarAttribute{ivar},'Role')
     413                        var_role=Field.VarAttribute{ivar}.Role;%'role' of the variable
     414                        if isequal(var_role,'errorflag')
     415                            disp_uvmat('ERROR','do not handle error flags in time series',checkrun)
     416                            return
    409417                        end
    410                         % check whether the variable ivar is a dimension variable
    411                         DimCell=Field.VarDimName{ivar};
    412                         if ischar(DimCell)
    413                             DimCell={DimCell};
     418                        if isequal(var_role,'warnflag')
     419                            testsum(ivar)=0;  % not recorded variable
     420                            eval(['DataOut=rmfield(DataOut,''' Field.ListVarName{ivar} ''');']);%remove variable
    414421                        end
    415                         if numel(DimCell)==1 && isequal(Field.ListVarName{ivar},DimCell{1})%detect dimension variables
    416                             testsum(ivar)=1;
     422                        if strcmp(var_role,'coord_x')||strcmp(var_role,'coord_y')||strcmp(var_role,'coord_z')||strcmp(var_role,'coord')
     423                            testsum(ivar)=1; %constant coordinates, record without time evolution
    417424                        end
    418425                    end
    419                 end
    420                 for ivar=1:nbvar
    421                     if testsum(ivar)==2
    422                         VarName=Field.ListVarName{ivar};
    423                         siz=size(Field.(VarName));
    424                         DataOut.(VarName)=zeros([nbfield siz]);
     426                    % check whether the variable ivar is a dimension variable
     427                    DimCell=Field.VarDimName{ivar};
     428                    if ischar(DimCell)
     429                        DimCell={DimCell};
    425430                    end
     431                    if numel(DimCell)==1 && isequal(Field.ListVarName{ivar},DimCell{1})%detect dimension variables
     432                        testsum(ivar)=1;
     433                    end
    426434                end
    427435            end
    428             DataOut.ListVarName=[{'Time'} DataOut.ListVarName];
    429         end
    430        
    431         % add data to the current field
    432         if checkhisto
    433             for ivar=1:length(Field.ListVarName)
    434                 VarName=Field.ListVarName{ivar};
    435                 if isfield(Data{1},VarName)
    436                     MaxValue=max(DataOut.(VarName));% current max of histogram absissa
    437                     MinValue=min(DataOut.(VarName));% current min of histogram absissa
    438                     MaxIndex=round(MaxValue/VarMesh);
    439                     MinIndex=round(MinValue/VarMesh);
    440                     MaxIndex_new=round(max(Field.(VarName)/VarMesh));% max of the current field
    441                     MinIndex_new=round(min(Field.(VarName)/VarMesh));
    442                     if MaxIndex_new>MaxIndex% the variable max for the current field exceeds the previous one
    443                         DataOut.(VarName)=[DataOut.(VarName) VarMesh*(MaxIndex+1:MaxIndex_new)];% append the new variable values
    444                         DataOut.([VarName 'Histo'])=[DataOut.([VarName 'Histo']) zeros(nbfield,MaxIndex_new-MaxIndex)]; % append the new histo values
    445                     end
    446                     if MinIndex_new <= MinIndex-1
    447                         DataOut.(VarName)=[VarMesh*(MinIndex_new:MinIndex-1) DataOut.(VarName)];% insert the new variable values
    448                         DataOut.([VarName 'Histo'])=[zeros(nbfield,MinIndex-MinIndex_new) DataOut.([VarName 'Histo'])];% insert the new histo values
    449                         ind_start=1;
    450                     else
    451                         ind_start=MinIndex_new-MinIndex+1;
    452                     end
    453                     DataOut.([VarName 'Histo'])(index,ind_start:ind_start+MaxIndex_new-MinIndex_new)=...
    454                         DataOut.([VarName 'Histo'])(index,ind_start:ind_start+MaxIndex_new-MinIndex_new)+Field.([VarName 'Histo']);
     436            for ivar=1:nbvar
     437                if testsum(ivar)==2
     438                    VarName=Field.ListVarName{ivar};
     439                    siz=size(Field.(VarName));
     440                    DataOut.(VarName)=zeros([nbfield siz]);
    455441                end
    456442            end
     443        end
     444        DataOut.ListVarName=[{'Time'} DataOut.ListVarName];
     445    end
     446    % end initialisation for index==1
     447   
     448    % append data from the current field
     449
     450    if checkhisto % case of histogram (projection mode=inside or outside)
     451        for ivar=1:length(Field.ListVarName)
     452            VarName=Field.ListVarName{ivar};
     453            if isfield(Data{1},VarName)
     454                MaxValue=max(DataOut.(VarName));% current max of histogram absissa
     455                MinValue=min(DataOut.(VarName));% current min of histogram absissa
     456                MaxIndex=round(MaxValue/VarMesh);
     457                MinIndex=round(MinValue/VarMesh);
     458                MaxIndex_new=round(max(Field.(VarName)/VarMesh));% max of the current field
     459                MinIndex_new=round(min(Field.(VarName)/VarMesh));
     460                if MaxIndex_new>MaxIndex% the variable max for the current field exceeds the previous one
     461                    DataOut.(VarName)=[DataOut.(VarName) VarMesh*(MaxIndex+1:MaxIndex_new)];% append the new variable values
     462                    DataOut.([VarName 'Histo'])=[DataOut.([VarName 'Histo']) zeros(nbfield,MaxIndex_new-MaxIndex)]; % append the new histo values
     463                end
     464                if MinIndex_new <= MinIndex-1
     465                    DataOut.(VarName)=[VarMesh*(MinIndex_new:MinIndex-1) DataOut.(VarName)];% insert the new variable values
     466                    DataOut.([VarName 'Histo'])=[zeros(nbfield,MinIndex-MinIndex_new) DataOut.([VarName 'Histo'])];% insert the new histo values
     467                    ind_start=1;
     468                else
     469                    ind_start=MinIndex_new-MinIndex+1;
     470                end
     471                DataOut.([VarName 'Histo'])(index,ind_start:ind_start+MaxIndex_new-MinIndex_new)=...
     472                    DataOut.([VarName 'Histo'])(index,ind_start:ind_start+MaxIndex_new-MinIndex_new)+Field.([VarName 'Histo']);
     473                VarVal=pdf2stat(Field.(VarName),Field.([VarName 'Histo']));% max of the current field
     474                for istat=1:numel(VarVal)
     475                    DataOut.([VarName StatName{istat}])(index)=VarVal(istat);
     476                end
     477            end
     478        end
     479    else % not histogram
     480        for ivar=1:length(Field.ListVarName)
     481            VarName=Field.ListVarName{ivar};
     482            VarVal=Field.(VarName);
     483            if testsum(ivar)==2% test for recorded variable
     484                if isempty(errormsg)
     485                    VarVal=shiftdim(VarVal,-1); %shift dimension
     486                    DataOut.(VarName)(index,:,:)=VarVal;%concanete the current field to the time series
     487                end
     488            elseif testsum(ivar)==1% variable representing fixed coordinates
     489                VarInit=DataOut.(VarName);
     490                if isempty(errormsg) && ~isequal(VarVal,VarInit)
     491                    disp_uvmat('ERROR',['time series requires constant coordinates ' VarName ': use projection mode interp'],checkrun)
     492                    return
     493                end
     494            end
     495        end
     496    end
     497   
     498    % record the time:
     499    if isempty(time)% time not set by xml filer(s)
     500        if isfield(Data{1},'Time')
     501            DataOut.Time(index,1)=Field.Time;
    457502        else
    458             for ivar=1:length(Field.ListVarName)
    459                 VarName=Field.ListVarName{ivar};
    460                 VarVal=Field.(VarName);
    461                 if testsum(ivar)==2% test for recorded variable
    462                     if isempty(errormsg)
    463                         VarVal=shiftdim(VarVal,-1); %shift dimension
    464                         DataOut.(VarName)(index,:,:)=VarVal;%concanete the current field to the time series
    465                     end
    466                 elseif testsum(ivar)==1% variable representing fixed coordinates
    467                     VarInit=DataOut.(VarName);
    468                     if isempty(errormsg) && ~isequal(VarVal,VarInit)
    469                         disp_uvmat('ERROR',['time series requires constant coordinates ' VarName ': use projection mode interp'],checkrun)
    470                         return
    471                     end
    472                 end
    473             end
    474         end
    475        
    476         % record the time:
    477         if isempty(time)% time not set by xml filer(s)
    478             if isfield(Data{1},'Time')
    479                 DataOut.Time(index,1)=Field.Time;
    480             else
    481                 DataOut.Time(index,1)=index;%default
    482             end
    483         else % time from ImaDoc prevails  TODO: correct
    484             DataOut.Time(index,1)=time(index);%
    485         end
    486        
    487         % record the number of missing input fields
    488         if ~isempty(errormsg)
    489             nbmissing=nbmissing+1;
    490             display(['index=' num2str(index) ':' errormsg])
    491         end
    492     end
    493 end
     503            DataOut.Time(index,1)=index;%default
     504        end
     505    else % time from ImaDoc prevails  TODO: correct
     506        DataOut.Time(index,1)=time(index);%
     507    end
     508   
     509    % record the number of missing input fields
     510    if ~isempty(errormsg)
     511        nbmissing=nbmissing+1;
     512        display(['index=' num2str(index) ':' errormsg])
     513    end
     514end
     515
    494516%%%%%%% END OF LOOP WITHIN A SLICE
    495517
     
    535557end
    536558
    537 %case of histograms
    538 if checkhisto
    539     for ivar=1:numel(Field.ListVarName)
    540         VarName=Field.ListVarName{ivar};
    541         if isfield(Data{1},VarName)
    542             DataOut.ListVarName=[DataOut.ListVarName {[VarName 'Histo']}];
    543             DataOut.VarDimName=[DataOut.VarDimName {{'Time',VarName}}];
    544         end
    545     end
    546 end
     559% %case of histograms
     560% if checkhisto
     561%     for ivar=1:numel(Field.ListVarName)
     562%         VarName=Field.ListVarName{ivar};
     563%         if isfield(Data{1},VarName)
     564%             DataOut.ListVarName=[DataOut.ListVarName {[VarName 'Histo']}];
     565%             DataOut.VarDimName=[DataOut.VarDimName {{'Time',VarName}}];
     566%         end
     567%     end
     568% end
    547569% display nbmissing
    548570if ~isequal(nbmissing,0)
     
    559581end
    560582
    561 %% plot the time series (the last one in case of multislices)
    562 if checkrun
    563     figure
    564     haxes=axes;
    565     plot_field(DataOut,haxes)
    566    
    567     %% display the result file using the GUI get_field
    568     hget_field=findobj(allchild(0),'name','get_field');
    569     if ~isempty(hget_field)
    570         delete(hget_field)
    571     end
    572     get_field(OutputFile,DataOut)
    573 end
    574 
    575 
     583%% plot the time series for  (the last one in case of multislices)
     584if checkrun %&& isfield(Param,'ProjObject') && strcmp(Param.ProjObject.Type,'points')
     585    %% open the result file with uvmat (in RUN mode)
     586    uvmat(OutputFile)% open the last result file with uvmat
     587end
     588%     figure
     589%     haxes=axes;
     590%     plot_field(DataOut,haxes)
     591%     
     592%     %% display the result file using the GUI get_field
     593%     hget_field=findobj(allchild(0),'name','get_field');
     594%     if ~isempty(hget_field)
     595%         delete(hget_field)
     596%     end
     597%     get_field(OutputFile,DataOut)
     598% end
     599
     600
Note: See TracChangeset for help on using the changeset viewer.