Changeset 876 for trunk


Ignore:
Timestamp:
Feb 22, 2015, 11:07:07 PM (10 years ago)
Author:
sommeria
Message:

deformation modified in civ_series and bug corrections

Location:
trunk/src
Files:
10 edited

Legend:

Unmodified
Added
Removed

  • trunk/src/get_field.m

    r874 r876  
    268268drawnow
    269269uiwait(handles.get_field);
    270 
    271 %------------------------------------------------------------------------
    272 % --- Executes when user attempts to close get_field.
    273 %------------------------------------------------------------------------
    274 %------------------------------------------------------------------------
    275 
    276 
    277 % Use UIRESUME instead of delete because the OutputFcn needs
    278 % to get the updated handles structure.
    279 
    280 % function get_field_CloseRequestFcn(hObject, eventdata)
    281 % handles.output=[];
    282 % guidata(hObject, handles);% Update handles structure
    283 % %delete(handles.get_field);
    284 % uiresume(handles.get_field);
    285 %drawnow
    286 % if isequal(get(handles.get_field, 'waitstatus'), 'waiting')
    287 %     % The GUI is still in UIWAIT, us UIRESUME
    288 %     uiresume(handles.get_field);
    289 % else
    290 %     % The GUI is no longer waiting, just close it
    291 %     delete(handles.get_field);
    292 % end
    293270
    294271% -----------------------------------------------------------------------
    … …  
    775752    end
    776753    if numel(var_component)<2
    777         if numel(test_coord)<2
     754        if numel(find(test_coord))<2
    778755            ListCoord={''};
    779756        else
    780             if numel(test_coord)>=3
     757            if numel(find(test_coord))>=3
    781758                set(handles.Coord_x,'Value',3)
    782759                set(handles.Coord_y,'Value',2)
    … …  
    10681045set(handles.vector_z,'Visible',status)
    10691046set(handles.W_title,'Visible',status)   
     1047if strcmp(status,'on')
    10701048   Field=get(handles.get_field,'UserData');
    10711049    if Field.MaxDim>=3% for 3D fields, propose to use the third variable as time
    … …  
    10771055        end
    10781056    end
    1079 
     1057end
    10801058%------------------------------------------------------------------------
    10811059% --- Executes on button press in OK.

  • trunk/src/series.m

    r872 r876  
    24192419        GetFieldData=get_field(FirstFileName,ParamIn);
    24202420        FieldList={};
     2421        if isfield(GetFieldData,'FieldOption')% if a field has been selected
    24212422        switch GetFieldData.FieldOption
    24222423            case 'vectors'
    … …  
    24822483        set(handles.Coord_x,'Visible','on')
    24832484        set(handles.Coord_y,'Visible','on')
     2485        end
    24842486    else
    24852487        msgbox_uvmat('ERROR',[FirstFileName ' does not exist'])

  • trunk/src/series/civ_input.m

    r874 r876  
    18601860     end
    18611861     if isfield(Param,'OutputSubDir')
    1862      Param=rmfield(Param,'OutputSubDir'); %remove output file option from civ_series
     1862         Param=rmfield(Param,'OutputSubDir'); %remove output file option from civ_series
    18631863     end
    18641864     ParamPatch1=Param.ActionInput.Patch1; %store the patch1 parameters
    18651865     Param.ActionInput=rmfield(Param.ActionInput,'Patch1');% does not execute Patch
     1866     Param.IndexRange.first_i=str2num(get(handles.ref_i,'String'));
     1867     Param.IndexRange.last_i=Param.IndexRange.first_i;
     1868     if strcmp(get(handles.ref_j,'Visible'),'on')
     1869         Param.IndexRange.first_j=str2num(get(handles.ref_j,'String'));
     1870         Param.IndexRange.last_j=Param.IndexRange.first_j;
     1871     else
     1872         Param.IndexRange.first_j=1;
     1873         Param.IndexRange.last_j=1;
     1874     end
    18661875     [Data,errormsg]=civ_series(Param);% get the civ1+fix1 results
    18671876     bckcolor=get(handles.civ_input,'Color');
    … …  
    19381947     set(handles.Fix1,'BackgroundColor',[1 1 0])
    19391948     set(handles.Patch1,'BackgroundColor',[1 1 0])
     1949     Param.IndexRange.first_i=str2num(get(handles.ref_i,'String'));
     1950     Param.IndexRange.last_i=Param.IndexRange.first_i;
     1951     if strcmp(get(handles.ref_j,'Visible'),'on')
     1952         Param.IndexRange.first_j=str2num(get(handles.ref_j,'String'));
     1953         Param.IndexRange.last_j=Param.IndexRange.first_j;
     1954     else
     1955         Param.IndexRange.first_j=1;
     1956         Param.IndexRange.last_j=1;
     1957     end
    19401958     [Data,errormsg]=civ_series(Param);% get the civ1+fix1 results
    19411959     
    … …  
    20312049     ParamPatch2=Param.ActionInput.Patch2; %store the patch1 parameters
    20322050     Param.ActionInput=rmfield(Param.ActionInput,'Patch2');% does not execute Patch
     2051     Param.IndexRange.first_i=str2num(get(handles.ref_i,'String'));
     2052     Param.IndexRange.last_i=Param.IndexRange.first_i;
     2053     if strcmp(get(handles.ref_j,'Visible'),'on')
     2054         Param.IndexRange.first_j=str2num(get(handles.ref_j,'String'));
     2055         Param.IndexRange.last_j=Param.IndexRange.first_j;
     2056     else
     2057         Param.IndexRange.first_j=1;
     2058         Param.IndexRange.last_j=1;
     2059     end
    20332060     [Data,errormsg]=civ_series(Param);% get the civ1+fix1 results
    20342061     bckcolor=get(handles.civ_input,'Color');

  • trunk/src/series/civ_series.m

    r873 r876  
    691691            end
    692692        end
    693 %         ibx2=ceil(par_civ2.CorrBoxSize(1)/2);
    694 %         iby2=ceil(par_civ2.CorrBoxSize(2)/2);
    695        % par_civ2.SearchBoxSize(1)=2*ibx2+9;% search ara +-4 pixels around the guess
    696         %par_civ2.SearchBoxSize(2)=2*iby2+9;
    697693        if CheckInputFile % else Dt given by par_civ2
    698694            if strcmp(Param.ActionInput.ListCompareMode,'displacement')
    … …  
    707703        par_civ2.Grid=[par_civ2.Grid(nbval>=1,1)-par_civ2.SearchBoxShift(:,1)/2 par_civ2.Grid(nbval>=1,2)-par_civ2.SearchBoxShift(:,2)/2];
    708704        if par_civ2.CheckDeformation
    709             par_civ2.DUDX=DUDX./nbval;
    710             par_civ2.DUDY=DUDY./nbval;
    711             par_civ2.DVDX=DVDX./nbval;
    712             par_civ2.DVDY=DVDY./nbval;
     705            par_civ2.DUDX=DUDX(nbval>=1)./nbval(nbval>=1);
     706            par_civ2.DUDY=DUDY(nbval>=1)./nbval(nbval>=1);
     707            par_civ2.DVDX=DVDX(nbval>=1)./nbval(nbval>=1);
     708            par_civ2.DVDY=DVDY(nbval>=1)./nbval(nbval>=1);
    713709        end
    714710       
    … …  
    910906
    911907%% prepare correlation and search boxes
    912 ibx2=ceil(par_civ.CorrBoxSize(1)/2);
    913 iby2=ceil(par_civ.CorrBoxSize(2)/2);
    914 isx2=ceil(par_civ.SearchBoxSize(1)/2);
    915 isy2=ceil(par_civ.SearchBoxSize(2)/2);
     908ibx2=floor(par_civ.CorrBoxSize(1)/2);
     909iby2=floor(par_civ.CorrBoxSize(2)/2);
     910isx2=floor(par_civ.SearchBoxSize(1)/2);
     911isy2=floor(par_civ.SearchBoxSize(2)/2);
    916912shiftx=round(par_civ.SearchBoxShift(:,1));%use the input shift estimate, rounded to the next integer value
    917913shifty=-round(par_civ.SearchBoxShift(:,2));% sign minus because image j index increases when y decreases
    … …  
    968964    end
    969965    check_undefined=(par_civ.Mask<200 & par_civ.Mask>=20 );
    970     par_civ.ImageA(check_undefined)=0;% put image A to zero (i.e. the min image value) in the undefined  area
    971     par_civ.ImageB(check_undefined)=0;% put image B to zero (i.e. the min image value) in the undefined  area
     966%     par_civ.ImageA(check_undefined)=0;% put image A to zero (i.e. the min image value) in the undefined  area
     967%     par_civ.ImageB(check_undefined)=0;% put image B to zero (i.e. the min image value) in the undefined  area
    972968end
    973969
    … …  
    993989        image2_crop=MinA*ones(numel(subrange2_y),numel(subrange2_x));% default value=min of image A
    994990        mask1_crop=ones(numel(subrange1_y),numel(subrange1_x));% default value=1 for mask
    995         mask2_crop=ones(numel(subrange2_y),numel(subrange2_x));% default value=min for mask
     991        mask2_crop=ones(numel(subrange2_y),numel(subrange2_x));% default value=1 for mask
    996992        check1_x=subrange1_x>=1 & subrange1_x<=par_civ.ImageWidth;% check which points in the subimage 1 are contained in the initial image 1
    997993        check1_y=subrange1_y>=1 & subrange1_y<=par_civ.ImageHeight;
    … …  
    10061002            F(ivec)=3; %
    10071003        else
     1004            image1_crop=image1_crop.*~mask1_crop;% put to zero the masked pixels (mask1_crop='true'=1)
     1005            image2_crop=image2_crop.*~mask2_crop;
    10081006            image1_mean=mean(mean(image1_crop))/(1-sizemask);
    10091007            image2_mean=mean(mean(image2_crop))/(1-sizemask);
    … …  
    10261024                XIant=XI-par_civ.DUDX(ivec)*XI-par_civ.DUDY(ivec)*YI+ceil(size(image1_crop,2)/2);
    10271025                YIant=YI-par_civ.DVDX(ivec)*XI-par_civ.DVDY(ivec)*YI+ceil(size(image1_crop,1)/2);
     1026                        image1_crop_1=image1_crop;
    10281027                image1_crop=interp2(image1_crop,XIant,YIant);
     1028                        image1_crop_2=image1_crop;
    10291029                image1_crop(isnan(image1_crop))=0;
     1030                  image1_crop_3=image1_crop;
    10301031                xi=(1:mesh:size(image2_crop,2));
    10311032                yi=(1:mesh:size(image2_crop,1))';
    10321033                image2_crop=interp2(image2_crop,xi,yi);
    10331034                image2_crop(isnan(image2_crop))=0;
     1035                par_civ.CorrSmooth=3;
    10341036                %%
    10351037            end
    … …  
    10411043            %Find the correlation max, at 255
    10421044            [y,x] = find(result_conv==255,1);
     1045            subimage2_crop=image2_crop(y:y+2*iby2/mesh,x:x+2*ibx2/mesh);%subimage of image 2 corresponding to the optimum displacement of first image
     1046            sum_square=sum_square*sum(sum(subimage2_crop.*subimage2_crop));% product of variances of image 1 and 2
     1047            sum_square=sqrt(sum_square);% srt of the variance product to normalise correlation
    10431048            if ~isempty(y) && ~isempty(x)
    10441049                try
    … …  
    10471052                    elseif par_civ.CorrSmooth==2
    10481053                        [vector,F(ivec)] = SUBPIX2DGAUSS (result_conv,x,y);
     1054                    else
     1055                        [vector,F(ivec)] = quadr_fit(result_conv,x,y);
    10491056                    end
    10501057                    utable(ivec)=vector(1)*mesh+shiftx(ivec);
    … …  
    10521059                    xtable(ivec)=iref+utable(ivec)/2-0.5;% convec flow (velocity taken at the point middle from imgae 1 and 2)
    10531060                    ytable(ivec)=jref+vtable(ivec)/2-0.5;% and position of pixel 1=0.5 (convention for image coordinates=0 at the edge)
    1054                     iref=round(xtable(ivec));% image index for the middle of the vector
     1061                    iref=round(xtable(ivec));% nearest image index for the middle of the vector
    10551062                    jref=round(ytable(ivec));
     1063                    % eliminate vectors located in the mask
    10561064                    if checkmask && par_civ.Mask(jref,iref)<200 && par_civ.Mask(jref,iref)>=100
    10571065                        utable(ivec)=0;
    … …  
    11501158end
    11511159vector=[peakx-floor(npx/2)-1 peaky-floor(npy/2)-1];
     1160
     1161%------------------------------------------------------------------------
     1162% --- Find the maximum of the correlation function after quadratic interpolation
     1163function [vector,F] = quadr_fit(result_conv,x,y)
     1164[npy,npx]=size(result_conv);
     1165if x<4 || y<4 || npx-x<4 ||npy-y <4
     1166    F=-2;
     1167    vector=[x y];
     1168else
     1169    F=0;
     1170    x_ind=x-4:x+4;
     1171    y_ind=y-4:y+4;
     1172    x_vec=0.25*(x_ind-x);
     1173    y_vec=0.25*(y_ind-y);
     1174    [X,Y]=meshgrid(x_vec,y_vec);
     1175    coord=[reshape(X,[],1) reshape(Y,[],1)];
     1176    result_conv=reshape(result_conv(y_ind,x_ind),[],1);
     1177   
     1178   
     1179    % n=numel(X);
     1180    % x=[X Y];
     1181    % X=X-0.5;
     1182    % Y=Y+0.5;
     1183    % y = (X.*X+2*Y.*Y+X.*Y+6) + 0.1*rand(n,1);
     1184    p = polyfitn(coord,result_conv,2);
     1185    A(1,1)=2*p.Coefficients(1);
     1186    A(1,2)=p.Coefficients(2);
     1187    A(2,1)=p.Coefficients(2);
     1188    A(2,2)=2*p.Coefficients(4);
     1189    vector=[x y]'-A\[p.Coefficients(3) p.Coefficients(5)]';
     1190    vector=vector'-[floor(npx/2) floor(npy/2)]-1 ;
     1191    % zg = polyvaln(p,coord);
     1192    % figure
     1193    % surf(x_vec,y_vec,reshape(zg,9,9))
     1194    % hold on
     1195    % plot3(X,Y,reshape(result_conv,9,9),'o')
     1196    % hold off
     1197end
    11521198
    11531199%'RUN_FIX': function for fixing velocity fields:

  • trunk/src/series/stereo_civ.m

    r864 r876  
    1 %'civ_series': PIV function activated by the general GUI series
     1%'stereo_civ': determination of topography by image correlation of two stereo views
    22% --- call the sub-functions:
    33%   civ: PIV function itself

  • trunk/src/series/stereo_input.m

    r865 r876  
    1 %'stereo_input': function associated with the GUI 'stereo_input.fig' to set the input parameters for civ_series
     1%'stereo_input': function associated with the GUI 'stereo_input.fig' to set the input parameters for stereo_civ
    22%------------------------------------------------------------------------
    33% function ParamOut = stereo_input(Param)

  • trunk/src/sub_field.m

    r866 r876  
    4343    SubData=Field;
    4444    return
     45end
     46if ~isfield(Field_1,'VarAttribute')
     47    Field_1.VarAttribute={};
    4548end
    4649
    … …  
    107110end
    108111
     112
    109113%% rename the dimensions of the second field if identical to those of the first
    110114for ilist=1:numel(Field_1.VarDimName)
    … …  
    127131ind_remove=false(size(Field_1.ListVarName));
    128132% loop on the variables of the second field Field_1
    129 for ilist=1:numel(Field_1.ListVarName)
     133for ilist=1:numel(Field_1.VarAttribute)
    130134    % case of variable with a single dimension
    131135    if ~isempty(Field_1.VarAttribute{ilist}) && ~isempty(regexp(Field_1.VarAttribute{ilist}.Role,'^coord'))% if variable with Role coord... is found.
    … …  
    173177ListVarNameNew=ListVarNameNew(~check_remove); % %list of renaimed varaibles corresponding to ListVarNameSub
    174178VarDimNameSub=Field_1.VarDimName(~check_remove);
     179if numel(Field_1.VarAttribute)<max(find(~check_remove))
     180    for ilist=numel(Field_1.VarAttribute)+1:max(find(~check_remove))
     181        Field_1.VarAttribute{ilist}={};
     182    end
     183end
    175184VarAttributeSub=Field_1.VarAttribute(~check_remove);
    176185check_rename=check_rename(~check_remove);

  • trunk/src/transform_field/diff_vel.m

    r870 r876  
    1 %'sub_field': combines two input fields
     1%'diff_vel': calculate the difference of two input velocity fields.
    22%
    3 % the two fields are subtstracted when of the same nature (scalar or
    4 % vector), if the coordinates do not coincide, the second field is
    5 % interpolated on the cooridintes of the first one
    6 %
    7 % when scalar and vectors are combined, the fields are just merged in a single matlab structure for common visualisation
     3% the second velocity field is linearly interpolated
     4% (after elimination of the vectors marked with an error flag) to the positions of
     5% the first one before subtraction. The ancilary data of the first field
     6% are preserved while those of the second one are lost.
     7
    88%-----------------------------------------------------------------------
    9 % function SubData=sub_field(Field,XmlData,Field_1)
     9% function SubData=diff_vel(Field,XmlData,Field_1)
    1010%
    1111% OUPUT:
    … …  
    1414% INPUT:
    1515% Field: matlab structure representing the first field
     16% XmlData: not used, needed for consistency with the call of transform fct.
    1617% Field_1:matlab structure representing the second field
    1718

  • trunk/src/transform_field/signal_spectrum.m

    r810 r876  
    1 % 'signal_spectrum': calculate and display spectrum of the current field
     1% 'signal_spectrum': calculate and display power spectrum of the current field
    22%  operate on a 1D signal or the first dimension of a higher dimensional matrix (then average over other dimensions)
    33%  this function aplies the Welch method and call the function of the matlab signal processing toolbox

  • trunk/src/transform_field/sub_field.m

    r810 r876  
    1 %'sub_field': combines two input fields
     1%'sub_field': combines two input fields, taking the difference if of the same nature
    22%
    3 % the two fields are subtstracted when of the same nature (scalar or
    4 % vector), if the coordinates do not coincide, the second field is
    5 % interpolated on the cooridintes of the first one
     3% the two fields are subtracted when of the same nature (scalar or
     4% vector),and defined at the same points
    65%
    76% when scalar and vectors are combined, the fields are just merged in a single matlab structure for common visualisation
Note: See TracChangeset for help on using the changeset viewer.