[654]  1  %'calc_field_interp': calculate fields (velocity, vort, div...) using linear interpolation if requested


[515]  2  %


[579]  3  % [VarVal,ListVarName,VarAttribute,errormsg]=calc_field_interp(Coord,Data,FieldName,XI,YI)


[515]  4  %


 5  % OUTPUT:


[575]  6  % VarVal: array giving the values of the calculated field


[579]  7  % ListVarName: corresponding list of variable names


 8  % VarAttribute: corresponding list of variable attributes, each term #ilist is of the form VarAttribute{ilist}.tag=value


[515]  9  %


 10  % INPUT:


[654]  11  % Coord(nbpoints,2): matrix of x,y coordinates of the input data points


[579]  12  % Data: inputfield structure


 13  % FieldName: string representing the field to calculate, or cell array of fields (as displayed in uvmat/FieldName)


[654]  14  % XI, YI: set of x and y coordinates where the fields need to be linearly interpolated,


 15  % if XI, YI are missing, there is no interpolation (case of colors in vector plots)


 16  %


[579]  17  function [VarVal,ListVarName,VarAttribute,errormsg]=calc_field_interp(Coord,Data,FieldName,XI,YI)


[515]  18 


[644]  19  %% initialization


[579]  20  VarVal={};


[515]  21  ListVarName={};


 22  VarAttribute={};


 23  errormsg='';


[517]  24  InputVarList={};


[579]  25  if ischar(FieldName),FieldName={FieldName};end


 26  check_skipped=zeros(size(FieldName));% default, =1 to mark the variables which can be calculated


[654]  27  check_interp=ones(size(FieldName));% default, =1 to mark the variables which can be interpolated (not ancillary)


[579]  28  Operator=cell(size(FieldName));


[644]  29 


 30  %% analyse the list of input fields: needed variables and requested operations


[579]  31  for ilist=1:numel(FieldName)


 32  Operator{ilist}='';%default empty operator (vec, norm,...)


 33  r=regexp(FieldName{ilist},'(?<Operator>(^vec^norm^curl^div^strain))\((?<UName>.+),(?<VName>.+)\)$','names');% analyse the field name


 34  if isempty(r) % the field name is a variable itself


[650]  35  ivar=find(strcmp(FieldName{ilist},Data.ListVarName));


[648]  36  if isempty(ivar)


[579]  37  check_skipped(ilist)=1; %variable not found


[654]  38  elseif isempty(find(strcmp(FieldName{ilist},InputVarList), 1));


[650]  39  if isfield(Data.VarAttribute{ivar},'Role') &&...


 40  (strcmp(Data.VarAttribute{ivar}.Role,'ancillary')strcmp(Data.VarAttribute{ivar}.Role,'warnflag')strcmp(Data.VarAttribute{ivar}.Role,'errorflag'))


[654]  41  check_interp(ilist)=0; % ancillary variable, not interpolated


[521]  42  end


[654]  43  InputVarList=[InputVarList FieldName{ilist}];% the variable is added to the list of input variables if it is not already in the list


[517]  44  end


 45  else


[579]  46  if ~isfield(Data,r.UName)~isfield(Data,r.VName)%needed input variable not found


[521]  47  check_skipped(ilist)=1;


[579]  48  elseif strcmp(r.Operator,'curl')strcmp(r.Operator,'div')strcmp(r.Operator,'strain')


 49  Operator{ilist}=r.Operator;


 50  switch r.Operator


 51  case 'curl'% case of CivX data format


[675]  52  if ~isfield(Data,'DjUi'), errormsg='field DjUi needed to get curl through linear interp: use ProjMode=interp_tps'; return; end


[579]  53  UName{ilist}='vort';


 54  Data.vort=Data.DjUi(:,1,2)Data.DjUi(:,2,1);


 55  case 'div'


[675]  56  if ~isfield(Data,'DjUi'), errormsg='field DjUi needed to get div through linear interp: use ProjMode=interp_tps'; return; end


[579]  57  UName{ilist}='div';


 58  Data.div=Data.DjUi(:,1,1)+Data.DjUi(:,2,2);


 59  case 'strain'


[675]  60  if ~isfield(Data,'DjUi'), errormsg='field DjUi needed to get strain through linear interp: use ProjMode=interp_tps'; return; end


[579]  61  UName{ilist}='strain';


 62  Data.strain=Data.DjUi(:,1,2)+Data.DjUi(:,2,1);


 63  end


 64  InputVarList=[InputVarList UName{ilist}]; %the variable is added to the list if it is not already in the list


[675]  65  else % case 'norm' for instance


[521]  66  UName{ilist}=r.UName;


 67  VName{ilist}=r.VName;


 68  if isempty(find(strcmp(r.UName,InputVarList)));


 69  InputVarList=[InputVarList UName{ilist}]; %the variable is added to the list if it is not already in the list


 70  end


[579]  71  if isempty(find(strcmp(r.VName,InputVarList), 1));


[521]  72  InputVarList=[InputVarList VName{ilist}]; %the variable is added to the list if it is not already in the list


 73  end


 74  Operator{ilist}=r.Operator;


[517]  75  end


[515]  76  end


 77  end


[644]  78 


 79  %% create interpolator for each variable to interpolate


[517]  80  if exist('XI','var')


 81  for ilist=1:numel(InputVarList)


 82  F.(InputVarList{ilist})=TriScatteredInterp(Coord,Data.(InputVarList{ilist}),'linear');


 83  end


[515]  84  end


[644]  85 


 86  %% perform the linear interpolation for the requested variables


[579]  87  for ilist=1:numel(FieldName)


[521]  88  if ~check_skipped(ilist)


[533]  89  nbvar=numel(ListVarName);


 90  switch Operator{ilist}


 91  case 'vec'


[517]  92  if exist('XI','var')


[667]  93  if check_interp(ilist)


[533]  94  VarVal{nbvar+1}=F.(UName{ilist})(XI,YI);


 95  VarVal{nbvar+2}=F.(VName{ilist})(XI,YI);


[654]  96  end


[517]  97  else


[533]  98  VarVal{nbvar+1}=Data.(UName{ilist});


 99  VarVal{nbvar+2}=Data.(VName{ilist});


[517]  100  end


[533]  101  ListVarName{nbvar+1}=UName{ilist};


 102  ListVarName{nbvar+2}=VName{ilist};


 103  VarAttribute{nbvar+1}.Role='vector_x';


 104  VarAttribute{nbvar+2}.Role='vector_y';


 105  case 'norm'


 106  if exist('XI','var')


[667]  107  if check_interp(ilist)


[533]  108  U2=F.(UName{ilist})(XI,YI).*F.(UName{ilist})(XI,YI);


 109  V2=F.(VName{ilist})(XI,YI).*F.(VName{ilist})(XI,YI);


[654]  110  end


[533]  111  else


 112  U2=Data.(UName{ilist}).*Data.(UName{ilist});


 113  V2=Data.(VName{ilist}).*Data.(VName{ilist});


 114  end


 115  VarVal{nbvar+1}=sqrt(U2+V2);


 116  ListVarName{nbvar+1}='norm';


[517]  117  VarAttribute{nbvar+1}.Role='scalar';


[579]  118  case {'curl','div','strain'}


 119  if exist('XI','var')


[667]  120  if check_interp(ilist)


[579]  121  VarVal{nbvar+1}=F.(UName{ilist})(XI,YI);


[654]  122  end


[579]  123  else


 124  VarVal{nbvar+1}=Data.(UName{ilist});


 125  end


 126  ListVarName{nbvar+1}=UName{ilist};


 127  VarAttribute{nbvar+1}.Role='scalar';


[533]  128  otherwise


[579]  129  if ~isempty(FieldName{ilist})


[533]  130  if exist('XI','var')


[667]  131  if check_interp(ilist)


[579]  132  VarVal{nbvar+1}=F.(FieldName{ilist})(XI,YI);


[654]  133  end


[533]  134  else


[579]  135  VarVal{nbvar+1}= Data.(FieldName{ilist});


[533]  136  end


[579]  137  ListVarName{nbvar+1}=FieldName{ilist};


[533]  138  VarAttribute{nbvar+1}.Role='scalar';


 139  end


 140  end


[515]  141  end


 142  end


[644]  143 


 144  %% put an error flag to indicate NaN data


[579]  145  if exist('XI','var')&&~isempty(VarVal)


[667]  146  nbvar=numel(VarVal);


[533]  147  ListVarName{nbvar+1}='FF';


 148  VarVal{nbvar+1}=isnan(VarVal{nbvar});


 149  VarAttribute{nbvar+1}.Role='errorflag';


[517]  150  end


[515]  151 


 152 


 153 


 154 


 155 

