Ignore:
Timestamp:
Aug 10, 2020, 3:20:45 PM (4 years ago)
Author:
sommeria
Message:

various updates

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/transform_field/phys_polar.m

    r1081 r1084  
    4949%% request input parameters
    5050if isfield(DataIn,'Action') && isfield(DataIn.Action,'RUN') && isequal(DataIn.Action.RUN,0)
    51     prompt = {'origin [x y] of polar coordinates';'reference radius';'reference angle(degrees)'};
     51    prompt = {'origin [x y] of polar coordinates';'reference radius';'reference angle(degrees)';'angle direction and switch x y(+/-)'};
    5252    dlg_title = 'set the parameters for the polar coordinates';
    5353    num_lines= 2;
    54     def     = { '[0 0]';'0';'0'};
     54    def     = { '[0 0]';'0';'0';'+'};
    5555    if isfield(XmlData,'TransformInput')
    5656        if isfield(XmlData.TransformInput,'PolarCentre')
     
    6262        if isfield(XmlData.TransformInput,'PolarReferenceAngle')
    6363            def{3}=num2str(XmlData.TransformInput.PolarReferenceAngle);
     64        end
     65        if isfield(XmlData.TransformInput,'PolarAngleDirection')
     66            def{4}=XmlData.TransformInput.PolarAngleDirection;
    6467        end
    6568    end
     
    6871    Data.TransformInput.PolarReferenceRadius=str2num(answer{2});
    6972    Data.TransformInput.PolarReferenceAngle=str2num(answer{3});
    70 %     if isfield(XmlData,'GeometryCalib')&& isfield(XmlData.GeometryCalib,'CoordUnit')
    71 %         Data.CoordUnit=XmlData.GeometryCalib.CoordUnit;% states that the output is in unit defined by GeometryCalib, then erased all projection objects with different units
    72 %     end
     73    Data.TransformInput.PolarAngleDirection=answer{4};
    7374    return
    7475end
     
    136137        angle_offset=(pi/180)*XmlData.TransformInput.PolarReferenceAngle; %offset angle (in unit of the final angle, degrees or arc length along the reference radius))
    137138    end
     139    check_reverse=isfield(XmlData.TransformInput,'PolarAngleDirection')&& strcmp(XmlData.TransformInput.PolarAngleDirection,'-');
    138140end
    139141
    140142%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    141143%% get fields
    142 check_scalar=0;
    143 check_vector=0;
     144
    144145nbvar=0;%counter for the number of output variables
    145146nbcoord=0;%counter for the number of variables for radial coordiantes (case of multiple field inputs)
     
    163164        ZIndex=0;
    164165    end
     166    check_scalar=zeros(1,numel(CellInfo));
     167    check_vector=zeros(1,numel(CellInfo));
    165168    for icell=1:numel(CellInfo)
    166169        if NbDim(icell)==2
     
    252255                    Data.VarDimName=[Data.VarDimName {radius_name} {theta_name}];
    253256                    nbvar=nbvar+2;
     257                    if check_reverse
     258                                            Data.VarAttribute{nbvar-1}.Role='coord_y';
     259                    Data.VarAttribute{nbvar}.Role='coord_x';
     260                    else
    254261                    Data.VarAttribute{nbvar-1}.Role='coord_x';
    255262                    Data.VarAttribute{nbvar}.Role='coord_y';
     263                    end
    256264                    check_unit=1;
    257265                    if isfield(DataCell{ifield},'CoordUnit')
     
    315323    [A,Data.radius,Data.theta]=phys_Ima_polar(A,coord_x,coord_y,Calib_new,ZInd,origin_xy,radius_offset,angle_offset,angle_scale);
    316324    for icell=1:numel(A)
    317         if icell<=numel(A)-1 && check_vector(icell)==1 && check_vector(icell+1)==1   %transform u,v into polar coordiantes
     325        if icell<=numel(A)-1 && check_vector(icell)==1 && check_vector(icell+1)==1   %transform u,v into polar coordinates
    318326            theta=Data.theta/angle_scale-angle_offset;
    319327            [~,Theta]=meshgrid(Data.radius,theta);%grid in physical coordinates
    320328            U_r_name= rename_indexing(U_r_name,Data.ListVarName);
    321             U_theta_name= rename_indexing(U_theta_name,Data.ListVarName);
    322             Data.ListVarName=[Data.ListVarName {U_r_name,U_theta_name}];
    323             Data.VarDimName=[Data.VarDimName {{theta_name,radius_name}} {{theta_name,radius_name}}];
    324             Data.(U_r_name)=A{icell}.*cos(Theta)+A{icell+1}.*sin(Theta);%radial velocity
    325             Data.(U_theta_name)=(-A{icell}.*sin(Theta)+A{icell+1}.*cos(Theta));%./(Data.X)%+radius_ref);% azimuthal velocity component
     329            U_theta_name= rename_indexing(U_theta_name,Data.ListVarName);       
     330                Data.(U_r_name)=A{icell}.*cos(Theta)+A{icell+1}.*sin(Theta);%radial velocity
     331                Data.(U_theta_name)=(-A{icell}.*sin(Theta)+A{icell+1}.*cos(Theta));% azimuthal velocity component
     332            if check_reverse
     333                Data.(U_theta_name)=(Data.(U_theta_name))';
     334                Data.(U_r_name)=Data.(U_r_name)';
     335                Data.ListVarName=[Data.ListVarName {U_theta_name,U_r_name}];
     336                Data.VarDimName=[Data.VarDimName {{radius_name,theta_name}} {{radius_name,theta_name}}];
     337            else
     338                Data.ListVarName=[Data.ListVarName {U_r_name,U_theta_name}];
     339                Data.VarDimName=[Data.VarDimName {{theta_name,radius_name}} {{theta_name,radius_name}}];
     340            end
    326341        elseif ~check_vector(icell)% for scalar fields
    327342            FieldName{icell}= rename_indexing(FieldName{icell},Data.ListVarName);
    328             Data.ListVarName=[Data.ListVarName {FieldName{icell}}];
    329             Data.VarDimName=[Data.VarDimName {{theta_name,radius_name}}];
    330             Data.(FieldName{icell})=A{icell};
    331         end
    332     end
     343            Data.ListVarName=[Data.ListVarName FieldName(icell)];       
     344            if check_reverse
     345                Data.(FieldName{icell})=A{icell}';
     346                Data.VarDimName=[Data.VarDimName {{radius_name,theta_name}}];
     347            else
     348                Data.VarDimName=[Data.VarDimName {{theta_name,radius_name}}];
     349                Data.(FieldName{icell})=A{icell};
     350            end
     351        end
     352    end
     353end
     354if check_reverse
     355    Data.(theta_name)=-Data.(theta_name);
    333356end
    334357
     
    406429radius=Min_r:Dr:Max_r;% polar coordinates for projections
    407430theta=Min_theta:Dtheta:Max_theta;
     431%theta=Max_theta:-Dtheta:Min_theta;
    408432[Radius,Theta]=meshgrid(radius,theta/angle_scale);%grid in polar coordinates (angles in radians)
    409433%transform X, Y in cartesian
     
    436460    Dy=(coord_y{icell}(end)-coord_y{icell}(1))/(npy(icell)-1);
    437461    indx_ima=1+round((XIMA-coord_x{icell}(1))/Dx);%indices of the initial matrix close to the points of the new grid
    438     indy_ima=1+round((YIMA-coord_y{icell}(1))/Dy);
    439     Delta_x=1+(XIMA-coord_x{icell}(1))/Dx-indx_ima;%
    440     Delta_y=1+(YIMA-coord_y{icell}(1))/Dy-indy_ima;
     462    %indy_ima=1+round((YIMA-coord_y{icell}(1))/Dy);
     463    indy_ima=1+round((coord_y{icell}(end)-YIMA)/Dy);
     464     Delta_x=1+(XIMA-coord_x{icell}(1))/Dx-indx_ima;%error in the index discretisation
     465     Delta_y=1+(coord_y{icell}(end)-YIMA)/Dy-indy_ima;
    441466    XIMA=reshape(indx_ima,1,[]);%indices reorganized in 'line'
    442467    YIMA=reshape(indy_ima,1,[]);%indices reorganized in 'line'
     
    450475        ind_in=find(flagin);
    451476        ind_out=find(~flagin);
    452         ICOMB=((XIMA-1)*npy(icell)+(npy(icell)+1-YIMA));
     477        ICOMB=((XIMA-1)*npy(icell)+(npy(icell)+1-YIMA));% indices in vec_A
    453478        ICOMB=ICOMB(flagin);%index corresponding to XIMA and YIMA in the aligned original image vec_A
    454479        vec_B(ind_in)=vec_A(ICOMB);
    455480        vec_B(ind_out)=zeros(size(ind_out));
    456481        A_out{icell}=reshape(vec_B,np_theta,np_r);%new image in real coordinates
    457         DA_y=circshift(A_out{icell},-1,1)-A_out{icell};
     482        DA_y=circshift(A_out{icell},-1,1)-A_out{icell};% derivative
    458483        DA_y(end,:)=0;
    459484        DA_x=circshift(A_out{icell},-1,2)-A_out{icell};
Note: See TracChangeset for help on using the changeset viewer.