Changeset 1131 for trunk


Ignore:
Timestamp:
Apr 11, 2024, 7:30:52 PM (8 months ago)
Author:
sommeria
Message:

updated to read pivdata from fluidimage

Location:
trunk/src
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/fileparts_uvmat.m

    r1127 r1131  
    9494            delim2=r.delim2;
    9595            num3=r.num3;
     96            ipair=1;
     97            jpair=1;
    9698            switch delim1
    9799                case '_'
     
    101103                            i1=str2double(num3);
    102104                            i2=str2double(num2);
     105                            if i2<=i1
     106                                ipair=0;
     107                            end
    103108                    end
    104109                case '-'
    105110                    j1=str2double(num2);
    106111                    j2=str2double(num1);
     112                    if j2<=j1
     113                                jpair=0;
     114                    end
    107115                    switch delim2
    108116                        case '_'
     
    110118                    end
    111119            end
     120            if ipair && jpair
    112121            NomType=[get_type(num3) delim2 get_type(num2) delim1 get_type(num1)];
    113122            RootFile=regexprep(FileName,[num3 delim2 num2 delim1 num1 '$'],'');
     123            else
     124                if ipair==0
     125                    i1=i2;
     126                    i2=[];
     127                    NomType=[get_type(num2) delim1 get_type(num1)];
     128                    RootFile=regexprep(FileName,[ num2 delim1 num1 '$'],'');
     129                elseif jpair==0
     130                    j1=j2;
     131                    j2=[];
     132                    NomType=get_type(num1);
     133                    RootFile=regexprep(FileName,[num1 '$'],'');
     134                end
     135            end
     136       
    114137        else
     138            ipair=1;
    115139            switch delim1
    116140                case '_'
     
    120144                    i1=str2double(num2);
    121145                    i2=str2double(num1);
    122             end
     146                    if i2<=i1
     147                        ipair=0;
     148                    end
     149            end
     150            if ipair
    123151            NomType=[get_type(num2) delim1 get_type(num1)];
    124152            RootFile=regexprep(FileName,[num2 delim1 num1 '$'],'');
     153            else
     154                i1=i2;
     155                i2=[];% no pair detected if i2<=i1
     156                NomType= get_type(num1);
     157            RootFile=regexprep(FileName,[ num1 '$'],'');
     158            end
    125159        end
    126160        NomType=regexprep(NomType,'-1','-2'); %set 1-2 instead of 1-1
     
    152186
    153187%% suppress '_' at the end of RootFile, put it on NomType
    154 % if strcmp(RootFile(end),'_')
    155 %     RootFile(end)=[];
    156 detect=regexp(RootFile,'_$'); %detect '_' at the end of RootFILE
     188detect=regexp(RootFile,'_$','once'); %detect '_' at the end of RootFILE
    157189if ~isempty(detect)
    158190    RootFile=regexprep(RootFile,'_$','');
     
    161193
    162194%% extract subdirectory for pairs i1-i2 or j1-j2 (or ab, AB)
    163 % if ~isempty(i2) || ~isempty(j2)
    164195    r=regexp(RootPath,'\<(?<newrootpath>.+)(\\|/)(?<subdir>[^\\^/]+)(\\|/)*\>','names');
    165196    if ~isempty(r)
     
    167198        RootPath=r.newrootpath;
    168199    end
    169 % end
    170200
    171201
  • trunk/src/plot_field.m

    r1127 r1131  
    981981    else %usual images (no contour)
    982982        % set  colormap for  image display
    983         if strcmp(ColorMap,'grayscale')
    984             vec=linspace(0,1,255);%define a linear greyscale colormap
    985             map=[vec' vec' vec'];
    986             colormap(map);  %grey scale color map
    987             if siz==3% true color images visualized in BW
    988                 A=uint16(sum(A,3));%sum the three color components for color images displayed with BW option
    989             end
    990         elseif strcmp(ColorMap,'BuYlRd')
    991             hh=load('BuYlRd.mat');
    992             colormap(hh.BuYlRd);
    993         else
    994             if siz==3 && CheckFixScalar % true color images rescaled by MaxA
    995                   A=uint8(255*double(A)/double(MaxA));
    996             end
    997             colormap(ColorMap); % standard false colors for div, vort , scalar fields
    998         end
     983        switch ColorMap
     984            case 'grayscale'
     985                vec=linspace(0,1,255);%define a linear greyscale colormap
     986                map=[vec' vec' vec'];
     987                colormap(map);  %grey scale color map
     988                if siz==3% true color images visualized in BW
     989                    A=uint16(sum(A,3));%sum the three color components for color images displayed with BW option
     990                end
     991            case 'BuYlRd'
     992                hh=load('BuYlRd.mat');
     993                colormap(hh.BuYlRd);
     994            case 'truecolor'
     995                if siz==3 && CheckFixScalar % true color images rescaled by MaxA
     996                    A=uint8(255*double(A)/double(MaxA));
     997                end
     998        end
     999 
    9991000       
    10001001        % interpolate field to increase resolution of image display
  • trunk/src/read_pivdata_fluidimage.m

    r1127 r1131  
    9999end
    100100
    101 %% reading data
    102 
    103 
    104 % Data=nc2struct(FileName,'ListGlobalAttribute','CivStage');
    105 % if isfield(Data,'Txt')
    106 %      erromsg=['error in read_civdata: ' Data.Txt];
    107 %     return
    108 % end
    109 % % set the list of variables to read and their role
    110 %[varlist,role,VelTypeOut]=varcivx_generator(ProjModeRequest,VelType,Datasets.CivStage);
    111 % if isempty(varlist)
    112 %     erromsg=['error in read_civdata: unknow velocity type ' VelType];
    113 %     return
    114 % else
    115 Datasets=hdf5load(FileName);%read the variables in the netcdf file
    116 %[varlist,role,VelTypeOut]=varcivx_generator(ProjModeRequest,VelType,Datasets.CivStage);
    117 role={'coord_x','coord_y','vector_x','vector_y','ancillary','warnflag','errorflag'};
    118 vardetect=ones(size(role));
    119 Field.ListGlobalAttribute= {'Dt','Time','CivStage'};
     101%% read the hdf file
     102Datasets=hdf5load(FileName);%read the hdf file
     103
     104%% Global attributes
     105Field.ListGlobalAttribute={'Dt','Time','CivStage','NbCoord','NbDim','TimeUnit','CoordUnit'};
    120106Field.Dt=1;
    121107Field.Time=0;
     
    125111    Field.CivStage=3;
    126112end
     113Field.NbCoord=2;
     114Field.NbDim=2;
     115Field.TimeUnit='s';
     116Field.CoordUnit='pixel';
     117
     118%% reading data
     119Field.ListVarName={'X'  'Y'  'U'  'V'  'C'  'F'  'FF'};
     120Field.VarDimName={'nb_vec' 'nb_vec' 'nb_vec' 'nb_vec' 'nb_vec' 'nb_vec' 'nb_vec'};
    127121VelTypeOut=VelType;
    128122switch VelType
    129123    case {'civ1','filter1'}
    130124        Data=Datasets.piv0;
    131         VelTypeOut='filter1';
     125        %VelTypeOut='filter1';
    132126    case {'civ2','filter2'}
    133127        Data=Datasets.piv1;
    134         VelTypeOut='filter2';
    135     case ''
     128        %VelTypeOut='filter2';
     129    case ''% no field specified as input, choose the most appropriate
    136130        if isfield(Datasets, 'piv1')
    137131            Data=Datasets.piv1;
     
    142136        end
    143137end
     138npy=double(Datasets.couple.shape_images(1)); %number of pixels along y for the image sources
    144139switch VelType
    145140    case {'civ1','civ2'}
    146141        Field.X= double(Data.xs);
    147         Field.Y= double(Data.ys);
     142        Field.Y= npy-double(Data.ys);
    148143        Field.U= double(Data.deltaxs);
    149         Field.U(isnan(Field.U)) = 0;
    150         Field.V= double(Data.deltays);
    151         Field.V(isnan(Field.V)) = 0;
    152     case {'filter1','filter2',''}
     144        Field.V= -double(Data.deltays);
     145        checkcolor=1;%color representation of the correlation and errors
     146    case 'filter1'
     147        Field.X= double(Data.ixvecs_approx);
     148        Field.Y= npy-double(Data.iyvecs_approx);
     149        Field.U= double(Data.deltaxs_approx);
     150        Field.V= -double(Data.deltays_approx);
     151        checkcolor=0;%no color representation of the correlation and errors
     152    case {'filter2',''}
    153153        Field.X= double(Data.ixvecs_final);
    154         Field.Y= double(Data.iyvecs_final);
     154        Field.Y= npy-double(Data.iyvecs_final);
    155155        Field.U= double(Data.deltaxs_final);
    156         Field.U(isnan(Field.U)) = 0;
    157         Field.V= double(Data.deltays_final);
    158         Field.V(isnan(Field.V)) = 0;
    159 end
    160 Field.ListVarName={'X'  'Y'  'U'  'V'  'C'  'F'  'FF'};
    161 Field.VarDimName={'nb_vec' 'nb_vec' 'nb_vec' 'nb_vec' 'nb_vec' 'nb_vec' 'nb_vec'};
    162 Field.C = double(Data.correls_max);
    163 Field.F=zeros(size(Field.U)); Field.F(Data.errors.keys + 1)=1; % !!! convention matlab vs python
    164 Field.FF=zeros(size(Field.U)); Field.FF(Data.errors.keys + 1)=1;
    165 % if vardetect(1)==0
    166 %      errormsg=[ 'requested field not available in ' FileName '/' VelType ': need to run patch'];
    167 %      return
    168 % end
    169 switch VelTypeOut
    170     case {'civ1','filter1'}
    171         if isfield(Field,'Patch1_SubDomain')
    172             Field.SubDomain=Field.Patch1_SubDomain;
    173             Field.ListGlobalAttribute=[Field.ListGlobalAttribute {'SubDomain'}];
    174         end
    175         if isfield(Field,'Civ1_Dt')
    176             Field.Dt=Field.Civ1_Dt;
    177         end
    178         if isfield(Field,'Civ1_Time')
    179             Field.Time=Field.Civ1_Time;
    180         end
    181     case {'civ2','filter2'}
    182         if isfield(Field,'Patch2_SubDomain')
    183             Field.SubDomain=Field.Patch2_SubDomain;
    184             Field.ListGlobalAttribute=[Field.ListGlobalAttribute {'SubDomain'}];
    185         end
    186         if isfield(Field,'Civ2_Dt')
    187         Field.Dt=Field.Civ2_Dt;
    188         end
    189         if isfield(Field,'Civ2_Time')
    190         Field.Time=Field.Civ2_Time;
    191         end
    192 end
    193 %Field.ListGlobalAttribute=[Field.ListGlobalAttribute {'Dt','Time'}];
     156        Field.V= -double(Data.deltays_final);
     157        checkcolor=1;%color representation of the correlation and errors
     158end
     159Field.U(isnan(Field.U)) = 0;
     160Field.V(isnan(Field.V)) = 0;
     161Field.C=ones(size(Field.U));%default
     162Field.F=zeros(size(Field.U));
     163Field.FF=zeros(size(Field.U));
     164if checkcolor
     165    Field.C = double(Data.correls_max);
     166    Field.F(Data.errors.keys + 1)=1; % !!! convention matlab vs python
     167    Field.FF(Data.errors.keys + 1)=1;
     168end
     169
     170%% set variable attributes
    194171ivar_U_tps=[];
    195172ivar_V_tps=[];
     173role={'coord_x','coord_y','vector_x','vector_y','ancillary','warnflag','errorflag'};
     174vardetect=ones(size(role));
    196175var_ind=find(vardetect);
    197176for ivar=1:numel(var_ind)
    198177    Field.VarAttribute{ivar}.Role=role{var_ind(ivar)};
    199     %Field.VarAttribute{ivar}.Mesh=0.025;%typical mesh for histograms O.025 pixel (used in series)
     178    Field.VarAttribute{ivar}.Mesh=0.025;%typical mesh for histograms O.025 pixel (used in series)
    200179    Field.VarAttribute{ivar}.ProjModeRequest=ProjModeRequest;
    201180    if strcmp(role{var_ind(ivar)},'vector_x')
     
    223202end
    224203
    225 %% update list of global attributes
    226 Field.ListGlobalAttribute=[Field.ListGlobalAttribute {'NbCoord','NbDim','TimeUnit','CoordUnit'}];
    227 Field.NbCoord=2;
    228 Field.NbDim=2;
    229 Field.TimeUnit='s';
    230 Field.CoordUnit='pixel';
     204
  • trunk/src/script_readnc.m

    r1119 r1131  
    1 DataFolder='.fsnet/project/meige/2019/19TORE';
     1DataFolder='.fsnet/project/coriolis/2024/24PLUME/VectorsEyemotion';
    22fileinput={'Tore_a8b40_T22s_N047_A1cm.png.civ.mproj.tfilter_1'}
    33fileinput=uigetfile_uvmat('pick an input file',DataFolder);
  • trunk/src/series/civ_series.m

    r1127 r1131  
    2121%           if absent no file is produced, result in the output structure Data (test mode)
    2222%     Param.ActionInput: substructure with the parameters provided by the GUI civ_input
    23 %                      .Civ1: parameters for civ1
     23%                      .Civ1: parameters for civ1cc
    2424%                      .Fix1: parameters for fix1
    2525%                      .Patch1:
     
    286286end
    287287for ifield=1:NbField
    288     tic
     288    tstart=tic;
     289    time_civ1=0;
     290  time_patch1=0;
     291  time_civ2=0;
     292  time_patch2=0;
    289293    if ~isempty(RUNHandle)% update the waitbar in interactive mode with GUI series  (checkrun=1)
    290294        update_waitbar(WaitbarHandle,ifield/NbField)
     
    382386                        continue
    383387                    end
     388                    tsart_input=tic;
    384389                    [par_civ1.ImageA,VideoObject_A] = read_image(ImageName_A,FileType_A,VideoObject_A,FrameIndex_A_Civ1(ifield));
     390                    time_input=toc(tsart_input);
    385391                end
    386392                ImageName_B=fullfile_uvmat(RootPath_B,SubDir_B,RootFile_B,FileExt_B,NomType_B,i2_series_Civ1(ifield),[],j2_series_Civ1(ifield));
     
    501507        else %usual PIV
    502508            % caluclate velocity data (y and v in indices, reverse to y component)
     509            tstart_civ1=tic;
    503510            [xtable, ytable, utable, vtable, ctable, F, result_conv, errormsg] = civ (par_civ1);
    504511            if ~isempty(errormsg)
     
    512519            Data.Civ1_C=reshape(ctable,[],1);
    513520            Data.Civ1_F=reshape(F,[],1);
     521            time_civ1=toc(tstart_civ1);
    514522        end
    515523    else% we use existing Civ1 data
     
    561569    if isfield (Param.ActionInput,'Patch1')
    562570        disp('patch1 started')
     571         tstart_patch1=tic;
    563572        if check_civx
    564573            errormsg='Civ Matlab input needed for patch';
     
    600609
    601610        % perform Patch calculation using the UVMAT fct 'filter_tps'
     611       
    602612        [Data.Civ1_SubRange,Data.Civ1_NbCentres,Data.Civ1_Coord_tps,Data.Civ1_U_tps,Data.Civ1_V_tps,tild,Ures, Vres,tild,FFres]=...
    603613            filter_tps([Data.Civ1_X(ind_good) Data.Civ1_Y(ind_good)],Data.Civ1_U(ind_good),Data.Civ1_V(ind_good),[],Data.Patch1_SubDomainSize,Data.Patch1_FieldSmooth,Data.Patch1_MaxDiff);
     
    605615        Data.Civ1_V_smooth(ind_good)=Vres;
    606616        Data.Civ1_FF(ind_good)=int8(FFres);
     617        time_patch1=toc(tstart_patch1);
    607618        disp('patch1 performed')
    608619    end
     
    611622    if isfield (Param.ActionInput,'Civ2')
    612623        disp('civ2 started')
     624        tstart_civ2=tic;
    613625        par_civ2=Param.ActionInput.Civ2;
    614626        if CheckInputFile % read input images (except in mode Test where it is introduced directly in Param.ActionInput.Civ1.ImageNameA and B)
     
    765777       
    766778        % calculate velocity data (y and v in image indices, reverse to y component)
     779       
    767780        [xtable, ytable, utable, vtable, ctable, F,result_conv,errormsg] = civ (par_civ2);
    768781       
     
    806819        Data.Civ2_F=reshape(F,[],1);
    807820        disp('civ2 performed')
     821        time_civ2=toc(tstart_civ2);
    808822    elseif ~isfield(Data,'ListVarName') % we start there, using existing Civ2 data
    809823        if exist('ncfile','var')
     
    836850        end
    837851        Data.ListGlobalAttribute=[Data.ListGlobalAttribute Fix2_param];
    838         %
    839         %         ListFixParam=fieldnames(Param.ActionInput.Fix2);
    840         %         for ilist=1:length(ListFixParam)
    841         %             ParamName=ListFixParam{ilist};
    842         %             ListName=['Fix2_' ParamName];
    843         %             eval(['Data.ListGlobalAttribute=[Data.ListGlobalAttribute ''' ParamName '''];'])
    844         %             eval(['Data.' ListName '=Param.ActionInput.Fix2.' ParamName ';'])
    845         %         end
    846852        if check_civx
    847853            if ~isfield(Data,'fix2')
     
    865871    if isfield (Param.ActionInput,'Patch2')
    866872        disp('patch2 started')
     873        tstart_patch2=tic;
    867874        list_param=fieldnames(Param.ActionInput.Patch2)';
    868875        list_param(strcmp('TestPatch2',list_param))=[];% remove the parameter TestCiv1 from the list
     
    894901                        disp_uvmat('ERROR','all vectors of civ2 are bad, check input parameters' ,checkrun)
    895902                        return
    896         end
     903                end
     904             
    897905        [Data.Civ2_SubRange,Data.Civ2_NbCentres,Data.Civ2_Coord_tps,Data.Civ2_U_tps,Data.Civ2_V_tps,tild,Ures,Vres,tild,FFres]=...
    898906            filter_tps([Data.Civ2_X(ind_good) Data.Civ2_Y(ind_good)],Data.Civ2_U(ind_good),Data.Civ2_V(ind_good),[],Data.Patch2_SubDomainSize,Data.Patch2_FieldSmooth,Data.Patch2_MaxDiff);
     
    901909        Data.Civ2_FF(ind_good)=FFres;
    902910        Data.CivStage=Data.CivStage+1;
     911        time_patch2=toc(tstart_patch2);
    903912        disp('patch2 performed')
    904913    end
     
    913922            disp(errormsg)
    914923        end
    915         disp(['ellapsed time ' num2str(toc/60,2) ' minutes'])
     924        time_total=toc(tstart);
     925        disp(['ellapsed time ' num2str(time_total/60,2) ' minutes'])
     926        disp(['time image reading ' num2str(time_input/60,2) ' minutes'])
     927        disp(['time civ1 ' num2str(time_civ1/60,2) ' minutes'])
     928        disp(['time patch1 ' num2str(time_patch1/60,2) ' minutes'])
     929        disp(['time civ2 ' num2str(time_civ2/60,2) ' minutes'])
     930        disp(['time patch2 ' num2str(time_patch2/60,2) ' minutes'])
     931        disp(['time other ' num2str((time_total-time_input-time_civ1-time_patch1-time_civ2-time_patch2)/60,2) ' minutes'])
    916932    end
    917933end
Note: See TracChangeset for help on using the changeset viewer.