source: trunk/src/imadoc2struct.m @ 606

Last change on this file since 606 was 591, checked in by sommeria, 11 years ago

various updates, in particular modification of series to do calculations in the cluster

File size: 16.3 KB
Line 
1%'imadoc2struct': reads the xml file for image documentation
2%------------------------------------------------------------------------
3% function [s,errormsg]=imadoc2struct(ImaDoc,option)
4%
5% OUTPUT:
6% s: structure representing ImaDoc
7%   s.Heading: information about the data hierarchical structure
8%   s.Time: matrix of times
9%   s.TimeUnit
10%  s.GeometryCalib: substructure containing the parameters for geometric calibration
11% errormsg: error message
12%
13% INPUT:
14% ImaDoc: full name of the xml input file with head key ImaDoc
15% varargin: optional list of strings to restrict the reading to a selection of subtrees, for instance 'GeometryCalib' (save time)
16
17function [s,errormsg]=imadoc2struct(ImaDoc,varargin)
18%% default input and output
19errormsg='';%default
20s=[];
21% s.Heading=[];%default
22% s.Time=[]; %default
23% s.TimeUnit=[]; %default
24% s.GeometryCalib=[];
25% tsai=[];%default
26
27%% opening the xml file
28[tild,tild,FileExt]=fileparts(ImaDoc);
29%% case of .civ files (obsolete)
30if strcmp(FileExt,'.civ')
31    [errormsg,time,TimeUnit,mode,npx,npy,s.GeometryCalib]=read_imatext(ImaDoc);
32    return
33end
34
35%% case of xml files
36if nargin ==1
37    [s,Heading]=xml2struct(ImaDoc);% convert the whole xml file in a structure s
38elseif nargin ==2
39    [s,Heading]=xml2struct(ImaDoc,varargin{1});% convert the xml file in a structure s, keeping only the subtree defined in input
40else %TODO: deal with more than two subtrees?
41    [s,Heading]=xml2struct(ImaDoc,varargin{1},varargin{2});% convert the xml file in a structure s, keeping only the subtree defined in input
42end
43if ~strcmp(Heading,'ImaDoc')
44    errormsg='the input xml file is not ImaDoc';
45    return
46end
47%% reading timing
48if isfield(s,'Camera')
49    Timing=s.Camera.BurstTiming;
50    if ~iscell(Timing)
51        Timing={Timing};
52    end
53    s.Time=[];
54    for k=1:length(Timing)
55        Frequency=1;
56        if isfield(Timing{k},'FrameFrequency')
57            Frequency=Timing{k}.FrameFrequency;
58        end
59        Dtj=[];
60        if isfield(Timing{k},'Dtj')
61            Dtj=Timing{k}.Dtj/Frequency;%Dtj converted from frame unit to TimeUnit (e.g. 's');
62        end
63        NbDtj=1;
64        if isfield(Timing{k},'NbDtj')&&~isempty(Timing{k}.NbDtj)
65            NbDtj=Timing{k}.NbDtj;
66        end
67        Dti=[];
68        if isfield(Timing{k},'Dti')
69            Dti=Timing{k}.Dti/Frequency;%Dti converted from frame unit to TimeUnit (e.g. 's');
70        end
71        NbDti=1;
72        if isfield(Timing{k},'NbDti')&&~isempty(Timing{k}.NbDti)
73            NbDti=Timing{k}.NbDti;
74        end
75        Time_val=Timing{k}.Time;%time in TimeUnit
76        if ~isempty(Dti)
77            Dti=reshape(Dti'*ones(1,NbDti),NbDti*numel(Dti),1); %concatene Dti vector NbDti times
78            Time_val=[Time_val;Time_val(end)+cumsum(Dti)];%append the times defined by the intervals  Dti
79        end
80        if ~isempty(Dtj)
81            Dtj=reshape(Dtj'*ones(1,NbDtj),1,NbDtj*numel(Dtj)); %concatene Dtj vector NbDtj times
82            Dtj=[0 Dtj];
83            Time_val=Time_val*ones(1,numel(Dtj))+ones(numel(Time_val),1)*cumsum(Dtj);% produce a time matrix with Dtj
84        end
85        % reading Dtk
86        Dtk=[];%default
87        NbDtk=1;%default
88        if isfield(Timing{k},'Dtk')
89            Dtk=Timing{k}.Dtk;
90        end
91        if isfield(Timing{k},'NbDtk')&&~isempty(Timing{k}.NbDtk)
92            NbDtk=Timing{k}.NbDtk;
93        end
94        if isempty(Dtk)
95            s.Time=[s.Time;Time_val];
96        else
97            for kblock=1:NbDtk+1
98                Time_val_k=Time_val+(kblock-1)*Dtk;
99                s.Time=[s.Time;Time_val_k];
100            end
101        end
102    end
103end
104
105% try
106%     t=xmltree(ImaDoc);
107% catch ME
108%     errormsg={['error reading ' ImaDoc ': ']; ME.message};
109%     display(errormsg);
110%     return
111% end
112% uid_root=find(t,'/ImaDoc');
113% if isempty(uid_root), return; end;%not an ImaDoc .xml file
114
115% %% Heading
116% uid_Heading=find(t,'/ImaDoc/Heading');
117% if ~isempty(uid_Heading),
118%     uid_Campaign=find(t,'/ImaDoc/Heading/Campaign');
119%     uid_Exp=find(t,'/ImaDoc/Heading/Experiment');
120%     uid_Device=find(t,'/ImaDoc/Heading/Device');
121%     uid_Record=find(t,'/ImaDoc/Heading/Record');
122%     uid_FirstImage=find(t,'/ImaDoc/Heading/ImageName');
123%     s.Heading.Campaign=get(t,children(t,uid_Campaign),'value');
124%     s.Heading.Experiment=get(t,children(t,uid_Exp),'value');
125%     s.Heading.Device=get(t,children(t,uid_Device),'value');
126%     if ~isempty(uid_Record)
127%         s.Heading.Record=get(t,children(t,uid_Record),'value');
128%     end
129%     s.Heading.ImageName=get(t,children(t,uid_FirstImage),'value');
130% end
131
132%% Camera  and timing
133% if strcmp(option,'*') || strcmp(option,'Camera')
134%     uid_Camera=find(t,'/ImaDoc/Camera');
135%     if ~isempty(uid_Camera)
136%         uid_ImageSize=find(t,'/ImaDoc/Camera/ImageSize');
137%         if ~isempty(uid_ImageSize);
138%             ImageSize=get(t,children(t,uid_ImageSize),'value');
139%             xindex=findstr(ImageSize,'x');
140%             if length(xindex)>=2
141%                 s.Npx=str2double(ImageSize(1:xindex(1)-1));
142%                 s.Npy=str2double(ImageSize(xindex(1)+1:xindex(2)-1));
143%             end
144%         end
145%         uid_TimeUnit=find(t,'/ImaDoc/Camera/TimeUnit');
146%         if ~isempty(uid_TimeUnit)
147%             s.TimeUnit=get(t,children(t,uid_TimeUnit),'value');
148%         end
149%         uid_BurstTiming=find(t,'/ImaDoc/Camera/BurstTiming');
150%         if ~isempty(uid_BurstTiming)
151%             for k=1:length(uid_BurstTiming)
152%                 subt=branch(t,uid_BurstTiming(k));%subtree under BurstTiming
153%                 % reading Dtk
154%                 Frequency=get_value(subt,'/BurstTiming/FrameFrequency',1);
155%                 Dtj=get_value(subt,'/BurstTiming/Dtj',[]);
156%                 Dtj=Dtj/Frequency;%Dtj converted from frame unit to TimeUnit (e.g. 's')
157%                 NbDtj=get_value(subt,'/BurstTiming/NbDtj',1);
158%                 Dti=get_value(subt,'/BurstTiming/Dti',[]);
159%                 Dti=Dti/Frequency;%Dtj converted from frame unit to TimeUnit (e.g. 's')
160%                 NbDti=get_value(subt,'/BurstTiming/NbDti',1);
161%                 Time_val=get_value(subt,'/BurstTiming/Time',0);%time in TimeUnit
162%                 if ~isempty(Dti)
163%                     Dti=reshape(Dti'*ones(1,NbDti),NbDti*numel(Dti),1); %concatene Dti vector NbDti times
164%                     Time_val=[Time_val;Time_val(end)+cumsum(Dti)];%append the times defined by the intervals  Dti
165%                 end
166%                 if ~isempty(Dtj)
167%                     Dtj=reshape(Dtj'*ones(1,NbDtj),1,NbDtj*numel(Dtj)); %concatene Dtj vector NbDtj times
168%                     Dtj=[0 Dtj];
169%                     Time_val=Time_val*ones(1,numel(Dtj))+ones(numel(Time_val),1)*cumsum(Dtj);% produce a time matrix with Dtj
170%                 end
171%                 % reading Dtk
172%                 Dtk=get_value(subt,'/BurstTiming/Dtk',[]);
173%                 NbDtk=get_value(subt,'/BurstTiming/NbDtk',1);
174%                 if isempty(Dtk)
175%                     s.Time=[s.Time;Time_val];
176%                 else
177%                     for kblock=1:NbDtk+1
178%                         Time_val_k=Time_val+(kblock-1)*Dtk;
179%                         s.Time=[s.Time;Time_val_k];
180%                     end
181%                 end
182%             end
183%         end
184%     end
185% end
186
187%% motor
188% if strcmp(option,'*') || strcmp(option,'GeometryCalib')
189%     uid_subtree=find(t,'/ImaDoc/TranslationMotor');
190%     if length(uid_subtree)==1
191%         subt=branch(t,uid_subtree);%subtree under GeometryCalib
192%        [s.TranslationMotor,errormsg]=read_subtree(subt,{'Nbslice','ZStart','ZEnd'},[1 1 1],[1 1 1]);
193%     end
194% end
195%%  geometric calibration
196% if strcmp(option,'*') || strcmp(option,'GeometryCalib')
197%     uid_GeometryCalib=find(t,'/ImaDoc/GeometryCalib');
198%     if ~isempty(uid_GeometryCalib)
199%         if length(uid_GeometryCalib)>1
200%             errormsg=['More than one GeometryCalib in ' filecivxml];
201%             return
202%         end
203%         subt=branch(t,uid_GeometryCalib);%subtree under GeometryCalib
204%         cont=get(subt,1,'contents');
205%         if ~isempty(cont)
206%             uid_CalibrationType=find(subt,'/GeometryCalib/CalibrationType');
207%             if isequal(length(uid_CalibrationType),1)
208%                 tsai.CalibrationType=get(subt,children(subt,uid_CalibrationType),'value');
209%             end
210%             uid_CoordUnit=find(subt,'/GeometryCalib/CoordUnit');
211%             if isequal(length(uid_CoordUnit),1)
212%                 tsai.CoordUnit=get(subt,children(subt,uid_CoordUnit),'value');
213%             end
214%             uid_fx_fy=find(subt,'/GeometryCalib/fx_fy');
215%             focal=[];%default fro old convention (Reg Wilson)
216%             if isequal(length(uid_fx_fy),1)
217%                 tsai.fx_fy=str2num(get(subt,children(subt,uid_fx_fy),'value'));
218%             else %old convention (Reg Wilson)
219%                 uid_focal=find(subt,'/GeometryCalib/focal');
220%                 uid_dpx_dpy=find(subt,'/GeometryCalib/dpx_dpy');
221%                 uid_sx=find(subt,'/GeometryCalib/sx');
222%                 if ~isempty(uid_focal) && ~isempty(uid_dpx_dpy) && ~isempty(uid_sx)
223%                     dpx_dpy=str2num(get(subt,children(subt,uid_dpx_dpy),'value'));
224%                     sx=str2num(get(subt,children(subt,uid_sx),'value'));
225%                     focal=str2num(get(subt,children(subt,uid_focal),'value'));
226%                     tsai.fx_fy(1)=sx*focal/dpx_dpy(1);
227%                     tsai.fx_fy(2)=focal/dpx_dpy(2);
228%                 end
229%             end
230%             uid_Cx_Cy=find(subt,'/GeometryCalib/Cx_Cy');
231%             if ~isempty(uid_Cx_Cy)
232%                 tsai.Cx_Cy=str2num(get(subt,children(subt,uid_Cx_Cy),'value'));
233%             end
234%             uid_kc=find(subt,'/GeometryCalib/kc');
235%             if ~isempty(uid_kc)
236%                 tsai.kc=str2double(get(subt,children(subt,uid_kc),'value'));
237%             else %old convention (Reg Wilson)
238%                 uid_kappa1=find(subt,'/GeometryCalib/kappa1');
239%                 if ~isempty(uid_kappa1)&& ~isempty(focal)
240%                     kappa1=str2double(get(subt,children(subt,uid_kappa1),'value'));
241%                     tsai.kc=-kappa1*focal*focal;
242%                 end
243%             end
244%             uid_Tx_Ty_Tz=find(subt,'/GeometryCalib/Tx_Ty_Tz');
245%             if ~isempty(uid_Tx_Ty_Tz)
246%                 tsai.Tx_Ty_Tz=str2num(get(subt,children(subt,uid_Tx_Ty_Tz),'value'));
247%             end
248%             uid_R=find(subt,'/GeometryCalib/R');
249%             if ~isempty(uid_R)
250%                 RR=get(subt,children(subt,uid_R),'value');
251%                 if length(RR)==3
252%                     tsai.R=[str2num(RR{1});str2num(RR{2});str2num(RR{3})];
253%                 end
254%             end
255%             
256%             %look for laser plane definitions
257%             uid_Angle=find(subt,'/GeometryCalib/PlaneAngle');
258%             uid_Pos=find(subt,'/GeometryCalib/SliceCoord');
259%             if isempty(uid_Pos)
260%                 uid_Pos=find(subt,'/GeometryCalib/PlanePos');%old convention
261%             end
262%             if ~isempty(uid_Angle)
263%                 tsai.PlaneAngle=str2num(get(subt,children(subt,uid_Angle),'value'));
264%             end
265%             if ~isempty(uid_Pos)
266%                 for j=1:length(uid_Pos)
267%                     tsai.SliceCoord(j,:)=str2num(get(subt,children(subt,uid_Pos(j)),'value'));
268%                 end
269%                 uid_DZ=find(subt,'/GeometryCalib/SliceDZ');
270%                 uid_NbSlice=find(subt,'/GeometryCalib/NbSlice');
271%                 if ~isempty(uid_DZ) && ~isempty(uid_NbSlice)
272%                     DZ=str2double(get(subt,children(subt,uid_DZ),'value'));
273%                     NbSlice=get(subt,children(subt,uid_NbSlice),'value');
274%                     if isequal(NbSlice,'volume')
275%                         tsai.NbSlice='volume';
276%                         NbSlice=NbDtj+1;
277%                     else
278%                         tsai.NbSlice=str2double(NbSlice);
279%                     end
280%                     tsai.SliceCoord=ones(NbSlice,1)*tsai.SliceCoord+DZ*(0:NbSlice-1)'*[0 0 1];
281%                 end
282%             end   
283%             tsai.SliceAngle=get_value(subt,'/GeometryCalib/SliceAngle',[0 0 0]);
284%             tsai.VolumeScan=get_value(subt,'/GeometryCalib/VolumeScan','n');
285%             tsai.InterfaceCoord=get_value(subt,'/GeometryCalib/InterfaceCoord',[0 0 0]);
286%             tsai.RefractionIndex=get_value(subt,'/GeometryCalib/RefractionIndex',1);
287%             
288%             if strcmp(option,'GeometryCalib')
289%                 tsai.PointCoord=get_value(subt,'/GeometryCalib/SourceCalib/PointCoord',[0 0 0 0 0]);
290%             end
291%             s.GeometryCalib=tsai;
292%         end
293%     end
294% end
295
296%--------------------------------------------------
297%  read a subtree
298% INPUT:
299% t: xltree
300% head_element: head elelemnt of the subtree
301% Data, structure containing
302%    .Key: element name
303%    .Type: type of element ('charg', 'float'....)
304%    .NbOccur: nbre of occurrence, NaN for un specified number
305function [s,errormsg]=read_subtree(subt,Data,NbOccur,NumTest)
306%--------------------------------------------------
307s=[];%default
308errormsg='';
309head_element=get(subt,1,'name');
310    cont=get(subt,1,'contents');
311    if ~isempty(cont)
312        for ilist=1:length(Data)
313            uid_key=find(subt,[head_element '/' Data{ilist}]);
314            if ~isequal(length(uid_key),NbOccur(ilist))
315                errormsg=['wrong number of occurence for ' Data{ilist}];
316                return
317            end
318            for ival=1:length(uid_key)
319                val=get(subt,children(subt,uid_key(ival)),'value');
320                if ~NumTest(ilist)
321                    eval(['s.' Data{ilist} '=val;']);
322                else
323                    eval(['s.' Data{ilist} '=str2double(val);'])
324                end
325            end
326        end
327    end
328
329
330%--------------------------------------------------
331%  read an xml element
332function val=get_value(t,label,default)
333%--------------------------------------------------
334val=default;
335uid=find(t,label);%find the element iud(s)
336if ~isempty(uid) %if the element named label exists
337   uid_child=children(t,uid);%find the children
338   if ~isempty(uid_child)
339       data=get(t,uid_child,'type');%get the type of child
340       if iscell(data)% case of multiple element
341           for icell=1:numel(data)
342               val_read=str2num(get(t,uid_child(icell),'value'));
343               if ~isempty(val_read)
344                   val(icell,:)=val_read;
345               end
346           end
347%           val=val';
348       else % case of unique element value
349           val_read=str2num(get(t,uid_child,'value'));
350           if ~isempty(val_read)
351               val=val_read;
352           else
353              val=get(t,uid_child,'value');%char string data
354           end
355       end
356   end
357end
358
359%------------------------------------------------------------------------
360%'read_imatext': reads the .civ file for image documentation (obsolete)
361% fileinput: name of the documentation file
362% time: matrix of times for the set of images
363%pxcmx: scale along x in pixels/cm
364%pxcmy: scale along y in pixels/cm
365function [error,time,TimeUnit,mode,npx,npy,GeometryCalib]=read_imatext(fileinput)
366%------------------------------------------------------------------------
367error='';%default
368time=[]; %default
369TimeUnit='s';
370mode='pairs';
371npx=[]; %default
372npy=[]; %default
373GeometryCalib=[];
374if ~exist(fileinput,'file'), error=['image doc file ' fileinput ' does not exist']; return;end;%input file does not exist
375dotciv=textread(fileinput);
376sizdot=size(dotciv);
377if ~isequal(sizdot(1)-8,dotciv(1,1));
378    error=1; %inconsistent number of bursts
379end
380nbfield=sizdot(1)-8;
381npx=(dotciv(2,1));
382npy=(dotciv(2,2));
383pxcmx=(dotciv(6,1));% pixels/cm in the .civ file
384pxcmy=(dotciv(6,2));
385% nburst=dotciv(3,1); % nbre of bursts
386abs_time1=dotciv([9:nbfield+8],2);
387dtime=dotciv(5,1)*(dotciv([9:nbfield+8],[3:end-1])+1);
388timeshift=[abs_time1 dtime];
389time=cumsum(timeshift,2);
390GeometryCalib.CalibrationType='rescale';
391GeometryCalib.R=[pxcmx 0 0; 0 pxcmy 0;0 0 0];
392GeometryCalib.Tx=0;
393GeometryCalib.Ty=0;
394GeometryCalib.Tz=1;
395GeometryCalib.dpx=1;
396GeometryCalib.dpy=1;
397GeometryCalib.sx=1;
398GeometryCalib.Cx=0;
399GeometryCalib.Cy=0;
400GeometryCalib.f=1;
401GeometryCalib.kappa1=0;
402GeometryCalib.CoordUnit='cm';
Note: See TracBrowser for help on using the repository browser.