source: trunk/src/find_field_cells.m @ 576

Last change on this file since 576 was 576, checked in by sommeria, 11 years ago

grid improved for civ: computation done closer to the edge. set_grid improved.

File size: 16.4 KB
Line 
1%'find_field_cells': test field structure for input in proj_field and plot_field
2%    group the variables  into 'fields' with common dimensions
3%------------------------------------------------------------------------
4% function  [CellInfo,NbDim,errormsg]=find_field_cells(Data)
5%
6% OUTPUT:
7% CellInfo: cell of structures describing field cells
8%     .CoordType:  type of coordinates for each field cell = 'scattered','grid','tps';
9%     .CoordIndex: array of the indices of the variables representing the coordinates (in the order z,y,x)
10%     .CoordSize: array of the nbre of values for each  coordinate in a grid, nbre of points in the unstructured case
11%     .NbSite_tps:
12%     .SubRange_tps
13%     .VarIndex: arrays of the variable indices in the field cell
14%     .VarIndex_ancillary: indices of ancillary variables
15%              _color : color image, the last index, which is not a coordinate variable, represent the 3 color components rgb
16%              _discrete: like scalar, but set of data points without continuity, represented as dots in a usual plot, instead of continuous lines otherwise
17%              _errorflag: index of error flag
18%              _image   : B/W image, (behaves like scalar)
19%              _vector_x,_y,_z: indices of variables giving the vector components x, y, z
20%              _warnflag: index of warnflag   
21%      .FieldRequest= 'interp_lin', 'interp_tps' indicate whether lin interpolation  or derivatives (tps) is needed to calculate the requested field
22%      .FieldName = operation to be performed to finalise the field cell after projection
23%      .SubCheck=0 /1 indicate that the field must be substracted (second  entry in uvmat)
24% NbDim: array with the length of CellVarIndex, giving the space dimension of each field cell
25% errormsg: error message
26%   
27% INPUT:
28% Data: structure representing fields, output of check_field_structure
29%            .ListGlobalAttributes
30%            .ListVarName: cell array listing the names (cahr strings) of the variables
31%            .VarDimName: cell array of cells containing the set of dimension names for each variable of .ListVarName
32%            .VarAttribute: cell array of structures containing the variable attributes:
33%                     .VarAttribute{ilist}.key=value, where ilist is the variable
34%                     index, key is the name of the attribute, value its value (char string or number)
35%            .Attr1, .Attr2....
36%        case of actual data:
37%            .Var1, .Var2...
38%        case of data structure, without the  data themselves
39%            .LisDimName: list of dimension names
40%            .DimValue: list of corresponding dimension values
41% HELP:
42% to get the dimensions of arrays common to the field #icell
43%         VarIndex=CellVarIndex{icell}; % list of variable indices
44%         DimIndex=Data.VarDimIndex{VarIndex(1)} % list of dimensions for each variable in the cell #icell
45%
46%AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
47%  Copyright Joel Sommeria, 2008, LEGI / CNRS-UJF-INPG, sommeria@coriolis-legi.org.
48%AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
49%     This file is part of the toolbox UVMAT.
50%
51%     UVMAT is free software; you can redistribute it and/or modify
52%     it under the terms of the GNU General Public License as published by
53%     the Free Software Foundation; either version 2 of the License, or
54%     (at your option) any later version.
55%
56%     UVMAT is distributed in the hope that it will be useful,
57%     but WITHOUT ANY WARRANTY; without even the implied warranty of
58%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
59%     GNU General Public License (file UVMAT/COPYING.txt) for more details.for input in proj_field and plot_field
60%    group the variables  into 'fields' with common dimensions
61
62function [CellInfo,NbDim,errormsg]=find_field_cells(Data)
63
64NbDim=0;
65errormsg='';
66if ~isfield(Data,'ListVarName'), errormsg='the list of variables .ListVarName is missing';return;end
67if ~isfield(Data,'VarDimName'), errormsg='the list of dimensions .VarDimName is missing';return;end
68nbvar=numel(Data.ListVarName);%number of field variables
69if ~isequal(numel(Data.VarDimName),nbvar), errormsg='.ListVarName and .VarDimName have unequal length';return;end
70% check the existence of variable data
71check_var=1;
72for ilist=1:numel(Data.ListVarName)
73    if ~isfield(Data,Data.ListVarName{ilist})
74        check_var=0;% dimensions of data defined, data not needed for this function
75        break
76    end
77end
78if ~check_var &&  ~(isfield(Data,'ListDimName')&& isfield(Data,'DimValue')&&isequal(numel(Data.ListDimName),numel(Data.DimValue)))
79    errormsg=['missing variable or values of dimensions' Data.ListVarName{ilist}];
80    return
81end
82
83
84%% role of variables and list of requested operations
85%ListRole={'coord_x','coord_y','coord_z','vector_x','vector_y','vector_z','vector_x_tps','vector_y_tps','warnflag','errorflag',...
86%   'ancillary','image','color','discrete','scalar','coord_tps'};% rmq vector_x_tps and vector_y_tps to be replaced by vector_x and vector_y
87Role=num2cell(blanks(nbvar));%initialize a cell array of nbvar blanks
88FieldRequest=regexprep(Role,' ',''); % fieldRequest set to '' by default
89FieldName=cell(size(Role)); % fieldRequest set to {} by default
90CheckSub=zeros(size(Role));% =1 for fields to substract
91Role=regexprep(Role,' ','scalar'); % Role set to 'scalar' by default
92if isfield(Data,'VarAttribute')
93    for ivar=1:numel(Data.VarAttribute)
94        if isfield(Data.VarAttribute{ivar},'Role')
95            Role{ivar}=Data.VarAttribute{ivar}.Role;
96        end
97        if isfield(Data.VarAttribute{ivar},'FieldRequest')
98            FieldRequest{ivar}=Data.VarAttribute{ivar}.FieldRequest;
99        end
100        if isfield(Data.VarAttribute{ivar},'FieldName')
101            FieldName{ivar}=Data.VarAttribute{ivar}.FieldName;
102        end
103        if isfield(Data.VarAttribute{ivar},'CheckSub')
104            CheckSub(ivar)=Data.VarAttribute{ivar}.CheckSub;
105        end
106    end
107end
108
109%% find scattered (unstructured) coordinates
110ivar_coord_x=find(strcmp('coord_x',Role));
111% VarDimCell=cell(numel(ivar_coord_x));
112check_select=false(1,nbvar);
113check_coord=false(1,nbvar);
114CellInfo=cell(1,numel(ivar_coord_x));
115NbDim=zeros(1,numel(ivar_coord_x));
116for icell=1:numel(ivar_coord_x)
117    DimCell=Data.VarDimName{ivar_coord_x(icell)};
118    if ischar(DimCell),DimCell={DimCell};end
119    check_cell=zeros(numel(DimCell),nbvar);
120    for idim=1:numel(DimCell)
121        for ivar=1:nbvar
122            check_cell(idim,ivar)=max(strcmp(DimCell{idim},Data.VarDimName{ivar}));
123        end
124    end
125    check_cell=sum(check_cell,1)==numel(DimCell);%logical array=1 for variables belonging to the current cell
126    VarIndex=find(check_cell);
127    if ~(numel(VarIndex)==1 && numel(DimCell)==1)% exclude case of isolated coord_x variable (treated later)
128        if ~(numel(VarIndex)==1 && numel(DimCell)>1)% a variable is associated to coordinate
129            CellInfo{icell}.CoordIndex=ivar_coord_x(icell);
130            % size of coordinate var
131            if check_var
132                CellInfo{icell}.CoordSize=numel(Data.(Data.ListVarName{ivar_coord_x(icell)}));
133            else
134                for idim=1:numel(DimCell)
135                 check_index= strcmp(DimCell{idim},Data.ListDimName);
136                 CellInfo{icell}.CoordSize(idim)=Data.DimValue(check_index);
137                end
138                CellInfo{icell}.CoordSize=prod(CellInfo{icell}.CoordSize);
139            end
140            ind_y=find(strcmp('coord_y',Role(VarIndex)));
141            if numel(VarIndex)==2||isempty(ind_y)% no variable, except possibly y
142                NbDim(icell)=1;
143            else
144                CellInfo{icell}.CoordType='scattered';
145                ind_z=find(strcmp('coord_z',Role(VarIndex)));
146                if numel(VarIndex)==3||isempty(ind_z)% no z variable, except possibly as a fct z(x,y)
147                    CellInfo{icell}.CoordIndex=[VarIndex(ind_y) CellInfo{icell}.CoordIndex];
148                    NbDim(icell)=2;
149                else
150                    CellInfo{icell}.CoordIndex=[VarIndex(ind_z) CellInfo{icell}.CoordIndex];
151                    NbDim(icell)=3;
152                end
153            end
154        end
155        CellInfo{icell}.VarIndex=VarIndex;
156        check_select=check_select|check_cell;
157    end
158end
159
160%% look for tps coordinates
161ivar_remain=find(~check_select);% indices of remaining variables (not already selected)
162check_coord_tps= strcmp('coord_tps',Role(~check_select));
163ivar_tps=ivar_remain(check_coord_tps);% variable indices corresponding to tps coordinates
164for icell_tps=1:numel(ivar_tps)
165    DimCell=Data.VarDimName{ivar_tps(icell_tps)};% dimension names for the current tps coordinate variable
166    icell=numel(CellInfo)+icell_tps; % new field cell index
167    CellInfo{icell}.CoordIndex=ivar_tps(icell_tps);% index of the  tps coordinate variable
168    %CellInfo{icell}.VarIndex_subrange_tps=[];
169    %CellInfo{icell}.VarIndex_nbsites_tps=[];
170    if numel(DimCell)==3
171        VarDimName=Data.VarDimName(~check_select);
172        for ivardim=1:numel(VarDimName)
173            if strcmp(VarDimName{ivardim},DimCell{3})
174                CellInfo{icell}.NbSite_tps= ivar_remain(ivardim);
175                check_cell(ivar_remain(ivardim))=1;% nbre of sites for each tps subdomain
176            elseif strcmp(VarDimName{ivardim}{1},DimCell{2}) && strcmp(VarDimName{ivardim}{3},DimCell{3})
177                CellInfo{icell}.SubRange_tps=ivar_remain(ivardim);
178                check_cell(ivar_remain(ivardim))=1;% subrange definiton for tps
179            elseif strcmp(VarDimName{ivardim}{1},DimCell{1}) && strcmp(VarDimName{ivardim}{2},DimCell{3})
180                check_cell(ivar_remain(ivardim))=1;% variable
181            end
182        end
183    end
184    CellInfo{icell}.CoordType='tps';
185    CellInfo{icell}.VarIndex=find(check_cell);
186    if check_var
187        NbDim(icell)=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),2);
188        CellInfo{icell}.CoordSize=size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),1)*size(Data.(Data.ListVarName{CellInfo{icell}.CoordIndex}),3);
189    else
190        check_index_1= strcmp(DimCell{1},Data.ListDimName);
191        check_index_2= strcmp(DimCell{2},Data.ListDimName);
192        check_index_3= strcmp(DimCell{3},Data.ListDimName);
193        NbDim(icell)=Data.DimValue(check_index_2);
194        CellInfo{icell}.CoordSize=Data.DimValue(check_index_1)*Data.DimValue(check_index_3);
195    end
196    check_select=check_select|check_cell;
197end
198
199%% look for coordinate variables and corresponding gridded data:
200% coordinate variables are variables associated with a single dimension, defining the coordinate values
201% two cases: 1)the coordiante variable represents the set of coordiante values
202%            2)the coordinate variable contains only two elements, representing the coordinate bounds for the dimension with the same name as the cordinate
203ivar_remain=find(~check_select);% indices of remaining variables, not already taken into account
204ListVarName=Data.ListVarName(~check_select);%list of remaining variables
205VarDimName=Data.VarDimName(~check_select);%dimensions of remaining variables
206check_coord_select= cellfun(@numel,VarDimName)==1|cellfun(@ischar,VarDimName)==1;% find remaining variables with a single dimension
207check_coord(~check_select)=check_coord_select;
208ListCoordIndex=ivar_remain(check_coord);% indices of remaining variables with a single dimension
209ListCoordName=Data.ListVarName(ListCoordIndex);% corresponding names of remaining variables with a single dimension
210ListDimName=Data.VarDimName(ListCoordIndex);% dimension names of remaining variables with a single dimension
211
212%remove redondant variables -> keep only one variable per dimension
213check_keep=logical(ones(size(ListDimName)));
214for idim=1:numel(ListDimName)
215    prev_ind=strcmp(ListDimName{idim},ListDimName(1:idim-1));% check whether the dimension is already taken into account
216    if ~isempty(prev_ind)
217        if strcmp(ListCoordName{idim},ListDimName{idim}) %variable with the same name as the coordinate taken in priority
218            check_keep(prev_ind)=0;
219        else
220           check_keep(idim)=0;
221        end
222    end
223end
224ListCoordIndex=ListCoordIndex(check_keep);% list of coordinate variable indices
225ListCoordName=ListCoordName(check_keep);% list of coordinate variable names
226ListDimName=ListDimName(check_keep);% list of coordinate dimension names
227
228% determine dimension sizes
229CoordSize=zeros(size(ListCoordIndex));
230for ilist=1:numel(ListCoordIndex)
231    if iscell(ListDimName{ilist})
232        ListDimName(ilist)=ListDimName{ilist};%transform cell to string
233    end
234    if check_var% if the list of dimensions has been directly defined, no variable data available
235        CoordSize(ilist)=numel(Data.(ListCoordName{ilist}));% number of elements in the variable corresponding to the dimension #ilist
236    else
237        check_index= strcmp(ListDimName{ilist},Data.ListDimName);% find the  index in the list of dimensions
238        CoordSize(ilist)=Data.DimValue(check_index);% find the  corresponding dimension value
239    end
240    if CoordSize(ilist)==2% case of uniform grid coordinate defined by lower and upper bounds only
241        ListDimName{ilist}=ListCoordName{ilist};% replace the dimension name by the coordinate variable name
242    end
243end
244
245% group the remaining variables in cells sharing the same coordinate variables
246NewCellInfo={};
247NewCellDimIndex={};
248NewNbDim=[];
249for ivardim=1:numel(VarDimName) % loop at the list of remaining variables
250    DimCell=VarDimName{ivardim};% dimension names of the current variable
251    if ischar(DimCell), DimCell={DimCell}; end %transform char to cell if needed
252    DimIndices=[];
253    for idim=1:numel(DimCell)
254        ind_dim=find(strcmp(DimCell{idim},ListDimName));%find the dim index in the list of coordinate variables
255        if ~isempty(ind_dim)
256            DimIndices=[DimIndices ind_dim]; %update the list of coord dimensions included in DimCell
257            if check_var && CoordSize(ind_dim)==2 % determine the size of the coordinate in case of coordiante variable limited to lower and upper bounds
258                if isvector(Data.(ListVarName{ivardim}))
259                    if numel(Data.(ListVarName{ivardim}))>2
260                        CoordSize(ind_dim)=numel(Data.(ListVarName{ivardim}));
261                    end
262                else
263                    CoordSize(ind_dim)=size(Data.(ListVarName{ivardim}),idim);
264                end
265            end
266        end
267    end
268    % look for cells of variables with the same coordinate variables
269    check_previous=0;
270    for iprev=1:numel(NewCellInfo)
271        if isequal(DimIndices,NewCellDimIndex{iprev})
272            check_previous=1;
273            NewCellInfo{iprev}.VarIndex=[NewCellInfo{iprev}.VarIndex ivar_remain(ivardim)];%append the current variable index to the found field cell
274            break
275        end
276    end
277    % create a new cell if no previous one contains the coordinate variables
278    if ~check_previous
279        nbcell=numel(NewCellInfo)+1;
280        NewCellDimIndex{nbcell}=DimIndices;
281        NewCellInfo{nbcell}.VarIndex=ivar_remain(ivardim);% create a new field cell with the current variable index
282        NewNbDim(nbcell)=numel(DimIndices);
283        NewCellInfo{nbcell}.CoordType='grid';
284        NewCellInfo{nbcell}.CoordSize=CoordSize;
285        NewCellInfo{nbcell}.CoordIndex=ListCoordIndex(DimIndices);
286    end
287end
288NbDim=[NbDim NewNbDim];
289CellInfo=[CellInfo NewCellInfo];
290%
291% %% suppress empty cells
292% check_empty=cellfun(@isempty,CellInfo);
293% CellInfo(check_empty)=[];
294% NbDim(check_empty)=[];
295
296%% suppress empty cells or cells with a single coordinate variable
297check_remove=false(size(CellInfo));
298for icell=1:numel(check_remove)
299    if isempty(CellInfo{icell})||(numel(CellInfo{icell}.VarIndex)==1 && check_coord(icell))
300        check_remove(icell)=1;
301    end
302end
303CellInfo(check_remove)=[];
304NbDim(check_remove)=[];
305
306%% document roles of non-coordinate variables
307% ListRole={'vector_x','vector_y','vector_z','vector_x_tps','vector_y_tps','warnflag','errorflag',...
308%    'ancillary','image','color','discrete','scalar'};% except coord,coord_x,_y,_z,Coord_tps already taken, into account
309for 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(FieldRequest{ivar})
318            CellInfo{icell}.FieldRequest=FieldRequest{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
327end
Note: See TracBrowser for help on using the repository browser.