Changeset 1144 for trunk/src/series
- Timestamp:
- May 13, 2024, 9:49:09 PM (7 months ago)
- Location:
- trunk/src/series
- Files:
-
- 3 edited
-
civ2vel_3C.m (modified) (13 diffs)
-
civ_series.m (modified) (22 diffs)
-
extract_multitif.m (modified) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/src/series/civ2vel_3C.m
r1141 r1144 58 58 ParamOut.WholeIndexRange='off';% prescribes the file index ranges from min to max (options 'off'/'on', 'off' by default) 59 59 ParamOut.NbSlice='off'; %nbre of slices ('off' by default) 60 ParamOut.VelType='o ne';% menu for selecting the velocity type (options 'off'/'one'/'two', 'off' by default)60 ParamOut.VelType='off';% menu for selecting the velocity type (options 'off'/'one'/'two', 'off' by default) 61 61 ParamOut.FieldName='off';% menu for selecting the field (s) in the input file(options 'off'/'one'/'two', 'off' by default) 62 62 ParamOut.FieldTransform = 'off';%use the phys transform function without choice … … 108 108 hdisp=disp_uvmat('WAITING...','checking the file series',checkrun); 109 109 [filecell,i1_series,i2_series,j1_series,j2_series]=get_file_series(Param); 110 if ~isempty(hdisp),delete(hdisp),end ;110 if ~isempty(hdisp),delete(hdisp),end 111 111 %%%%%%%%%%%% 112 112 % The cell array filecell is the list of input file names, while … … 124 124 %% define the directory for result file (with path=RootPath{1}) 125 125 OutputDir=[Param.OutputSubDir Param.OutputDirExt];% subdirectory for output files 126 %127 % if ~isfield(Param,'InputFields')128 % Param.InputFields.FieldName='';129 % end130 126 131 127 %% calibration data and timing: read the ImaDoc files … … 175 171 CheckOverwrite=Param.CheckOverwrite; 176 172 end 173 177 174 for index=1:NbField 178 175 179 180 181 %% generating the name of the merged field 176 %% generating the name of the merged field 182 177 i1=i1_series{1}(index); 183 178 if ~isempty(i2_series{end}) … … 198 193 OutputFile=fullfile_uvmat(RootPath{1},OutputDir,RootFile{1},'.nc','_1-2',i1,i2,j1,j2); 199 194 200 %%201 202 203 195 if ~isempty(RUNHandle) && ~strcmp(get(RUNHandle,'BusyAction'),'queue') 204 196 disp('program stopped by user') … … 206 198 end 207 199 208 if (~CheckOverwrite && exist(OutputFile,'file'))209 210 211 end212 200 if (~CheckOverwrite && exist(OutputFile,'file')) 201 disp('existing output file already exists, skip to next field') 202 continue% skip iteration if the mode overwrite is desactivated and the result file already exists 203 end 204 213 205 %%%%%%%%%%%%%%%% loop on views (input lines) %%%%%%%%%%%%%%%% 214 206 Data=cell(1,NbView);%initiate the set Data … … 217 209 %get Xphys,Yphys,Zphys from 1 or 2 stereo folders. Positions are taken 218 210 %at the middle between to time step 219 clear ZItemp 220 ZItemp=zeros(size(XI,1),size(XI,2),2); 221 222 if index==1 211 clear ZItemp 212 ZItemp=zeros(size(XI,1),size(XI,2),2); 213 CheckZ=0; 214 215 if index==1 223 216 first_img=i1_series{1,1}(1,1); %id of the first image of the series 224 end 225 226 idtemp=0; 227 for indextemp=index:index+1; 228 idtemp=idtemp+1; 229 if NbView==3 % if there is only 1 stereo folder, extract directly Xphys,Yphys and Zphys 230 217 end 218 219 idtemp=0; 220 for indextemp=index:index+1 221 idtemp=idtemp+1; 222 if NbView==3 % if there is only 1 stereo folder, extract directly Xphys,Yphys and Zphys 223 224 [Data{3},tild,errormsg] = nc2struct([Param.InputTable{3,1},'/',Param.InputTable{3,2},'/',Param.InputTable{3,3},'_',int2str(first_img+indextemp-1),'.nc']); 225 226 if exist('Data{3}.Civ3_FF','var') % FF is present, remove wrong vector 227 temp=find(Data{3}.Civ3_FF==0); 228 Zphys=Data{3}.Zphys(temp); 229 Yphys=Data{3}.Yphys(temp); 230 Xphys=Data{3}.Xphys(temp); 231 else 232 Zphys=Data{3}.Zphys; 233 Yphys=Data{3}.Yphys; 234 Xphys=Data{3}.Xphys; 235 end 236 237 elseif NbView==4 % is there is 2 stereo folders, get global U and V and compute Zphys 238 239 %test if the seconde camera is the same for both folder 240 for i=3:4 241 indpt(i)=strfind(Param.InputTable{i,2},'.'); % indice of the "." is the folder name 1 242 indline(i)=strfind(Param.InputTable{i,2},'-'); % indice of the "-" is the folder name1 243 camname{i}=Param.InputTable{i,2}(indline(i)+1:indpt(i)-1);% extract the second camera name 244 end 245 246 if strcmp(camname{3},camname{4})==0 247 disp_uvmat('ERROR','The 2 stereo folders should have the same camera for the second position',checkrun) 248 return 249 end 250 251 [Data{3},tild,errormsg] = nc2struct([Param.InputTable{3,1},'/',Param.InputTable{3,2},'/',Param.InputTable{3,3},'_',int2str(first_img+indextemp-1),'.nc']); 252 253 if exist('Data{3}.Civ3_FF','var') % if FF is present, remove wrong vector 254 temp=find(Data{3}.Civ3_FF==0); 255 Xmid3=Data{3}.Xmid(temp); 256 Ymid3=Data{3}.Ymid(temp); 257 U3=Data{3}.Uphys(temp); 258 V3=Data{3}.Vphys(temp); 259 else 260 Xmid3=Data{3}.Xmid; 261 Ymid3=Data{3}.Ymid; 262 U3=Data{3}.Uphys; 263 V3=Data{3}.Vphys; 264 end 265 %temporary gridd of merging the 2 stereos datas 266 [xq,yq] = meshgrid(min(Xmid3+(U3)/2):(max(Xmid3+(U3)/2)-min(Xmid3+(U3)/2))/128:max(Xmid3+(U3)/2),min(Ymid3+(V3)/2):(max(Ymid3+(V3)/2)-min(Ymid3+(V3)/2))/128:max(Ymid3+(V3)/2)); 267 268 %1st folder : interpolate the first camera (Dalsa1) points on the second (common) camera 269 %(Dalsa 3) 270 x3Q=griddata(Xmid3+(U3)/2,Ymid3+(V3)/2,Xmid3-(U3)/2,xq,yq); 271 y3Q=griddata(Xmid3+(U3)/2,Ymid3+(V3)/2,Ymid3-(V3)/2,xq,yq); 272 273 [Data{4},tild,errormsg] = nc2struct([Param.InputTable{4,1},'/',Param.InputTable{4,2},'/',Param.InputTable{4,3},'_',int2str(first_img+indextemp-1),'.nc']); 274 if exist('Data{4}.Civ3_FF','var') % if FF is present, remove wrong vector 275 temp=find(Data{4}.Civ3_FF==0); 276 Xmid4=Data{4}.Xmid(temp); 277 Ymid4=Data{4}.Ymid(temp); 278 U4=Data{4}.Uphys(temp); 279 V4=Data{4}.Vphys(temp); 280 else 281 Xmid4=Data{4}.Xmid; 282 Ymid4=Data{4}.Ymid; 283 U4=Data{4}.Uphys; 284 V4=Data{4}.Vphys; 285 end 286 287 %2nd folder :interpolate the first camera (Dalsa2) points on the second (common) camera 288 %(Dalsa 3) 289 x4Q=griddata(Xmid4+(U4)/2,Ymid4+(V4)/2,Xmid4-(U4)/2,xq,yq); 290 y4Q=griddata(Xmid4+(U4)/2,Ymid4+(V4)/2,Ymid4-(V4)/2,xq,yq); 291 292 xmid=reshape((x4Q+x3Q)/2,length(xq(:,1)).*length(xq(1,:)),1); 293 ymid=reshape((y4Q+y3Q)/2,length(yq(:,1)).*length(yq(1,:)),1); 294 u=reshape(x4Q-x3Q,length(xq(:,1)).*length(xq(1,:)),1); 295 v=reshape(y4Q-y3Q,length(yq(:,1)).*length(yq(1,:)),1); 296 297 298 [Zphys,Xphys,Yphys,error]=shift2z(xmid, ymid, u, v,XmlData); %get Xphy,Yphy and Zphys 299 %remove NaN 300 tempNaN=isnan(Zphys);tempind=find(tempNaN==1); 301 Zphys(tempind)=[]; 302 Xphys(tempind)=[]; 303 Yphys(tempind)=[]; 304 end 231 305 232 233 [Data{3},tild,errormsg] = nc2struct([Param.InputTable{3,1},'/',Param.InputTable{3,2},'/',Param.InputTable{3,3},'_',int2str(first_img+indextemp-1),'.nc']); 234 235 if exist('Data{3}.Civ3_FF','var') % FF is present, remove wrong vector 236 temp=find(Data{3}.Civ3_FF==0); 237 Zphys=Data{3}.Zphys(temp); 238 Yphys=Data{3}.Yphys(temp); 239 Xphys=Data{3}.Xphys(temp); 240 else 241 Zphys=Data{3}.Zphys; 242 Yphys=Data{3}.Yphys; 243 Xphys=Data{3}.Xphys; 244 end 245 246 247 248 elseif NbView==4 % is there is 2 stereo folders, get global U and V and compute Zphys 249 250 251 %test if the seconde camera is the same for both folder 252 for i=3:4 253 indpt(i)=strfind(Param.InputTable{i,2},'.'); % indice of the "." is the folder name 1 254 indline(i)=strfind(Param.InputTable{i,2},'-'); % indice of the "-" is the folder name1 255 camname{i}=Param.InputTable{i,2}(indline(i)+1:indpt(i)-1);% extract the second camera name 256 end 257 258 if strcmp(camname{3},camname{4})==0 259 disp_uvmat('ERROR','The 2 stereo folders should have the same camera for the second position',checkrun) 260 return 261 end 262 263 264 265 [Data{3},tild,errormsg] = nc2struct([Param.InputTable{3,1},'/',Param.InputTable{3,2},'/',Param.InputTable{3,3},'_',int2str(first_img+indextemp-1),'.nc']); 266 267 if exist('Data{3}.Civ3_FF','var') % if FF is present, remove wrong vector 268 temp=find(Data{3}.Civ3_FF==0); 269 Xmid3=Data{3}.Xmid(temp); 270 Ymid3=Data{3}.Ymid(temp); 271 U3=Data{3}.Uphys(temp); 272 V3=Data{3}.Vphys(temp); 273 else 274 Xmid3=Data{3}.Xmid; 275 Ymid3=Data{3}.Ymid; 276 U3=Data{3}.Uphys; 277 V3=Data{3}.Vphys; 278 end 279 %temporary gridd of merging the 2 stereos datas 280 [xq,yq] = meshgrid(min(Xmid3+(U3)/2):(max(Xmid3+(U3)/2)-min(Xmid3+(U3)/2))/128:max(Xmid3+(U3)/2),min(Ymid3+(V3)/2):(max(Ymid3+(V3)/2)-min(Ymid3+(V3)/2))/128:max(Ymid3+(V3)/2)); 281 282 %1st folder : interpolate the first camera (Dalsa1) points on the second (common) camera 283 %(Dalsa 3) 284 x3Q=griddata(Xmid3+(U3)/2,Ymid3+(V3)/2,Xmid3-(U3)/2,xq,yq); 285 y3Q=griddata(Xmid3+(U3)/2,Ymid3+(V3)/2,Ymid3-(V3)/2,xq,yq); 286 287 288 289 [Data{4},tild,errormsg] = nc2struct([Param.InputTable{4,1},'/',Param.InputTable{4,2},'/',Param.InputTable{4,3},'_',int2str(first_img+indextemp-1),'.nc']); 290 if exist('Data{4}.Civ3_FF','var') % if FF is present, remove wrong vector 291 temp=find(Data{4}.Civ3_FF==0); 292 Xmid4=Data{4}.Xmid(temp); 293 Ymid4=Data{4}.Ymid(temp); 294 U4=Data{4}.Uphys(temp); 295 V4=Data{4}.Vphys(temp); 296 else 297 Xmid4=Data{4}.Xmid; 298 Ymid4=Data{4}.Ymid; 299 U4=Data{4}.Uphys; 300 V4=Data{4}.Vphys; 301 end 302 303 %2nd folder :interpolate the first camera (Dalsa2) points on the second (common) camera 304 %(Dalsa 3) 305 x4Q=griddata(Xmid4+(U4)/2,Ymid4+(V4)/2,Xmid4-(U4)/2,xq,yq); 306 y4Q=griddata(Xmid4+(U4)/2,Ymid4+(V4)/2,Ymid4-(V4)/2,xq,yq); 307 308 xmid=reshape((x4Q+x3Q)/2,length(xq(:,1)).*length(xq(1,:)),1); 309 ymid=reshape((y4Q+y3Q)/2,length(yq(:,1)).*length(yq(1,:)),1); 310 u=reshape(x4Q-x3Q,length(xq(:,1)).*length(xq(1,:)),1); 311 v=reshape(y4Q-y3Q,length(yq(:,1)).*length(yq(1,:)),1); 312 313 314 [Zphys,Xphys,Yphys,error]=shift2z(xmid, ymid, u, v,XmlData); %get Xphy,Yphy and Zphys 315 %remove NaN 316 tempNaN=isnan(Zphys);tempind=find(tempNaN==1); 317 Zphys(tempind)=[]; 318 Xphys(tempind)=[]; 319 Yphys(tempind)=[]; 320 error(tempind)=[]; 321 322 end 323 324 if NbView>2 325 ZItemp(:,:,idtemp)=griddata(Xphys,Yphys,Zphys,XI,YI); %interpolation on the choosen gridd 326 end 327 328 end 329 ZI=mean(ZItemp,3); %mean between two the two time step 330 Vtest=ZItemp(:,:,2)-ZItemp(:,:,1); 306 if NbView>2 307 ZItemp(:,:,idtemp)=griddata(Xphys,Yphys,Zphys,XI,YI); %interpolation on the choosen gridd 308 CheckZ=1; 309 end 310 end 311 ZI=mean(ZItemp,3); %mean between the two time step 331 312 332 313 [Xa,Ya]=px_XYZ(XmlData{1}.GeometryCalib,[],XI,YI,ZI);% set of image coordinates on view a 333 314 [Xb,Yb]=px_XYZ(XmlData{2}.GeometryCalib,[],XI,YI,ZI);% set of image coordinates on view b 334 315 335 316 336 317 for iview=1:2 337 318 %% reading input file(s) … … 342 323 end 343 324 % get the time defined in the current file if not already defined from the xml file 344 if isfield(Data{iview},'Time')&& isequal(Data{iview}.Time,Data{1}.Time)325 if isfield(Data{iview},'Time')&& (Data{iview}.Time-Data{1}.Time)<0.0001 345 326 Time=Data{iview}.Time; 346 327 else … … 355 336 end 356 337 end 357 %remove wrong vector 338 %remove wrong vector 358 339 if isfield(Data{1},'FF') 359 340 temp=find(Data{1}.FF==0); … … 373 354 [A]=get_coeff(XmlData{1}.GeometryCalib,Xa,Ya,XI,YI,ZI); %get coef A~ 374 355 375 %remove wrong vector 356 %remove wrong vector 376 357 if isfield(Data{2},'FF') 377 358 temp=find(Data{2}.FF==0); … … 396 377 S=ones(size(XI,1),size(XI,2),3); 397 378 D=ones(size(XI,1),size(XI,2),3,3); 398 379 399 380 S(:,:,1)=A(:,:,1,1).*Ua+A(:,:,2,1).*Va+B(:,:,1,1).*Ub+B(:,:,2,1).*Vb; 400 381 S(:,:,2)=A(:,:,1,2).*Ua+A(:,:,2,2).*Va+B(:,:,1,2).*Ub+B(:,:,2,2).*Vb; … … 416 397 W(indj,indi)=dxyz(3); 417 398 end 418 end 399 end 419 400 Error=zeros(size(XI,1),size(XI,2),4); 420 401 Error(:,:,1)=A(:,:,1,1).*U+A(:,:,1,2).*V+A(:,:,1,3).*W-Ua; … … 424 405 425 406 426 427 428 429 430 407 %% recording the merged field 431 408 if index==1% initiate the structure at first index 432 MergeData.ListGlobalAttribute={'Conventions','Time','Dt' };409 MergeData.ListGlobalAttribute={'Conventions','Time','Dt','CoordUnit'}; 433 410 MergeData.Conventions='uvmat'; 434 MergeData.Time=Time; 435 MergeData.Dt=Dt; 436 MergeData.ListVarName={'coord_x','coord_y','Z','U','V','W','Error'}; 411 if isfield (XmlData{1}.GeometryCalib,'CoordUnit') && isfield (XmlData{2}.GeometryCalib,'CoordUnit') && strcmp(XmlData{1}.GeometryCalib.CoordUnit, XmlData{2}.GeometryCalib.CoordUnit) 412 MergeData.CoordUnit=XmlData{1}.GeometryCalib.CoordUnit; 413 else 414 disp_uvmat('ERROR','inconsistent coord units in the two input velocity series',checkrun) 415 return 416 end 417 MergeData.ListVarName={'coord_x','coord_y','U','V','W','Error'}; 437 418 MergeData.VarDimName={'coord_x','coord_y',{'coord_y','coord_x'},{'coord_y','coord_x'}... 438 {'coord_y','coord_x'},{'coord_y','coord_x'},{'coord_y','coord_x'}}; 419 {'coord_y','coord_x'},{'coord_y','coord_x'}}; 420 MergeData.VarAttribute{6}.unit='pixel'; %error estimate expressed in pixel 421 if CheckZ 422 MergeData.ListVarName=[MergeData.ListVarName {'Z'}]; 423 MergeData.VarDimName=[MergeData.ListVarName {'coord_y','coord_x'}]; 424 MergeData.Z=ZI; 425 end 439 426 MergeData.coord_x=xI; 440 427 MergeData.coord_y=yI; 441 428 end 429 MergeData.Time=Time; 430 MergeData.Dt=Dt; 442 431 MergeData.U=U/Dt; 443 432 MergeData.V=V/Dt; 444 433 MergeData.W=W/Dt; 445 MergeData.Z=ZI;446 447 %mfx=(XmlData{1}.GeometryCalib.fx_fy(1)+XmlData{2}.GeometryCalib.fx_fy(1))/2;448 %mfy=(XmlData{1}.GeometryCalib.fx_fy(2)+XmlData{2}.GeometryCalib.fx_fy(2))/2;449 MergeData.Error=0. 5*sqrt(sum(Error.^2,3));434 435 436 mfx=(XmlData{1}.GeometryCalib.fx_fy(1)+XmlData{2}.GeometryCalib.fx_fy(1))/2; 437 mfy=(XmlData{1}.GeometryCalib.fx_fy(2)+XmlData{2}.GeometryCalib.fx_fy(2))/2; 438 MergeData.Error=0.25*(mfx+mfy)*sqrt(sum(Error.^2,3)); 450 439 errormsg=struct2nc(OutputFile,MergeData);%save result file 451 440 if isempty(errormsg) -
trunk/src/series/civ_series.m
r1143 r1144 249 249 Data.Program='civ_series'; 250 250 Data.CivStage=0;%default 251 check_civx=0;%default252 251 253 252 %% get timing from the ImaDoc file or input video … … 332 331 if strcmp(Param.ActionInput.ListCompareMode,'PIV') 333 332 ncfile_out=fullfile_uvmat(OutputPath,OutputDir,RootFile_A,'.nc',NomTypeNc,i1_civ2,i2_civ2,j1_civ2,j2_civ2); 334 % ncfile_out=fullfile_uvmat(RootPath_A,OutputDir,RootFile_A,'.nc',NomTypeNc,i1_civ2,i2_civ2,j1_civ2,j2_civ2);335 333 else % displacement 336 334 ncfile_out=fullfile_uvmat(OutputPath,OutputDir,RootFile_A,'.nc',NomTypeNc,i2_civ2,[],j2_civ2); … … 448 446 end 449 447 % set the list of variables 450 Data.ListVarName={'Civ1_X','Civ1_Y','Civ1_U','Civ1_V','Civ1_ F','Civ1_C'};% cell array containing the names of the fields to record448 Data.ListVarName={'Civ1_X','Civ1_Y','Civ1_U','Civ1_V','Civ1_C','Civ1_FF'};% cell array containing the names of the fields to record 451 449 Data.VarDimName={'nb_vec_1','nb_vec_1','nb_vec_1','nb_vec_1','nb_vec_1','nb_vec_1'}; 452 450 Data.VarAttribute{1}.Role='coord_x'; … … 454 452 Data.VarAttribute{3}.Role='vector_x'; 455 453 Data.VarAttribute{4}.Role='vector_y'; 456 Data.VarAttribute{5}.Role='warnflag'; 454 Data.VarAttribute{5}.Role='ancillary'; 455 Data.VarAttribute{6}.Role='errorflag'; 457 456 % case of mask 458 457 if par_civ1.CheckMask&&~isempty(par_civ1.Mask) … … 487 486 Data.ListVarName=[Data.ListVarName 'Civ1_Z']; 488 487 Data.Civ1_X=[];Data.Civ1_Y=[];Data.Civ1_Z=[]; 489 Data.Civ1_U=[];Data.Civ1_V=[];Data.Civ1_C=[]; Data.Civ1_F=[];488 Data.Civ1_U=[];Data.Civ1_V=[];Data.Civ1_C=[]; 490 489 for ivol=1:NbSlice 491 490 % caluclate velocity data (y and v in indices, reverse to y component) … … 501 500 Data.Civ1_V=[Data.Civ1_V reshape(-vtable,[],1)]; 502 501 Data.Civ1_C=[Data.Civ1_C reshape(ctable,[],1)]; 503 Data.Civ1_F =[Data.Civ1_Creshape(F,[],1)];502 Data.Civ1_FF=[Data.Civ1_FF reshape(F,[],1)]; 504 503 end 505 504 else %usual PIV … … 516 515 Data.Civ1_V=reshape(-vtable,[],1); 517 516 Data.Civ1_C=reshape(ctable,[],1); 518 Data.Civ1_F =reshape(F,[],1);517 Data.Civ1_FF=reshape(F,[],1); 519 518 time_civ1=toc(tstart_civ1); 520 519 end … … 557 556 end 558 557 Data.ListGlobalAttribute=[Data.ListGlobalAttribute Fix1_param]; 559 Data.ListVarName=[Data.ListVarName {'Civ1_FF'}]; 560 Data.VarDimName=[Data.VarDimName {'nb_vec_1'}]; 561 nbvar=length(Data.ListVarName); 562 Data.VarAttribute{nbvar}.Role='errorflag'; 563 Data.Civ1_FF=int8(detect_false(Param.ActionInput.Fix1,Data.Civ1_F,Data.Civ1_C,Data.Civ1_U,Data.Civ1_V)); 558 Data.Civ1_FF=uint8(detect_false(Param.ActionInput.Fix1,Data.Civ1_C,Data.Civ1_U,Data.Civ1_V,Data.Civ1_FF)); 564 559 Data.CivStage=2; 565 560 end … … 568 563 disp('patch1 started') 569 564 tstart_patch1=tic; 570 if check_civx 571 errormsg='Civ Matlab input needed for patch'; 572 disp_uvmat('ERROR',errormsg,checkrun) 573 return 574 end 575 565 576 566 % record the processing parameters of Patch1 as global attributes in the result nc file 577 567 list_param=fieldnames(Param.ActionInput.Patch1)'; … … 612 602 Data.Civ1_U_smooth(ind_good)=Ures;% take the interpolated (smoothed) velocity values for good vectors, keep civ1 data for the other 613 603 Data.Civ1_V_smooth(ind_good)=Vres; 614 Data.Civ1_FF(ind_good)= int8(FFres);604 Data.Civ1_FF(ind_good)=uint8(FFres); 615 605 time_patch1=toc(tstart_patch1); 616 606 disp('patch1 performed') … … 802 792 % define the Civ2 variable (if Civ2 data are not replaced from previous calculation) 803 793 if isempty(find(strcmp('Civ2_X',Data.ListVarName),1)) 804 Data.ListVarName=[Data.ListVarName {'Civ2_X','Civ2_Y','Civ2_U','Civ2_V','Civ2_ F','Civ2_C'}];% cell array containing the names of the fields to record794 Data.ListVarName=[Data.ListVarName {'Civ2_X','Civ2_Y','Civ2_U','Civ2_V','Civ2_C','Civ2_FF'}];% cell array containing the names of the fields to record 805 795 Data.VarDimName=[Data.VarDimName {'nb_vec_2','nb_vec_2','nb_vec_2','nb_vec_2','nb_vec_2','nb_vec_2'}]; 806 796 Data.VarAttribute{nbvar+1}.Role='coord_x'; … … 808 798 Data.VarAttribute{nbvar+3}.Role='vector_x'; 809 799 Data.VarAttribute{nbvar+4}.Role='vector_y'; 810 Data.VarAttribute{nbvar+5}.Role='warnflag'; 800 Data.VarAttribute{nbvar+5}.Role='ancillary'; 801 Data.VarAttribute{nbvar+6}.Role='errorflag'; 811 802 end 812 803 Data.Civ2_X=reshape(xtable,[],1); … … 815 806 Data.Civ2_V=reshape(-vtable,[],1); 816 807 Data.Civ2_C=reshape(ctable,[],1); 817 Data.Civ2_F =reshape(F,[],1);808 Data.Civ2_FF=reshape(F,[],1); 818 809 disp('civ2 performed') 819 810 time_civ2=toc(tstart_civ2); … … 847 838 Data.(Fix2_param{ilist})=Param.ActionInput.Fix2.(list_param{ilist}); 848 839 end 849 Data.ListGlobalAttribute=[Data.ListGlobalAttribute Fix2_param]; 850 if check_civx 851 if ~isfield(Data,'fix2') 852 Data.ListGlobalAttribute=[Data.ListGlobalAttribute 'fix2']; 853 Data.fix2=1; 854 Data.ListVarName=[Data.ListVarName {'vec2_FixFlag'}]; 855 Data.VarDimName=[Data.VarDimName {'nb_vectors2'}]; 856 end 857 Data.vec_FixFlag=detect_false(Param.Fix2,Data.vec2_F,Data.vec2_C,Data.vec2_U,Data.vec2_V); 858 else 859 Data.ListVarName=[Data.ListVarName {'Civ2_FF'}]; 860 Data.VarDimName=[Data.VarDimName {'nb_vec_2'}]; 861 nbvar=length(Data.ListVarName); 862 Data.VarAttribute{nbvar}.Role='errorflag'; 863 Data.Civ2_FF=double(detect_false(Param.ActionInput.Fix2,Data.Civ2_F,Data.Civ2_C,Data.Civ2_U,Data.Civ2_V)); 864 Data.CivStage=Data.CivStage+1; 865 end 840 Data.ListGlobalAttribute=[Data.ListGlobalAttribute Fix2_param]; 841 Data.Civ2_FF=double(detect_false(Param.ActionInput.Fix2,Data.Civ2_C,Data.Civ2_U,Data.Civ2_V,Data.Civ2_FF)); 842 Data.CivStage=Data.CivStage+1; 866 843 end 867 844 … … 1082 1059 if sizemask > 1/2% eliminate point if more than half of the correlation box is masked 1083 1060 F(ivec)=3; % 1084 utable(ivec)= 0;1085 vtable(ivec)= 0;1061 utable(ivec)=NaN; 1062 vtable(ivec)=NaN; 1086 1063 else 1087 1064 image1_crop=image1_crop.*~mask1_crop;% put to zero the masked pixels (mask1_crop='true'=1) … … 1103 1080 end 1104 1081 if F(ivec)==3 1105 utable(ivec)= 0;1106 vtable(ivec)= 0;1082 utable(ivec)=NaN; 1083 vtable(ivec)=NaN; 1107 1084 else 1108 1085 %mask … … 1197 1174 peaky = peaky+ (f1-f2)/(2*f1-4*f0+2*f2); 1198 1175 else 1199 F= -2; % warning flag for vector truncated by the limited search box1176 F=1; % warning flag for vector truncated by the limited search box 1200 1177 end 1201 1178 peakx=x; … … 1206 1183 peakx = peakx+ (f1-f2)/(2*f1-4*f0+2*f2); 1207 1184 else 1208 F= -2; % warning flag for vector truncated by the limited search box1185 F=1; % warning flag for vector truncated by the limited search box 1209 1186 end 1210 1187 vector=[peakx-floor(npx/2)-1 peaky-floor(npy/2)-1]; … … 1215 1192 %------------------------------------------------------------------------ 1216 1193 % vector=[0 0]; %default 1217 F= -2;1194 F=1; 1218 1195 peaky=y; 1219 1196 peakx=x; … … 1258 1235 [npy,npx]=size(result_conv); 1259 1236 if x<4 || y<4 || npx-x<4 ||npy-y <4 1260 F= -2;1237 F=1; 1261 1238 vector=[x y]; 1262 1239 else … … 1292 1269 1293 1270 1294 function FF=detect_false(Param, F,C,U,V)1295 FF= zeros(size(F));%default1296 % FF= -2, for correlation max at edge1297 % FF= -1, for too small correlation1298 % FF= 1, for velocity outside bounds1299 % FF= 2 for exclusion by difference with the smoothed field1300 FF(F==-2)=-2; 1271 function FF=detect_false(Param,C,U,V,FFIn) 1272 FF=FFIn;%default, good vectors 1273 % FF=1, for correlation max at edge, not set in this function 1274 % FF=2, for too small correlation 1275 % FF=3, for velocity outside bounds 1276 % FF=4 for exclusion by difference with the smoothed field, not set in this function 1277 1301 1278 if isfield (Param,'MinCorr') 1302 FF(C<Param.MinCorr )=-1;1279 FF(C<Param.MinCorr & FFIn==0)=2; 1303 1280 end 1304 1281 if (isfield(Param,'MinVel')&&~isempty(Param.MinVel))||(isfield (Param,'MaxVel')&&~isempty(Param.MaxVel)) … … 1306 1283 if isfield (Param,'MinVel')&&~isempty(Param.MinVel) 1307 1284 U2Min=Param.MinVel*Param.MinVel; 1308 FF(Umod<U2Min )=1;1285 FF(Umod<U2Min & FFIn==0)=3; 1309 1286 end 1310 1287 if isfield (Param,'MaxVel')&&~isempty(Param.MaxVel) 1311 1288 U2Max=Param.MaxVel*Param.MaxVel; 1312 FF(Umod>U2Max)=1; 1313 end 1314 end 1315 1316 1317 %criterium on warn flags 1318 % FlagName={'CheckFmin2','CheckF2','CheckF3','CheckF4'}; 1319 % FlagVal=[-2 2 3 4]; 1320 % for iflag=1:numel(FlagName) 1321 % if isfield(Param,FlagName{iflag}) && Param.(FlagName{iflag}) 1322 % FF=(FF==1| F==FlagVal(iflag)); 1323 % end 1324 % end 1325 % %criterium on correlation values 1326 % if isfield (Param,'MinCorr') 1327 % FF=FF==1 | C<Param.MinCorr; 1328 % end 1329 % if (isfield(Param,'MinVel')&&~isempty(Param.MinVel))||(isfield (Param,'MaxVel')&&~isempty(Param.MaxVel)) 1330 % Umod= U.*U+V.*V; 1331 % if isfield (Param,'MinVel')&&~isempty(Param.MinVel) 1332 % FF=FF==1 | Umod<(Param.MinVel*Param.MinVel); 1333 % end 1334 % if isfield (Param,'MaxVel')&&~isempty(Param.MaxVel) 1335 % FF=FF==1 | Umod>(Param.MaxVel*Param.MaxVel); 1336 % end 1337 % end 1338 1289 FF(Umod>U2Max & FFIn==0)=3; 1290 end 1291 end 1339 1292 1340 1293 %------------------------------------------------------------------------ -
trunk/src/series/extract_multitif.m
r1129 r1144 77 77 ParamOut.AllowInputSort='off';% allow alphabetic sorting of the list of input file SubDir (options 'off'/'on', 'off' by default) 78 78 ParamOut.WholeIndexRange='on';% prescribes the file index ranges from min to max (options 'off'/'on', 'off' by default) 79 ParamOut.NbSlice='o ff'; % impose calculation in a single process (no parallel processing to avoid 'holes'))79 ParamOut.NbSlice='on'; % 80 80 ParamOut.VelType='off';% menu for selecting the velocity type (options 'off'/'one'/'two', 'off' by default) 81 81 ParamOut.FieldName='off';% menu for selecting the field (s) in the input file(options 'off'/'one'/'two', 'off' by default) … … 91 91 first_j=[];% note that the function will propose to cover the whole range of indices 92 92 if isfield(Param.IndexRange,'MinIndex_j'); first_j=Param.IndexRange.MinIndex_j; end 93 last_j=[];94 if isfield(Param.IndexRange,'MaxIndex_j'); last_j=Param.IndexRange.MaxIndex_j; end93 % % last_j=[]; 94 % % if isfield(Param.IndexRange,'MaxIndex_j'); last_j=Param.IndexRange.MaxIndex_j; end 95 95 PairString=''; 96 96 if isfield(Param.IndexRange,'PairString'); PairString=Param.IndexRange.PairString; end … … 101 101 msgbox_uvmat('WARNING',['the first input file ' FirstFileName ' does not exist']) 102 102 end 103 103 ParamOut.NbSlice=Param.IndexRange.MaxIndex_i; 104 104 %% check the validity of input file types 105 105 FileInfo=get_file_info(FirstFileName); 106 106 if ~strcmp(FileInfo.FileType,'multimage') 107 msgbox_uvmat('ERROR',['invalid file type input: ' FileInfo.FileType ' not a n image'])107 msgbox_uvmat('ERROR',['invalid file type input: ' FileInfo.FileType ' not a tiff image series']) 108 108 return 109 109 end
Note: See TracChangeset
for help on using the changeset viewer.
