source: trunk/src/geometry_calib.m @ 363

Last change on this file since 363 was 356, checked in by sommeria, 13 years ago

civ updated with new functions for opening files, consistently with uvmat
Bugs to be expected (use previous version then)

File size: 60.4 KB
Line 
1%'geometry_calib': associated to the GUI geometry_calib to perform geometric calibration from a set of reference points
2%------------------------------------------------------------------------
3% function hgeometry_calib = geometry_calib(inputfile,pos)
4%
5%OUTPUT:
6% hgeometry_calib=current handles of the GUI geometry_calib.fig
7%
8%INPUT:
9% inputfile: (optional) name of an xml file containing coordinates of reference points
10% pos: (optional) 4 element vector setting the 'Position' of the GUI
11%
12%A%AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
13%  Copyright Joel Sommeria, 2008, LEGI / CNRS-UJF-INPG, sommeria@coriolis-legi.org.
14%AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
15%     This file is part of the toolbox UVMAT.
16%
17%     UVMAT is free software; you can redistribute it and/or modify
18%     it under the terms of the GNU General Public License as published by
19%     the Free Software Foundation; either version 2 of the License, or
20%     (at your option) any later version.
21%
22%     UVMAT is distributed in the hope that it will be useful,
23%     but WITHOUT ANY WARRANTY; without even the implied warranty of
24%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
25%     GNU General Public License (file UVMAT/COPYING.txt) for more details.
26%AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
27
28function varargout = geometry_calib(varargin)
29% GEOMETRY_CALIB M-file for geometry_calib.fig
30%      GEOMETRY_CALIB, by itself, creates a MenuCoord GEOMETRY_CALIB or raises the existing
31%      singleton*.
32%
33%      H = GEOMETRY_CALIB returns the handle to a MenuCoord GEOMETRY_CALIB or the handle to
34%      the existing singleton*.
35%
36%      GEOMETRY_CALIB('CALLBACK',hObject,eventData,handles,...) calls the local
37%      function named CALLBACK in GEOMETRY_CALIB.M with the given input arguments.
38%
39%      GEOMETRY_CALIB('Property','Value',...) creates a MenuCoord GEOMETRY_CALIB or raises the
40%      existing singleton*.  Starting from the left, property value pairs are
41%      applied to the GUI before geometry_calib_OpeningFunction gets called.  An
42%      unrecognized property name or invalid value makes property application
43%      stop.  All inputs are passed to geometry_calib_OpeningFcn via varargin.
44%
45%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one
46%      instance to run (singleton)".
47%
48% See also: GUIDE, GUIDATA, GUIHANDLES
49
50% Edit the above text to modify the response to help geometry_calib
51
52% Last Modified by GUIDE v2.5 05-Oct-2010 13:47:00
53
54% Begin initialization code - DO NOT edit
55gui_Singleton = 1;
56gui_State = struct('gui_Name',       mfilename, ...
57                   'gui_Singleton',  gui_Singleton, ...
58                   'gui_OpeningFcn', @geometry_calib_OpeningFcn, ...
59                   'gui_OutputFcn',  @geometry_calib_OutputFcn, ...
60                   'gui_LayoutFcn',  [] , ...
61                   'gui_Callback',   []);
62if nargin && ischar(varargin{1}) && ~isempty(regexp(varargin{1},'_Callback','once'))
63    gui_State.gui_Callback = str2func(varargin{1});
64end
65
66if nargout
67    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
68else
69    gui_mainfcn(gui_State, varargin{:});
70end
71% End initialization code - DO NOT edit
72
73
74% --- Executes just before geometry_calib is made visible.
75%INPUT:
76%handles: handles of the geometry_calib interface elements
77% PlotHandles: set of handles of the elements contolling the plotting
78% parameters on the uvmat interface (obtained by 'get_plot_handle.m')
79%------------------------------------------------------------------------
80function geometry_calib_OpeningFcn(hObject, eventdata, handles,inputfile,pos)
81%------------------------------------------------------------------------
82% Choose default command line output for geometry_calib
83
84handles.output = hObject;
85
86% Update handles structure
87guidata(hObject, handles);
88set(hObject,'DeleteFcn',{@closefcn})%
89set(hObject,'WindowButtonDownFcn',{'mouse_alt_gui',handles}) % allows mouse action with right button (zoom for uicontrol display)
90
91%set the position of the interface
92if exist('pos','var')&& length(pos)>=4
93%     %pos_gui=get(hObject,'Position');
94%     pos_gui(1)=pos(1);
95%     pos_gui(2)=pos(2);
96    set(hObject,'Position',pos);
97end
98
99%set menu of calibration options
100set(handles.calib_type,'String',{'rescale';'linear';'3D_linear';'3D_quadr';'3D_extrinsic'})
101% inputxml='';
102if exist('inputfile','var')&& ~isempty(inputfile)
103    struct.XmlInputFile=inputfile;
104%     [Pathsub,RootFile,field_count,str2,str_a,str_b,ext]=name2display(inputfile);
105    [RootPath,~,RootFile,~,~,~,~,FileExt]=fileparts_uvmat(inputfile);
106    if ~strcmp(FileExt,'.xml')
107        inputfile=[fullfile(RootPath,RootFile) '.xml'];%xml file corresponding to the input file
108    end
109    set(handles.ListCoord,'String',{'......'})
110    if exist(inputfile,'file')
111        Heading=loadfile(handles,inputfile);% load the point coordiantes existing in the xml file
112        if isfield(Heading,'Campaign')&& ischar(Heading.Campaign)
113            struct.Campaign=Heading.Campaign;
114        end
115    end   
116    set(hObject,'UserData',struct)
117end
118
119set(handles.ListCoord,'KeyPressFcn',{@key_press_fcn,handles})%set keyboard action function
120
121
122%------------------------------------------------------------------------
123% --- Outputs from this function are returned to the command line.
124function varargout = geometry_calib_OutputFcn(hObject, eventdata, handles)
125%------------------------------------------------------------------------
126% Get default command line output from handles structure
127varargout{1} = handles.output;
128varargout{2}=handles;
129%
130%------------------------------------------------------------------------
131% executed when closing: set the parent interface button to value 0
132function closefcn(gcbo,eventdata)
133%------------------------------------------------------------------------
134huvmat=findobj(allchild(0),'Name','uvmat');
135if ~isempty(huvmat)
136    handles=guidata(huvmat);
137    hobject=findobj(handles.axes3,'tag','calib_points');
138    if ~isempty(hobject)
139        delete(hobject)
140    end
141    hobject=findobj(handles.axes3,'tag','calib_marker');
142    if ~isempty(hobject)
143        delete(hobject)
144    end   
145end
146
147%------------------------------------------------------------------------
148% --- Executes on button press in calibrate_lin.
149function APPLY_Callback(hObject, eventdata, handles)
150%------------------------------------------------------------------------
151%read the current calibration points
152Coord_cell=get(handles.ListCoord,'String');
153Object=read_geometry_calib(Coord_cell);
154Coord=Object.Coord;
155% apply the calibration, whose type is selected in  handles.calib_type
156if ~isempty(Coord)
157    calib_cell=get(handles.calib_type,'String');
158    val=get(handles.calib_type,'Value');
159    GeometryCalib=feval(['calib_' calib_cell{val}],Coord,handles);
160else
161    msgbox_uvmat('ERROR','No calibration points, abort')
162    return
163end
164Z_plane=[];
165if ~isempty(Coord)
166    %check error
167    X=Coord(:,1);
168    Y=Coord(:,2);
169    Z=Coord(:,3);
170    x_ima=Coord(:,4);
171    y_ima=Coord(:,5);
172    [Xpoints,Ypoints]=px_XYZ(GeometryCalib,X,Y,Z);
173    GeometryCalib.ErrorRms(1)=sqrt(mean((Xpoints-x_ima).*(Xpoints-x_ima)));
174    [GeometryCalib.ErrorMax(1),index(1)]=max(abs(Xpoints-x_ima));
175    GeometryCalib.ErrorRms(2)=sqrt(mean((Ypoints-y_ima).*(Ypoints-y_ima)));
176    [GeometryCalib.ErrorMax(2),index(2)]=max(abs(Ypoints-y_ima));
177    [~,ind_dim]=max(GeometryCalib.ErrorMax);
178    index=index(ind_dim);
179    %set the Z position of the reference plane used for calibration
180    if isequal(max(Z),min(Z))%Z constant
181        Z_plane=Z(1);
182        GeometryCalib.NbSlice=1;
183        GeometryCalib.SliceCoord=[0 0 Z_plane];
184    end
185end
186%set the coordinate unit
187unitlist=get(handles.CoordUnit,'String');
188unit=unitlist{get(handles.CoordUnit,'value')};
189GeometryCalib.CoordUnit=unit;
190%record the points
191GeometryCalib.SourceCalib.PointCoord=Coord;
192display_intrinsic(GeometryCalib,handles)%display calibration intrinsic parameters
193
194% Display extrinsinc parameters (rotation and translation of camera with  respect to the phys coordiantes)
195set(handles.Tx,'String',num2str(GeometryCalib.Tx_Ty_Tz(1),4))
196set(handles.Ty,'String',num2str(GeometryCalib.Tx_Ty_Tz(2),4))
197set(handles.Tz,'String',num2str(GeometryCalib.Tx_Ty_Tz(3),4))
198set(handles.Phi,'String',num2str(GeometryCalib.omc(1),4))
199set(handles.Theta,'String',num2str(GeometryCalib.omc(2),4))
200set(handles.Psi,'String',num2str(GeometryCalib.omc(3),4))
201
202% store the calibration data, by default in the xml file of the currently displayed image
203huvmat=findobj(allchild(0),'Name','uvmat');
204UvData=get(huvmat,'UserData');
205NbSlice_j=1;%default
206ZStart=Z_plane;
207ZEnd=Z_plane;
208volume_scan='n';
209if isfield(UvData,'XmlData')
210    UvData.XmlData
211    if isfield(UvData.XmlData,'TranslationMotor')
212        NbSlice_j=UvData.XmlData.TranslationMotor.Nbslice;
213        ZStart=UvData.XmlData.TranslationMotor.ZStart/10;
214        ZEnd=UvData.XmlData.TranslationMotor.ZEnd/10;
215        volume_scan='y';
216    end
217end
218hhuvmat=guidata(huvmat);%handles of elements in the GUI uvmat
219RootPath='';
220% RootFile='';
221if ~isempty(hhuvmat.RootPath)&& ~isempty(hhuvmat.RootFile)
222%     testhandle=1;
223    RootPath=get(hhuvmat.RootPath,'String');
224    RootFile=get(hhuvmat.RootFile,'String');
225    filebase=fullfile(RootPath,RootFile);
226    outputfile=[filebase '.xml'];%xml file associated with the currently displayed image
227else
228    question={'save the calibration data and point coordinates in'};
229    def={fullfile(RootPath,'ObjectCalib.xml')};
230    options.Resize='on';
231    answer=inputdlg(question,'save average in a new file',1,def,options);
232    outputfile=answer{1};
233end
234answer=msgbox_uvmat('INPUT_Y-N',{[outputfile ' updated with calibration data'];...
235    ['Error rms (along x,y)=' num2str(GeometryCalib.ErrorRms) ' pixels'];...
236    ['Error max (along x,y)=' num2str(GeometryCalib.ErrorMax) ' pixels']});
237
238%% record the calibration parameters and display the current image of uvmat in the new phys coordinates
239if strcmp(answer,'Yes')
240    if strcmp(calib_cell{val}(1:2),'3D')%set the plane position for 3D (projection) calibration
241       input_key={'Z (first position)','Z (last position)','Z (water surface)', 'refractive index','NbSlice','volume scan (y/n)','tilt angle y axis','tilt angle x axis'};
242       input_val=[{num2str(ZEnd)} {num2str(ZStart)} {num2str(ZStart)} {'1.333'} num2str(NbSlice_j) {volume_scan} {'0'} {'0'}];
243        answer=inputdlg(input_key,'slice position(s)',ones(1,8), input_val,'on');
244        %answer_1=msgbox_uvmat('INPUT_TXT',' Z= ',num2str(Z_plane));
245        GeometryCalib.NbSlice=str2double(answer{5});
246        GeometryCalib.VolumeScan=answer{6};
247        if isempty(answer)
248            Z_plane=0; %default
249        else
250            Z_plane=linspace(str2double(answer{1}),str2double(answer{2}),GeometryCalib.NbSlice);
251        end     
252        GeometryCalib.SliceCoord=Z_plane'*[0 0 1];
253        GeometryCalib.SliceAngle(:,3)=0;
254        GeometryCalib.SliceAngle(:,2)=str2double(answer{7})*ones(GeometryCalib.NbSlice,1);%rotation around y axis (to generalise)
255        GeometryCalib.SliceAngle(:,1)=str2double(answer{8})*ones(GeometryCalib.NbSlice,1);%rotation around x axis (to generalise)
256        GeometryCalib.InterfaceCoord=[0 0 str2double(answer{3})];
257        GeometryCalib.RefractionIndex=str2double(answer{4});     
258    end
259    errormsg=update_imadoc(GeometryCalib,outputfile);% introduce the calibration data in the xml file
260    if ~strcmp(errormsg,'')
261        msgbox_uvmat('ERROR',errormsg);
262    end
263   
264    %display image with new calibration in the currently opened uvmat interface
265    hhh=findobj(hhuvmat.axes3,'Tag','calib_marker');% delete calib points and markers
266    if ~isempty(hhh)
267        delete(hhh);     
268    end
269    hhh=findobj(hhuvmat.axes3,'Tag','calib_points');
270    if ~isempty(hhh)
271        delete(hhh);
272    end
273    set(hhuvmat.CheckFixLimits,'Value',0)% put FixedLimits option to 'off'
274    set(hhuvmat.CheckFixLimits,'BackgroundColor',[0.7 0.7 0.7])
275    UserData=get(handles.geometry_calib,'UserData');
276    UserData.XmlInputFile=outputfile;%save the current xml file name
277    set(handles.geometry_calib,'UserData',UserData)
278    uvmat('RootPath_Callback',hObject,eventdata,hhuvmat); %file input with xml reading  in uvmat, show the image in phys coordinates
279    MenuPlot_Callback(hObject, eventdata, handles)
280    set(handles.ListCoord,'Value',index)% indicate in the list the point with max deviation (possible mistake)
281    ListCoord_Callback(hObject, eventdata, handles)
282    figure(handles.geometry_calib)
283end
284
285%------------------------------------------------------------------
286% --- Executes on button press in calibrate_lin.
287
288function REPLICATE_Callback(hObject, eventdata, handles)
289%------------------------------------------------------------------------
290
291%% Apply calibration
292calib_cell=get(handles.calib_type,'String'); %#ok<NASGU>
293val=get(handles.calib_type,'Value'); %#ok<NASGU>
294
295%read the current calibration points
296Coord_cell=get(handles.ListCoord,'String');
297Object=read_geometry_calib(Coord_cell);
298Coord=Object.Coord;
299
300% apply the calibration, whose type is selected in  handles.calib_type
301if ~isempty(Coord)
302    calib_cell=get(handles.calib_type,'String');
303    val=get(handles.calib_type,'Value');
304    GeometryCalib=feval(['calib_' calib_cell{val}],Coord,handles);
305else
306    msgbox_uvmat('ERROR','No calibration points, abort')
307    return
308end
309
310if ~isempty(Coord)
311    %check error
312    X=Coord(:,1);
313    Y=Coord(:,2);
314    Z=Coord(:,3);
315    x_ima=Coord(:,4);
316    y_ima=Coord(:,5);
317    [Xpoints,Ypoints]=px_XYZ(GeometryCalib,X,Y,Z);
318    GeometryCalib.ErrorRms(1)=sqrt(mean((Xpoints-x_ima).*(Xpoints-x_ima)));
319    [GeometryCalib.ErrorMax(1)]=max(abs(Xpoints-x_ima));
320    GeometryCalib.ErrorRms(2)=sqrt(mean((Ypoints-y_ima).*(Ypoints-y_ima)));
321    [GeometryCalib.ErrorMax(2)]=max(abs(Ypoints-y_ima));
322%     [EM,ind_dim]=max(GeometryCalib.ErrorMax);
323    %set the Z position of the reference plane used for calibration
324    Z_plane=[];
325    if isequal(max(Z),min(Z))
326        Z_plane=Z(1);
327    end
328    answer_1=msgbox_uvmat('INPUT_TXT',' Z= ',num2str(Z_plane));
329    Z_plane=str2double(answer_1);
330    GeometryCalib.NbSlice=1;
331    GeometryCalib.SliceCoord=[0 0 Z_plane];
332    %set the coordinate unit
333    unitlist=get(handles.CoordUnit,'String');
334    unit=unitlist{get(handles.CoordUnit,'value')};
335    GeometryCalib.CoordUnit=unit;
336    %record the points
337    GeometryCalib.SourceCalib.PointCoord=Coord;
338end
339
340%% display calibration paprameters
341display_intrinsic(GeometryCalib,handles)%display calibration intrinsic parameters
342
343% Display extrinsinc parameters (rotation and translation of camera with  respect to the phys coordiantes)
344set(handles.Tx,'String',num2str(GeometryCalib.Tx_Ty_Tz(1),4))
345set(handles.Ty,'String',num2str(GeometryCalib.Tx_Ty_Tz(2),4))
346set(handles.Tz,'String',num2str(GeometryCalib.Tx_Ty_Tz(3),4))
347set(handles.Phi,'String',num2str(GeometryCalib.omc(1),4))
348set(handles.Theta,'String',num2str(GeometryCalib.omc(2),4))
349set(handles.Psi,'String',num2str(GeometryCalib.omc(3),4))
350
351%% open the GUI dataview
352h_dataview=findobj(allchild(0),'name','dataview');
353if ~isempty(h_dataview)
354    delete(h_dataview)
355end
356CalibData=get(handles.geometry_calib,'UserData');%read the calibration image source on the interface userdata
357% InputFile='';
358if isfield(CalibData,'XmlInputFile')
359    InputDir=fileparts(CalibData.XmlInputFile);
360    [InputDir,DirName]=fileparts(InputDir);
361end
362SubCampaignTest='n'; %default
363testup=0;
364if isfield(CalibData,'SubCampaign')
365    SubCampaignTest='y';
366    dir_ref=CalibData.SubCampaign;
367    testup=1;
368elseif isfield(CalibData,'Campaign')
369    dir_ref=CalibData.Campaign;
370    testup=1;
371end
372while testup
373    [InputDir,DirName]=fileparts(InputDir);
374    if strcmp(DirName,dir_ref)
375        break
376    end
377end
378InputDir=fullfile(InputDir,DirName);
379answer=msgbox_uvmat('INPUT_TXT','Campaign ?',InputDir);
380if strcmp(answer,'Cancel')
381    return
382end
383
384dataview(answer,SubCampaignTest,GeometryCalib);
385
386%------------------------------------------------------------------------
387% determine the parameters for a calibration by an affine function (rescaling and offset, no rotation)
388function GeometryCalib=calib_rescale(Coord,handles)
389%------------------------------------------------------------------------
390X=Coord(:,1);
391Y=Coord(:,2);% Z not used
392x_ima=Coord(:,4);
393y_ima=Coord(:,5);
394[px]=polyfit(X,x_ima,1);
395[py]=polyfit(Y,y_ima,1);
396% T_x=px(2);
397% T_y=py(2);
398GeometryCalib.CalibrationType='rescale';
399GeometryCalib.fx_fy=[px(1) py(1)];%.fx_fy corresponds to pxcm along x and y
400GeometryCalib.CoordUnit=[];% default value, to be updated by the calling function
401GeometryCalib.Tx_Ty_Tz=[px(2)/px(1) py(2)/py(1) 1];
402GeometryCalib.omc=[0 0 0];
403
404%------------------------------------------------------------------------
405% determine the parameters for a calibration by a linear transform matrix (rescale and rotation)
406
407
408function GeometryCalib=calib_linear(Coord,handles)
409%------------------------------------------------------------------------
410X=Coord(:,1);
411Y=Coord(:,2);% Z not used
412x_ima=Coord(:,4);
413y_ima=Coord(:,5);
414XY_mat=[ones(size(X)) X Y];
415a_X1=XY_mat\x_ima; %transformation matrix for X
416a_Y1=XY_mat\y_ima;%transformation matrix for X
417R=[a_X1(2),a_X1(3);a_Y1(2),a_Y1(3)];
418epsilon=sign(det(R));
419norm=abs(det(R));
420GeometryCalib.CalibrationType='linear';
421if (a_X1(2)/a_Y1(3))>0
422    GeometryCalib.fx_fy(1)=sqrt((a_X1(2)/a_Y1(3))*norm);
423else
424    GeometryCalib.fx_fy(1)=-sqrt(-(a_X1(2)/a_Y1(3))*norm);
425end
426GeometryCalib.fx_fy(2)=(a_Y1(3)/a_X1(2))*GeometryCalib.fx_fy(1);
427GeometryCalib.CoordUnit=[];% default value, to be updated by the calling function
428GeometryCalib.Tx_Ty_Tz=[a_X1(1)/GeometryCalib.fx_fy(1) a_Y1(1)/GeometryCalib.fx_fy(2) 1];
429R(1,:)=R(1,:)/GeometryCalib.fx_fy(1);
430R(2,:)=R(2,:)/GeometryCalib.fx_fy(2);
431R=[R;[0 0]];
432GeometryCalib.R=[R [0;0;-epsilon]];
433GeometryCalib.omc=(180/pi)*[acos(GeometryCalib.R(1,1)) 0 0];
434%------------------------------------------------------------------------
435% determine the tsai parameters for a view normal to the grid plane
436% NOT USED
437function GeometryCalib=calib_normal(Coord,handles)
438%------------------------------------------------------------------------
439Calib.f1=str2num(get(handles.fx,'String'));
440Calib.f2=str2num(get(handles.fy,'String'));
441Calib.k=str2num(get(handles.kc,'String'));
442Calib.Cx=str2num(get(handles.Cx,'String'));
443Calib.Cy=str2num(get(handles.Cy,'String'));
444%default
445if isempty(Calib.f1)
446    Calib.f1=25/0.012;
447end
448if isempty(Calib.f2)
449    Calib.f2=25/0.012;
450end
451if isempty(Calib.k)
452    Calib.k=0;
453end
454if isempty(Calib.Cx)||isempty(Calib.Cy)
455    huvmat=findobj(allchild(0),'Tag','uvmat');
456    hhuvmat=guidata(huvmat);
457    Calib.Cx=str2num(get(hhuvmat.num_Npx,'String'))/2;
458    Calib.Cx=str2num(get(hhuvmat.num_Npy,'String'))/2;
459end   
460%tsai parameters
461Calib.dpx=0.012;%arbitrary
462Calib.dpy=0.012;
463Calib.sx=Calib.f1*Calib.dpx/(Calib.f2*Calib.dpy);
464Calib.f=Calib.f2*Calib.dpy;
465Calib.kappa1=Calib.k/(Calib.f*Calib.f);
466
467%initial guess
468X=Coord(:,1);
469Y=Coord(:,2);
470Zmean=mean(Coord(:,3));
471x_ima=Coord(:,4)-Calib.Cx;
472y_ima=Coord(:,5)-Calib.Cy;
473XY_mat=[ones(size(X)) X Y];
474a_X1=XY_mat\x_ima; %transformation matrix for X
475a_Y1=XY_mat\y_ima;%transformation matrix for Y
476R=[a_X1(2),a_X1(3),0;a_Y1(2),a_Y1(3),0;0,0,-1];% rotation+ z axis reversal (upward)
477norm=sqrt(det(-R));
478calib_param(1)=0;% quadratic distortion
479calib_param(2)=a_X1(1);
480calib_param(3)=a_Y1(1);
481calib_param(4)=Calib.f/(norm*Calib.dpx)-R(3,3)*Zmean;
482calib_param(5)=angle(a_X1(2)+1i*a_X1(3));
483display(['initial guess=' num2str(calib_param)])
484
485%optimise the parameters: minimisation of error
486calib_param = fminsearch(@(calib_param) error_calib(calib_param,Calib,Coord),calib_param);
487
488GeometryCalib.CalibrationType='tsai_normal';
489GeometryCalib.focal=Calib.f;
490GeometryCalib.dpx_dpy=[Calib.dpx Calib.dpy];
491GeometryCalib.Cx_Cy=[Calib.Cx Calib.Cy];
492GeometryCalib.sx=Calib.sx;
493GeometryCalib.kappa1=calib_param(1);
494GeometryCalib.CoordUnit=[];% default value, to be updated by the calling function
495GeometryCalib.Tx_Ty_Tz=[calib_param(2) calib_param(3) calib_param(4)];
496alpha=calib_param(5);
497GeometryCalib.R=[cos(alpha) sin(alpha) 0;-sin(alpha) cos(alpha) 0;0 0 -1];
498
499%------------------------------------------------------------------------
500function GeometryCalib=calib_3D_linear(Coord,handles)
501%------------------------------------------------------------------
502path_uvmat=which('uvmat');% check the path detected for source file uvmat
503path_UVMAT=fileparts(path_uvmat); %path to UVMAT
504huvmat=findobj(allchild(0),'Tag','uvmat');
505hhuvmat=guidata(huvmat);
506coord_files=get(handles.coord_files,'String');
507if ischar(coord_files)
508    coord_files={coord_files};
509end
510if isempty(coord_files{1}) || isequal(coord_files,{''})
511    coord_files={};
512end
513%retrieve the calibration points stored in the files listed in the popup list coord_files
514x_1=Coord(:,4:5)';%px coordinates of the ref points
515nx=str2num(get(hhuvmat.num_Npx,'String'));
516ny=str2num(get(hhuvmat.num_Npy,'String'));
517x_1(2,:)=ny-x_1(2,:);%reverse the y image coordinates
518X_1=Coord(:,1:3)';%phys coordinates of the ref points
519n_ima=numel(coord_files)+1;
520if ~isempty(coord_files)
521    msgbox_uvmat('CONFIRMATION',['The xy coordinates of the calibration points in ' num2str(n_ima) ' planes will be used'])
522    for ifile=1:numel(coord_files)
523    t=xmltree(coord_files{ifile});
524    s=convert(t);%convert to matlab structure
525        if isfield(s,'GeometryCalib')
526            if isfield(s.GeometryCalib,'SourceCalib')
527                if isfield(s.GeometryCalib.SourceCalib,'PointCoord')
528                PointCoord=s.GeometryCalib.SourceCalib.PointCoord;
529                Coord_file=zeros(length(PointCoord),5);%default
530                for i=1:length(PointCoord)
531                    line=str2num(PointCoord{i});
532                    Coord_file(i,4:5)=line(4:5);%px x
533                    Coord_file(i,1:3)=line(1:3);%phys x
534                end
535                eval(['x_' num2str(ifile+1) '=Coord_file(:,4:5)'';']);
536                eval(['x_' num2str(ifile+1) '(2,:)=ny-x_' num2str(ifile+1) '(2,:);' ]);
537                eval(['X_' num2str(ifile+1) '=Coord_file(:,1:3)'';']);
538                end
539            end
540        end
541    end
542end
543n_ima=numel(coord_files)+1;
544est_dist=[0;0;0;0;0];
545est_aspect_ratio=0;
546est_fc=[1;1];
547%fc=[25;25]/0.012;
548center_optim=0;
549run(fullfile(path_UVMAT,'toolbox_calib','go_calib_optim'));
550GeometryCalib.CalibrationType='3D_linear';
551GeometryCalib.fx_fy=fc';
552%GeometryCalib.focal=fc(2);
553%GeometryCalib.dpx_dpy=[1 1];
554GeometryCalib.Cx_Cy=cc';
555%GeometryCalib.sx=fc(1)/fc(2);
556GeometryCalib.kc=kc(1);
557%GeometryCalib.kappa1=-kc(1)/fc(2)^2;
558GeometryCalib.CoordUnit=[];% default value, to be updated by the calling function
559GeometryCalib.Tx_Ty_Tz=Tc_1';
560GeometryCalib.R=Rc_1;
561GeometryCalib.R(2,1:3)=-GeometryCalib.R(2,1:3);%inversion of the y image coordinate
562GeometryCalib.Tx_Ty_Tz(2)=-GeometryCalib.Tx_Ty_Tz(2);%inversion of the y image coordinate
563GeometryCalib.Cx_Cy(2)=ny-GeometryCalib.Cx_Cy(2);%inversion of the y image coordinate
564GeometryCalib.omc=(180/pi)*omc_1;%angles in degrees
565GeometryCalib.ErrorRMS=[];
566GeometryCalib.ErrorMax=[];
567
568%------------------------------------------------------------------------
569function GeometryCalib=calib_3D_quadr(Coord,handles)
570%------------------------------------------------------------------
571
572path_uvmat=which('uvmat');% check the path detected for source file uvmat
573path_UVMAT=fileparts(path_uvmat); %path to UVMAT
574huvmat=findobj(allchild(0),'Tag','uvmat');
575hhuvmat=guidata(huvmat);
576% check_cond=0;
577coord_files=get(handles.coord_files,'String');
578if ischar(coord_files)
579    coord_files={coord_files};
580end
581if isempty(coord_files{1}) || isequal(coord_files,{''})
582    coord_files={};
583end
584
585%retrieve the calibration points stored in the files listed in the popup list coord_files
586x_1=Coord(:,4:5)';%px coordinates of the ref points
587nx=str2num(get(hhuvmat.num_Npx,'String'));
588ny=str2num(get(hhuvmat.num_Npy,'String'));
589x_1(2,:)=ny-x_1(2,:);%reverse the y image coordinates
590X_1=Coord(:,1:3)';%phys coordinates of the ref points
591n_ima=numel(coord_files)+1;
592if ~isempty(coord_files)
593    msgbox_uvmat('CONFIRMATION',['The xy coordinates of the calibration points in ' num2str(n_ima) ' planes will be used'])
594    for ifile=1:numel(coord_files)
595    t=xmltree(coord_files{ifile});
596    s=convert(t);%convert to matlab structure
597        if isfield(s,'GeometryCalib')
598            if isfield(s.GeometryCalib,'SourceCalib')
599                if isfield(s.GeometryCalib.SourceCalib,'PointCoord')
600                PointCoord=s.GeometryCalib.SourceCalib.PointCoord;
601                Coord_file=zeros(length(PointCoord),5);%default
602                for i=1:length(PointCoord)
603                    line=str2num(PointCoord{i});
604                    Coord_file(i,4:5)=line(4:5);%px x
605                    Coord_file(i,1:3)=line(1:3);%phys x
606                end
607                eval(['x_' num2str(ifile+1) '=Coord_file(:,4:5)'';']);
608                eval(['x_' num2str(ifile+1) '(2,:)=ny-x_' num2str(ifile+1) '(2,:);' ]);
609                eval(['X_' num2str(ifile+1) '=Coord_file(:,1:3)'';']);
610                end
611            end
612        end
613    end
614end
615n_ima=numel(coord_files)+1;
616est_dist=[1;0;0;0;0];
617est_aspect_ratio=1;
618%est_fc=[0;0];
619%fc=[25;25]/0.012;
620center_optim=0;
621run(fullfile(path_UVMAT,'toolbox_calib','go_calib_optim'));
622
623GeometryCalib.CalibrationType='3D_quadr';
624GeometryCalib.fx_fy=fc';
625%GeometryCalib.focal=fc(2);
626%GeometryCalib.dpx_dpy=[1 1];
627GeometryCalib.Cx_Cy=cc';
628%GeometryCalib.sx=fc(1)/fc(2);
629GeometryCalib.kc=kc(1);
630%GeometryCalib.kappa1=-kc(1)/fc(2)^2;
631GeometryCalib.CoordUnit=[];% default value, to be updated by the calling function
632GeometryCalib.Tx_Ty_Tz=Tc_1';
633GeometryCalib.R=Rc_1;
634GeometryCalib.R(2,1:3)=-GeometryCalib.R(2,1:3);%inversion of the y image coordinate
635GeometryCalib.Tx_Ty_Tz(2)=-GeometryCalib.Tx_Ty_Tz(2);%inversion of the y image coordinate
636GeometryCalib.Cx_Cy(2)=ny-GeometryCalib.Cx_Cy(2);%inversion of the y image coordinate
637GeometryCalib.omc=(180/pi)*omc_1;%angles in degrees
638GeometryCalib.ErrorRMS=[];
639GeometryCalib.ErrorMax=[];
640
641
642%------------------------------------------------------------------------
643function GeometryCalib=calib_3D_extrinsic(Coord,handles)
644%------------------------------------------------------------------
645path_uvmat=which('geometry_calib');% check the path detected for source file uvmat
646path_UVMAT=fileparts(path_uvmat); %path to UVMAT
647x_1=double(Coord(:,4:5)');%image coordiantes
648X_1=double(Coord(:,1:3)');% phys coordinates
649huvmat=findobj(allchild(0),'Tag','uvmat');
650hhuvmat=guidata(huvmat);
651ny=str2double(get(hhuvmat.num_Npy,'String'));
652x_1(2,:)=ny-x_1(2,:);%reverse the y image coordinates
653n_ima=1;
654GeometryCalib.CalibrationType='3D_extrinsic';
655GeometryCalib.fx_fy(1)=str2num(get(handles.fx,'String'));
656GeometryCalib.fx_fy(2)=str2num(get(handles.fy,'String'));
657GeometryCalib.Cx_Cy(1)=str2num(get(handles.Cx,'String'));
658GeometryCalib.Cx_Cy(2)=str2num(get(handles.Cy,'String'));
659GeometryCalib.kc=str2num(get(handles.kc,'String'));
660fct_path=fullfile(path_UVMAT,'toolbox_calib');
661addpath(fct_path)
662GeometryCalib.Cx_Cy(2)=ny-GeometryCalib.Cx_Cy(2);%reverse Cx_Cy(2) for calibration (inversion of px ordinate)
663% [omc1,Tc1,Rc1,H,x,ex,JJ] = compute_extrinsic(x_1,X_1,...
664%     [Calib.f Calib.f*Calib.sx]',...
665%     [Calib.Cx Calib.Cy]',...
666%     [-Calib.kappa1*Calib.f^2 0 0 0 0]);
667[omc,Tc1,Rc1,H,x,ex,JJ] = compute_extrinsic(x_1,X_1,...
668    (GeometryCalib.fx_fy)',GeometryCalib.Cx_Cy',[GeometryCalib.kc 0 0 0 0]);
669rmpath(fct_path);
670GeometryCalib.CoordUnit=[];% default value, to be updated by the calling function
671GeometryCalib.Tx_Ty_Tz=Tc1';
672%inversion of z axis
673GeometryCalib.R=Rc1;
674GeometryCalib.R(2,1:3)=-GeometryCalib.R(2,1:3);%inversion of the y image coordinate
675GeometryCalib.Tx_Ty_Tz(2)=-GeometryCalib.Tx_Ty_Tz(2);%inversion of the y image coordinate
676GeometryCalib.Cx_Cy(2)=ny-GeometryCalib.Cx_Cy(2);%inversion of the y image coordinate
677GeometryCalib.omc=(180/pi)*omc';
678%GeometryCalib.R(3,1:3)=-GeometryCalib.R(3,1:3);%inversion for z upward
679
680
681
682%------------------------------------------------------------------------
683%function GeometryCalib=calib_tsai_heikkila(Coord)
684% TEST: NOT IMPLEMENTED
685%------------------------------------------------------------------
686% path_uvmat=which('uvmat');% check the path detected for source file uvmat
687% path_UVMAT=fileparts(path_uvmat); %path to UVMAT
688% path_calib=fullfile(path_UVMAT,'toolbox_calib_heikkila');
689% addpath(path_calib)
690% npoints=size(Coord,1);
691% Coord(:,1:3)=10*Coord(:,1:3);
692% Coord=[Coord zeros(npoints,2) -ones(npoints,1)];
693% [par,pos,iter,res,er,C]=cacal('dalsa',Coord);
694% GeometryCalib.CalibrationType='tsai';
695% GeometryCalib.focal=par(2);
696
697
698%------------------------------------------------------------------------
699% --- determine the rms of calibration error
700function ErrorRms=error_calib(calib_param,Calib,Coord)
701%calib_param: vector of free calibration parameters (to optimise)
702%Calib: structure of the given calibration parameters
703%Coord: list of phys coordinates (columns 1-3, and pixel coordinates (columns 4-5)
704Calib.f=25;
705Calib.dpx=0.012;
706Calib.dpy=0.012;
707Calib.sx=1;
708Calib.Cx=512;
709Calib.Cy=512;
710Calib.kappa1=calib_param(1);
711Calib.Tx=calib_param(2);
712Calib.Ty=calib_param(3);
713Calib.Tz=calib_param(4);
714alpha=calib_param(5);
715Calib.R=[cos(alpha) sin(alpha) 0;-sin(alpha) cos(alpha) 0;0 0 -1];
716
717X=Coord(:,1);
718Y=Coord(:,2);
719Z=Coord(:,3);
720x_ima=Coord(:,4);
721y_ima=Coord(:,5);
722[Xpoints,Ypoints]=px_XYZ(Calib,X,Y,Z);
723ErrorRms(1)=sqrt(mean((Xpoints-x_ima).*(Xpoints-x_ima)));
724ErrorRms(2)=sqrt(mean((Ypoints-y_ima).*(Ypoints-y_ima)));
725ErrorRms=mean(ErrorRms);
726
727%------------------------------------------------------------------------
728function XImage_Callback(hObject, eventdata, handles)
729%------------------------------------------------------------------------
730update_list(hObject, eventdata,handles)
731
732%------------------------------------------------------------------------
733function YImage_Callback(hObject, eventdata, handles)
734%------------------------------------------------------------------------
735update_list(hObject, eventdata,handles)
736
737%------------------------------------------------------------------------
738% --- Executes on button press in STORE.
739function STORE_Callback(hObject, eventdata, handles)
740Coord_cell=get(handles.ListCoord,'String');
741Object=read_geometry_calib(Coord_cell);
742unitlist=get(handles.CoordUnit,'String');
743unit=unitlist{get(handles.CoordUnit,'value')};
744GeometryCalib.CoordUnit=unit;
745GeometryCalib.SourceCalib.PointCoord=Object.Coord;
746huvmat=findobj(allchild(0),'Name','uvmat');
747hhuvmat=guidata(huvmat);%handles of elements in the GUI uvmat
748% RootPath='';
749% RootFile='';
750if ~isempty(hhuvmat.RootPath)&& ~isempty(hhuvmat.RootFile)
751%     testhandle=1;
752    RootPath=get(hhuvmat.RootPath,'String');
753    RootFile=get(hhuvmat.RootFile,'String');
754    filebase=fullfile(RootPath,RootFile);
755    while exist([filebase '.xml'],'file')
756        filebase=[filebase '~'];
757    end
758    outputfile=[filebase '.xml'];
759    errormsg=update_imadoc(GeometryCalib,outputfile);
760    if ~strcmp(errormsg,'')
761        msgbox_uvmat('ERROR',errormsg);
762    end
763    listfile=get(handles.coord_files,'string');
764    if isequal(listfile,{''})
765        listfile={outputfile};
766    else
767        listfile=[listfile;{outputfile}];%update the list of coord files
768    end
769    set(handles.coord_files,'string',listfile);
770end
771set(handles.ListCoord,'Value',1)% refresh the display of coordinates
772set(handles.ListCoord,'String',{'......'})
773
774% --------------------------------------------------------------------
775% --- Executes on button press in CLEAR_PTS: clear the list of calibration points
776function CLEAR_PTS_Callback(hObject, eventdata, handles)
777% --------------------------------------------------------------------
778set(handles.ListCoord,'Value',1)% refresh the display of coordinates
779set(handles.ListCoord,'String',{'......'})
780MenuPlot_Callback(hObject, eventdata, handles)
781
782%------------------------------------------------------------------------
783% --- Executes on button press in CLEAR.
784function CLEAR_Callback(hObject, eventdata, handles)
785%------------------------------------------------------------------------
786set(handles.coord_files,'Value',1)
787set(handles.coord_files,'String',{''})
788
789%------------------------------------------------------------------------
790function XObject_Callback(hObject, eventdata, handles)
791%------------------------------------------------------------------------
792update_list(hObject, eventdata,handles)
793
794%------------------------------------------------------------------------
795function YObject_Callback(hObject, eventdata, handles)
796%------------------------------------------------------------------------
797update_list(hObject, eventdata,handles)
798
799%------------------------------------------------------------------------
800function ZObject_Callback(hObject, eventdata, handles)
801%------------------------------------------------------------------------
802update_list(hObject, eventdata,handles)
803
804%------------------------------------------------------------------------
805function update_list(hObject, eventdata, handles)
806%------------------------------------------------------------------------
807newval(4)=str2double(get(handles.XImage,'String'));
808newval(5)=str2double(get(handles.YImage,'String'));
809newval(1)=str2double(get(handles.XObject,'String'));
810newval(2)=str2double(get(handles.YObject,'String'));
811newval(3)=str2double(get(handles.ZObject,'String'));
812if isnan(newval(3))
813    newval(3)=0;%put z to 0 by default
814end
815Coord=get(handles.ListCoord,'String');
816Coord(end)=[]; %remove last string '.....'
817val=get(handles.ListCoord,'Value');
818data=read_geometry_calib(Coord);
819data.Coord(val,:)=newval;
820for i=1:size(data.Coord,1)
821    for j=1:5
822          Coord_cell{i,j}=num2str(data.Coord(i,j),4);%display coordiantes with 4 digits
823    end
824end
825
826Tabchar=cell2tab(Coord_cell,' | ');
827Tabchar=[Tabchar ;{'......'}];
828set(handles.ListCoord,'String',Tabchar)
829
830%update the plot
831ListCoord_Callback(hObject, eventdata, handles)
832MenuPlot_Callback(hObject, eventdata, handles)
833
834%------------------------------------------------------------------------
835% --- Executes on selection change in ListCoord.
836function ListCoord_Callback(hObject, eventdata, handles)
837%------------------------------------------------------------------------
838huvmat=findobj(allchild(0),'Name','uvmat');%find the current uvmat interface handle
839hplot=findobj(huvmat,'Tag','axes3');%main plotting axis of uvmat
840hhh=findobj(hplot,'Tag','calib_marker');
841Coord_cell=get(handles.ListCoord,'String');
842val=get(handles.ListCoord,'Value');
843if numel(val)>1
844    return %no action if several lines have been selected
845end
846coord_str=Coord_cell{val};
847k=findstr(' | ',coord_str);
848if isempty(k)%last line '.....' selected
849    if ~isempty(hhh)
850        delete(hhh)%delete the circle marker
851    end
852    return
853end
854%fill the edit boxex
855set(handles.XObject,'String',coord_str(1:k(1)-1))
856set(handles.YObject,'String',coord_str(k(1)+3:k(2)-1))
857set(handles.ZObject,'String',coord_str(k(2)+3:k(3)-1))
858set(handles.XImage,'String',coord_str(k(3)+3:k(4)-1))
859set(handles.YImage,'String',coord_str(k(4)+3:end))
860h_menu_coord=findobj(huvmat,'Tag','transform_fct');
861menu=get(h_menu_coord,'String');
862choice=get(h_menu_coord,'Value');
863if iscell(menu)
864    option=menu{choice};
865else
866    option='px'; %default
867end
868if isequal(option,'phys')
869    XCoord=str2double(coord_str(1:k(1)-1));
870    YCoord=str2double(coord_str(k(1)+3:k(2)-1));
871elseif isequal(option,'px')|| isequal(option,'')
872    XCoord=str2double(coord_str(k(3)+3:k(4)-1));
873    YCoord=str2double(coord_str(k(4)+3:end));
874else
875    msgbox_uvmat('ERROR','the choice in menu_coord of uvmat must be px or phys ')
876end
877if isempty(XCoord)||isempty(YCoord)
878     if ~isempty(hhh)
879        delete(hhh)%delete the circle marker
880    end
881    return
882end
883xlim=get(hplot,'XLim');
884ylim=get(hplot,'YLim');
885ind_range=max(abs(xlim(2)-xlim(1)),abs(ylim(end)-ylim(1)))/20;%defines the size of the circle marker
886if isempty(hhh)
887    set(0,'CurrentFig',huvmat)
888    set(huvmat,'CurrentAxes',hplot)
889    rectangle('Curvature',[1 1],...
890              'Position',[XCoord-ind_range/2 YCoord-ind_range/2 ind_range ind_range],'EdgeColor','m',...
891              'LineStyle','-','Tag','calib_marker');
892else
893    set(hhh,'Position',[XCoord-ind_range/2 YCoord-ind_range/2 ind_range ind_range])
894end
895
896%------------------------------------------------------------------------
897% --- Executes on selection change in edit_append.
898function edit_append_Callback(hObject, eventdata, handles)
899%------------------------------------------------------------------------
900choice=get(handles.edit_append,'Value');
901if choice
902    set(handles.edit_append,'BackgroundColor',[1 1 0])
903    huvmat=findobj(allchild(0),'tag','uvmat');
904    if ishandle(huvmat)
905        hhuvmat=guidata(huvmat);
906        set(hhuvmat.edit_object,'Value',0)
907        set(hhuvmat.edit_object,'BackgroundColor',[0.7 0.7 0.7])
908    end
909else
910    set(handles.edit_append,'BackgroundColor',[0.7 0.7 0.7])
911end
912   
913function NEW_Callback(hObject, eventdata, handles)
914%A METTRE SOUS UN BOUTON
915huvmat=findobj(allchild(0),'Name','uvmat');
916hchild=get(huvmat,'children');
917hcoord=findobj(hchild,'Tag','menu_coord');
918coordtype=get(hcoord,'Value');
919haxes=findobj(hchild,'Tag','axes3');
920AxeData=get(haxes,'UserData');
921if ~isequal(hcoord,2)
922    set(hcoord,'Value',2)
923    huvmat=uvmat(AxeData);
924    'relancer uvmat';
925end
926if ~isfield(AxeData,'ZoomAxes')
927    msgbox_uvmat('ERROR','first draw a window around a grid marker')
928    return
929end
930XLim=get(AxeData.ZoomAxes,'XLim');
931YLim=get(AxeData.ZoomAxes,'YLim');
932np=size(AxeData.A);
933ind_sub_x=round(XLim);
934ind_sub_y=np(1)-round(YLim);
935Mfiltre=AxeData.A(ind_sub_y(2):ind_sub_y(1) ,ind_sub_x,:);
936Mfiltre_norm=double(Mfiltre);
937Mfiltre_norm=Mfiltre_norm/sum(sum(Mfiltre_norm));
938Mfiltre_norm=100*(Mfiltre_norm-mean(mean(Mfiltre_norm)));
939Atype=class(AxeData.A);
940Data.NbDim=2;
941Data.A=filter2(Mfiltre_norm,double(AxeData.A));
942Data.A=feval(Atype,Data.A);
943Data.AName='image';
944Data.AX=AxeData.AX;
945Data.AY=AxeData.AY;
946Data.CoordType='px';
947plot_field(Data)
948
949%------------------------------------------------------------------------
950function MenuPlot_Callback(hObject, eventdata, handles)
951%------------------------------------------------------------------------
952huvmat=findobj(allchild(0),'Name','uvmat');%find the current uvmat interface handle
953%UvData=get(huvmat,'UserData');%Data associated to the current uvmat interface
954hhuvmat=guidata(huvmat); %handles of GUI elements in uvmat
955%hplot=findobj(huvmat,'Tag','axes3');%main plotting axis of uvmat
956h_menu_coord=findobj(huvmat,'Tag','transform_fct');
957menu=get(h_menu_coord,'String');
958choice=get(h_menu_coord,'Value');
959if iscell(menu)
960    option=menu{choice};
961else
962    option='px'; %default
963end
964Coord_cell=get(handles.ListCoord,'String');
965ObjectData=read_geometry_calib(Coord_cell);
966%ObjectData=read_geometry_calib(handles);%read the interface input parameters defining the object
967if ~isempty(ObjectData.Coord)
968    if isequal(option,'phys')
969        ObjectData.Coord=ObjectData.Coord(:,1:3);
970    elseif isequal(option,'px')||isequal(option,'')
971        ObjectData.Coord=ObjectData.Coord(:,4:5);
972    else
973        msgbox_uvmat('ERROR','the choice in menu_coord of uvmat must be '''', px or phys ')
974    end
975end
976%axes(hhuvmat.axes3)
977set(0,'CurrentFigure',huvmat)
978set(huvmat,'CurrentAxes',hhuvmat.axes3)
979hh=findobj('Tag','calib_points');
980if  ~isempty(ObjectData.Coord) && isempty(hh)
981    hh=line(ObjectData.Coord(:,1),ObjectData.Coord(:,2),'Color','m','Tag','calib_points','LineStyle','.','Marker','+');
982elseif isempty(ObjectData.Coord)%empty list of points, suppress the plot
983    delete(hh)
984else
985    set(hh,'XData',ObjectData.Coord(:,1))
986    set(hh,'YData',ObjectData.Coord(:,2))
987end
988pause(.1)
989figure(handles.geometry_calib)
990
991% --------------------------------------------------------------------
992function MenuHelp_Callback(hObject, eventdata, handles)
993path_to_uvmat=which('uvmat');% check the path of uvmat
994pathelp=fileparts(path_to_uvmat);
995helpfile=fullfile(pathelp,'uvmat_doc','uvmat_doc.html');
996if isempty(dir(helpfile)), msgbox_uvmat('ERROR','Please put the help file uvmat_doc.html in the sub-directory /uvmat_doc of the UVMAT package')
997else
998   addpath (fullfile(pathelp,'uvmat_doc'))
999   web([helpfile '#geometry_calib'])
1000end
1001
1002%------------------------------------------------------------------------
1003function MenuCreateGrid_Callback(hObject, eventdata, handles)
1004%------------------------------------------------------------------------
1005%hcalib=get(handles.calib_type,'parent');%handles of the GUI geometry_calib
1006CalibData=get(handles.geometry_calib,'UserData');
1007Tinput=[];%default
1008if isfield(CalibData,'grid')
1009    Tinput=CalibData.grid;
1010end
1011[T,CalibData.grid]=create_grid(Tinput);%display the GUI create_grid
1012set(handles.geometry_calib,'UserData',CalibData)
1013
1014%grid in phys space
1015Coord=get(handles.ListCoord,'String');
1016val=get(handles.ListCoord,'Value');
1017data=read_geometry_calib(Coord);
1018%nbpoints=size(data.Coord,1); %nbre of calibration points
1019data.Coord(val:val+size(T,1)-1,1:3)=T(end:-1:1,:);%update the existing list of phys coordinates from the GUI create_grid
1020% for i=1:nbpoints
1021%    for j=1:5
1022%           Coord{i,j}=num2str(data.Coord(i,j),4);%display coordiantes with 4 digits
1023%    end
1024% end
1025%update the phys coordinates starting from the selected point (down in the
1026Coord(end,:)=[]; %remove last string '.....'
1027for i=1:size(data.Coord,1)
1028    for j=1:5
1029          Coord{i,j}=num2str(data.Coord(i,j),4);%display coordiantes with 4 digits
1030    end
1031end
1032
1033%size(data.Coord,1)
1034Tabchar=cell2tab(Coord,' | ');
1035Tabchar=[Tabchar ;{'......'}];
1036set(handles.ListCoord,'String',Tabchar)
1037
1038% -----------------------------------------------------------------------
1039% --- automatic grid dectection from local maxima of the images
1040function MenuDetectGrid_Callback(hObject, eventdata, handles)
1041%------------------------------------------------------------------------
1042%% read the four last point coordinates in pixels
1043Coord_cell=get(handles.ListCoord,'String');%read list of coordinates on geometry_calib
1044data=read_geometry_calib(Coord_cell);
1045nbpoints=size(data.Coord,1); %nbre of calibration points
1046if nbpoints~=4
1047    msgbox_uvmat('ERROR','four points must have be selected by the mouse to delimitate the phys grid area; the Ox axis will be defined by the two first points')
1048    return
1049end
1050corners_X=(data.Coord(end:-1:end-3,4)); %pixel absissa of the four corners
1051corners_Y=(data.Coord(end:-1:end-3,5));
1052
1053%reorder the last two points (the two first in the list) if needed
1054angles=angle((corners_X-corners_X(1))+1i*(corners_Y-corners_Y(1)));
1055if abs(angles(4)-angles(2))>abs(angles(3)-angles(2))
1056      X_end=corners_X(4);
1057      Y_end=corners_Y(4);
1058      corners_X(4)=corners_X(3);
1059      corners_Y(4)=corners_Y(3);
1060      corners_X(3)=X_end;
1061      corners_Y(3)=Y_end;
1062end
1063
1064%% initiate the grid
1065CalibData=get(handles.geometry_calib,'UserData');%get information stored on the GUI geometry_calib
1066grid_input=[];%default
1067if isfield(CalibData,'grid')
1068    grid_input=CalibData.grid;%retrieve the previously used grid
1069end
1070[T,CalibData.grid,white_test]=create_grid(grid_input,'detect_grid');%display the GUI create_grid, read the set of phys coordinates T
1071set(handles.geometry_calib,'UserData',CalibData)%store the phys grid parameters for later use
1072
1073
1074
1075%% read the current image, displayed in the GUI uvmat
1076huvmat=findobj(allchild(0),'Name','uvmat');
1077UvData=get(huvmat,'UserData');
1078A=UvData.Field.A;
1079npxy=size(A);
1080X=[CalibData.grid.x_0 CalibData.grid.x_1 CalibData.grid.x_0 CalibData.grid.x_1]';%corner absissa in the phys coordinates (cm)
1081Y=[CalibData.grid.y_0 CalibData.grid.y_0 CalibData.grid.y_1 CalibData.grid.y_1]';%corner ordinates in the phys coordinates (cm)
1082
1083%calculate transform matrices for plane projection
1084% reference: http://alumni.media.mit.edu/~cwren/interpolator/ by Christopher R. Wren
1085B = [ X Y ones(size(X)) zeros(4,3)        -X.*corners_X -Y.*corners_X ...
1086      zeros(4,3)        X Y ones(size(X)) -X.*corners_Y -Y.*corners_Y ];
1087B = reshape (B', 8 , 8 )';
1088D = [ corners_X , corners_Y ];
1089D = reshape (D', 8 , 1 );
1090l = (B' * B)\B' * D;
1091Amat = reshape([l(1:6)' 0 0 1 ],3,3)';
1092C = [l(7:8)' 1];
1093
1094% transform grid image into 'phys' coordinates
1095GeometryCalib.fx_fy=[1 1];
1096GeometryCalib.Tx_Ty_Tz=[Amat(1,3) Amat(2,3) 1];
1097GeometryCalib.R=[Amat(1,1),Amat(1,2),0;Amat(2,1),Amat(2,2),0;C(1),C(2),0];
1098GeometryCalib.CoordUnit='cm';
1099path_uvmat=which('uvmat');% check the path detected for source file uvmat
1100path_UVMAT=fileparts(path_uvmat); %path to UVMAT
1101addpath(fullfile(path_UVMAT,'transform_field'))
1102Data.ListVarName={'AY','AX','A'};
1103Data.VarDimName={'AY','AX',{'AY','AX'}};
1104if ndims(A)==3
1105    A=mean(A,3);
1106end
1107Data.A=A-min(min(A));
1108Data.AY=[npxy(1)-0.5 0.5];
1109Data.AX=[0.5 npxy(2)];
1110Data.CoordUnit='pixel';
1111Calib.GeometryCalib=GeometryCalib;
1112DataOut=phys(Data,Calib);
1113rmpath(fullfile(path_UVMAT,'transform_field'))
1114Amod=DataOut.A;
1115Rangx=DataOut.AX;
1116Rangy=DataOut.AY;
1117if white_test
1118    Amod=double(Amod);%case of white grid markers: will look for image maxima
1119else
1120    Amod=-double(Amod);%case of black grid markers: will look for image minima
1121end
1122% figure(12) %display corrected image
1123% Amax=max(max(Amod));
1124% image(Rangx,Rangy,uint8(255*Amod/Amax))
1125
1126%% detection of local image extrema in each direction
1127Dx=(Rangx(2)-Rangx(1))/(npxy(2)-1); %x mesh in real space
1128Dy=(Rangy(2)-Rangy(1))/(npxy(1)-1); %y mesh in real space
1129ind_range_x=ceil(abs(GeometryCalib.R(1,1)*CalibData.grid.Dx/3));% range of search of image ma around each point obtained by linear interpolation from the marked points
1130ind_range_y=ceil(abs(GeometryCalib.R(2,2)*CalibData.grid.Dy/3));% range of search of image ma around each point obtained by linear interpolation from the marked points
1131nbpoints=size(T,1);
1132for ipoint=1:nbpoints
1133    i0=1+round((T(ipoint,1)-Rangx(1))/Dx);%round(Xpx(ipoint));
1134    j0=1+round((T(ipoint,2)-Rangy(1))/Dy);%round(Xpx(ipoint));
1135    j0min=max(j0-ind_range_y,1);
1136    j0max=min(j0+ind_range_y,size(Amod,1));
1137    i0min=max(i0-ind_range_x,1);
1138    i0max=min(i0+ind_range_x,size(Amod,2));
1139    Asub=Amod(j0min:j0max,i0min:i0max);
1140    x_profile=sum(Asub,1);
1141    y_profile=sum(Asub,2);
1142    [Amax,ind_x_max]=max(x_profile);
1143    [Amax,ind_y_max]=max(y_profile);
1144    %sub-pixel improvement using moments
1145    x_shift=0;
1146    y_shift=0;
1147    if ind_x_max+2<=numel(x_profile) && ind_x_max-2>=1
1148        Atop=x_profile(ind_x_max-2:ind_x_max+2);
1149        x_shift=sum(Atop.*[-2 -1 0 1 2])/sum(Atop);
1150    end
1151    if ind_y_max+2<=numel(y_profile) && ind_y_max-2>=1
1152        Atop=y_profile(ind_y_max-2:ind_y_max+2);
1153        y_shift=sum(Atop.*[-2 -1 0 1 2]')/sum(Atop);
1154    end
1155    Delta(ipoint,1)=(x_shift+ind_x_max+i0min-i0-1)*Dx;%shift from the initial guess
1156    Delta(ipoint,2)=(y_shift+ind_y_max+j0min-j0-1)*Dy;
1157end
1158Tmod=T(:,(1:2))+Delta;
1159[Xpx,Ypx]=px_XYZ(GeometryCalib,Tmod(:,1),Tmod(:,2));
1160for ipoint=1:nbpoints
1161     Coord{ipoint,1}=num2str(T(ipoint,1),4);%display coordiantes with 4 digits
1162     Coord{ipoint,2}=num2str(T(ipoint,2),4);%display coordiantes with 4 digits
1163     Coord{ipoint,3}=num2str(T(ipoint,3),4);%display coordiantes with 4 digits;
1164     Coord{ipoint,4}=num2str(Xpx(ipoint),4);%display coordiantes with 4 digits
1165     Coord{ipoint,5}=num2str(Ypx(ipoint),4);%display coordiantes with 4 digits
1166end
1167Tabchar=cell2tab(Coord(end:-1:1,:),' | ');
1168Tabchar=[Tabchar ;{'......'}];
1169set(handles.ListCoord,'Value',1)
1170set(handles.ListCoord,'String',Tabchar)
1171MenuPlot_Callback(hObject, eventdata, handles)
1172
1173%-----------------------------------------------------------------------
1174function MenuTranslatePoints_Callback(hObject, eventdata, handles)
1175%-----------------------------------------------------------------------
1176%hcalib=get(handles.calib_type,'parent');%handles of the GUI geometry_calib
1177CalibData=get(handles.geometry_calib,'UserData');
1178Tinput=[];%default
1179if isfield(CalibData,'translate')
1180    Tinput=CalibData.translate;
1181end
1182T=translate_points(Tinput);%display translate_points GUI and get shift parameters
1183CalibData.translate=T;
1184set(handles.geometry_calib,'UserData',CalibData)
1185%translation
1186Coord_cell=get(handles.ListCoord,'String');
1187data=read_geometry_calib(Coord_cell);
1188data.Coord(:,1)=T(1)+data.Coord(:,1);
1189data.Coord(:,2)=T(2)+data.Coord(:,2);
1190data.Coord(:,3)=T(3)+data.Coord(:,3);
1191data.Coord(:,[4 5])=data.Coord(:,[4 5]);
1192for i=1:size(data.Coord,1)
1193    for j=1:5
1194          Coord{i,j}=num2str(data.Coord(i,j),4);%phys x,y,z
1195   end
1196end
1197Tabchar=cell2tab(Coord,' | ');
1198Tabchar=[Tabchar; {'.....'}];
1199%set(handles.ListCoord,'Value',1)
1200set(handles.ListCoord,'String',Tabchar)
1201
1202
1203% --------------------------------------------------------------------
1204function MenuRotatePoints_Callback(hObject, eventdata, handles)
1205%hcalib=get(handles.calib_type,'parent');%handles of the GUI geometry_calib
1206CalibData=get(handles.geometry_calib,'UserData');
1207Tinput=[];%default
1208if isfield(CalibData,'rotate')
1209    Tinput=CalibData.rotate;
1210end
1211T=rotate_points(Tinput);%display rotate_points GUI to introduce rotation parameters
1212CalibData.rotate=T;
1213set(handles.geometry_calib,'UserData',CalibData)
1214%-----------------------------------------------------
1215%rotation
1216Phi=T(1);
1217O_x=0;%default
1218O_y=0;%default
1219if numel(T)>=2
1220    O_x=T(2);%default
1221end
1222if numel(T)>=3
1223    O_y=T(3);%default
1224end
1225Coord_cell=get(handles.ListCoord,'String');
1226data=read_geometry_calib(Coord_cell);
1227r1=cos(pi*Phi/180);
1228r2=-sin(pi*Phi/180);
1229r3=sin(pi*Phi/180);
1230r4=cos(pi*Phi/180);
1231x=data.Coord(:,1)-O_x;
1232y=data.Coord(:,2)-O_y;
1233data.Coord(:,1)=r1*x+r2*y;
1234data.Coord(:,2)=r3*x+r4*y;
1235% data.Coord(:,[4 5])=data.Coord(:,[4 5]);
1236for i=1:size(data.Coord,1)
1237    for j=1:5
1238          Coord{i,j}=num2str(data.Coord(i,j),4);%phys x,y,z
1239   end
1240end
1241Tabchar=cell2tab(Coord,' | ');
1242Tabchar=[Tabchar;{'......'}];
1243set(handles.ListCoord,'Value',1)
1244set(handles.ListCoord,'String',Tabchar)
1245
1246
1247% %------------------------------------------------------------------------
1248% % --- Executes on button press in rotation.
1249% function rotation_Callback(hObject, eventdata, handles)
1250% %------------------------------------------------------------------------
1251% angle_rot=(pi/180)*str2num(get(handles.Phi,'String'));
1252% Coord_cell=get(handles.ListCoord,'String');
1253% data=read_geometry_calib(Coord_cell);
1254% data.Coord(:,1)=cos(angle_rot)*data.Coord(:,1)+sin(angle_rot)*data.Coord(:,2);
1255% data.Coord(:,1)=-sin(angle_rot)*data.Coord(:,1)+cos(angle_rot)*data.Coord(:,2);
1256% set(handles.XObject,'String',num2str(data.Coord(:,1),4));
1257% set(handles.YObject,'String',num2str(data.Coord(:,2),4));
1258
1259
1260%------------------------------------------------------------------------
1261% image transform from px to phys
1262%INPUT:
1263%Zindex: index of plane
1264% function [A_out,Rangx,Rangy]=phys_Ima(A,Calib,ZIndex)
1265% %------------------------------------------------------------------------
1266% xcorner=[];
1267% ycorner=[];
1268% npx=[];
1269% npy=[];
1270% siz=size(A)
1271% npx=[npx siz(2)];
1272% npy=[npy siz(1)]
1273% xima=[0.5 siz(2)-0.5 0.5 siz(2)-0.5];%image coordinates of corners
1274% yima=[0.5 0.5 siz(1)-0.5 siz(1)-0.5];
1275% [xcorner,ycorner]=phys_XYZ(Calib,xima,yima,ZIndex);%corresponding physical coordinates
1276% Rangx(1)=min(xcorner);
1277% Rangx(2)=max(xcorner);
1278% Rangy(2)=min(ycorner);
1279% Rangy(1)=max(ycorner);
1280% test_multi=(max(npx)~=min(npx)) | (max(npy)~=min(npy));
1281% npx=max(npx);
1282% npy=max(npy);
1283% x=linspace(Rangx(1),Rangx(2),npx);
1284% y=linspace(Rangy(1),Rangy(2),npy);
1285% [X,Y]=meshgrid(x,y);%grid in physical coordiantes
1286% vec_B=[];
1287%
1288% zphys=0; %default
1289% if isfield(Calib,'SliceCoord') %.Z= index of plane
1290%    SliceCoord=Calib.SliceCoord(ZIndex,:);
1291%    zphys=SliceCoord(3); %to generalize for non-parallel planes
1292% end
1293% [XIMA,YIMA]=px_XYZ(Calib,X,Y,zphys);%corresponding image indices for each point in the real space grid
1294% XIMA=reshape(round(XIMA),1,npx*npy);%indices reorganized in 'line'
1295% YIMA=reshape(round(YIMA),1,npx*npy);
1296% flagin=XIMA>=1 & XIMA<=npx & YIMA >=1 & YIMA<=npy;%flagin=1 inside the original image
1297% testuint8=isa(A,'uint8');
1298% testuint16=isa(A,'uint16');
1299% if numel(siz)==2 %(B/W images)
1300%     vec_A=reshape(A,1,npx*npy);%put the original image in line
1301%     ind_in=find(flagin);
1302%     ind_out=find(~flagin);
1303%     ICOMB=((XIMA-1)*npy+(npy+1-YIMA));
1304%     ICOMB=ICOMB(flagin);%index corresponding to XIMA and YIMA in the aligned original image vec_A
1305%     vec_B(ind_in)=vec_A(ICOMB);
1306%     vec_B(ind_out)=zeros(size(ind_out));
1307%     A_out=reshape(vec_B,npy,npx);%new image in real coordinates
1308% elseif numel(siz)==3     
1309%     for icolor=1:siz(3)
1310%         vec_A=reshape(A{icell}(:,:,icolor),1,npx*npy);%put the original image in line
1311%         ind_in=find(flagin);
1312%         ind_out=find(~flagin);
1313%         ICOMB=((XIMA-1)*npy+(npy+1-YIMA));
1314%         ICOMB=ICOMB(flagin);%index corresponding to XIMA and YIMA in the aligned original image vec_A
1315%         vec_B(ind_in)=vec_A(ICOMB);
1316%         vec_B(ind_out)=zeros(size(ind_out));
1317%         A_out(:,:,icolor)=reshape(vec_B,npy,npx);%new image in real coordinates
1318%     end
1319% end
1320% if testuint8
1321%     A_out=uint8(A_out);
1322% end
1323% if testuint16
1324%     A_out=uint16(A_out);
1325% end
1326
1327%------------------------------------------------------------------------
1328% pointwise transform from px to phys
1329%INPUT:
1330%Z: index of plane
1331% function [Xphys,Yphys,Zphys]=phys_XYZ(Calib,X,Y,Z)
1332% %------------------------------------------------------------------------
1333% if exist('Z','var')& isequal(Z,round(Z))& Z>0 & isfield(Calib,'SliceCoord')&length(Calib.SliceCoord)>=Z
1334%     Zindex=Z;
1335%     Zphys=Calib.SliceCoord(Zindex,3);%GENERALISER AUX CAS AVEC ANGLE
1336% else
1337%     Zphys=0;
1338% end
1339% if ~exist('X','var')||~exist('Y','var')
1340%     Xphys=[];
1341%     Yphys=[];%default
1342%     return
1343% end
1344% Xphys=X;%default
1345% Yphys=Y;
1346% %image transform
1347% if isfield(Calib,'R')
1348%     R=(Calib.R)';
1349%     Dx=R(5)*R(7)-R(4)*R(8);
1350%     Dy=R(1)*R(8)-R(2)*R(7);
1351%     D0=Calib.f*(R(2)*R(4)-R(1)*R(5));
1352%     Z11=R(6)*R(8)-R(5)*R(9);
1353%     Z12=R(2)*R(9)-R(3)*R(8); 
1354%     Z21=R(4)*R(9)-R(6)*R(7);
1355%     Z22=R(3)*R(7)-R(1)*R(9);
1356%     Zx0=R(3)*R(5)-R(2)*R(6);
1357%     Zy0=R(1)*R(6)-R(3)*R(4);
1358%     A11=R(8)*Calib.Ty-R(5)*Calib.Tz+Z11*Zphys;
1359%     A12=R(2)*Calib.Tz-R(8)*Calib.Tx+Z12*Zphys;
1360%     A21=-R(7)*Calib.Ty+R(4)*Calib.Tz+Z21*Zphys;
1361%     A22=-R(1)*Calib.Tz+R(7)*Calib.Tx+Z11*Zphys;
1362%     X0=Calib.f*(R(5)*Calib.Tx-R(2)*Calib.Ty+Zx0*Zphys);
1363%     Y0=Calib.f*(-R(4)*Calib.Tx+R(1)*Calib.Ty+Zy0*Zphys);
1364%         %px to camera:
1365%     Xd=(Calib.dpx/Calib.sx)*(X-Calib.Cx); % sensor coordinates
1366%     Yd=Calib.dpy*(Y-Calib.Cy);
1367%     dist_fact=1+Calib.kappa1*(Xd.*Xd+Yd.*Yd); %distortion factor
1368%     Xu=dist_fact.*Xd;%undistorted sensor coordinates
1369%     Yu=dist_fact.*Yd;
1370%     denom=Dx*Xu+Dy*Yu+D0;
1371%     % denom2=denom.*denom;
1372%     Xphys=(A11.*Xu+A12.*Yu+X0)./denom;%world coordinates
1373%     Yphys=(A21.*Xu+A22.*Yu+Y0)./denom;
1374% end
1375
1376
1377% --------------------------------------------------------------------
1378function MenuImportPoints_Callback(hObject, eventdata, handles)
1379fileinput=browse_xml(hObject, eventdata, handles);
1380if isempty(fileinput)
1381    return
1382end
1383[s,errormsg]=imadoc2struct(fileinput,'GeometryCalib');
1384GeometryCalib=s.GeometryCalib;
1385%GeometryCalib=load_calib(hObject, eventdata, handles)
1386calib=reshape(GeometryCalib.PointCoord,[],1);
1387for ilist=1:numel(calib)
1388    CoordCell{ilist}=num2str(calib(ilist));
1389end
1390CoordCell=reshape(CoordCell,[],5);
1391Tabchar=cell2tab(CoordCell,' | ');%transform cells into table ready for display
1392Tabchar=[Tabchar;{'......'}];
1393set(handles.ListCoord,'Value',1)
1394set(handles.ListCoord,'String',Tabchar)
1395MenuPlot_Callback(handles.geometry_calib, [], handles)
1396
1397% -----------------------------------------------------------------------
1398function MenuImportIntrinsic_Callback(hObject, eventdata, handles)
1399%------------------------------------------------------------------------
1400fileinput=browse_xml(hObject, eventdata, handles);
1401if isempty(fileinput)
1402    return
1403end
1404[s,errormsg]=imadoc2struct(fileinput,'GeometryCalib');
1405GeometryCalib=s.GeometryCalib;
1406display_intrinsic(GeometryCalib,handles)
1407
1408% -----------------------------------------------------------------------
1409function MenuImportAll_Callback(hObject, eventdata, handles)
1410%------------------------------------------------------------------------
1411fileinput=browse_xml(hObject, eventdata, handles);
1412if ~isempty(fileinput)
1413    loadfile(handles,fileinput)
1414end
1415
1416% -----------------------------------------------------------------------
1417% --- Executes on menubar option Import/Grid file: introduce previous grid files
1418function MenuGridFile_Callback(hObject, eventdata, handles)
1419% -----------------------------------------------------------------------
1420inputfile=browse_xml(hObject, eventdata, handles);
1421listfile=get(handles.coord_files,'string');
1422if isequal(listfile,{''})
1423    listfile={inputfile};
1424else
1425    listfile=[listfile;{inputfile}];%update the list of coord files
1426end
1427set(handles.coord_files,'string',listfile);
1428
1429%------------------------------------------------------------------------
1430% --- 'key_press_fcn:' function activated when a key is pressed on the keyboard
1431function key_press_fcn(hObject,eventdata,handles)
1432%------------------------------------------------------------------------
1433xx=double(get(handles.geometry_calib,'CurrentCharacter')); %get the keyboard character
1434if ismember(xx,[8 127])%backspace or delete
1435    Coord_cell=get(handles.ListCoord,'String');
1436    val=get(handles.ListCoord,'Value');
1437     if max(val)<numel(Coord_cell) % the last element '...' has not been selected
1438        Coord_cell(val)=[];%remove the selected line
1439        set(handles.ListCoord,'Value',min(val))
1440        set(handles.ListCoord,'String',Coord_cell)         
1441        ListCoord_Callback(hObject, eventdata, handles)
1442        MenuPlot_Callback(hObject,eventdata,handles)
1443     end
1444end
1445
1446%------------------------------------------------------------------------
1447function fileinput=browse_xml(hObject, eventdata, handles)
1448%------------------------------------------------------------------------
1449fileinput=[];%default
1450oldfile=''; %default
1451UserData=get(handles.geometry_calib,'UserData');
1452if isfield(UserData,'XmlInputFile')
1453    oldfile=UserData.XmlInputFile;
1454end
1455[FileName, PathName, filterindex] = uigetfile( ...
1456       {'*.xml;*.mat', ' (*.xml,*.mat)';
1457       '*.xml',  '.xml files '; ...
1458        '*.mat',  '.mat matlab files '}, ...
1459        'Pick a file',oldfile);
1460fileinput=[PathName FileName];%complete file name
1461testblank=findstr(fileinput,' ');%look for blanks
1462if ~isempty(testblank)
1463    msgbox_uvmat('ERROR','forbidden input file name or path: no blank character allowed')
1464    return
1465end
1466sizf=size(fileinput);
1467if (~ischar(fileinput)||~isequal(sizf(1),1)),return;end
1468UserData.XmlInputFile=fileinput;
1469set(handles.geometry_calib,'UserData',UserData)%record current file foer further use of browser
1470
1471% -----------------------------------------------------------------------
1472function Heading=loadfile(handles,fileinput)
1473%------------------------------------------------------------------------
1474Heading=[];%default
1475[s,errormsg]=imadoc2struct(fileinput,'GeometryCalib');
1476if ~isempty(errormsg)
1477    msgbox_uvmat('ERROR',['Error for reading ' fileinput ': '  errormsg])
1478    return
1479end
1480if ~isempty(s.Heading)
1481    Heading=s.Heading;
1482end
1483   
1484GeometryCalib=s.GeometryCalib;
1485fx=1;fy=1;Cx=0;Cy=0;kc=0; %default
1486%     Tabchar={};
1487CoordCell={};
1488%     kc=0;%default
1489%     f1=1000;
1490%     f2=1000;
1491%     hhuvmat=guidata(findobj(allchild(0),'Name','uvmat'));
1492%     Cx=str2num(get(hhuvmat.num_Npx,'String'))/2;
1493%     Cy=str2num(get(hhuvmat.num_Npy,'String'))/2;
1494Tabchar={};%default
1495val_cal=1;%default
1496if ~isempty(GeometryCalib)
1497    % choose the calibration option
1498    if isfield(GeometryCalib,'CalibrationType')
1499       calib_list=get(handles.calib_type,'String');
1500       for ilist=1:numel(calib_list)
1501           if strcmp(calib_list{ilist},GeometryCalib.CalibrationType)
1502               val_cal=ilist;
1503               break
1504           end
1505       end
1506    end
1507    display_intrinsic(GeometryCalib,handles)%intrinsic param
1508    %extrinsic param
1509    if isfield(GeometryCalib,'Tx_Ty_Tz')
1510        Tx_Ty_Tz=GeometryCalib.Tx_Ty_Tz;
1511        set(handles.Tx,'String',num2str(GeometryCalib.Tx_Ty_Tz(1),4))
1512        set(handles.Ty,'String',num2str(GeometryCalib.Tx_Ty_Tz(2),4))
1513        set(handles.Tz,'String',num2str(GeometryCalib.Tx_Ty_Tz(3),4))
1514    end
1515    if isfield(GeometryCalib,'omc')
1516        set(handles.Phi,'String',num2str(GeometryCalib.omc(1),4))
1517        set(handles.Theta,'String',num2str(GeometryCalib.omc(2),4))
1518        set(handles.Psi,'String',num2str(GeometryCalib.omc(3),4))
1519    end
1520    calib=reshape(GeometryCalib.PointCoord,[],1);
1521    for ilist=1:numel(calib)
1522        CoordCell{ilist}=num2str(calib(ilist));
1523    end
1524    CoordCell=reshape(CoordCell,[],5);
1525    Tabchar=cell2tab(CoordCell,' | ');%transform cells into table ready for display
1526    MenuPlot_Callback(handles.geometry_calib, [], handles)
1527end
1528set(handles.calib_type,'Value',val_cal)
1529Tabchar=[Tabchar;{'......'}];
1530set(handles.ListCoord,'Value',1)
1531set(handles.ListCoord,'String',Tabchar)
1532
1533if isempty(CoordCell)% allow mouse action by default in the absence of input points
1534    set(handles.edit_append,'Value',1)
1535    set(handles.edit_append,'BackgroundColor',[1 1 0])
1536else % does not allow mouse action by default in the presence of input points
1537    set(handles.edit_append,'Value',0)
1538    set(handles.edit_append,'BackgroundColor',[0.7 0.7 0.7])
1539end
1540
1541%------------------------------------------------------------------------
1542%---display calibration intrinsic parameters
1543function display_intrinsic(GeometryCalib,handles)
1544%------------------------------------------------------------------------
1545fx=[];
1546fy=[];
1547if isfield(GeometryCalib,'fx_fy')
1548    fx=GeometryCalib.fx_fy(1);
1549    fy=GeometryCalib.fx_fy(2);
1550end
1551Cx_Cy=[0 0];%default
1552if isfield(GeometryCalib,'Cx_Cy')
1553    Cx_Cy=GeometryCalib.Cx_Cy;
1554end
1555kc=0;
1556if isfield(GeometryCalib,'kc')
1557    kc=GeometryCalib.kc; %* GeometryCalib.focal*GeometryCalib.focal;
1558end
1559set(handles.fx,'String',num2str(fx,5))
1560set(handles.fy,'String',num2str(fy,5))
1561set(handles.Cx,'String',num2str(Cx_Cy(1),'%1.1f'))
1562set(handles.Cy,'String',num2str(Cx_Cy(2),'%1.1f'))
1563set(handles.kc,'String',num2str(kc,'%1.4f'))
1564
1565
Note: See TracBrowser for help on using the repository browser.