source: trunk/src/geometry_calib.m @ 491

Last change on this file since 491 was 444, checked in by sommeria, 12 years ago

masks and calibration updated to fit with the new conventions on file organisation.
bug corrected in mouse_up (object creation)
bug corrected in civ in mode TESTciv

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