Changeset 864 for trunk/src/series/civ2vel_3C.m
 Timestamp:
 Feb 6, 2015, 7:55:10 PM (9 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/src/series/civ2vel_3C.m
r863 r864 169 169 %% MAIN LOOP ON FIELDS 170 170 warning off 171 for index=1:NbField 171 if NbField<2 && length(filecell{:,1})>2 172 disp_uvmat('ERROR','you need at least 2 images to compute the mean position for the stereo.',checkrun) 173 return 174 end 175 for index=1:NbField1 172 176 update_waitbar(WaitbarHandle,index/NbField) 173 177 if ~isempty(RUNHandle) && ~strcmp(get(RUNHandle,'BusyAction'),'queue') … … 179 183 Data=cell(1,NbView);%initiate the set Data 180 184 timeread=zeros(1,NbView); 181 [Data{3},tild,errormsg] = nc2struct(filecell{3,index}); 182 ZI=griddata(Data{3}.Xphys,Data{3}.Yphys,Data{3}.Zphys,XI,YI); 185 186 %get Xphys,Yphys,Zphys from 1 or 2 stereo folders. Positions are taken 187 %at the middle between to time step 188 clear ZItemp 189 ZItemp=zeros(size(XI,1),size(XI,2),2); 190 for indextemp=index:index+1; 191 if NbView==3 % if there is only 1 stereo folder, extract directly Xphys,Yphys and Zphys 192 193 [Data{3},tild,errormsg] = nc2struct(filecell{3,indextemp}); 194 195 if exist('Data{3}.Civ3_FF','var') % FF is present, remove wrong vector 196 temp=find(Data{3}.Civ3_FF==0); 197 Zphys=Data{3}.Zphys(temp); 198 Yphys=Data{3}.Yphys(temp); 199 Xphys=Data{3}.Xphys(temp); 200 else 201 Zphys=Data{3}.Zphys; 202 Yphys=Data{3}.Yphys; 203 Xphys=Data{3}.Xphys; 204 end 205 206 207 208 elseif NbView==4 % is there is 2 stereo folders, get global U and V and compute Zphys 209 210 211 %test if the seconde camera is the same for both folder 212 for i=3:4 213 indpt(i)=strfind(Param.InputTable{i,2},'.'); % indice of the "." is the folder name 1 214 indline(i)=strfind(Param.InputTable{i,2},''); % indice of the "" is the folder name1 215 camname{i}=Param.InputTable{i,2}(indline(i)+1:indpt(i)1);% extract the second camera name 216 end 217 218 if strcmp(camname{3},camname{4})==0 219 disp_uvmat('ERROR','The 2 stereo folders should have the same camera for the second position',checkrun) 220 return 221 end 222 223 [Data{3},tild,errormsg] = nc2struct(filecell{3,indextemp}); 224 if exist('Data{3}.Civ3_FF','var') % if FF is present, remove wrong vector 225 temp=find(Data{3}.Civ3_FF==0); 226 Xmid3=Data{3}.Xmid(temp); 227 Ymid3=Data{3}.Ymid(temp); 228 U3=Data{3}.Uphys(temp); 229 V3=Data{3}.Vphys(temp); 230 else 231 Xmid3=Data{3}.Xmid; 232 Ymid3=Data{3}.Ymid; 233 U3=Data{3}.Uphys; 234 V3=Data{3}.Vphys; 235 end 236 %temporary gridd of merging the 2 stereos datas 237 [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)); 238 239 %1st folder : interpolate the first camera (Dalsa1) points on the second (common) camera 240 %(Dalsa 3) 241 x3Q=griddata(Xmid3+(U3)/2,Ymid3+(V3)/2,Xmid3(U3)/2,xq,yq); 242 y3Q=griddata(Xmid3+(U3)/2,Ymid3+(V3)/2,Ymid3(V3)/2,xq,yq); 243 244 245 [Data{4},tild,errormsg] = nc2struct(filecell{4,indextemp}); 246 if exist('Data{4}.Civ3_FF','var') % if FF is present, remove wrong vector 247 temp=find(Data{4}.Civ3_FF==0); 248 Xmid4=Data{4}.Xmid(temp); 249 Ymid4=Data{4}.Ymid(temp); 250 U4=Data{4}.Uphys(temp); 251 V4=Data{4}.Vphys(temp); 252 else 253 Xmid4=Data{4}.Xmid; 254 Ymid4=Data{4}.Ymid; 255 U4=Data{4}.Uphys; 256 V4=Data{4}.Vphys; 257 end 258 259 %2nd folder :interpolate the first camera (Dalsa2) points on the second (common) camera 260 %(Dalsa 3) 261 x4Q=griddata(Xmid4+(U4)/2,Ymid4+(V4)/2,Xmid4(U4)/2,xq,yq); 262 y4Q=griddata(Xmid4+(U4)/2,Ymid4+(V4)/2,Ymid4(V4)/2,xq,yq); 263 264 xmid=reshape((x4Q+x3Q)/2,length(xq(:,1)).*length(xq(1,:)),1); 265 ymid=reshape((y4Q+y3Q)/2,length(yq(:,1)).*length(yq(1,:)),1); 266 u=reshape(x4Qx3Q,length(xq(:,1)).*length(xq(1,:)),1); 267 v=reshape(y4Qy3Q,length(yq(:,1)).*length(yq(1,:)),1); 268 269 270 [Zphys,Xphys,Yphys,error]=shift2z(xmid, ymid, u, v,XmlData); %get Xphy,Yphy and Zphys 271 %remove NaN 272 tempNaN=isnan(Zphys);tempind=find(tempNaN==1); 273 Zphys(tempind)=[]; 274 Xphys(tempind)=[]; 275 Yphys(tempind)=[]; 276 error(tempind)=[]; 277 278 end 279 280 281 282 283 ZItemp(:,:,indextemp)=griddata(Xphys,Yphys,Zphys,XI,YI); %interpolation on the choosen gridd 284 285 end 286 ZI=mean(ZItemp,3); %mean between two the two time step 287 183 288 [Xa,Ya]=px_XYZ(XmlData{1}.GeometryCalib,XI,YI,ZI);% set of image coordinates on view a 184 289 [Xb,Yb]=px_XYZ(XmlData{2}.GeometryCalib,XI,YI,ZI);% set of image coordinates on view b … … 208 313 end 209 314 end 210 Ua=griddata(Data{1}.X,Data{1}.Y,Data{1}.U,Xa,Ya); 211 Va=griddata(Data{1}.X,Data{1}.Y,Data{1}.V,Xa,Ya); 315 %remove wrong vector 316 temp=find(Data{1}.FF==0); 317 X1=Data{1}.X(temp); 318 Y1=Data{1}.Y(temp); 319 U1=Data{1}.U(temp); 320 V1=Data{1}.V(temp); 321 322 Ua=griddata(X1,Y1,U1,Xa,Ya); 323 Va=griddata(X1,Y1,V1,Xa,Ya); 212 324 A=get_coeff(XmlData{1}.GeometryCalib,Xa,Ya,YI,YI,ZI); 213 325 214 Ub=griddata(Data{2}.X,Data{2}.Y,Data{1}.U,Xb,Yb); 215 Vb=griddata(Data{2}.X,Data{2}.Y,Data{2}.V,Xb,Yb); 326 327 temp=find(Data{2}.FF==0); 328 X2=Data{2}.X(temp); 329 Y2=Data{2}.Y(temp); 330 U2=Data{2}.U(temp); 331 V2=Data{2}.V(temp); 332 Ub=griddata(X2,Y2,U2,Xb,Yb); 333 Vb=griddata(X2,Y2,V2,Xb,Yb); 216 334 B=get_coeff(XmlData{2}.GeometryCalib,Xb,Yb,YI,YI,ZI); 217 335 S=ones(size(XI,1),size(XI,2),3); … … 302 420 303 421 304 305 306 422 function [z,Xphy,Yphy,error]=shift2z(xmid, ymid, u, v,XmlData) 423 z=0; 424 error=0; 425 426 427 %% first image 428 Calib_A=XmlData{1}.GeometryCalib; 429 R=(Calib_A.R)'; 430 x_a=xmid u/2; 431 y_a=ymid u/2; 432 z_a=R(7)*x_a+R(8)*y_a+Calib_A.Tx_Ty_Tz(1,3); 433 Xa=(R(1)*x_a+R(2)*y_a+Calib_A.Tx_Ty_Tz(1,1))./z_a; 434 Ya=(R(4)*x_a+R(5)*y_a+Calib_A.Tx_Ty_Tz(1,2))./z_a; 435 436 A_1_1=R(1)R(7)*Xa; 437 A_1_2=R(2)R(8)*Xa; 438 A_1_3=R(3)R(9)*Xa; 439 A_2_1=R(4)R(7)*Ya; 440 A_2_2=R(5)R(8)*Ya; 441 A_2_3=R(6)R(9)*Ya; 442 Det=A_1_1.*A_2_2A_1_2.*A_2_1; 443 Dxa=(A_1_2.*A_2_3A_2_2.*A_1_3)./Det; 444 Dya=(A_2_1.*A_1_3A_1_1.*A_2_3)./Det; 445 446 %% second image 447 Calib_B=XmlData{2}.GeometryCalib; 448 R=(Calib_B.R)'; 449 x_b=xmid+ u/2; 450 y_b=ymid+ v/2; 451 z_b=R(7)*x_b+R(8)*y_b+Calib_B.Tx_Ty_Tz(1,3); 452 Xb=(R(1)*x_b+R(2)*y_b+Calib_B.Tx_Ty_Tz(1,1))./z_b; 453 Yb=(R(4)*x_b+R(5)*y_b+Calib_B.Tx_Ty_Tz(1,2))./z_b; 454 B_1_1=R(1)R(7)*Xb; 455 B_1_2=R(2)R(8)*Xb; 456 B_1_3=R(3)R(9)*Xb; 457 B_2_1=R(4)R(7)*Yb; 458 B_2_2=R(5)R(8)*Yb; 459 B_2_3=R(6)R(9)*Yb; 460 Det=B_1_1.*B_2_2B_1_2.*B_2_1; 461 Dxb=(B_1_2.*B_2_3B_2_2.*B_1_3)./Det; 462 Dyb=(B_2_1.*B_1_3B_1_1.*B_2_3)./Det; 463 464 %% result 465 Den=(DxbDxa).*(DxbDxa)+(DybDya).*(DybDya); 466 error=((DybDya).*u(DxbDxa).*v)./Den; 467 468 xnew(1,:)=Dxa.*z+x_a; 469 xnew(2,:)=Dxb.*z+x_b; 470 ynew(1,:)=Dya.*z+y_a; 471 ynew(2,:)=Dyb.*z+y_b; 472 473 Xphy=mean(xnew,1); %on moyenne les 2 valeurs 474 Yphy=mean(ynew,1); 475 z=((DxbDxa).*u(DybDya).*v)./Den; 476 477
Note: See TracChangeset
for help on using the changeset viewer.