[124] | 1 | %'calc_field': defines fields (velocity, vort, div...) from civx data and calculate them
|
---|
[8] | 2 | %---------------------------------------------------------------------
|
---|
[399] | 3 | % [DataOut,errormsg]=calc_field(FieldList,DataIn,Coord_interp)
|
---|
[8] | 4 | %
|
---|
[124] | 5 | % OUTPUT:
|
---|
| 6 | % Scal: matlab vector representing the scalar values (length nbvec defined by var_read)
|
---|
[123] | 7 | % if no input, Scal=list of programmed scalar names (to put in menus)
|
---|
| 8 | % if only the field name is put as input, vec_A=type of scalar, which can be:
|
---|
| 9 | % 'discrete': related to the individual velocity vectors, not interpolated by patch
|
---|
[8] | 10 | % 'vel': scalar calculated solely from velocity components
|
---|
[124] | 11 | % 'der': needs spatial derivatives
|
---|
[8] | 12 | % 'var': the scalar name directly corresponds to a field name in the netcdf files
|
---|
| 13 | % error: error flag
|
---|
[89] | 14 | % error = 0; OK
|
---|
| 15 | % error = 1; the prescribed scalar cannot be read or calculated from available fields
|
---|
| 16 | %
|
---|
[8] | 17 | % INPUT:
|
---|
[399] | 18 | % FieldList: cell array of strings representing the name(s) of the field(s) to calculate
|
---|
[8] | 19 | % DataIn: structure representing the field, as defined in check_field_srtructure.m
|
---|
[399] | 20 | % Coord_interp(:,nb_coord) optional set of coordinates to interpolate the field (use with thin plate shell)
|
---|
[89] | 21 | %
|
---|
[8] | 22 | % FUNCTION related
|
---|
[124] | 23 | % varname_generator.m: determines the field names to read in the netcdf
|
---|
| 24 | % file, depending on the scalar
|
---|
[517] | 25 | function [FieldOptions,ColorList]=calc_field
|
---|
| 26 | %function [DataOut,errormsg]=calc_field(FieldList,DataIn,Coord_interp)
|
---|
[8] | 27 |
|
---|
| 28 | %list of defined scalars to display in menus (in addition to 'ima_cor').
|
---|
[124] | 29 | % a type is associated to each scalar:
|
---|
[8] | 30 | % 'discrete': related to the individual velocity vectors, not interpolated by patch
|
---|
[124] | 31 | % 'vel': calculated from velocity components, continuous field (interpolated with velocity)
|
---|
| 32 | % 'der': needs spatial derivatives
|
---|
[8] | 33 | % 'var': the scalar name corresponds to a field name in the netcdf files
|
---|
| 34 | % a specific variable name for civ1 and civ2 fields are also associated, if
|
---|
[123] | 35 | % the scalar is calculated from other fields, as explicited below
|
---|
| 36 |
|
---|
[404] | 37 | %% list of field options implemented
|
---|
[515] | 38 | FieldOptions={'vec(U,V)';...%image correlation corresponding to a vel vector
|
---|
| 39 | 'C';...%image correlation corresponding to a vel vector
|
---|
| 40 | 'norm(U,V)';...%norm of the velocity
|
---|
| 41 | 'curl(U,V)';...%vorticity
|
---|
| 42 | 'div(U,V)';...%divergence
|
---|
| 43 | 'strain(U,V)';...%rate of strain
|
---|
| 44 | 'U';... %u velocity component
|
---|
| 45 | 'V';... %v velocity component
|
---|
[124] | 46 | 'w';... %w velocity component
|
---|
| 47 | 'w_normal';... %w velocity component normal to the plane
|
---|
| 48 | 'error'}; %error associated to a vector (for stereo or patch)
|
---|
[517] | 49 | ColorList={'C';...%image correlation corresponding to a vel vector
|
---|
| 50 | 'norm(U,V)';...%norm of the velocity
|
---|
| 51 | 'U';... %u velocity component
|
---|
| 52 | 'V';... %v velocity component
|
---|
| 53 | }
|
---|
| 54 | % errormsg=[]; %default error message
|
---|
| 55 | % if ~exist('FieldList','var')
|
---|
| 56 | % DataOut=FieldOptions;% gives the list of possible field inputs in the absence of input
|
---|
| 57 | % return
|
---|
| 58 | % end
|
---|
| 59 | return
|
---|
[404] | 60 | %% check input
|
---|
| 61 | if ~exist('DataIn','var')
|
---|
| 62 | DataIn=[];
|
---|
| 63 | end
|
---|
| 64 | if ischar(FieldList)
|
---|
| 65 | FieldList={FieldList};%convert a string input to a cell with one string element
|
---|
| 66 | end
|
---|
| 67 | check_grid=0;
|
---|
| 68 | check_der=0;
|
---|
| 69 | check_calc=ones(size(FieldList));
|
---|
| 70 | for ilist=1:length(FieldList)
|
---|
| 71 | switch FieldList{ilist}
|
---|
[491] | 72 | case {'u','v','velocity','norm_vel','ima_cor'}
|
---|
[404] | 73 | check_grid=1;% needs a regular grid
|
---|
| 74 | case{'vort','div','strain'}% needs spatial derivatives spatial derivatives
|
---|
| 75 | check_der=1;
|
---|
[491] | 76 | % case {'velocity','norm_vel','ima_cor'};
|
---|
[404] | 77 | otherwise
|
---|
| 78 | check_calc(ilist)=0;
|
---|
| 79 | end
|
---|
| 80 | end
|
---|
| 81 | FieldList=FieldList(check_calc==1);
|
---|
| 82 | if isfield(DataIn,'Z')&& isequal(size(DataIn.Z),size(DataIn.X))
|
---|
| 83 | nbcoord=3;
|
---|
[8] | 84 | else
|
---|
[404] | 85 | nbcoord=2;
|
---|
| 86 | end
|
---|
| 87 | ListVarName={};
|
---|
| 88 | ValueList={};
|
---|
| 89 | RoleList={};
|
---|
| 90 | units_cell={};
|
---|
| 91 |
|
---|
[515] | 92 | %% reproduce global attributes
|
---|
| 93 | DataOut.ListGlobalAttribute=DataIn.ListGlobalAttribute;
|
---|
| 94 | for ilist=1:numel(DataOut.ListGlobalAttribute)
|
---|
| 95 | DataOut.(DataOut.ListGlobalAttribute{ilist})=DataIn.(DataIn.ListGlobalAttribute{ilist});
|
---|
| 96 | end
|
---|
| 97 |
|
---|
[404] | 98 | %% interpolation with new civ data
|
---|
[515] | 99 | [CellVarIndex,NbDimVec,VarTypeCell,errormsg]=find_field_cells(DataIn);
|
---|
| 100 | icell_tps=[];
|
---|
| 101 | for icell=1:numel(CellVarIndex)
|
---|
| 102 | VarType=VarTypeCell{icell};
|
---|
| 103 | if NbDimVec(icell)>=2 && ~isempty(VarType.subrange_tps)
|
---|
| 104 | icell_tps=[icell_tps icell];
|
---|
[8] | 105 | end
|
---|
[515] | 106 | end
|
---|
[406] | 107 |
|
---|
[515] | 108 | %if isfield(DataIn,'SubRange') && isfield(DataIn,'Coord_tps') && (exist('Coord_interp','var') || check_grid ||check_der)
|
---|
| 109 | if ~isempty(icell_tps)
|
---|
[405] | 110 | %create a default grid if needed
|
---|
| 111 | if ~exist('Coord_interp','var')
|
---|
[515] | 112 | for ind=1:numel(icell_tps)
|
---|
| 113 | SubRange=DataIn.(DataIn.ListVarName{VarType{icell_tps(ind)}.subrange_tps});
|
---|
| 114 | XMax(ind)=max(max(SubRange(1,:,:)));% extrema of the coordinates
|
---|
| 115 | YMax(ind)=max(max(SubRange(2,:,:)));
|
---|
| 116 | XMin(ind)=min(min(SubRange(1,:,:)));
|
---|
| 117 | YMin(ind)=min(min(SubRange(2,:,:)));
|
---|
| 118 | end
|
---|
| 119 | XMax=max(XMax);
|
---|
| 120 | YMax=max(YMax);
|
---|
| 121 | XMin=min(XMin);
|
---|
| 122 | YMin=min(YMin);
|
---|
[405] | 123 | if ~isfield(DataIn,'Mesh')
|
---|
| 124 | DataIn.Mesh=sqrt(2*(XMax-XMin)*(YMax-YMin)/numel(DataIn.Coord_tps));
|
---|
| 125 | % adjust the mesh to a value 1, 2 , 5 *10^n
|
---|
| 126 | ord=10^(floor(log10(DataIn.Mesh)));%order of magnitude
|
---|
| 127 | if DataIn.Mesh/ord>=5
|
---|
| 128 | DataIn.Mesh=5*ord;
|
---|
| 129 | elseif DataIn.Mesh/ord>=2
|
---|
| 130 | DataIn.Mesh=2*ord;
|
---|
| 131 | else
|
---|
| 132 | DataIn.Mesh=ord;
|
---|
| 133 | end
|
---|
| 134 | end
|
---|
[501] | 135 | coord_x=XMin:DataIn.Mesh:XMax;% increase the recommanded mesh to visualisation
|
---|
[405] | 136 | coord_y=YMin:DataIn.Mesh:YMax;
|
---|
[406] | 137 | DataOut.coord_x=[XMin XMax];
|
---|
| 138 | DataOut.coord_y=[YMin YMax];
|
---|
[405] | 139 | [XI,YI]=meshgrid(coord_x,coord_y);
|
---|
[491] | 140 | Coord_interp=cat(3,XI,YI);%[XI YI];
|
---|
[405] | 141 | end
|
---|
[491] | 142 | npx=size(Coord_interp,2);
|
---|
| 143 | npy=size(Coord_interp,1);
|
---|
| 144 | Coord_interp=reshape(Coord_interp,npx*npy,size(Coord_interp,3));
|
---|
[515] | 145 |
|
---|
[406] | 146 | %initialise output
|
---|
[404] | 147 | nb_sites=size(Coord_interp,1);
|
---|
| 148 | nb_coord=size(Coord_interp,2);
|
---|
[406] | 149 | nbval=zeros(nb_sites,1);
|
---|
| 150 | NbSubDomain=size(DataIn.SubRange,3);
|
---|
| 151 | DataOut.ListVarName={'coord_y','coord_x'};
|
---|
| 152 | DataOut.VarDimName={{'coord_y'},{'coord_x'}};
|
---|
| 153 | DataOut.VarAttribute{1}=[];
|
---|
| 154 | DataOut.VarAttribute{2}=[];
|
---|
[404] | 155 | for ilist=1:length(FieldList)
|
---|
| 156 | switch FieldList{ilist}
|
---|
| 157 | case 'velocity'
|
---|
| 158 | DataOut.U=zeros(nb_sites,1);
|
---|
| 159 | DataOut.V=zeros(nb_sites,1);
|
---|
| 160 | otherwise
|
---|
| 161 | DataOut.(FieldList{ilist})=zeros(nb_sites,1);
|
---|
[389] | 162 | end
|
---|
[404] | 163 | end
|
---|
[406] | 164 |
|
---|
| 165 | % interpolate data in each subdomain
|
---|
[515] | 166 | for icell=icell_tps
|
---|
| 167 | for isub=1:NbSubDomain
|
---|
| 168 | nbvec_sub=DataIn.(DataIn.ListVarName{VarType{icell}.nbsites_tps});
|
---|
| 169 | SubRange=DataIn.(DataIn.ListVarName{VarType{icell}.subrange_tps});
|
---|
| 170 | check_range=(Coord_interp >=ones(nb_sites,1)*SubRange(:,1,isub)' & Coord_interp<=ones(nb_sites,1)*SubRange(:,2,isub)');
|
---|
| 171 | ind_sel=find(sum(check_range,2)==nb_coord);
|
---|
| 172 | nbval(ind_sel)=nbval(ind_sel)+1;% records the number of values for eacn interpolation point (in case of subdomain overlap)
|
---|
| 173 | Coord_tps=DataIn.(DataIn.ListVarName{VarType{icell}.coord_tps});
|
---|
| 174 | if check_grid
|
---|
| 175 | EM = tps_eval(Coord_interp(ind_sel,:),Coord_tps(1:nbvec_sub,:,isub));%kernels for calculating the velocity from tps 'sources'
|
---|
| 176 | end
|
---|
| 177 | if check_der
|
---|
| 178 | [EMDX,EMDY] = tps_eval_dxy(Coord_interp(ind_sel,:),Coord_tps(1:nbvec_sub,:,isub));%kernels for calculating the spatial derivatives from tps 'sources'
|
---|
| 179 | end
|
---|
| 180 | ListVar={};
|
---|
| 181 | U_tps=DataIn.(DataIn.ListVarName{VarType{icell}.vector_x_tps});
|
---|
| 182 | V_tps=DataIn.(DataIn.ListVarName{VarType{icell}.vector_y_tps});
|
---|
| 183 | for ilist=1:length(FieldList)
|
---|
| 184 | var_count=numel(ListVar);
|
---|
| 185 | switch FieldList{ilist}
|
---|
| 186 | case 'velocity'
|
---|
| 187 | ListVar=[ListVar {'U', 'V'}];
|
---|
| 188 | VarAttribute{var_count+1}.Role='vector_x';
|
---|
| 189 | VarAttribute{var_count+2}.Role='vector_y';
|
---|
| 190 | DataOut.U(ind_sel)=DataOut.U(ind_sel)+EM *U_tps(1:nbvec_sub+3,isub);
|
---|
| 191 | DataOut.V(ind_sel)=DataOut.V(ind_sel)+EM *V_tps(1:nbvec_sub+3,isub);
|
---|
| 192 | case 'u'
|
---|
| 193 | ListVar=[ListVar {'u'}];
|
---|
| 194 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 195 | DataOut.u(ind_sel)=DataOut.u(ind_sel)+EM *U_tps(1:nbvec_sub+3,isub);
|
---|
| 196 | case 'v'
|
---|
| 197 | ListVar=[ListVar {'v'}];
|
---|
| 198 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 199 | DataOut.v(ind_sel)=DataOut.v(ind_sel)+EM *V_tps(1:nbvec_sub+3,isub);
|
---|
| 200 | case 'norm_vel'
|
---|
| 201 | ListVar=[ListVar {'norm_vel'}];
|
---|
| 202 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 203 | U=DataOut.U(ind_sel)+EM *U_tps(1:nbvec_sub+3,isub);
|
---|
| 204 | V=DataOut.V(ind_sel)+EM *V_tps(1:nbvec_sub+3,isub);
|
---|
| 205 | DataOut.norm_vel(ind_sel)=sqrt(U.*U+V.*V);
|
---|
| 206 | case 'vort'
|
---|
| 207 | ListVar=[ListVar {'vort'}];
|
---|
| 208 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 209 | DataOut.vort(ind_sel)=DataOut.vort(ind_sel)-EMDY *DataIn.U_tps(1:nbvec_sub+3,isub)+EMDX *DataIn.V_tps(1:nbvec_sub+3,isub);
|
---|
| 210 | case 'div'
|
---|
| 211 | ListVar=[ListVar {'div'}];
|
---|
| 212 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 213 | DataOut.div(ind_sel)=DataOut.div(ind_sel)+EMDX*DataIn.U_tps(1:nbvec_sub+3,isub)+EMDY *DataIn.V_tps(1:nbvec_sub+3,isub);
|
---|
| 214 | case 'strain'
|
---|
| 215 | ListVar=[ListVar {'strain'}];
|
---|
| 216 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 217 | DataOut.strain(ind_sel)=DataOut.strain(ind_sel)+EMDY*DataIn.U_tps(1:nbvec_sub+3,isub)+EMDX *DataIn.V_tps(1:nbvec_sub+3,isub);
|
---|
| 218 | end
|
---|
| 219 | end
|
---|
[246] | 220 | end
|
---|
[515] | 221 | DataOut.FF=nbval==0; %put errorflag to 1 for points outside the interpolation rang
|
---|
| 222 | nbval(nbval==0)=1;% to avoid division by zero for averaging
|
---|
| 223 | if isempty(find(strcmp('FF',DataOut.ListVarName),1))% if FF is not already listed
|
---|
| 224 | DataOut.ListVarName=[DataOut.ListVarName {'FF'}];
|
---|
| 225 | DataOut.VarDimName=[DataOut.VarDimName {{'coord_y','coord_x'}}];
|
---|
| 226 | DataOut.VarAttribute{length(DataOut.ListVarName)}.Role='errorflag';
|
---|
[404] | 227 | end
|
---|
[515] | 228 | DataOut.ListVarName=[DataOut.ListVarName ListVar];
|
---|
| 229 | for ifield=1:numel(ListVar)
|
---|
| 230 | VarDimName{ifield}={'coord_y','coord_x'};
|
---|
| 231 | DataOut.(ListVar{ifield})=DataOut.(ListVar{ifield})./nbval;
|
---|
| 232 | DataOut.(ListVar{ifield})=reshape(DataOut.(ListVar{ifield}),npy,npx);
|
---|
[246] | 233 | end
|
---|
[515] | 234 | DataOut.FF=reshape(DataOut.FF,npy,npx);
|
---|
| 235 | DataOut.VarDimName=[DataOut.VarDimName VarDimName];
|
---|
| 236 | DataOut.VarAttribute=[DataOut.VarAttribute VarAttribute];
|
---|
[404] | 237 | end
|
---|
| 238 | else
|
---|
| 239 |
|
---|
| 240 | %% civx data
|
---|
| 241 | DataOut=DataIn;
|
---|
| 242 | for ilist=1:length(FieldList)
|
---|
| 243 | if ~isempty(FieldList{ilist})
|
---|
| 244 | [VarName,Value,Role,units]=feval(FieldList{ilist},DataIn);%calculate field with appropriate function named FieldList{ilist}
|
---|
| 245 | ListVarName=[ListVarName VarName];
|
---|
| 246 | ValueList=[ValueList Value];
|
---|
| 247 | RoleList=[RoleList Role];
|
---|
| 248 | units_cell=[units_cell units];
|
---|
[246] | 249 | end
|
---|
[124] | 250 | end
|
---|
[404] | 251 | %erase previous data (except coordinates)
|
---|
| 252 | for ivar=nbcoord+1:length(DataOut.ListVarName)
|
---|
| 253 | VarName=DataOut.ListVarName{ivar};
|
---|
| 254 | DataOut=rmfield(DataOut,VarName);
|
---|
| 255 | end
|
---|
| 256 | DataOut.ListVarName=DataOut.ListVarName(1:nbcoord);
|
---|
| 257 | if isfield(DataOut,'VarDimName')
|
---|
| 258 | DataOut.VarDimName=DataOut.VarDimName(1:nbcoord);
|
---|
| 259 | else
|
---|
| 260 | errormsg='element .VarDimName missing in input data';
|
---|
| 261 | return
|
---|
| 262 | end
|
---|
| 263 | DataOut.VarAttribute=DataOut.VarAttribute(1:nbcoord);
|
---|
| 264 | %append new data
|
---|
| 265 | DataOut.ListVarName=[DataOut.ListVarName ListVarName];
|
---|
| 266 | for ivar=1:length(ListVarName)
|
---|
| 267 | DataOut.VarDimName{nbcoord+ivar}=DataOut.VarDimName{1};
|
---|
| 268 | DataOut.VarAttribute{nbcoord+ivar}.Role=RoleList{ivar};
|
---|
| 269 | DataOut.VarAttribute{nbcoord+ivar}.units=units_cell{ivar};
|
---|
| 270 | DataOut.(ListVarName{ivar})=ValueList{ivar};
|
---|
| 271 | end
|
---|
[8] | 272 | end
|
---|
| 273 |
|
---|
[246] | 274 |
|
---|
[404] | 275 |
|
---|
[8] | 276 | %%%%%%%%%%%%% velocity fieldn%%%%%%%%%%%%%%%%%%%%
|
---|
| 277 | function [VarName,ValCell,Role,units_cell]=velocity(DataIn)
|
---|
| 278 | VarName={};
|
---|
| 279 | ValCell={};
|
---|
| 280 | Role={};
|
---|
| 281 | units_cell={};
|
---|
| 282 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 283 | units=[DataIn.CoordUnit '/' DataIn.TimeUnit];
|
---|
| 284 | else
|
---|
| 285 | units='pixel';
|
---|
| 286 | end
|
---|
| 287 | if isfield(DataIn,'U')
|
---|
| 288 | VarName=[VarName {'U'}];
|
---|
| 289 | ValCell=[ValCell {DataIn.U}];
|
---|
| 290 | Role=[Role {'vector_x'}];
|
---|
| 291 | units_cell=[units_cell {units}];
|
---|
| 292 | end
|
---|
| 293 | if isfield(DataIn,'V')
|
---|
| 294 | VarName=[VarName {'V'}];
|
---|
| 295 | ValCell=[ValCell {DataIn.V}];
|
---|
[124] | 296 | Role=[Role {'vector_y'}];
|
---|
[8] | 297 | units_cell=[units_cell {units}];
|
---|
| 298 | end
|
---|
| 299 | if isfield(DataIn,'W')
|
---|
| 300 | VarName=[VarName {'W'}];
|
---|
| 301 | ValCell=[ValCell {DataIn.W}];
|
---|
| 302 | Role=[Role {'vector_z'}];
|
---|
| 303 | units_cell=[units_cell {units}];
|
---|
| 304 | end
|
---|
| 305 | if isfield(DataIn,'F')
|
---|
| 306 | VarName=[VarName {'F'}];
|
---|
| 307 | ValCell=[ValCell {DataIn.F}];
|
---|
| 308 | Role=[Role {'warnflag'}];
|
---|
| 309 | units_cell=[units_cell {[]}];
|
---|
| 310 | end
|
---|
| 311 | if isfield(DataIn,'FF')
|
---|
| 312 | VarName=[VarName,{'FF'}];
|
---|
| 313 | ValCell=[ValCell {DataIn.FF}];
|
---|
| 314 | Role=[Role {'errorflag'}];
|
---|
| 315 | units_cell=[units_cell {[]}];
|
---|
| 316 | end
|
---|
| 317 |
|
---|
| 318 | %%%%%%%%%%%%% ima cor%%%%%%%%%%%%%%%%%%%%
|
---|
| 319 | function [VarName,ValCell,Role,units]=ima_cor(DataIn)
|
---|
| 320 | VarName={};
|
---|
| 321 | ValCell={};
|
---|
| 322 | Role={};
|
---|
| 323 | units={};
|
---|
| 324 | if isfield(DataIn,'C')
|
---|
| 325 | VarName{1}='C';
|
---|
| 326 | ValCell{1}=DataIn.C;
|
---|
| 327 | Role={'ancillary'};
|
---|
| 328 | units={[]};
|
---|
| 329 | end
|
---|
| 330 |
|
---|
| 331 | %%%%%%%%%%%%% norm_vec %%%%%%%%%%%%%%%%%%%%
|
---|
| 332 | function [VarName,ValCell,Role,units]=norm_vel(DataIn)
|
---|
| 333 | VarName={};
|
---|
| 334 | ValCell={};
|
---|
| 335 | Role={};
|
---|
| 336 | units={};
|
---|
| 337 | if isfield(DataIn,'U') && isfield(DataIn,'V')
|
---|
| 338 | VarName{1}='norm_vel';
|
---|
[124] | 339 | ValCell{1}=DataIn.U.*DataIn.U+ DataIn.V.*DataIn.V;
|
---|
| 340 | if isfield(DataIn,'W') && isequal(size(DataIn.W),size(DataIn.U))
|
---|
| 341 | ValCell{1}=ValCell{1}+DataIn.W.*DataIn.W;
|
---|
| 342 | end
|
---|
| 343 | ValCell{1}=sqrt(ValCell{1});
|
---|
| 344 | Role{1}='scalar';
|
---|
| 345 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
[8] | 346 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
[124] | 347 | else
|
---|
[8] | 348 | units={'pixel'};
|
---|
[124] | 349 | end
|
---|
| 350 | end
|
---|
[8] | 351 |
|
---|
| 352 |
|
---|
| 353 |
|
---|
| 354 | %%%%%%%%%%%%% vorticity%%%%%%%%%%%%%%%%%%%%
|
---|
| 355 | function [VarName,ValCell,Role,units]=vort(DataIn)
|
---|
| 356 | VarName={};
|
---|
| 357 | ValCell={};
|
---|
| 358 | Role={};
|
---|
| 359 | units={};
|
---|
| 360 | if isfield(DataIn,'DjUi')
|
---|
| 361 | VarName{1}='vort';
|
---|
| 362 | ValCell{1}=DataIn.DjUi(:,1,2)-DataIn.DjUi(:,2,1); %vorticity
|
---|
| 363 | siz=size(ValCell{1});
|
---|
| 364 | ValCell{1}=reshape(ValCell{1},siz(1),1);
|
---|
| 365 | Role{1}='scalar';
|
---|
| 366 | if isfield(DataIn,'TimeUnit')
|
---|
| 367 | units={[DataIn.TimeUnit '-1']};
|
---|
| 368 | else
|
---|
| 369 | units={[]};
|
---|
| 370 | end
|
---|
[124] | 371 | end
|
---|
[8] | 372 |
|
---|
| 373 | %%%%%%%%%%%%% divergence%%%%%%%%%%%%%%%%%%%%
|
---|
| 374 | function [VarName,ValCell,Role,units]=div(DataIn)
|
---|
| 375 | VarName={};
|
---|
| 376 | ValCell={};
|
---|
| 377 | Role={};
|
---|
| 378 | units={};
|
---|
| 379 | if isfield(DataIn,'DjUi')
|
---|
| 380 | VarName{1}='div';
|
---|
| 381 | ValCell{1}=DataIn.DjUi(:,1,1)+DataIn.DjUi(:,2,2); %DUDX+DVDY
|
---|
| 382 | siz=size(ValCell{1});
|
---|
| 383 | ValCell{1}=reshape(ValCell{1},siz(1),1);
|
---|
| 384 | Role{1}='scalar';
|
---|
| 385 | if isfield(DataIn,'TimeUnit')
|
---|
| 386 | units={[DataIn.TimeUnit '-1']};
|
---|
| 387 | else
|
---|
| 388 | units={[]};
|
---|
| 389 | end
|
---|
[124] | 390 | end
|
---|
[8] | 391 |
|
---|
| 392 | %%%%%%%%%%%%% strain %%%%%%%%%%%%%%%%%%%%
|
---|
| 393 | function [VarName,ValCell,Role,units]=strain(DataIn)
|
---|
| 394 | VarName={};
|
---|
| 395 | ValCell={};
|
---|
| 396 | Role={};
|
---|
| 397 | units={};
|
---|
| 398 | if isfield(DataIn,'DjUi')
|
---|
[124] | 399 | VarName{1}='strain';
|
---|
| 400 | ValCell{1}=DataIn.DjUi(:,1,2)+DataIn.DjUi(:,2,1);%DVDX+DUDY
|
---|
| 401 | siz=size(ValCell{1});
|
---|
| 402 | ValCell{1}=reshape(ValCell{1},siz(1),1);
|
---|
[156] | 403 | Role{1}='scalar';
|
---|
[124] | 404 | if isfield(DataIn,'TimeUnit')
|
---|
[8] | 405 | units={[DataIn.TimeUnit '-1']};
|
---|
[124] | 406 | else
|
---|
[8] | 407 | units={[]};
|
---|
[124] | 408 | end
|
---|
| 409 | end
|
---|
[8] | 410 |
|
---|
| 411 | %%%%%%%%%%%%% u %%%%%%%%%%%%%%%%%%%%
|
---|
| 412 | function [VarName,ValCell,Role,units]=u(DataIn)
|
---|
| 413 | VarName={};
|
---|
| 414 | ValCell={};
|
---|
| 415 | Role={};
|
---|
| 416 | units={};
|
---|
| 417 | if isfield(DataIn,'U')
|
---|
| 418 | VarName{1}='U';
|
---|
| 419 | ValCell{1}=DataIn.U;
|
---|
| 420 | Role{1}='scalar';
|
---|
| 421 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 422 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 423 | else
|
---|
| 424 | units={'pixel'};
|
---|
| 425 | end
|
---|
| 426 | end
|
---|
| 427 |
|
---|
| 428 | %%%%%%%%%%%%% v %%%%%%%%%%%%%%%%%%%%
|
---|
| 429 | function [VarName,ValCell,Role,units]=v(DataIn)
|
---|
| 430 | VarName={};
|
---|
| 431 | ValCell={};
|
---|
| 432 | Role={};
|
---|
| 433 | units={};
|
---|
| 434 | if isfield(DataIn,'V')
|
---|
| 435 | VarName{1}='V';
|
---|
| 436 | ValCell{1}=DataIn.V;
|
---|
| 437 | Role{1}='scalar';
|
---|
| 438 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 439 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 440 | else
|
---|
| 441 | units={'pixel'};
|
---|
| 442 | end
|
---|
| 443 | end
|
---|
| 444 |
|
---|
| 445 | %%%%%%%%%%%%% w %%%%%%%%%%%%%%%%%%%%
|
---|
| 446 | function [VarName,ValCell,Role,units]=w(DataIn)
|
---|
| 447 | VarName={};
|
---|
| 448 | ValCell={};
|
---|
| 449 | Role={};
|
---|
| 450 | units={};
|
---|
| 451 | if isfield(DataIn,'W')
|
---|
| 452 | VarName{1}='W';
|
---|
| 453 | ValCell{1}=DataIn.W;
|
---|
| 454 | Role{1}='scalar';%will remain unchanged by projection
|
---|
| 455 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 456 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 457 | else
|
---|
| 458 | units={'pixel'};
|
---|
| 459 | end
|
---|
| 460 | end
|
---|
| 461 |
|
---|
| 462 | %%%%%%%%%%%%% w_normal %%%%%%%%%%%%%%%%%%%%
|
---|
| 463 | function [VarName,ValCell,Role,units]=w_normal(DataIn)
|
---|
| 464 | VarName={};
|
---|
| 465 | ValCell={};
|
---|
| 466 | Role={};
|
---|
| 467 | units={};
|
---|
| 468 | if isfield(DataIn,'W')
|
---|
| 469 | VarName{1}='W';
|
---|
| 470 | ValCell{1}=DataIn.W;
|
---|
[124] | 471 | Role{1}='vector_z';%will behave like a vector component by projection
|
---|
[8] | 472 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 473 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 474 | else
|
---|
| 475 | units={'pixel'};
|
---|
| 476 | end
|
---|
| 477 | end
|
---|
| 478 |
|
---|
| 479 | %%%%%%%%%%%%% error %%%%%%%%%%%%%%%%%%%%
|
---|
| 480 | function [VarName,ValCell,Role,units]=error(DataIn)
|
---|
| 481 | VarName={};
|
---|
| 482 | ValCell={};
|
---|
| 483 | Role={};
|
---|
| 484 | units={};
|
---|
| 485 | if isfield(DataIn,'E')
|
---|
| 486 | VarName{1}='E';
|
---|
| 487 | ValCell{1}=DataIn.E;
|
---|
| 488 | Role{1}='ancillary'; %TODO CHECK units in actual fields
|
---|
| 489 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 490 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 491 | else
|
---|
| 492 | units={'pixel'};
|
---|
| 493 | end
|
---|
| 494 | end
|
---|
| 495 |
|
---|