Changeset 862


Ignore:
Timestamp:
Jan 30, 2015, 8:37:03 PM (10 years ago)
Author:
sommeria
Message:

bugs repaired

Location:
trunk/src
Files:
2 deleted
7 edited

Legend:

Unmodified
Added
Removed

  • trunk/src/mouse_motion.m

    r859 r862  
    232232            if test_piv
    233233               [dd,ind_pt]=min(abs(Field.X-xy(1,1))+abs(Field.Y-xy(1,2)));
    234                if (isfield(hhciv,'TestCiv2')&& get(hhciv.TestCiv2,'Value'))% if TestCiv2 is activated
     234               if isfield(hhciv,'TestCiv2') && strcmp(get(hhciv.Civ2,'Visible'),'on')&& ...
     235                                     strcmp(get(hhciv.TestCiv2,'Visible'),'on')&& get(hhciv.TestCiv2,'Value')% if TestCiv2 is activated
    235236                   CivOption='Civ2';
    236237                   Param.CheckCiv1=0;
    237238                   par_civ=read_GUI(hhciv.Civ2);%read the Civ2 panel in civ_input
    238                    par_civ.Civ1_SubRange=Field.Civ1_SubRange;
    239                    par_civ.Civ1_NbCentres=Field.Civ1_NbCentres;
    240                    par_civ.Civ1_Coord_tps=Field.Civ1_Coord_tps;
    241                    par_civ.Civ1_U_tps=Field.Civ1_U_tps;
    242                    par_civ.Civ1_V_tps=Field.Civ1_V_tps;
     239                   par_civ.Civ1_SubRange=Field.Civ1_SubRange;% get the subranges used for patch
     240                   par_civ.Civ1_NbCentres=Field.Civ1_NbCentres;% get the nbre of civ1 data points in each subrange
     241                   par_civ.Civ1_Coord_tps=Field.Civ1_Coord_tps;% get the positions of the centres (Civ1 data)
     242                   par_civ.Civ1_U_tps=Field.Civ1_U_tps;% get the tps coefficients for spline, U component
     243                   par_civ.Civ1_V_tps=Field.Civ1_V_tps;% get the tps coefficients for spline, V component
    243244                   par_civ.Civ1_Dt=Field.Civ1_Dt;
    244245                   par_civ.Civ2_Dt=Field.Civ2_Dt;
    245                    shiftx=Field.ShiftX(ind_pt);
     246                   shiftx=Field.ShiftX(ind_pt);% get field info stored in the cuirrent figure
    246247                   shifty=Field.ShiftY(ind_pt);
    247248               else
    … …  
    306307                        if ~isempty(corrfig)
    307308                            set(0,'CurrentFigure',corrfig(1))
    308                             AxeData.CurrentCorrImage=imagesc(rangx,-rangy,result_conv,[0 1]);
     309                            AxeData.CurrentCorrImage=imagesc(rangx,-rangy,result_conv,[0 1]);% plot the correlation map, with full colormap from 0 to 1
    309310                            AxeData.CurrentVector=line([0 Data.([CivOption '_U'])],[0 Data.([CivOption '_V'])],'Tag','vector');
    310                             AxeData.TitleHandle=title(num2str(par_civ.Grid));
     311                            AxeData.TitleHandle=title(['[x,y]= ' num2str(par_civ.Grid) ' px']);
    311312                            colorbar
    312313                            set(CurrentAxes,'UserData',AxeData)
    313                             set(get(AxeData.CurrentCorrImage,'parent'),'YDir','normal')
     314                            ParentAxe=get(AxeData.CurrentCorrImage,'parent');
     315                            set(ParentAxe,'YDir','normal')
     316                             set(ParentAxe,'DataAspectRatio',[1 1 1]); %equal axis scale
    314317                        end
    315318                    else
    … …  
    317320                        set(AxeData.CurrentCorrImage,'XData',rangx)
    318321                        set(AxeData.CurrentCorrImage,'YData',-rangy)
     322                        ParentAxe=get(AxeData.CurrentCorrImage,'parent');
     323                        rangx(1)=min(rangx(1),0);
     324                        rangx(2)=max(rangx(2),0);
     325                        rangy(2)=max(rangy(2),0);
     326                        rangy(1)=min(rangy(1),0);
     327                         set(ParentAxe,'XLim',rangx)
     328                         set(ParentAxe,'YLim',[-rangy(2) -rangy(1)])
    319329                        set(AxeData.CurrentVector,'XData',[0 Data.([CivOption '_U'])],'YData',[0 Data.([CivOption '_V'])]);
    320                         set(AxeData.TitleHandle,'String',num2str(par_civ.Grid))
     330                        set(AxeData.TitleHandle,'String',['[x,y]= ' num2str(par_civ.Grid) ' px'])
    321331                    end
    322332                end

  • trunk/src/series/civ2vel_3C.m

    r840 r862  
    1 %'civ2vel_3C': combine velocity fields from two camerasto get three velocity components
     1%'civ2vel_3C': combine velocity fields from two cameras to get three velocity components
    22%------------------------------------------------------------------------
    33% function ParamOut=civ2vel_3C(Param)
    … …  
    116116     tsaiA=XmlData{1}.GeometryCalib;
    117117 else
    118      msgbox_uvmat('ERROR','no geometric calibration available for image A')
     118     disp_uvmat('ERROR','no geometric calibration available for image A',checkrun)
    119119     return
    120120 end
    … …  
    122122     tsaiB=XmlData{2}.GeometryCalib;
    123123 else
    124      msgbox_uvmat('ERROR','no geometric calibration available for image B')
     124     disp_uvmat('ERROR','no geometric calibration available for image B',checkrun)
    125125     return
    126126 end
    127127[filecell,i1_series,i2_series,j1_series,j2_series]=get_file_series(Param);
    128128
    129 %% grid of physical positions
     129%% grid of physical positions (given by projection plane)
     130if ~Param.CheckObject
     131         disp_uvmat('ERROR','a projection plane with interpolation is needed',checkrun)
     132     return
     133end
     134ObjectData=Param.ProjObject;
     135
    130136[x,y]=meshgrid(ObjectData.RangeX(1):ObjectData.DX:ObjectData.RangeX(2),ObjectData.RangeY(1):ObjectData.DY:ObjectData.RangeY(2));
     137z=zeros(size(x));
     138[Field,ParamOut,errormsg] = read_field(FileName,FileType,ParamIn,num)
    131139%camera coordinates:initialisation 2 cameras
    132140NbCamera=2;
    … …  
    178186        %% transform the input field (e.g; phys) if requested (no transform involving two input fields)
    179187            %camera coordinates
    180         xc=R(1)*Xphys+R(2)*Yphys+R(3)*Zphys+Calib.Tx_Ty_Tz(1);
    181         yc=R(4)*Xphys+R(5)*Yphys+R(6)*Zphys+Calib.Tx_Ty_Tz(2);
    182         zc=R(7)*Xphys+R(8)*Yphys+R(9)*Zphys+Calib.Tx_Ty_Tz(3);
    183 
    184         %undistorted image coordinates
    185         Xu=xc./zc;
    186         Yu=yc./zc;
    187188        Data{iview}=phys(Data{iview},XmlData{iview});
    188189       

  • trunk/src/series/civ_input.m

    r861 r862  
    397397%% Civ2 parameters
    398398Param.Civ2.CorrBoxSize=[21 21];
    399 Param.Civ2.SearchBoxSize=[31 31];
     399Param.Civ2.SearchBoxSize=[27 27];
    400400Param.Civ2.CorrSmooth=1;
    401401Param.Civ2.Dx=10;
    … …  
    539539    set(handles.ListPairCiv2,'Visible','off')
    540540end
    541 if max(checkbox(1:3))==0 && get(handles.CheckCiv2,'UserData')==6,% no operation asked before Civ2 and input file ready for civ3
     541hseries=findobj(allchild(0),'Tag','series');% find the parent GUI 'series'
     542hhseries=guidata(hseries); %handles of the elements in 'series'
     543InputTable=get(hhseries.InputTable,'Data');
     544if size(InputTable,1)>=1 && strcmp(InputTable{1,5},'.nc') && max(checkbox(1:3))==0 %&& get(handles.CheckCiv2,'UserData')==6,% no operation asked before Civ2 and input file ready for civ3
    542545    set(handles.CheckCiv3,'Visible','on')
    543546else
    … …  
    573576    checkeven=(mod(ActionInput.Civ2.CorrBoxSize,2)==0);
    574577    ActionInput.Civ2.CorrBoxSize(checkeven)=ActionInput.Civ2.CorrBoxSize(checkeven)+1;% set correlation box sizes to odd values
    575     ActionInput.Civ2.SearchBoxSize=max(ActionInput.Civ2.SearchBoxSize,ActionInput.Civ2.CorrBoxSize+6);
     578    ActionInput.Civ2.SearchBoxSize=max(ActionInput.Civ2.SearchBoxSize,ActionInput.Civ2.CorrBoxSize+4);
    576579    checkeven=(mod(ActionInput.Civ2.SearchBoxSize,2)==0);
    577580    ActionInput.Civ2.SearchBoxSize(checkeven)=ActionInput.Civ2.SearchBoxSize(checkeven)+1;% set search box sizes to odd values
    … …  
    11321135% --- Executes on button press in CheckDeformation.
    11331136function CheckDeformation_Callback(hObject, eventdata, handles)
    1134 set(handles.configSource,'String','NEW')
     1137set(handles.ConfigSource,'String','NEW')
    11351138set(handles.OK,'BackgroundColor',[1 0 1])
    11361139
    … …  
    17511754     Param.Action.RUN=1;
    17521755     Param.ActionInput=read_GUI(handles.civ_input);
     1756     if isfield(Param.ActionInput,'Fix1')
     1757         Param.ActionInput=rmfield(Param.ActionInput,'Fix1');
     1758     end
     1759     if isfield(Param.ActionInput,'Patch1')
     1760         Param.ActionInput=rmfield(Param.ActionInput,'Patch1');
     1761     end
     1762     if isfield(Param.ActionInput,'Civ2')%remove options that may be selected beyond Patch1
     1763         Param.ActionInput=rmfield(Param.ActionInput,'Civ2');
     1764     end
     1765     if isfield(Param.ActionInput,'Fix2')
     1766         Param.ActionInput=rmfield(Param.ActionInput,'Fix2');
     1767     end
     1768     if isfield(Param.ActionInput,'Patch2')
     1769         Param.ActionInput=rmfield(Param.ActionInput,'Patch2');
     1770     end
    17531771     if isfield(Param,'OutputSubDir')
    17541772     Param=rmfield(Param,'OutputSubDir'); %remove output file option from civ_series
    … …  
    18191837     Param.Action.RUN=1;
    18201838     Param.ActionInput=read_GUI(handles.civ_input);
     1839     if isfield(Param.ActionInput,'Civ2')%remove options that may be selected beyond Patch1
     1840         Param.ActionInput=rmfield(Param.ActionInput,'Civ2');
     1841     end
     1842     if isfield(Param.ActionInput,'Fix2')
     1843         Param.ActionInput=rmfield(Param.ActionInput,'Fix2');
     1844     end
     1845     if isfield(Param.ActionInput,'Patch2')
     1846         Param.ActionInput=rmfield(Param.ActionInput,'Patch2');
     1847     end
    18211848     if isfield(Param,'OutputSubDir')
    18221849     Param=rmfield(Param,'OutputSubDir'); %remove output file option from civ_series
    … …  
    19281955    ViewData.PlotAxes.X=Data.Civ2_X';
    19291956    ViewData.PlotAxes.Y=Data.Civ2_Y';
    1930     ViewData.PlotAxes.ShiftX=Data.Civ2_U';% shift at each point (from patch1) estimated by running civ2
     1957    ViewData.PlotAxes.ShiftX=Data.Civ2_U';% shift at each point (from patch1) estimated by the preliminary run of civ2
    19311958    ViewData.PlotAxes.ShiftY=Data.Civ2_V';
    19321959    ViewData.PlotAxes.Civ1_SubRange=Data.Civ1_SubRange;
    … …  
    22792306    drawnow
    22802307end
    2281 
    2282 
    2283 
    2284 
    2285 
    2286 

  • trunk/src/series/civ_series.m

    r861 r862  
    614614            end
    615615        end
    616         if par_civ2.CheckDeformation
    617             DUDX=zeros(size(par_civ2.Grid,1),1);
    618             DUDY=zeros(size(par_civ2.Grid,1),1);
    619             DVDX=zeros(size(par_civ2.Grid,1),1);
    620             DVDY=zeros(size(par_civ2.Grid,1),1);
    621         end
    622616       
    623617        % get the guess from patch1 or patch2 (case 'CheckCiv3')
    … …  
    661655        Shifty=zeros(size(par_civ2.Grid,1),1);
    662656        nbval=zeros(size(par_civ2.Grid,1),1);% nbre of interpolated values at each grid point (from the different patch subdomains)
     657        if par_civ2.CheckDeformation
     658            DUDX=zeros(size(par_civ2.Grid,1),1);
     659            DUDY=zeros(size(par_civ2.Grid,1),1);
     660            DVDX=zeros(size(par_civ2.Grid,1),1);
     661            DVDY=zeros(size(par_civ2.Grid,1),1);
     662        end
    663663        NbSubDomain=size(SubRange,3);
    664664        for isub=1:NbSubDomain% for each sub-domain of Patch1
    665665            nbvec_sub=NbCentres(isub);% nbre of Civ vectors in the subdomain
    666666            ind_sel=find(par_civ2.Grid(:,1)>=SubRange(1,1,isub) & par_civ2.Grid(:,1)<=SubRange(1,2,isub) &...
    667                 par_civ2.Grid(:,2)>=SubRange(2,1,isub) & par_civ2.Grid(:,2)<=SubRange(2,2,isub));
     667                par_civ2.Grid(:,2)>=SubRange(2,1,isub) & par_civ2.Grid(:,2)<=SubRange(2,2,isub));% grid points in the subdomain
    668668            if ~isempty(ind_sel)
    669                 epoints = par_civ2.Grid(ind_sel,:);% coordinates of interpolation sites
     669                epoints = par_civ2.Grid(ind_sel,:);% coordinates of interpolation sites (measurement grids)
    670670                ctrs=Coord_tps(1:nbvec_sub,:,isub) ;%(=initial points) ctrs
    671671                nbval(ind_sel)=nbval(ind_sel)+1;% records the number of values for each interpolation point (in case of subdomain overlap)
    672                 EM = tps_eval(epoints,ctrs);
    673                 Shiftx(ind_sel)=Shiftx(ind_sel)+EM*U_tps(1:nbvec_sub+3,isub);
     672                EM = tps_eval(epoints,ctrs);% thin plate spline (tps) coefficient
     673                Shiftx(ind_sel)=Shiftx(ind_sel)+EM*U_tps(1:nbvec_sub+3,isub);%velocity shift estimated by tps from civ1
    674674                Shifty(ind_sel)=Shifty(ind_sel)+EM*V_tps(1:nbvec_sub+3,isub);
    675675                if par_civ2.CheckDeformation
    … …  
    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;
     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;
    697697        if CheckInputFile % else Dt given by par_civ2
    698698            if strcmp(Param.ActionInput.ListCompareMode,'displacement')
    … …  
    704704        end
    705705        par_civ2.SearchBoxShift=(Civ2_Dt/Civ1_Dt)*[Shiftx(nbval>=1)./nbval(nbval>=1) Shifty(nbval>=1)./nbval(nbval>=1)];
    706         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];% grid taken at the extrapolated origin of the displacement vectors
     706        % shift the grid points by half the expected shift to provide the correlation box position in image A
     707        par_civ2.Grid=[par_civ2.Grid(nbval>=1,1)-par_civ2.SearchBoxShift(nbval>=1,1)/2 par_civ2.Grid(nbval>=1,2)-par_civ2.SearchBoxShift(nbval>=1,2)/2];
    707708        if par_civ2.CheckDeformation
    708709            par_civ2.DUDX=DUDX./nbval;
    … …  
    712713        end
    713714       
    714         % calculate velocity data (y and v in indices, reverse to y component)
     715        % calculate velocity data (y and v in image indices, reverse to y component)
    715716        [xtable, ytable, utable, vtable, ctable, F,result_conv,errormsg] = civ (par_civ2);
    716717       
    … …  
    880881%  .ImageWidth: nb of pixels of the image in x
    881882%  .Dx, Dy: mesh for the PIV calculation
    882 %  .Grid: grid giving the PIV calculation points (alternative to .Dx .Dy)
     883%  .Grid: grid giving the PIV calculation points (alternative to .Dx .Dy): centres of the correlation boxes in Image A
    883884%  .Mask: name of a mask file or mask image matrix itself
    884885%  .MinIma: thresholds for image luminosity
    … …  
    893894
    894895%% prepare measurement grid
    895 if isfield(par_civ,'Grid')% grid points set as input
    896     if ischar(par_civ.Grid)%read the grid file if the input is a file name
     896if isfield(par_civ,'Grid')% grid points set as input, central positions of the sub-images in image A
     897    if ischar(par_civ.Grid)%read the grid file if the input is a file name (grid in x, y image coordinates)
    897898        par_civ.Grid=dlmread(par_civ.Grid);
    898899        par_civ.Grid(1,:)=[];%the first line must be removed (heading in the grid file)
    899900    end
    900     % else par_civ.Grid is already an array
    901 else% automatic grid
     901    % else par_civ.Grid is already an array, no action here
     902else% automatic grid in x, y image coordinates
    902903    minix=floor(par_civ.Dx/2)-0.5;
    903904    maxix=minix+par_civ.Dx*floor((par_civ.ImageWidth-1)/par_civ.Dx);
    904     miniy=floor(par_civ.Dy/2)-0.5;
     905    miniy=floor(par_civ.Dy/2)-0.5;% first automatic grid point at half the mesh Dy
    905906    maxiy=minix+par_civ.Dy*floor((par_civ.ImageHeight-1)/par_civ.Dy);
    906907    [GridX,GridY]=meshgrid(minix:par_civ.Dx:maxix,miniy:par_civ.Dy:maxiy);
    907908    par_civ.Grid(:,1)=reshape(GridX,[],1);
    908     par_civ.Grid(:,2)=reshape(GridY,[],1);
     909    par_civ.Grid(:,2)=reshape(GridY,[],1);% increases with array index
    909910end
    910911nbvec=size(par_civ.Grid,1);
    … …  
    915916isx2=ceil(par_civ.SearchBoxSize(1)/2);
    916917isy2=ceil(par_civ.SearchBoxSize(2)/2);
    917 shiftx=round(par_civ.SearchBoxShift(:,1));
     918shiftx=round(par_civ.SearchBoxShift(:,1));%use the input shift estimate, rounded to the next integer value
    918919shifty=-round(par_civ.SearchBoxShift(:,2));% sign minus because image j index increases when y decreases
    919920if numel(shiftx)==1% case of a unique shift for the whole field( civ1)
    … …  
    922923end
    923924
    924 %% Default output
    925 xtable=par_civ.Grid(:,1);
    926 ytable=par_civ.Grid(:,2);
    927 utable=zeros(nbvec,1);
    928 vtable=zeros(nbvec,1);
     925%% Array initialisation and default output  if par_civ.CorrSmooth=0 (just the grid calculated, no civ computation)
     926xtable=round(par_civ.Grid(:,1)+0.5)-0.5;
     927ytable=round(par_civ.ImageHeight-par_civ.Grid(:,2)+0.5)-0.5;% y index corresponding to the position in image coordiantes
     928utable=shiftx;%zeros(nbvec,1);
     929vtable=shifty;%zeros(nbvec,1);
    929930ctable=zeros(nbvec,1);
    930931F=zeros(nbvec,1);
    … …  
    943944check_MaxIma=isfield(par_civ,'MaxIma') && ~isempty(par_civ.MaxIma);
    944945
    945 % %% prepare images
    946 % if isfield(par_civ,'reverse_pair')
    947 %     if par_civ.reverse_pair
    948 %         if ischar(par_civ.ImageB)
    949 %             temp=par_civ.ImageA;
    950 %             par_civ.ImageA=imread(par_civ.ImageB);
    951 %         end
    952 %         if ischar(temp)
    953 %             par_civ.ImageB=imread(temp);
    954 %         end
    955 %     end
    956 % else
    957 %     if ischar(par_civ.ImageA)
    958 %         par_civ.ImageA=imread(par_civ.ImageA);
    959 %     end
    960 %     if ischar(par_civ.ImageB)
    961 %         par_civ.ImageB=imread(par_civ.ImageB);
    962 %     end
    963 % end
    964946par_civ.ImageA=sum(double(par_civ.ImageA),3);%sum over rgb component for color images
    965947par_civ.ImageB=sum(double(par_civ.ImageB),3);
    … …  
    1000982    mesh=0.25;%mesh in pixels for subpixel image interpolation (x 4 in each direction)
    1001983end
    1002 if par_civ.CorrSmooth~=0 % par_civ.CorrSmooth=0 implies no civ computation (just input image and grid points viven)
     984
     985if par_civ.CorrSmooth~=0 % par_civ.CorrSmooth=0 implies no civ computation (just input image and grid points given)
    1003986    for ivec=1:nbvec
    1004987        iref=round(par_civ.Grid(ivec,1)+0.5);% xindex on the image A for the middle of the correlation box
    1005         jref=round(par_civ.ImageHeight-par_civ.Grid(ivec,2)+0.5);% yindex on the image B for the middle of the correlation box
     988        jref=round(par_civ.ImageHeight-par_civ.Grid(ivec,2)+0.5);%  j index  for the middle of the correlation box in the image A
    1006989        F(ivec)=0;
    1007990        subrange1_x=iref-ibx2:iref+ibx2;% x indices defining the first subimage
    … …  
    10761059            else
    10771060                F(ivec)=3;
    1078                 [y,x]
    10791061            end
    10801062        end
    … …  
    10851067%------------------------------------------------------------------------
    10861068% --- Find the maximum of the correlation function after interpolation
     1069% OUPUT:
     1070% vector = optimum displacement vector with subpixel correction
     1071% F =flag: =0 OK
     1072%           =-2 , warning: max too close to the edge of the search box (1 pixel margin)
     1073% INPUT:
     1074% x,y: position of the maximum correlation at integer values
     1075
    10871076function [vector,F] = SUBPIXGAUSS (result_conv,x,y)
    10881077%------------------------------------------------------------------------
    1089 vector=[0 0]; %default
     1078% vector=[0 0]; %default
    10901079F=0;
    10911080[npy,npx]=size(result_conv);
    … …  
    11171106function [vector,F] = SUBPIX2DGAUSS (result_conv,x,y)
    11181107%------------------------------------------------------------------------
    1119 vector=[0 0]; %default
     1108% vector=[0 0]; %default
    11201109F=-2;
    11211110peaky=y;

  • trunk/src/series/stereo_civ.m

    r854 r862  
    5050    path_series=fileparts(which('series'));
    5151    addpath(fullfile(path_series,'series'))
    52     Data=civ_input(Param);% introduce the civ parameters using the GUI civ_input
     52   % Data=civ_input(Param);% introduce the civ parameters using the GUI civ_input
     53   Data=stereo_input(Param)
    5354    if isempty(Data)
    5455        Data=Param;% if  civ_input has been cancelled, keep previous parameters

  • trunk/src/sub_field.m

    r811 r862  
    150150
    151151%% append the other variables of the second field, modifying their name if needed
    152 for ilist=1:numel(Field_1.ListVarName)
     152ListVarNameNew=Field_1.ListVarName;
     153check_rename=zeros(size(ListVarNameNew));
     154for ilist=1:numel(ListVarNameNew)
    153155    VarName=Field_1.ListVarName{ilist};
    154     ind_prev=find(strcmp(VarName,Field.ListVarName));
    155     if isempty(ind_prev)% variable name does not exist in Field
    156         VarNameNew=VarName;
    157     else  % variable name exists in Field     
    158             VarNameNew=[VarName '_1'];   
    159             if isfield(Field_1.VarAttribute{ilist},'FieldName')
    160                 Field_1.VarAttribute{ilist}.FieldName=regexprep_r(Field_1.VarAttribute{ilist}.FieldName,VarName,VarNameNew);
    161             end
    162     end
    163         SubData.ListVarName=[SubData.ListVarName {VarNameNew}];
    164         SubData.VarDimName=[SubData.VarDimName Field_1.VarDimName(ilist)];
    165         SubData.(VarNameNew)=Field_1.(VarName);
    166         SubData.VarAttribute=[SubData.VarAttribute Field_1.VarAttribute(ilist)];
    167         SubData.VarAttribute{end}.CheckSub=1;% mark that the field needs to be substracted
     156    ind_prev=find(strcmp(ListVarNameNew{ilist},Field.ListVarName),1);% look for duplicated variable name
     157    if ~isempty(ind_prev)% variable name exists in Field
     158        check_rename(ilist)=1;
     159        ListVarNameNew{ilist}=[ListVarNameNew{ilist} '_1'];
     160    end
     161    SubData.ListVarName=[SubData.ListVarName ListVarNameNew{ilist}];
     162    SubData.VarDimName=[SubData.VarDimName Field_1.VarDimName(ilist)];
     163    SubData.(ListVarNameNew{ilist})=Field_1.(VarName);% teke the values of the old variable for the newly named one
     164    %SubData.VarAttribute=[SubData.VarAttribute Field_1.VarDimName(ilist)];
     165end
     166
     167%% replace variable name in field expression FieldName, e.g. 'norm(U,V)'-> 'norm(U_1,V_1)'
     168for ilist=1:numel(ListVarNameNew)
     169    if check_rename(ilist)&&  isfield(Field_1.VarAttribute{ilist},'FieldName')
     170        for ivar=1:numel(find(check_rename))
     171            Field_1.VarAttribute{ilist}.FieldName=regexprep_r(Field_1.VarAttribute{ilist}.FieldName,...
     172                Field_1.ListVarName{ivar},ListVarNameNew{ivar});
     173        end
     174    end
     175    SubData.VarAttribute=[SubData.VarAttribute Field_1.VarAttribute(ilist)];
     176    SubData.VarAttribute{end}.CheckSub=1;% mark that the field needs to be substracted as an attribute
    168177end
    169178
    … …  
    198207SubData.VarDimName(find(ind_remove))=[];
    199208SubData.VarAttribute(find(ind_remove))=[];
    200 'end'
    201209
    202210function OutputCell=regexprep_r(InputCell,search_string,new_string)
Note: See TracChangeset for help on using the changeset viewer.