Ignore:
Timestamp:
May 21, 2018, 7:06:45 PM (7 years ago)
Author:
sommeria
Message:

find field cells improved

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/find_field_cells.m

    r1041 r1045  
    1818%              _vector_x,_y,_z: indices of variables giving the vector components x, y, z
    1919%              _warnflag: index of warnflag   
     20%       .DimIndex
    2021%      .ProjModeRequest= 'interp_lin', 'interp_tps' indicate whether lin interpolation  or derivatives (tps) is needed to calculate the requested field
    2122%      .FieldName = operation to be performed to finalise the field cell after projection
     
    6465
    6566function [CellInfo,NbDim,errormsg]=find_field_cells(Data)
    66 CellInfo={};
     67CellInfo={};%default output
    6768NbDim=0;
    6869errormsg='';
     
    8788FieldName=cell(size(Role)); % fieldRequest set to {} by default
    8889CheckSub=zeros(size(Role));% =1 for fields to substract
    89 Role=regexprep(Role,' ','scalar'); % Role set to 'scalar' by default
     90%Role=regexprep(Role,' ','scalar'); % Role set to 'scalar' by default
    9091if isfield(Data,'VarAttribute')
    9192    for ivar=1:numel(Data.VarAttribute)
     
    105106end
    106107
    107 %% find scattered (unstructured) coordinates
    108 ivar_coord_x=find(strcmp('coord_x',Role));%find variables with Role='coord_x'
    109 check_select=false(1,nbvar);
    110 check_coord=false(1,nbvar);
    111 CellInfo=cell(1,numel(ivar_coord_x));
    112 NbDim=zeros(1,numel(ivar_coord_x));
    113 % loop on unstructured coordinate x -> different field cells
    114 for icell=1:numel(ivar_coord_x)
    115     DimCell=Data.VarDimName{ivar_coord_x(icell)};% cell of dimension names for ivar_coord_x(icell)
    116     if ischar(DimCell),DimCell={DimCell};end % transform char to cell for a single dimension
    117     % look for variables sharing dimension(s) with ivar_coord_x(icell)
    118     check_cell=zeros(numel(DimCell),nbvar);
    119     for idim=1:numel(DimCell);% for each variable with role coord_x, look at which other variables contain the same dimension
    120         for ivar=1:nbvar
    121             check_cell(idim,ivar)=max(strcmp(DimCell{idim},Data.VarDimName{ivar}));
    122         end
    123     end
    124     check_cell=sum(check_cell,1)==numel(DimCell);%logical array=1 for variables belonging to the current cell
    125     VarIndex=find(check_cell);% list of detected variable indices
    126     if ~(numel(VarIndex)==1 && numel(DimCell)==1)% exclude case of isolated coord_x variable (treated later)
    127         if ~(numel(VarIndex)==1 && numel(DimCell)>1)% a variable is associated to coordinate
    128             CellInfo{icell}.CoordIndex=ivar_coord_x(icell);
    129             % size of coordinate var
     108%% detect  fields with different roles
     109ind_scalar=find(strcmp('scalar',Role));
     110ind_errorflag=find(strcmp('errorflag',Role));
     111ind_image=find(strcmp('image',Role));
     112ind_vector_x=[find(strcmp('vector_x',Role)) find(strcmp('vector_x_tps',Role))];
     113if ~isempty(ind_vector_x)
     114    ind_vector_y=[find(strcmp('vector_y',Role)) find(strcmp('vector_y_tps',Role))];
     115    ind_vector_z=find(strcmp('vector_z',Role));
     116    ind_warnflag=find(strcmp('warnflag',Role));
     117    ind_ancillary=find(strcmp('ancillary',Role));
     118end
     119ind_discrete=find(strcmp('discrete',Role));
     120ind_coord_x=find(strcmp('coord_x',Role));
     121ind_coord_y=find(strcmp('coord_y',Role));
     122ind_coord_z=find(strcmp('coord_z',Role));
     123ind_coord_tps=find(strcmp('coord_tps',Role));
     124check_string=cellfun(@ischar,Data.VarDimName)==1;
     125index_string=find(check_string);
     126for ivar=index_string
     127    Data.VarDimName{ivar}={Data.VarDimName{ivar}};%transform char strings into cells
     128end
     129check_coord_names= cellfun(@numel,Data.VarDimName)==1;
     130check_coord_raster=false(size(check_coord_names));% check variables describing regular mesh (raster coordinates), from two values, min and max.
     131if check_var
     132    for ivar=find(check_coord_names)
     133        if numel(Data.(Data.ListVarName{ivar}))==2
     134            check_coord_raster(ivar)=true;
     135        end
     136    end
     137else
     138    for ivar=find(check_coord_names)
     139        DimIndex=find(strcmp(Data.VarDimName{ivar},Data.ListDimName));
     140        if Data.DimValue(DimIndex)==2
     141            check_coord_raster(ivar)=true;
     142        end
     143    end
     144end
     145
     146
     147%% initate cells around each scalar field
     148index_remove=[];
     149cell_nbre=numel(ind_scalar)+numel(ind_vector_x);
     150flag_remove=false(1,cell_nbre);
     151NbDim=zeros(1,cell_nbre);
     152index_coord_x=zeros(size(ind_coord_x));
     153for icell=1:numel(ind_scalar)
     154    CellInfo{icell}.VarType='scalar';
     155    CellInfo{icell}.VarIndex_scalar=ind_scalar(icell);
     156    CellInfo{icell}.VarIndex=ind_scalar(icell);
     157    DimCell_var=Data.VarDimName{ind_scalar(icell)};% cell of dimension names for ivar_coord_x(icell)
     158    %look for errorflag
     159    for ivar=ind_errorflag
     160        DimCell=Data.VarDimName{ivar};
     161        if isequal(DimCell,DimCell_var)
     162            CellInfo{icell}.VarIndex(2)=ivar;
     163            CellInfo{icell}.VarIndex_errorflag=ivar;
     164            break
     165        end
     166    end
     167end
     168
     169%% initate cells around each vector field
     170for icell=numel(ind_scalar)+1:cell_nbre
     171    CellInfo{icell}.VarType='vector';
     172    CellInfo{icell}.VarIndex(1)=ind_vector_x(icell-numel(ind_scalar));
     173    CellInfo{icell}.VarIndex_vector_x=ind_vector_x(icell-numel(ind_scalar));
     174    DimCell_var=Data.VarDimName{ind_vector_x(icell-numel(ind_scalar))};% cell of dimension names for ivar_coord_x(icell)
     175    % look for the associated y vector component
     176    nbvar=1;
     177    for ivar=ind_vector_y
     178        DimCell=Data.VarDimName{ivar};
     179        if isequal(DimCell,DimCell_var)
     180            CellInfo{icell}.VarIndex(2)=ivar;
     181            nbvar=2;
     182            CellInfo{icell}.VarIndex_vector_y=ivar;
     183            break
     184        end
     185    end
     186    if ~isfield(CellInfo{icell},'VarIndex_vector_y')
     187        flag_remove(icell)=true;% no vector_y found , mark cell to remove
     188    end
     189    % look for the associated z vector component
     190    for ivar=ind_vector_z
     191        DimCell=Data.VarDimName{ivar};
     192        if isequal(DimCell,DimCell_var)
     193            CellInfo{icell}.VarIndex(3)=ivar;
     194            nbvar=3;
     195            break
     196        end
     197    end
     198    %look for the vector color scalar (ancillary)
     199    for ivar=ind_ancillary
     200        DimCell=Data.VarDimName{ivar};
     201        if isequal(DimCell,DimCell_var)
     202            nbvar=nbvar+1;
     203            CellInfo{icell}.VarIndex(nbvar)=ivar;
     204            CellInfo{icell}.VarIndex_ancillary=ivar;
     205            break
     206        end
     207    end
     208    %look for warnflag
     209    for ivar=ind_warnflag
     210        DimCell=Data.VarDimName{ivar};
     211        if isequal(DimCell,DimCell_var)
     212            nbvar=nbvar+1;
     213            CellInfo{icell}.VarIndex(nbvar)=ivar;
     214            CellInfo{icell}.VarIndex_warnflag=ivar;
     215            break
     216        end
     217    end
     218    %look for errorflag
     219    for ivar=ind_errorflag
     220        DimCell=Data.VarDimName{ivar};
     221        if isequal(DimCell,DimCell_var)
     222            nbvar=nbvar+1;
     223            CellInfo{icell}.VarIndex(nbvar)=ivar;
     224            CellInfo{icell}.VarIndex_errorflag=ivar;
     225            break
     226        end
     227    end
     228end
     229
     230%% find coordinates for each cell around field variables, scalars or vectors
     231for icell=1:cell_nbre
     232    CellInfo{icell}.CoordType='';
     233    ind_var=CellInfo{icell}.VarIndex(1);
     234    DimCell_var=Data.VarDimName{ind_var};% cell of dimension names for ivar_coord_x(icell)
     235    if ~check_var
     236        for idim=1:numel(DimCell_var)
     237            CellInfo{icell}.DimIndex(idim)=find(strcmp(DimCell_var{idim},Data.ListDimName));
     238        end
     239    end
     240    %look for z scattered coordinates
     241    if isempty(ind_coord_z)
     242        NbDim(icell)=2;
     243        CellInfo{icell}.CoordIndex=[0 0];
     244    else
     245        NbDim(icell)=3;
     246        CellInfo{icell}.CoordIndex=[0 0 0];
     247        for ivar=ind_coord_z
     248            DimCell=Data.VarDimName{ivar};
     249            if isequal(DimCell,DimCell_var)
     250                CellInfo{icell}.CoordType='scattered';
     251                CellInfo{icell}.CoordIndex(1)=ivar;
     252                CellInfo{icell}.ZName=Data.ListVarName{ivar};
     253                CellInfo{icell}.ZIndex=ivar;
     254                break
     255            end
     256        end
     257    end
     258    % look for y coordinate
     259    for ivar=ind_coord_y
     260        % detect scattered y coordinates, variable with the same dimension(s) as the field variable considered
     261        DimCell=Data.VarDimName{ivar};
     262        if isequal(DimCell,DimCell_var)
     263            CellInfo{icell}.CoordType='scattered';
     264            CellInfo{icell}.CoordIndex(NbDim(icell)-1)=ivar;
     265            CellInfo{icell}.YName=Data.ListVarName{ivar};
     266            CellInfo{icell}.YIndex=ivar;
     267            break
     268        end
     269    end
     270   
     271    %look for x coordinates
     272    if strcmp(CellInfo{icell}.CoordType,'scattered')
     273        for ivar=ind_coord_x
     274            DimCell=Data.VarDimName{ivar};
     275            if isequal(DimCell,DimCell_var)
     276                CellInfo{icell}.CoordIndex(NbDim(icell))=ivar;
     277                CellInfo{icell}.XName=Data.ListVarName{ivar};
     278                CellInfo{icell}.XIndex=ivar;
     279                break
     280            end
     281        end
     282    end
     283    if isfield(CellInfo{icell},'ZName')
     284        if isfield(CellInfo{icell},'YName')&& isfield(CellInfo{icell},'XName')
     285            continue %scattered coordinates OK
     286        end
     287    else
     288        if isfield(CellInfo{icell},'YName')
     289            if isfield(CellInfo{icell},'XName')
     290                NbDim(icell)=2;
     291                continue %scattered coordinates OK
     292            end
     293        else
     294            if isfield(CellInfo{icell},'XName'); % only one coordinate x, switch vector field to 1D plot
     295                for ind=1:numel(CellInfo{icell}.VarIndex)
     296                    Role{CellInfo{icell}.VarIndex(ind)}='coord_y';
     297                end
     298                continue
     299            end
     300        end
     301    end
     302   
     303    %look for grid  coordinates
     304    if isempty(CellInfo{icell}.CoordType)
     305        NbDim(icell)=numel(DimCell_var);
     306        CellInfo{icell}.DimOrder=[];
     307        if NbDim(icell)==3
     308            %coord z
     309            for ivar=ind_coord_z
     310                if check_coord_names(ivar)
     311                    DimRank=find(strcmp(Data.VarDimName{ivar},DimCell_var));
     312                check_coord=~isempty(DimRank);
     313            elseif check_coord_raster(ivar)
     314                DimRank=find(strcmp(Data.ListVarName{ivar},DimCell_var));
     315                check_coord=~isempty(DimRank);
     316            end
     317%                 check_coord= (check_coord_names(ivar) && strcmp(Data.VarDimName{ivar},DimCell_var{1}))||...% coord varbable
     318%                     (check_coord_raster(ivar) && strcmp(Data.ListVarName{ivar},DimCell_var{1})); % rasrewr coord defined by min and max
     319                if check_coord
     320                    CellInfo{icell}.CoordType='grid';
     321                    CellInfo{icell}.CoordIndex(1)=ivar;
     322                    CellInfo{icell}.ZName=Data.ListVarName{ivar};
     323                    CellInfo{icell}.ZIndex=ivar;
     324                    CellInfo{icell}.DimOrder=DimRank;
     325                    break
     326                end
     327            end
     328        end
     329        for ivar=ind_coord_y
     330              if check_coord_names(ivar)
     331                    DimRank=find(strcmp(Data.VarDimName{ivar},DimCell_var));
     332                check_coord=~isempty(DimRank);
     333            elseif check_coord_raster(ivar)
     334                DimRank=find(strcmp(Data.ListVarName{ivar},DimCell_var));
     335                check_coord=~isempty(DimRank);
     336            end
     337%             check_coord= (check_coord_names(ivar) && strcmp(Data.VarDimName{ivar},DimCell_var{NbDim(icell)-1}))||...% coord variable
     338%                 (check_coord_raster(ivar) && strcmp(Data.ListVarName{ivar},DimCell_var{NbDim(icell)-1})); % rasrewr coord defined by min and max
     339            if check_coord
     340                CellInfo{icell}.CoordType='grid';
     341                CellInfo{icell}.CoordIndex(NbDim(icell)-1)=ivar;
     342                CellInfo{icell}.YName=Data.ListVarName{ivar};
     343                CellInfo{icell}.YIndex=ivar;
     344                CellInfo{icell}.DimOrder=[CellInfo{icell}.DimOrder DimRank];
     345                break
     346            end
     347        end
     348        for ivar=ind_coord_x
     349            if check_coord_names(ivar)
     350                    DimRank=find(strcmp(Data.VarDimName{ivar},DimCell_var));
     351                check_coord=~isempty(DimRank);
     352            elseif check_coord_raster(ivar)
     353                DimRank=find(strcmp(Data.ListVarName{ivar},DimCell_var));
     354                check_coord=~isempty(DimRank);
     355            end
     356%             check_coord= (check_coord_names(ivar) && strcmp(Data.VarDimName{ivar},DimCell_var{NbDim(icell)}))||...% coord variable
     357%                 (check_coord_raster(ivar) && strcmp(Data.ListVarName{ivar},DimCell_var{NbDim(icell)})); % raster coord defined by min and max
     358            if check_coord
     359                CellInfo{icell}.CoordIndex(NbDim(icell))=ivar;
     360                CellInfo{icell}.XName=Data.ListVarName{ivar};
     361                CellInfo{icell}.XIndex=ivar;
     362                CellInfo{icell}.DimOrder=[CellInfo{icell}.DimOrder DimRank];
     363                break
     364            end
     365        end
     366    end
     367    %look for tps coordinates
     368    for ivar=ind_coord_tps
     369        DimCell=Data.VarDimName{ivar};
     370        if  numel(DimCell)==3 && strcmp(DimCell{1},DimCell_var{1})
     371            CellInfo{icell}.CoordType='tps';
     372            CellInfo{icell}.CoordIndex=ivar;
    130373            if check_var
    131                 CellInfo{icell}.CoordSize=numel(Data.(Data.ListVarName{ivar_coord_x(icell)}));
     374                NbDim(icell)=size(Data.(Data.ListVarName{ivar}),2);
    132375            else
    133                 for idim=1:numel(DimCell)
    134                     check_index= strcmp(DimCell{idim},Data.ListDimName);
    135                     CellInfo{icell}.CoordSize(idim)=Data.DimValue(check_index);
     376                DimIndex=find(strcmp(Data.VarDimName{ivar},Data.ListDimName));
     377                NbDim(icell)= Data.DimValue(DimIndex);
     378            end
     379            for ivardim=1:numel(Data.VarDimName)
     380                if strcmp(Data.VarDimName{ivardim},DimCell{3})
     381                    CellInfo{icell}.NbCentres_tps= ivardim;% nbre of sites for each tps subdomain
     382                elseif strcmp(Data.VarDimName{ivardim}{1},DimCell{2}) && numel(Data.VarDimName{ivardim})>=3 && strcmp(Data.VarDimName{ivardim}{3},DimCell{3})
     383                    CellInfo{icell}.SubRange_tps=ivardim;% subrange definiton for tps
    136384                end
    137                 CellInfo{icell}.CoordSize=prod(CellInfo{icell}.CoordSize);
    138             end
    139             ind_scalar=find(strcmp('scalar',Role(VarIndex)));
    140             ind_vector_x=find(strcmp('vector_x',Role(VarIndex)));
    141             ind_vector_y=find(strcmp('vector_y',Role(VarIndex)));
    142             ind_y=find(strcmp('coord_y',Role(VarIndex)));
    143             if numel([ind_scalar ind_vector_x ind_vector_y])==0           
    144 %             if numel(VarIndex)==2||isempty(ind_y)% no variable, except possibly y
    145                 NbDim(icell)=1;
    146             else
    147                 CellInfo{icell}.CoordType='scattered';
    148                 ind_z=find(strcmp('coord_z',Role(VarIndex)));
    149                 if numel(VarIndex)==3||isempty(ind_z)% no z variable, except possibly as a fct z(x,y)
    150                     CellInfo{icell}.CoordIndex=[VarIndex(ind_y) CellInfo{icell}.CoordIndex];
    151                     NbDim(icell)=2;
     385            end
     386        end
     387        break
     388    end
     389end
     390
     391%% get number of coordinate points for each cell
     392if check_var
     393    for icell=1:numel(CellInfo)
     394        switch CellInfo{icell}.CoordType
     395            case 'scattered'
     396                CellInfo{icell}.CoordSize=numel(Data.(CellInfo{icell}.XName));
     397            case 'grid'
     398                if NbDim(icell)==3
     399                    CellInfo{icell}.CoordSize=[numel(Data.(CellInfo{icell}.XName)) numel(Data.(CellInfo{icell}.YName)) numel(Data.(CellInfo{icell}.YName))];
    152400                else
    153                     CellInfo{icell}.CoordIndex=[VarIndex(ind_z) CellInfo{icell}.CoordIndex];
    154                     NbDim(icell)=3;
     401                    CellInfo{icell}.CoordSize=[numel(Data.(CellInfo{icell}.XName)) numel(Data.(CellInfo{icell}.YName))];
    155402                end
    156             end
    157         end
    158         CellInfo{icell}.VarIndex=VarIndex;
    159         check_select=check_select|check_cell;
    160     end
    161 end
    162 
    163 %% look for tps coordinates
    164 ivar_remain=find(~check_select);% indices of remaining variables (not already selected)
    165 check_coord_tps= strcmp('coord_tps',Role(~check_select));
    166 ivar_tps=ivar_remain(check_coord_tps);% variable indices corresponding to tps coordinates
    167 
    168 % loop on the tps coordinate sets
    169 for icell_tps=1:numel(ivar_tps)
    170     check_cell=zeros(1,nbvar);% =1 for the variables selected in the current cell
    171     check_cell(ivar_tps(icell_tps))=1;% mark the coordinate variable as selected
    172     DimCell=Data.VarDimName{ivar_tps(icell_tps)};% dimension names for the current tps coordinate variable
    173     icell=numel(CellInfo)+icell_tps; % new field cell index
    174     CellInfo{icell}.CoordIndex=ivar_tps(icell_tps);% index of the  tps coordinate variable
    175     if numel(DimCell)==3
    176         VarDimName=Data.VarDimName(~check_select);
    177         for ivardim=1:numel(VarDimName)
    178             if strcmp(VarDimName{ivardim},DimCell{3})
    179                 CellInfo{icell}.NbCentres_tps= ivar_remain(ivardim);% nbre of sites for each tps subdomain
    180                 check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
    181             elseif strcmp(VarDimName{ivardim}{1},DimCell{2}) && numel(VarDimName{ivardim})>=3 && strcmp(VarDimName{ivardim}{3},DimCell{3})
    182                 CellInfo{icell}.SubRange_tps=ivar_remain(ivardim);% subrange definiton for tps
    183                 check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
    184             elseif strcmp(VarDimName{ivardim}{1},DimCell{1}) && strcmp(VarDimName{ivardim}{2},DimCell{3})% variable
    185                 check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
    186             end
    187         end
    188     end
    189     CellInfo{icell}.CoordType='tps';
    190     CellInfo{icell}.VarIndex=find(check_cell);
    191     if check_var
    192         NbDim(icell)=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),2);
    193         CellInfo{icell}.CoordSize=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),1)*size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),3);
    194     else
    195         check_index_1= strcmp(DimCell{1},Data.ListDimName);
    196         check_index_2= strcmp(DimCell{2},Data.ListDimName);
    197         NbDim(icell)=Data.DimValue(check_index_2);
    198         if numel(DimCell)>=3
    199         check_index_3= strcmp(DimCell{3},Data.ListDimName); 
    200         CellInfo{icell}.CoordSize=Data.DimValue(check_index_1)*Data.DimValue(check_index_3);
    201         end
    202     end
    203     check_select=check_select|check_cell;
    204 end
    205 
    206 %% look for coordinate variables and corresponding gridded data:
    207 % coordinate variables are variables associated with a single dimension, defining the coordinate values
    208 % two cases: 1)the coordiante variable represents the set of coordiante values
    209 %            2)the coordinate variable contains only two elements, representing the coordinate bounds for the dimension with the same name as the cordinate
    210 ivar_remain=find(~check_select);% indices of remaining variables, not already taken into account
    211 ListVarName=Data.ListVarName(~check_select);%list of names of remaining variables
    212 VarDimName=Data.VarDimName(~check_select);%dimensions of remaining variables
    213 check_coord_select= cellfun(@numel,VarDimName)==1|cellfun(@ischar,VarDimName)==1;% find remaining variables with a single dimension
    214 check_coord_select=check_coord_select & ~strcmp('ancillary',Role(~check_select));% do not select ancillary variables as coordinates
    215 ListCoordIndex=ivar_remain(check_coord_select);% indices of remaining variables with a single dimension
    216 ListCoordName=ListVarName(check_coord_select);% corresponding names of remaining variables with a single dimension
    217 ListDimName=VarDimName(check_coord_select);% dimension names of remaining variables with a single dimension
    218 
    219 %remove redondant variables -> keep only one variable per dimension
    220 check_keep=true(size(ListDimName));
    221 for idim=1:numel(ListDimName)
    222     prev_ind=find(strcmp(ListDimName{idim},ListDimName(1:idim-1)));% check whether the dimension is already taken into account
    223     if ~isempty(prev_ind)% in case of multiple coord variable
    224         if strcmp(ListCoordName{idim},ListDimName{idim}) %variable with the same name as the coordinate taken in priority
    225             check_keep(prev_ind)=0;% choose a variable with the same name as coordinate in priority
    226         else
    227            check_keep(idim)=0; %keep the first coordiante variable found
    228         end
    229     end
    230 end
    231 ListCoordIndex=ListCoordIndex(check_keep);% list of coordinate variable indices
    232 ListCoordName=ListCoordName(check_keep);% list of coordinate variable names
    233 ListDimName=ListDimName(check_keep);% list of coordinate dimension names
    234 
    235 % determine dimension sizes
    236 CoordSize=zeros(size(ListCoordIndex));
    237 for ilist=1:numel(ListCoordIndex)
    238     if iscell(ListDimName{ilist})
    239         ListDimName(ilist)=ListDimName{ilist};%transform cell to string
    240     end
    241     if check_var% if the list of dimensions has been directly defined, no variable data available
    242         CoordSize(ilist)=numel(Data.(ListCoordName{ilist}));% number of elements in the variable corresponding to the dimension #ilist
    243     else
    244         check_index= strcmp(ListDimName{ilist},Data.ListDimName);% find the  index in the list of dimensions
    245         CoordSize(ilist)=Data.DimValue(check_index);% find the  corresponding dimension value
    246     end
    247     if CoordSize(ilist)==2% case of uniform grid coordinate defined by lower and upper bounds only
    248         ListDimName{ilist}=ListCoordName{ilist};% replace the dimension name by the coordinate variable name
    249     end
    250 end
    251 
    252 %% group the remaining variables in cells sharing the same coordinate variables
    253 NewCellInfo={};
    254 NewCellDimIndex={};
    255 NewNbDim=[];
    256 for ivardim=1:numel(VarDimName) % loop at the list of dimensions for the remaining variables
    257     DimCell=VarDimName{ivardim};% dimension names of the current variable
    258     if ischar(DimCell), DimCell={DimCell}; end %transform char to cell if needed
    259     DimIndices=[];
    260     for idim=1:numel(DimCell)
    261         ind_dim=find(strcmp(DimCell{idim},ListDimName));%find the dim index in the list of dimensions ListDimName
    262         if ~isempty(ind_dim)
    263             DimIndices=[DimIndices ind_dim]; %update the list of dim indices included in DimCell
    264             if check_var & CoordSize(ind_dim)==2 % determine the size of the coordinate in case of coordinate definition limited to lower and upper bounds
    265                 if isvector(Data.(ListVarName{ivardim}))
    266                     if numel(Data.(ListVarName{ivardim}))>2
    267                         CoordSize(ind_dim)=numel(Data.(ListVarName{ivardim}));
    268                     end
    269                 else
    270                     CoordSize(ind_dim)=size(Data.(ListVarName{ivardim}),idim);
    271                 end
    272             end
    273         end
    274     end
    275     % look for cells of variables with the same coordinate variables
    276     check_previous=0;
    277     for iprev=1:numel(NewCellInfo)
    278         if isequal(DimIndices,NewCellDimIndex{iprev})
    279             check_previous=1;
    280             NewCellInfo{iprev}.VarIndex=[NewCellInfo{iprev}.VarIndex ivar_remain(ivardim)];%append the current variable index to the found field cell
    281             break
    282         end
    283     end
    284     % create a new cell if no previous one contains the coordinate variables
    285     if ~check_previous
    286         nbcell=numel(NewCellInfo)+1;
    287         NewCellDimIndex{nbcell}=DimIndices;
    288         NewCellInfo{nbcell}.VarIndex=ivar_remain(ivardim);% create a new field cell with the current variable index
    289         NewNbDim(nbcell)=numel(DimIndices);
    290         NewCellInfo{nbcell}.CoordType='grid';
    291         NewCellInfo{nbcell}.CoordSize=CoordSize(DimIndices);
    292         NewCellInfo{nbcell}.CoordIndex=ListCoordIndex(DimIndices);
    293     end
    294 end
    295 NbDim=[NbDim NewNbDim];
    296 CellInfo=[CellInfo NewCellInfo];
    297 
    298 %% suppress empty cells or cells with a single coordinate variable
    299 check_remove=false(size(CellInfo));
    300 for icell=1:numel(check_remove)
    301     if isempty(CellInfo{icell})||(numel(CellInfo{icell}.VarIndex)==1 && numel(check_coord)>=icell && check_coord(icell))
    302         check_remove(icell)=1;
    303     end
    304 end
    305 CellInfo(check_remove)=[];
    306 NbDim(check_remove)=[];
     403            case 'tps'
     404                NbDim(icell)=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),2);
     405                CellInfo{icell}.CoordSize=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),1);
     406        end
     407    end
     408else
     409    for icell=1:numel(CellInfo)
     410        CellInfo{icell}.CoordSize=size(Data.DimValue(CellInfo{icell}.DimIndex));
     411    end
     412end
     413%
     414% %% loop on the tps coordinate sets
     415%
     416%     for icell_tps=1:numel(ind_coord_tps)
     417%         check_cell=zeros(1,nbvar);% =1 for the variables selected in the current cell
     418%         check_cell(ivar_tps(icell_tps))=1;% mark the coordinate variable as selected
     419%         DimCell=Data.VarDimName{ivar_tps(icell_tps)};% dimension names for the current tps coordinate variable
     420%         icell=numel(CellInfo)+icell_tps; % new field cell index
     421%         CellInfo{icell}.CoordIndex=ivar_tps(icell_tps);% index of the  tps coordinate variable
     422%         if numel(DimCell)==3
     423%             VarDimName=Data.VarDimName(~check_select);
     424%             for ivardim=1:numel(VarDimName)
     425%                 if strcmp(VarDimName{ivardim},DimCell{3})
     426%                     CellInfo{icell}.NbCentres_tps= ivar_remain(ivardim);% nbre of sites for each tps subdomain
     427%                     check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
     428%                 elseif strcmp(VarDimName{ivardim}{1},DimCell{2}) && numel(VarDimName{ivardim})>=3 && strcmp(VarDimName{ivardim}{3},DimCell{3})
     429%                     CellInfo{icell}.SubRange_tps=ivar_remain(ivardim);% subrange definiton for tps
     430%                     check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
     431%                 elseif strcmp(VarDimName{ivardim}{1},DimCell{1}) && strcmp(VarDimName{ivardim}{2},DimCell{3})% variable
     432%                     check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
     433%                 end
     434%             end
     435%         end
     436%         CellInfo{icell}.CoordType='tps';
     437%         CellInfo{icell}.VarIndex=find(check_cell);
     438%         if check_var
     439%             NbDim(icell)=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),2);
     440%             CellInfo{icell}.CoordSize=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),1)*size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),3);
     441%         else
     442%             check_index_1= strcmp(DimCell{1},Data.ListDimName);
     443%             check_index_2= strcmp(DimCell{2},Data.ListDimName);
     444%             NbDim(icell)=Data.DimValue(check_index_2);
     445%             if numel(DimCell)>=3
     446%                 check_index_3= strcmp(DimCell{3},Data.ListDimName);
     447%                 CellInfo{icell}.CoordSize=Data.DimValue(check_index_1)*Data.DimValue(check_index_3);
     448%             end
     449%         end
     450%         check_select=check_select|check_cell;
     451%     end
     452
     453
     454%% cell for ordinary plots
     455iremove=false(1,numel(ind_coord_y));
     456for ilist=1:numel(ind_coord_y)% remove the y coordinates which have been used yet in scalar or vector fields
     457    for icell=1:numel(CellInfo)
     458        if isfield(CellInfo{icell},'YIndex')&& isequal(CellInfo{icell}.YIndex,ind_coord_y(ilist))
     459            iremove(ilist)=true;
     460            continue
     461        end
     462    end
     463end
     464ind_coord_y(iremove)=[];
     465if ~isempty(ind_coord_x)
     466    y_nbre=zeros(1,numel(ind_coord_x));
     467    for icell=1:numel(ind_coord_x)
     468        Cell1DPlot{icell}.VarType='1DPlot';
     469        Cell1DPlot{icell}.XIndex=ind_coord_x(icell);
     470        Cell1DPlot{icell}.XName=Data.ListVarName{ind_coord_x(icell)};
     471        DimCell_x=Data.VarDimName{ind_coord_x(icell)};
     472        for ivar=[ind_coord_y ind_discrete]
     473            DimCell=Data.VarDimName{ivar};
     474            if  numel(DimCell)==1 && strcmp(DimCell_x{1},DimCell{1})
     475                y_nbre(icell)=y_nbre(icell)+1;
     476                Cell1DPlot{icell}.YIndex(y_nbre(icell))=ivar;
     477                break
     478            end
     479        end
     480    end
     481    Cell1DPlot(find(y_nbre==0))=[];
     482    CellInfo=[CellInfo Cell1DPlot];
     483    NbDim=[NbDim ones(1,numel(Cell1DPlot))];
     484end
    307485
    308486%% document roles of non-coordinate variables
    309487for icell=1:numel(CellInfo)
    310     VarIndex=CellInfo{icell}.VarIndex;
    311     for ivar=VarIndex
    312         if isfield(CellInfo{icell},['VarIndex_' Role{ivar}])
    313             CellInfo{icell}.(['VarIndex_' Role{ivar}])=[CellInfo{icell}.(['VarIndex_' Role{ivar}]) ivar];
    314         else
    315             CellInfo{icell}.(['VarIndex_' Role{ivar}])= ivar;
    316         end
    317         if ~isempty(ProjModeRequest{ivar})
    318             CellInfo{icell}.ProjModeRequest=ProjModeRequest{ivar};
    319         end
    320         if ~isempty(FieldName{ivar})
    321             CellInfo{icell}.FieldName=FieldName{ivar};
    322         end
    323         if CheckSub(ivar)==1
    324             CellInfo{icell}.CheckSub=1;
    325         end
    326     end
    327 end
     488    if isfield(CellInfo{icell},'VarIndex')
     489        VarIndex=CellInfo{icell}.VarIndex;
     490        for ivar=VarIndex
     491            %         if isfield(CellInfo{icell},['VarIndex_' Role{ivar}])
     492            %             CellInfo{icell}.(['VarIndex_' Role{ivar}])=[CellInfo{icell}.(['VarIndex_' Role{ivar}]) ivar];
     493            %         else
     494            %             CellInfo{icell}.(['VarIndex_' Role{ivar}])= ivar;
     495            %         end
     496            if ~isempty(ProjModeRequest{ivar})
     497                CellInfo{icell}.ProjModeRequest=ProjModeRequest{ivar};
     498            end
     499            if ~isempty(FieldName{ivar})
     500                CellInfo{icell}.FieldName=FieldName{ivar};
     501            end
     502            if CheckSub(ivar)==1
     503                CellInfo{icell}.CheckSub=1;
     504            end
     505        end
     506    end
     507end
     508% for icell=ind_coord_tps
     509%     VarIndex=CellInfo{icell}.VarIndex;
     510%     for ivar=VarIndex
     511%         if isfield(CellInfo{icell},['VarIndex_' Role{ivar}])
     512%             CellInfo{icell}.(['VarIndex_' Role{ivar}])=[CellInfo{icell}.(['VarIndex_' Role{ivar}]) ivar];
     513%         else
     514%             CellInfo{icell}.(['VarIndex_' Role{ivar}])= ivar;
     515%         end
     516%         if ~isempty(ProjModeRequest{ivar})
     517%             CellInfo{icell}.ProjModeRequest=ProjModeRequest{ivar};
     518%         end
     519%         if ~isempty(FieldName{ivar})
     520%             CellInfo{icell}.FieldName=FieldName{ivar};
     521%         end
     522%         if CheckSub(ivar)==1
     523%             CellInfo{icell}.CheckSub=1;
     524%         end
     525%     end
     526% end
     527
     528
     529
     530%
     531% %% analyse vector fields
     532% if ~isempty(ind_vector_x) && ~isempty(ind_vector_y)
     533%     if numel(ind_vector_x)>1
     534%         errormsg='multiply defined vector x component'
     535%         return
     536%     end
     537%     DimCell_vec=Data.VarDimName{ind_vector_x};% cell of dimension names for ivar_coord_x(icell)
     538%     if ischar(DimCell),DimCell={DimCell};end % transform char to cell for a single dimension
     539%     DimCell_y=Data.VarDimName{ind_vector_y};% cell of dimension names for ivar_coord_x(icell)
     540%     if ischar(DimCell_y),DimCell_y={DimCell_y};end % transform char to cell for a single dimension
     541%     if ~isequal(DimCell,DimCell_y)
     542%         errormsg='inconsistent x and y vector components';
     543%         return
     544%     end
     545%     %look for coordinates
     546%     for ivar=ind_coord_y
     547%         DimCell=Data.VarDimName{ivar};
     548%         if ischar(DimCell),DimCell={DimCell};end % transform char to cell for a single dimension
     549%         if isequal(DimCell,DimCell_vec)
     550%             CoordType='scattered';
     551%             coordy=ivar;
     552%         else
     553%             if isempty(ind_coord_z) && strcmp(DimCell{1},DimCell_vec{1})
     554%                 CoordType='grid';
     555%                 coordy=ivar;
     556%             elseif ~isempty(ind_coord_z) && strcmp(DimCell{1},DimCell_vec{2})
     557%                 CoordType='grid';
     558%                 coordy=ivar;
     559%                 coordz=ind_coord_z;
     560%             end
     561%         end
     562%         
     563%         %% find scattered (unstructured) coordinates
     564%         ivar_coord_x=find(strcmp('coord_x',Role));%find variables with Role='coord_x'
     565%         check_select=false(1,nbvar);
     566%         check_coord=false(1,nbvar);
     567%         CellInfo=cell(1,numel(ivar_coord_x));
     568%         NbDim=zeros(1,numel(ivar_coord_x));
     569%         % loop on unstructured coordinate x -> different field cells
     570%         for icell=1:numel(ivar_coord_x)
     571%             DimCell=Data.VarDimName{ivar_coord_x(icell)};% cell of dimension names for ivar_coord_x(icell)
     572%             if ischar(DimCell),DimCell={DimCell};end % transform char to cell for a single dimension
     573%             % look for variables sharing dimension(s) with ivar_coord_x(icell)
     574%             check_cell=zeros(numel(DimCell),nbvar);
     575%             for idim=1:numel(DimCell);% for each variable with role coord_x, look at which other variables contain the same dimension
     576%                 for ivar=1:nbvar
     577%                     check_cell(idim,ivar)=max(strcmp(DimCell{idim},Data.VarDimName{ivar}));
     578%                 end
     579%             end
     580%             check_cell=sum(check_cell,1)==numel(DimCell);%logical array=1 for variables belonging to the current cell
     581%             VarIndex=find(check_cell);% list of detected variable indices
     582%             if ~(numel(VarIndex)==1 && numel(DimCell)==1)% exclude case of isolated coord_x variable (treated later)
     583%                 if ~(numel(VarIndex)==1 && numel(DimCell)>1)% a variable is associated to coordinate
     584%                     CellInfo{icell}.CoordIndex=ivar_coord_x(icell);
     585%                     % size of coordinate var
     586%                     if check_var
     587%                         CellInfo{icell}.CoordSize=numel(Data.(Data.ListVarName{ivar_coord_x(icell)}));
     588%                     else
     589%                         for idim=1:numel(DimCell)
     590%                             check_index= strcmp(DimCell{idim},Data.ListDimName);
     591%                             CellInfo{icell}.CoordSize(idim)=Data.DimValue(check_index);
     592%                         end
     593%                         CellInfo{icell}.CoordSize=prod(CellInfo{icell}.CoordSize);
     594%                     end
     595%                     %             ind_scalar=find(strcmp('scalar',Role(VarIndex)));
     596%                     %             ind_vector_x=find(strcmp('vector_x',Role(VarIndex)));
     597%                     %             ind_vector_y=find(strcmp('vector_y',Role(VarIndex)));
     598%                     ind_y=find(strcmp('coord_y',Role(VarIndex)));
     599%                     if numel([ind_scalar ind_vector_x ind_vector_y])==0
     600%                         %             if numel(VarIndex)==2||isempty(ind_y)% no variable, except possibly y
     601%                         NbDim(icell)=1;
     602%                     else
     603%                         CellInfo{icell}.CoordType='scattered';
     604%                         ind_z=find(strcmp('coord_z',Role(VarIndex)));
     605%                         if numel(VarIndex)==3||isempty(ind_z)% no z variable, except possibly as a fct z(x,y)
     606%                             CellInfo{icell}.CoordIndex=[VarIndex(ind_y) CellInfo{icell}.CoordIndex];
     607%                             NbDim(icell)=2;
     608%                         else
     609%                             CellInfo{icell}.CoordIndex=[VarIndex(ind_z) CellInfo{icell}.CoordIndex];
     610%                             NbDim(icell)=3;
     611%                         end
     612%                     end
     613%                 end
     614%                 CellInfo{icell}.VarIndex=VarIndex;
     615%                 check_select=check_select|check_cell;
     616%             end
     617%         end
     618%         
     619%         %% look for tps coordinates
     620%         ivar_remain=find(~check_select);% indices of remaining variables (not already selected)
     621%         check_coord_tps= strcmp('coord_tps',Role(~check_select));
     622%         ivar_tps=ivar_remain(check_coord_tps);% variable indices corresponding to tps coordinates
     623%         
     624%         % loop on the tps coordinate sets
     625%         for icell_tps=1:numel(ivar_tps)
     626%             check_cell=zeros(1,nbvar);% =1 for the variables selected in the current cell
     627%             check_cell(ivar_tps(icell_tps))=1;% mark the coordinate variable as selected
     628%             DimCell=Data.VarDimName{ivar_tps(icell_tps)};% dimension names for the current tps coordinate variable
     629%             icell=numel(CellInfo)+icell_tps; % new field cell index
     630%             CellInfo{icell}.CoordIndex=ivar_tps(icell_tps);% index of the  tps coordinate variable
     631%             if numel(DimCell)==3
     632%                 VarDimName=Data.VarDimName(~check_select);
     633%                 for ivardim=1:numel(VarDimName)
     634%                     if strcmp(VarDimName{ivardim},DimCell{3})
     635%                         CellInfo{icell}.NbCentres_tps= ivar_remain(ivardim);% nbre of sites for each tps subdomain
     636%                         check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
     637%                     elseif strcmp(VarDimName{ivardim}{1},DimCell{2}) && numel(VarDimName{ivardim})>=3 && strcmp(VarDimName{ivardim}{3},DimCell{3})
     638%                         CellInfo{icell}.SubRange_tps=ivar_remain(ivardim);% subrange definiton for tps
     639%                         check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
     640%                     elseif strcmp(VarDimName{ivardim}{1},DimCell{1}) && strcmp(VarDimName{ivardim}{2},DimCell{3})% variable
     641%                         check_cell(ivar_remain(ivardim))=1;% mark the variable as selected
     642%                     end
     643%                 end
     644%             end
     645%             CellInfo{icell}.CoordType='tps';
     646%             CellInfo{icell}.VarIndex=find(check_cell);
     647%             if check_var
     648%                 NbDim(icell)=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),2);
     649%                 CellInfo{icell}.CoordSize=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),1)*size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),3);
     650%             else
     651%                 check_index_1= strcmp(DimCell{1},Data.ListDimName);
     652%                 check_index_2= strcmp(DimCell{2},Data.ListDimName);
     653%                 NbDim(icell)=Data.DimValue(check_index_2);
     654%                 if numel(DimCell)>=3
     655%                     check_index_3= strcmp(DimCell{3},Data.ListDimName);
     656%                     CellInfo{icell}.CoordSize=Data.DimValue(check_index_1)*Data.DimValue(check_index_3);
     657%                 end
     658%             end
     659%             check_select=check_select|check_cell;
     660%         end
     661%         
     662%     
     663%         
     664%         % determine dimension sizes
     665%         CoordSize=zeros(size(ListCoordIndex));
     666%         for ilist=1:numel(ListCoordIndex)
     667%             if iscell(ListDimName{ilist})
     668%                 ListDimName(ilist)=ListDimName{ilist};%transform cell to string
     669%             end
     670%             if check_var% if the list of dimensions has been directly defined, no variable data available
     671%                 CoordSize(ilist)=numel(Data.(ListCoordName{ilist}));% number of elements in the variable corresponding to the dimension #ilist
     672%             else
     673%                 check_index= strcmp(ListDimName{ilist},Data.ListDimName);% find the  index in the list of dimensions
     674%                 CoordSize(ilist)=Data.DimValue(check_index);% find the  corresponding dimension value
     675%             end
     676%             if CoordSize(ilist)==2% case of uniform grid coordinate defined by lower and upper bounds only
     677%                 ListDimName{ilist}=ListCoordName{ilist};% replace the dimension name by the coordinate variable name
     678%             end
     679%         end
     680%     end
     681% end
     682%
     683%
     684%
     685% %% suppress empty cells or cells with a single coordinate variable
     686% check_remove=false(size(CellInfo));
     687% for icell=1:numel(check_remove)
     688%     if isempty(CellInfo{icell})||(numel(CellInfo{icell}.VarIndex)==1 && numel(check_coord)>=icell && check_coord(icell))
     689%         check_remove(icell)=1;
     690%     end
     691% end
     692% CellInfo(check_remove)=[];
     693% NbDim(check_remove)=[];
     694%
     695%
Note: See TracChangeset for help on using the changeset viewer.