[575] | 1 | %'calc_field': defines fields (velocity, vort, div...) from civx data (old conventions) 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);
|
---|
[575] | 123 | if ~isfield(DataIn,'CoordMesh')
|
---|
| 124 | DataIn.Coord
|
---|
| 125 | Mesh=sqrt(2*(XMax-XMin)*(YMax-YMin)/numel(DataIn.Coord_tps));
|
---|
[405] | 126 | % adjust the mesh to a value 1, 2 , 5 *10^n
|
---|
[575] | 127 | ord=10^(floor(log10(DataIn.CoordMesh)));%order of magnitude
|
---|
| 128 | if DataIn.CoordMesh/ord>=5
|
---|
| 129 | DataIn.CoordMesh=5*ord;
|
---|
| 130 | elseif DataIn.CoordMesh/ord>=2
|
---|
| 131 | DataIn.CoordMesh=2*ord;
|
---|
[405] | 132 | else
|
---|
[575] | 133 | DataIn.CoordMesh=ord;
|
---|
[405] | 134 | end
|
---|
| 135 | end
|
---|
[575] | 136 | coord_x=XMin:DataIn.CoordMesh:XMax;% increase the recommanded mesh to visualisation
|
---|
| 137 | coord_y=YMin:DataIn.CoordMesh:YMax;
|
---|
[406] | 138 | DataOut.coord_x=[XMin XMax];
|
---|
| 139 | DataOut.coord_y=[YMin YMax];
|
---|
[405] | 140 | [XI,YI]=meshgrid(coord_x,coord_y);
|
---|
[491] | 141 | Coord_interp=cat(3,XI,YI);%[XI YI];
|
---|
[405] | 142 | end
|
---|
[491] | 143 | npx=size(Coord_interp,2);
|
---|
| 144 | npy=size(Coord_interp,1);
|
---|
| 145 | Coord_interp=reshape(Coord_interp,npx*npy,size(Coord_interp,3));
|
---|
[515] | 146 |
|
---|
[406] | 147 | %initialise output
|
---|
[404] | 148 | nb_sites=size(Coord_interp,1);
|
---|
| 149 | nb_coord=size(Coord_interp,2);
|
---|
[406] | 150 | nbval=zeros(nb_sites,1);
|
---|
| 151 | NbSubDomain=size(DataIn.SubRange,3);
|
---|
| 152 | DataOut.ListVarName={'coord_y','coord_x'};
|
---|
| 153 | DataOut.VarDimName={{'coord_y'},{'coord_x'}};
|
---|
| 154 | DataOut.VarAttribute{1}=[];
|
---|
| 155 | DataOut.VarAttribute{2}=[];
|
---|
[404] | 156 | for ilist=1:length(FieldList)
|
---|
| 157 | switch FieldList{ilist}
|
---|
| 158 | case 'velocity'
|
---|
| 159 | DataOut.U=zeros(nb_sites,1);
|
---|
| 160 | DataOut.V=zeros(nb_sites,1);
|
---|
| 161 | otherwise
|
---|
| 162 | DataOut.(FieldList{ilist})=zeros(nb_sites,1);
|
---|
[389] | 163 | end
|
---|
[404] | 164 | end
|
---|
[406] | 165 |
|
---|
| 166 | % interpolate data in each subdomain
|
---|
[515] | 167 | for icell=icell_tps
|
---|
| 168 | for isub=1:NbSubDomain
|
---|
| 169 | nbvec_sub=DataIn.(DataIn.ListVarName{VarType{icell}.nbsites_tps});
|
---|
| 170 | SubRange=DataIn.(DataIn.ListVarName{VarType{icell}.subrange_tps});
|
---|
| 171 | check_range=(Coord_interp >=ones(nb_sites,1)*SubRange(:,1,isub)' & Coord_interp<=ones(nb_sites,1)*SubRange(:,2,isub)');
|
---|
| 172 | ind_sel=find(sum(check_range,2)==nb_coord);
|
---|
| 173 | nbval(ind_sel)=nbval(ind_sel)+1;% records the number of values for eacn interpolation point (in case of subdomain overlap)
|
---|
| 174 | Coord_tps=DataIn.(DataIn.ListVarName{VarType{icell}.coord_tps});
|
---|
| 175 | if check_grid
|
---|
| 176 | EM = tps_eval(Coord_interp(ind_sel,:),Coord_tps(1:nbvec_sub,:,isub));%kernels for calculating the velocity from tps 'sources'
|
---|
| 177 | end
|
---|
| 178 | if check_der
|
---|
| 179 | [EMDX,EMDY] = tps_eval_dxy(Coord_interp(ind_sel,:),Coord_tps(1:nbvec_sub,:,isub));%kernels for calculating the spatial derivatives from tps 'sources'
|
---|
| 180 | end
|
---|
| 181 | ListVar={};
|
---|
| 182 | U_tps=DataIn.(DataIn.ListVarName{VarType{icell}.vector_x_tps});
|
---|
| 183 | V_tps=DataIn.(DataIn.ListVarName{VarType{icell}.vector_y_tps});
|
---|
| 184 | for ilist=1:length(FieldList)
|
---|
| 185 | var_count=numel(ListVar);
|
---|
| 186 | switch FieldList{ilist}
|
---|
| 187 | case 'velocity'
|
---|
| 188 | ListVar=[ListVar {'U', 'V'}];
|
---|
| 189 | VarAttribute{var_count+1}.Role='vector_x';
|
---|
| 190 | VarAttribute{var_count+2}.Role='vector_y';
|
---|
| 191 | DataOut.U(ind_sel)=DataOut.U(ind_sel)+EM *U_tps(1:nbvec_sub+3,isub);
|
---|
| 192 | DataOut.V(ind_sel)=DataOut.V(ind_sel)+EM *V_tps(1:nbvec_sub+3,isub);
|
---|
| 193 | case 'u'
|
---|
| 194 | ListVar=[ListVar {'u'}];
|
---|
| 195 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 196 | DataOut.u(ind_sel)=DataOut.u(ind_sel)+EM *U_tps(1:nbvec_sub+3,isub);
|
---|
| 197 | case 'v'
|
---|
| 198 | ListVar=[ListVar {'v'}];
|
---|
| 199 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 200 | DataOut.v(ind_sel)=DataOut.v(ind_sel)+EM *V_tps(1:nbvec_sub+3,isub);
|
---|
| 201 | case 'norm_vel'
|
---|
| 202 | ListVar=[ListVar {'norm_vel'}];
|
---|
| 203 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 204 | U=DataOut.U(ind_sel)+EM *U_tps(1:nbvec_sub+3,isub);
|
---|
| 205 | V=DataOut.V(ind_sel)+EM *V_tps(1:nbvec_sub+3,isub);
|
---|
| 206 | DataOut.norm_vel(ind_sel)=sqrt(U.*U+V.*V);
|
---|
| 207 | case 'vort'
|
---|
| 208 | ListVar=[ListVar {'vort'}];
|
---|
| 209 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 210 | 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);
|
---|
| 211 | case 'div'
|
---|
| 212 | ListVar=[ListVar {'div'}];
|
---|
| 213 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 214 | 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);
|
---|
| 215 | case 'strain'
|
---|
| 216 | ListVar=[ListVar {'strain'}];
|
---|
| 217 | VarAttribute{var_count+1}.Role='scalar';
|
---|
| 218 | 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);
|
---|
| 219 | end
|
---|
| 220 | end
|
---|
[246] | 221 | end
|
---|
[515] | 222 | DataOut.FF=nbval==0; %put errorflag to 1 for points outside the interpolation rang
|
---|
| 223 | nbval(nbval==0)=1;% to avoid division by zero for averaging
|
---|
| 224 | if isempty(find(strcmp('FF',DataOut.ListVarName),1))% if FF is not already listed
|
---|
| 225 | DataOut.ListVarName=[DataOut.ListVarName {'FF'}];
|
---|
| 226 | DataOut.VarDimName=[DataOut.VarDimName {{'coord_y','coord_x'}}];
|
---|
| 227 | DataOut.VarAttribute{length(DataOut.ListVarName)}.Role='errorflag';
|
---|
[404] | 228 | end
|
---|
[515] | 229 | DataOut.ListVarName=[DataOut.ListVarName ListVar];
|
---|
| 230 | for ifield=1:numel(ListVar)
|
---|
| 231 | VarDimName{ifield}={'coord_y','coord_x'};
|
---|
| 232 | DataOut.(ListVar{ifield})=DataOut.(ListVar{ifield})./nbval;
|
---|
| 233 | DataOut.(ListVar{ifield})=reshape(DataOut.(ListVar{ifield}),npy,npx);
|
---|
[246] | 234 | end
|
---|
[515] | 235 | DataOut.FF=reshape(DataOut.FF,npy,npx);
|
---|
| 236 | DataOut.VarDimName=[DataOut.VarDimName VarDimName];
|
---|
| 237 | DataOut.VarAttribute=[DataOut.VarAttribute VarAttribute];
|
---|
[404] | 238 | end
|
---|
| 239 | else
|
---|
| 240 |
|
---|
| 241 | %% civx data
|
---|
| 242 | DataOut=DataIn;
|
---|
| 243 | for ilist=1:length(FieldList)
|
---|
| 244 | if ~isempty(FieldList{ilist})
|
---|
| 245 | [VarName,Value,Role,units]=feval(FieldList{ilist},DataIn);%calculate field with appropriate function named FieldList{ilist}
|
---|
| 246 | ListVarName=[ListVarName VarName];
|
---|
| 247 | ValueList=[ValueList Value];
|
---|
| 248 | RoleList=[RoleList Role];
|
---|
| 249 | units_cell=[units_cell units];
|
---|
[246] | 250 | end
|
---|
[124] | 251 | end
|
---|
[404] | 252 | %erase previous data (except coordinates)
|
---|
| 253 | for ivar=nbcoord+1:length(DataOut.ListVarName)
|
---|
| 254 | VarName=DataOut.ListVarName{ivar};
|
---|
| 255 | DataOut=rmfield(DataOut,VarName);
|
---|
| 256 | end
|
---|
| 257 | DataOut.ListVarName=DataOut.ListVarName(1:nbcoord);
|
---|
| 258 | if isfield(DataOut,'VarDimName')
|
---|
| 259 | DataOut.VarDimName=DataOut.VarDimName(1:nbcoord);
|
---|
| 260 | else
|
---|
| 261 | errormsg='element .VarDimName missing in input data';
|
---|
| 262 | return
|
---|
| 263 | end
|
---|
| 264 | DataOut.VarAttribute=DataOut.VarAttribute(1:nbcoord);
|
---|
| 265 | %append new data
|
---|
| 266 | DataOut.ListVarName=[DataOut.ListVarName ListVarName];
|
---|
| 267 | for ivar=1:length(ListVarName)
|
---|
| 268 | DataOut.VarDimName{nbcoord+ivar}=DataOut.VarDimName{1};
|
---|
| 269 | DataOut.VarAttribute{nbcoord+ivar}.Role=RoleList{ivar};
|
---|
| 270 | DataOut.VarAttribute{nbcoord+ivar}.units=units_cell{ivar};
|
---|
| 271 | DataOut.(ListVarName{ivar})=ValueList{ivar};
|
---|
| 272 | end
|
---|
[8] | 273 | end
|
---|
| 274 |
|
---|
[246] | 275 |
|
---|
[404] | 276 |
|
---|
[8] | 277 | %%%%%%%%%%%%% velocity fieldn%%%%%%%%%%%%%%%%%%%%
|
---|
| 278 | function [VarName,ValCell,Role,units_cell]=velocity(DataIn)
|
---|
| 279 | VarName={};
|
---|
| 280 | ValCell={};
|
---|
| 281 | Role={};
|
---|
| 282 | units_cell={};
|
---|
| 283 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 284 | units=[DataIn.CoordUnit '/' DataIn.TimeUnit];
|
---|
| 285 | else
|
---|
| 286 | units='pixel';
|
---|
| 287 | end
|
---|
| 288 | if isfield(DataIn,'U')
|
---|
| 289 | VarName=[VarName {'U'}];
|
---|
| 290 | ValCell=[ValCell {DataIn.U}];
|
---|
| 291 | Role=[Role {'vector_x'}];
|
---|
| 292 | units_cell=[units_cell {units}];
|
---|
| 293 | end
|
---|
| 294 | if isfield(DataIn,'V')
|
---|
| 295 | VarName=[VarName {'V'}];
|
---|
| 296 | ValCell=[ValCell {DataIn.V}];
|
---|
[124] | 297 | Role=[Role {'vector_y'}];
|
---|
[8] | 298 | units_cell=[units_cell {units}];
|
---|
| 299 | end
|
---|
| 300 | if isfield(DataIn,'W')
|
---|
| 301 | VarName=[VarName {'W'}];
|
---|
| 302 | ValCell=[ValCell {DataIn.W}];
|
---|
| 303 | Role=[Role {'vector_z'}];
|
---|
| 304 | units_cell=[units_cell {units}];
|
---|
| 305 | end
|
---|
| 306 | if isfield(DataIn,'F')
|
---|
| 307 | VarName=[VarName {'F'}];
|
---|
| 308 | ValCell=[ValCell {DataIn.F}];
|
---|
| 309 | Role=[Role {'warnflag'}];
|
---|
| 310 | units_cell=[units_cell {[]}];
|
---|
| 311 | end
|
---|
| 312 | if isfield(DataIn,'FF')
|
---|
| 313 | VarName=[VarName,{'FF'}];
|
---|
| 314 | ValCell=[ValCell {DataIn.FF}];
|
---|
| 315 | Role=[Role {'errorflag'}];
|
---|
| 316 | units_cell=[units_cell {[]}];
|
---|
| 317 | end
|
---|
| 318 |
|
---|
| 319 | %%%%%%%%%%%%% ima cor%%%%%%%%%%%%%%%%%%%%
|
---|
| 320 | function [VarName,ValCell,Role,units]=ima_cor(DataIn)
|
---|
| 321 | VarName={};
|
---|
| 322 | ValCell={};
|
---|
| 323 | Role={};
|
---|
| 324 | units={};
|
---|
| 325 | if isfield(DataIn,'C')
|
---|
| 326 | VarName{1}='C';
|
---|
| 327 | ValCell{1}=DataIn.C;
|
---|
| 328 | Role={'ancillary'};
|
---|
| 329 | units={[]};
|
---|
| 330 | end
|
---|
| 331 |
|
---|
| 332 | %%%%%%%%%%%%% norm_vec %%%%%%%%%%%%%%%%%%%%
|
---|
| 333 | function [VarName,ValCell,Role,units]=norm_vel(DataIn)
|
---|
| 334 | VarName={};
|
---|
| 335 | ValCell={};
|
---|
| 336 | Role={};
|
---|
| 337 | units={};
|
---|
| 338 | if isfield(DataIn,'U') && isfield(DataIn,'V')
|
---|
| 339 | VarName{1}='norm_vel';
|
---|
[124] | 340 | ValCell{1}=DataIn.U.*DataIn.U+ DataIn.V.*DataIn.V;
|
---|
| 341 | if isfield(DataIn,'W') && isequal(size(DataIn.W),size(DataIn.U))
|
---|
| 342 | ValCell{1}=ValCell{1}+DataIn.W.*DataIn.W;
|
---|
| 343 | end
|
---|
| 344 | ValCell{1}=sqrt(ValCell{1});
|
---|
| 345 | Role{1}='scalar';
|
---|
| 346 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
[8] | 347 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
[124] | 348 | else
|
---|
[8] | 349 | units={'pixel'};
|
---|
[124] | 350 | end
|
---|
| 351 | end
|
---|
[8] | 352 |
|
---|
| 353 |
|
---|
| 354 |
|
---|
| 355 | %%%%%%%%%%%%% vorticity%%%%%%%%%%%%%%%%%%%%
|
---|
| 356 | function [VarName,ValCell,Role,units]=vort(DataIn)
|
---|
| 357 | VarName={};
|
---|
| 358 | ValCell={};
|
---|
| 359 | Role={};
|
---|
| 360 | units={};
|
---|
| 361 | if isfield(DataIn,'DjUi')
|
---|
| 362 | VarName{1}='vort';
|
---|
| 363 | ValCell{1}=DataIn.DjUi(:,1,2)-DataIn.DjUi(:,2,1); %vorticity
|
---|
| 364 | siz=size(ValCell{1});
|
---|
| 365 | ValCell{1}=reshape(ValCell{1},siz(1),1);
|
---|
| 366 | Role{1}='scalar';
|
---|
| 367 | if isfield(DataIn,'TimeUnit')
|
---|
| 368 | units={[DataIn.TimeUnit '-1']};
|
---|
| 369 | else
|
---|
| 370 | units={[]};
|
---|
| 371 | end
|
---|
[124] | 372 | end
|
---|
[8] | 373 |
|
---|
| 374 | %%%%%%%%%%%%% divergence%%%%%%%%%%%%%%%%%%%%
|
---|
| 375 | function [VarName,ValCell,Role,units]=div(DataIn)
|
---|
| 376 | VarName={};
|
---|
| 377 | ValCell={};
|
---|
| 378 | Role={};
|
---|
| 379 | units={};
|
---|
| 380 | if isfield(DataIn,'DjUi')
|
---|
| 381 | VarName{1}='div';
|
---|
| 382 | ValCell{1}=DataIn.DjUi(:,1,1)+DataIn.DjUi(:,2,2); %DUDX+DVDY
|
---|
| 383 | siz=size(ValCell{1});
|
---|
| 384 | ValCell{1}=reshape(ValCell{1},siz(1),1);
|
---|
| 385 | Role{1}='scalar';
|
---|
| 386 | if isfield(DataIn,'TimeUnit')
|
---|
| 387 | units={[DataIn.TimeUnit '-1']};
|
---|
| 388 | else
|
---|
| 389 | units={[]};
|
---|
| 390 | end
|
---|
[124] | 391 | end
|
---|
[8] | 392 |
|
---|
| 393 | %%%%%%%%%%%%% strain %%%%%%%%%%%%%%%%%%%%
|
---|
| 394 | function [VarName,ValCell,Role,units]=strain(DataIn)
|
---|
| 395 | VarName={};
|
---|
| 396 | ValCell={};
|
---|
| 397 | Role={};
|
---|
| 398 | units={};
|
---|
| 399 | if isfield(DataIn,'DjUi')
|
---|
[124] | 400 | VarName{1}='strain';
|
---|
| 401 | ValCell{1}=DataIn.DjUi(:,1,2)+DataIn.DjUi(:,2,1);%DVDX+DUDY
|
---|
| 402 | siz=size(ValCell{1});
|
---|
| 403 | ValCell{1}=reshape(ValCell{1},siz(1),1);
|
---|
[156] | 404 | Role{1}='scalar';
|
---|
[124] | 405 | if isfield(DataIn,'TimeUnit')
|
---|
[8] | 406 | units={[DataIn.TimeUnit '-1']};
|
---|
[124] | 407 | else
|
---|
[8] | 408 | units={[]};
|
---|
[124] | 409 | end
|
---|
| 410 | end
|
---|
[8] | 411 |
|
---|
| 412 | %%%%%%%%%%%%% u %%%%%%%%%%%%%%%%%%%%
|
---|
| 413 | function [VarName,ValCell,Role,units]=u(DataIn)
|
---|
| 414 | VarName={};
|
---|
| 415 | ValCell={};
|
---|
| 416 | Role={};
|
---|
| 417 | units={};
|
---|
| 418 | if isfield(DataIn,'U')
|
---|
| 419 | VarName{1}='U';
|
---|
| 420 | ValCell{1}=DataIn.U;
|
---|
| 421 | Role{1}='scalar';
|
---|
| 422 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 423 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 424 | else
|
---|
| 425 | units={'pixel'};
|
---|
| 426 | end
|
---|
| 427 | end
|
---|
| 428 |
|
---|
| 429 | %%%%%%%%%%%%% v %%%%%%%%%%%%%%%%%%%%
|
---|
| 430 | function [VarName,ValCell,Role,units]=v(DataIn)
|
---|
| 431 | VarName={};
|
---|
| 432 | ValCell={};
|
---|
| 433 | Role={};
|
---|
| 434 | units={};
|
---|
| 435 | if isfield(DataIn,'V')
|
---|
| 436 | VarName{1}='V';
|
---|
| 437 | ValCell{1}=DataIn.V;
|
---|
| 438 | Role{1}='scalar';
|
---|
| 439 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 440 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 441 | else
|
---|
| 442 | units={'pixel'};
|
---|
| 443 | end
|
---|
| 444 | end
|
---|
| 445 |
|
---|
| 446 | %%%%%%%%%%%%% w %%%%%%%%%%%%%%%%%%%%
|
---|
| 447 | function [VarName,ValCell,Role,units]=w(DataIn)
|
---|
| 448 | VarName={};
|
---|
| 449 | ValCell={};
|
---|
| 450 | Role={};
|
---|
| 451 | units={};
|
---|
| 452 | if isfield(DataIn,'W')
|
---|
| 453 | VarName{1}='W';
|
---|
| 454 | ValCell{1}=DataIn.W;
|
---|
| 455 | Role{1}='scalar';%will remain unchanged by projection
|
---|
| 456 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 457 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 458 | else
|
---|
| 459 | units={'pixel'};
|
---|
| 460 | end
|
---|
| 461 | end
|
---|
| 462 |
|
---|
| 463 | %%%%%%%%%%%%% w_normal %%%%%%%%%%%%%%%%%%%%
|
---|
| 464 | function [VarName,ValCell,Role,units]=w_normal(DataIn)
|
---|
| 465 | VarName={};
|
---|
| 466 | ValCell={};
|
---|
| 467 | Role={};
|
---|
| 468 | units={};
|
---|
| 469 | if isfield(DataIn,'W')
|
---|
| 470 | VarName{1}='W';
|
---|
| 471 | ValCell{1}=DataIn.W;
|
---|
[124] | 472 | Role{1}='vector_z';%will behave like a vector component by projection
|
---|
[8] | 473 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 474 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 475 | else
|
---|
| 476 | units={'pixel'};
|
---|
| 477 | end
|
---|
| 478 | end
|
---|
| 479 |
|
---|
| 480 | %%%%%%%%%%%%% error %%%%%%%%%%%%%%%%%%%%
|
---|
| 481 | function [VarName,ValCell,Role,units]=error(DataIn)
|
---|
| 482 | VarName={};
|
---|
| 483 | ValCell={};
|
---|
| 484 | Role={};
|
---|
| 485 | units={};
|
---|
| 486 | if isfield(DataIn,'E')
|
---|
| 487 | VarName{1}='E';
|
---|
| 488 | ValCell{1}=DataIn.E;
|
---|
| 489 | Role{1}='ancillary'; %TODO CHECK units in actual fields
|
---|
| 490 | if isfield(DataIn,'CoordUnit') && isfield(DataIn,'TimeUnit')
|
---|
| 491 | units={[DataIn.CoordUnit '/' DataIn.TimeUnit]};
|
---|
| 492 | else
|
---|
| 493 | units={'pixel'};
|
---|
| 494 | end
|
---|
| 495 | end
|
---|
| 496 |
|
---|