source: trunk/src/geometry_calib.m @ 177

Last change on this file since 177 was 177, checked in by sommeria, 10 years ago

bug with mouse object editing resolved. Display feature 'satus' for PIV task advancement introduced. Various bug repair and cleaning

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