source: trunk/src/calc_field_interp.m @ 644

Last change on this file since 644 was 644, checked in by sommeria, 7 years ago

various improvements: resize GUI uvmat, projection on lines

File size: 5.9 KB
Line 
1%'calc_field_interp': defines fields (velocity, vort, div...) from civ data and calculate them
2% for projection with linear interpolation
3%---------------------------------------------------------------------
4% [VarVal,ListVarName,VarAttribute,errormsg]=calc_field_interp(Coord,Data,FieldName,XI,YI)
5%
6% OUTPUT:
7% VarVal: array giving the values of the calculated field
8% ListVarName: corresponding list of variable names
9% VarAttribute: corresponding list of variable attributes, each term #ilist is of the form VarAttribute{ilist}.tag=value
10%
11% INPUT:
12% Coord(nbpoints,2): matrix of x,y coordinates of theinput data points
13% Data: inputfield structure
14% FieldName: string representing the field to calculate, or cell array of fields (as displayed in uvmat/FieldName)
15% XI, YI: set of x and y coordiantes where the fields need to be linearly interpolated
16
17function [VarVal,ListVarName,VarAttribute,errormsg]=calc_field_interp(Coord,Data,FieldName,XI,YI)
18
19%% initialization
20VarVal={};
21ListVarName={};
22VarAttribute={};
23errormsg='';
24InputVarList={};
25if ischar(FieldName),FieldName={FieldName};end
26check_skipped=zeros(size(FieldName));% default, =1 to mark the variables which can be calculated
27Operator=cell(size(FieldName));
28
29%% analyse the list of input fields: needed variables and requested operations
30for ilist=1:numel(FieldName)
31    Operator{ilist}='';%default empty operator (vec, norm,...)
32    r=regexp(FieldName{ilist},'(?<Operator>(^vec|^norm|^curl|^div|^strain))\((?<UName>.+),(?<VName>.+)\)$','names');% analyse the field name
33    if isempty(r) % the field name is a variable itself
34        if ~isfield(Data,FieldName{ilist})
35            check_skipped(ilist)=1; %variable not found
36        else
37            if isempty(find(strcmp(FieldName{ilist},InputVarList), 1));
38                InputVarList=[InputVarList FieldName{ilist}];% the variable is added to the list of input variables if it is not already in the list
39            end
40        end
41    else
42        if ~isfield(Data,r.UName)||~isfield(Data,r.VName)%needed input variable not found
43            check_skipped(ilist)=1;
44        elseif strcmp(r.Operator,'curl')||strcmp(r.Operator,'div')||strcmp(r.Operator,'strain')
45            Operator{ilist}=r.Operator;
46            switch r.Operator
47                case 'curl'% case of CivX data format
48                    UName{ilist}='vort';
49                    Data.vort=Data.DjUi(:,1,2)-Data.DjUi(:,2,1);
50                case 'div'
51                    UName{ilist}='div';
52                    Data.div=Data.DjUi(:,1,1)+Data.DjUi(:,2,2);
53                case 'strain'
54                    UName{ilist}='strain';
55                    Data.strain=Data.DjUi(:,1,2)+Data.DjUi(:,2,1);
56            end
57            InputVarList=[InputVarList UName{ilist}]; %the variable is added to the list if it is not already in the list
58        else
59            UName{ilist}=r.UName;
60            VName{ilist}=r.VName;
61            if isempty(find(strcmp(r.UName,InputVarList)));
62                InputVarList=[InputVarList UName{ilist}]; %the variable is added to the list if it is not already in the list
63            end
64            if isempty(find(strcmp(r.VName,InputVarList), 1));
65                InputVarList=[InputVarList VName{ilist}]; %the variable is added to the list if it is not already in the list
66            end
67            Operator{ilist}=r.Operator;
68        end
69    end
70end
71
72%% create interpolator for each variable to interpolate
73if exist('XI','var')
74    for ilist=1:numel(InputVarList)
75        F.(InputVarList{ilist})=TriScatteredInterp(Coord,Data.(InputVarList{ilist}),'linear');
76    end
77end
78
79%% perform the linear interpolation for the requested variables
80for ilist=1:numel(FieldName)
81    if ~check_skipped(ilist)
82        nbvar=numel(ListVarName);
83        switch Operator{ilist}
84            case 'vec'
85                if exist('XI','var')
86                    VarVal{nbvar+1}=F.(UName{ilist})(XI,YI);
87                    VarVal{nbvar+2}=F.(VName{ilist})(XI,YI);
88                else
89                    VarVal{nbvar+1}=Data.(UName{ilist});
90                    VarVal{nbvar+2}=Data.(VName{ilist});
91                end
92                ListVarName{nbvar+1}=UName{ilist};
93                ListVarName{nbvar+2}=VName{ilist};
94                VarAttribute{nbvar+1}.Role='vector_x';
95                VarAttribute{nbvar+2}.Role='vector_y';
96            case 'norm'
97                if exist('XI','var')
98                    U2=F.(UName{ilist})(XI,YI).*F.(UName{ilist})(XI,YI);
99                    V2=F.(VName{ilist})(XI,YI).*F.(VName{ilist})(XI,YI);
100                else
101                    U2=Data.(UName{ilist}).*Data.(UName{ilist});
102                    V2=Data.(VName{ilist}).*Data.(VName{ilist});
103                end
104                VarVal{nbvar+1}=sqrt(U2+V2);
105                ListVarName{nbvar+1}='norm';
106                VarAttribute{nbvar+1}.Role='scalar';
107            case {'curl','div','strain'}
108                if exist('XI','var')
109                    VarVal{nbvar+1}=F.(UName{ilist})(XI,YI);
110                else
111                    VarVal{nbvar+1}=Data.(UName{ilist});
112                end
113                ListVarName{nbvar+1}=UName{ilist};
114                VarAttribute{nbvar+1}.Role='scalar';
115            otherwise
116                if ~isempty(FieldName{ilist})
117                    if exist('XI','var')
118                        VarVal{nbvar+1}=F.(FieldName{ilist})(XI,YI);
119                    else
120                        VarVal{nbvar+1}= Data.(FieldName{ilist});
121                    end
122                    ListVarName{nbvar+1}=FieldName{ilist};
123                    VarAttribute{nbvar+1}.Role='scalar';
124                end
125        end
126    end
127end
128
129%% put an error flag to indicate NaN data
130if exist('XI','var')&&~isempty(VarVal)
131    nbvar=numel(ListVarName);
132    ListVarName{nbvar+1}='FF';
133    VarVal{nbvar+1}=isnan(VarVal{nbvar});
134    VarAttribute{nbvar+1}.Role='errorflag';
135end
136
137
138
139
140
Note: See TracBrowser for help on using the repository browser.