source: trunk/src/read_civxdata.m @ 206

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

bug fixes to deal with volumes, storage of ACTION menu in series fixed

File size: 10.6 KB
Line 
1%'read_civxdata': reads civx data from netcdf files
2%------------------------------------------------------------------
3%
4% function [Field,VelTypeOut]=read_civxdata(filename,FieldNames,VelType)
5%
6% OUTPUT:
7% Field: structure representing the selected field, containing
8%            .Txt: (char string) error message if any
9%            .ListGlobalAttribute: list of global attributes containing:
10%                    .NbCoord: number of vector components
11%                    .NbDim: number of dimensions (=2 or 3)
12%                    .dt: time interval for the corresponding image pair
13%                    .Time: absolute time (average of the initial image pair)
14%                    .CivStage: =0,
15%                       =1, civ1 has been performed only
16%                       =2, fix1 has been performed
17%                       =3, pacth1 has been performed
18%                       =4, civ2 has been performed
19%                       =5, fix2 has been performed
20%                       =6, pacth2 has been performed
21%                     .CoordUnit: 'pixel'
22%             .ListVarName: {'X'  'Y'  'U'  'V'  'F'  'FF'}
23%                   .X, .Y, .Z: set of vector coordinates
24%                    .U,.V,.W: corresponding set of vector components
25%                    .F: warning flags
26%                    .FF: error flag, =0 for good vectors
27%                     .C: scalar associated with velocity (used for vector colors)
28%                    .DijU; matrix of spatial derivatives (DijU(1,1,:)=DUDX,
29%                        DijU(1,2,:)=DUDY, Dij(2,1,:)=DVDX, DijU(2,2,:)=DVDY
30%
31% VelTypeOut: velocity type corresponding to the selected field: ='civ1','interp1','interp2','civ2'....
32%
33% INPUT:
34% filename: file name (string).
35% FieldNames =cell of field names to get, which can contain the strings:
36%             'ima_cor': image correlation, vec_c or vec2_C
37%             'vort','div','strain': requires velocity derivatives DUDX...
38%             'error': error estimate (vec_E or vec2_E)
39%             
40% VelType : character string indicating the types of velocity fields to read ('civ1','civ2'...)
41%            if vel_type=[] or'*', a  priority choice, given by vel_type_out{1,2}, is done depending
42%            if vel_type='filter'; a structured field is sought (filter2 in priority, then filter1)
43
44
45% FUNCTIONS called:
46% 'varcivx_generator':, sets the names of vaiables to read in the netcdf file
47% 'nc2struct': reads a netcdf file
48
49function [Field,VelTypeOut]=read_civxdata(filename,FieldNames,VelType)
50%% default input
51if ~exist('VelType','var')
52    VelType=[];
53end
54if isequal(VelType,'*')
55    VelType=[];
56end
57if ~exist('FieldNames','var')
58    FieldNames=[]; %default
59end
60
61%% reading data
62VelTypeOut=VelType;%default
63[var,role,units,vel_type_out_cell]=varcivx_generator(FieldNames,VelType);%determine the names of constants and variables to read
64[Field,vardetect,ichoice]=nc2struct(filename,var);%read the variables in the netcdf file
65if isfield(Field,'Txt')
66    return
67end
68if vardetect(1)==0
69     Field.Txt=[ 'requested field not available in ' filename '/' VelType];
70     return
71end
72var_ind=find(vardetect);
73for ivar=1:length(var_ind)
74    Field.VarAttribute{ivar}.Role=role{var_ind(ivar)};
75    Field.VarAttribute{ivar}.Mesh=0.1;%typical mesh for histograms O.1 pixel
76end
77VelTypeOut=VelType;
78if ~isempty(ichoice)
79    VelTypeOut=vel_type_out_cell{ichoice};
80end
81
82%% adjust for Djui:
83if isfield(Field,'DjUi')
84    Field.ListVarName(end-2:end)=[];
85    Field.ListVarName{end}='DjUi';
86    Field.VarDimName(end-2:end)=[];
87    Field.VarAttribute(end-2:end)=[];
88end
89
90%% renaming for standard conventions
91Field.NbCoord=Field.nb_coord;
92Field.NbDim=Field.nb_dim;
93
94%% CivStage
95if isfield(Field,'patch2')&& isequal(Field.patch2,1)
96    Field.CivStage=6;
97elseif isfield(Field,'fix2')&& isequal(Field.fix2,1)
98    Field.CivStage=5;
99elseif isfield(Field,'civ2')&& isequal(Field.civ2,1)
100    Field.CivStage=4;
101elseif isfield(Field,'patch')&& isequal(Field.patch,1)
102    Field.CivStage=3;
103elseif isfield(Field,'fix')&& isequal(Field.fix,1)
104    Field.CivStage=2;
105else
106    Field.CivStage=1;
107end
108
109%determine the appropriate constant for time and dt for the PIV pair
110test_civ1=isequal(VelTypeOut,'civ1')||isequal(VelTypeOut,'interp1')||isequal(VelTypeOut,'filter1');
111test_civ2=isequal(VelTypeOut,'civ2')||isequal(VelTypeOut,'interp2')||isequal(VelTypeOut,'filter2');
112Field.Time=0; %default
113if test_civ1
114    if isfield(Field,'absolut_time_T0')
115        Field.Time=double(Field.absolut_time_T0);
116        Field.dt=double(Field.dt);
117    else
118       Field.Txt='the input file is not civx';
119       Field.CivStage=0;
120       Field.dt=0;
121    end
122elseif test_civ2
123    Field.Time=double(Field.absolut_time_T0_2);
124    Field.dt=double(Field.dt2);
125else
126    Field.Txt='the input file is not civx';
127    Field.CivStage=0;
128    Field.dt=0;
129end
130
131
132
133%% rescale fields to pixel coordinates
134if isfield(Field,'pixcmx')
135    Field.pixcmx=double(Field.pixcmx);
136    Field.pixcmy=double(Field.pixcmy);
137    Field.U=Field.U*Field.pixcmx;
138    Field.V=Field.V*Field.pixcmy;
139    Field.X=Field.X*Field.pixcmx;
140    Field.Y=Field.Y*Field.pixcmy;
141end
142if ~isequal(Field.dt,0)
143    Field.U=Field.U*Field.dt;%translate in px displacement
144    Field.V=Field.V*Field.dt;
145    if isfield(Field,'DjUi')
146       Field.DjUi(:,1,1)=Field.dt*Field.DjUi(:,1,1);
147       Field.DjUi(:,2,2)=Field.dt*Field.DjUi(:,2,2);
148       Field.DjUi(:,1,2)=(Field.pixcmy/Field.pixcmx)*Field.dt*Field.DjUi(:,1,2);
149       Field.DjUi(:,2,1)=(Field.pixcmx/Field.pixcmy)*Field.dt*Field.DjUi(:,2,1);
150    end
151end
152
153%% update list of global attributes
154List=Field.ListGlobalAttribute;
155ind_remove=[];
156for ilist=1:length(List)
157    switch(List{ilist})
158        case {'patch2','fix2','civ2','patch','fix','dt2','absolut_time_T0','absolut_time_T0_2','nb_coord','nb_dim','pixcmx','pixcmy'}
159            ind_remove=[ind_remove ilist];
160            Field=rmfield(Field,List{ilist});
161    end
162end
163List(ind_remove)=[];
164Field.ListGlobalAttribute=[{'NbCoord'},{'NbDim'} List {'Time','CivStage','CoordUnit'}];
165Field.CoordUnit='pixel';
166
167
168%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
169% [var,role,units,vel_type_out]=varcivx_generator(FieldNames,vel_type)
170%INPUT:
171% FieldNames =cell of field names to get, which can contain the strings:
172%             'ima_cor': image correlation, vec_c or vec2_C
173%             'vort','div','strain': requires velocity derivatives DUDX...
174%             'error': error estimate (vec_E or vec2_E)
175%             
176% vel_type: character string indicating the types of velocity fields to read ('civ1','civ2'...)
177%            if vel_type=[] or'*', a  priority choice, given by vel_type_out{1,2}, is done depending
178%            if vel_type='filter'; a structured field is sought (filter2 in priority, then filter1)
179
180function [var,role,units,vel_type_out]=varcivx_generator(FieldNames,vel_type)
181
182%% default input values
183if ~exist('vel_type','var'),vel_type=[];end;
184if iscell(vel_type),vel_type=vel_type{1}; end;%transform cell to string if needed
185if ~exist('FieldNames','var'),FieldNames={'ima_cor'};end;%default scalar
186if ischar(FieldNames), FieldNames={FieldNames}; end;
187
188%% select the priority order for automatic vel_type selection
189testder=0;
190for ilist=1:length(FieldNames)
191    if ~isempty(FieldNames{ilist})
192    switch FieldNames{ilist}
193        case {'vort','div','strain'}
194            testder=1;
195    end
196    end
197end     
198if isempty(vel_type) || isequal(vel_type,'*') %undefined velocity type (civ1,civ2...)
199    if testder
200         vel_type_out{1}='filter2'; %priority to filter2 for scalar reading, filter1 as second
201        vel_type_out{2}='filter1';
202    else
203        vel_type_out{1}='civ2'; %priority to civ2 for vector reading, civ1 as second priority     
204        vel_type_out{2}='civ1';
205    end
206elseif isequal(vel_type,'filter')
207        vel_type_out{1}='filter2'; %priority to filter2 for scalar reading, filter1 as second
208        vel_type_out{2}='filter1';
209        if ~testder
210            vel_type_out{3}='civ1';%civ1 as third priority if derivatives are not needed
211        end
212elseif testder
213    test_civ1=isequal(vel_type,'civ1')||isequal(vel_type,'interp1')||isequal(vel_type,'filter1');
214    if test_civ1
215        vel_type_out{1}='filter1'; %switch to filter for reading spatial derivatives
216    else
217        vel_type_out{1}='filter2';
218    end
219else   
220    vel_type_out{1}=vel_type;%imposed velocity field
221end
222vel_type_out=vel_type_out';
223
224%% determine names of netcdf variables to read
225var={'X','Y','Z','U','V','W','C','F','FF'};
226role={'coord_x','coord_y','coord_z','vector_x','vector_y','vector_z','ancillary','warnflag','errorflag'};
227units={'pixel','pixel','pixel','pixel','pixel','pixel',[],[],[]};
228if testder
229    var=[var {'DjUi(:,1,1)','DjUi(:,1,2)','DjUi(:,2,1)','DjUi(:,2,2)'}];
230    role=[role {'tensor','tensor','tensor','tensor'}];
231    units=[units {'pixel','pixel','pixel','pixel'}];
232end
233for ilist=1:length(vel_type_out)
234    var=[var;varname1(vel_type_out{ilist},FieldNames)];
235end
236
237%------------------------------------------------------------------------ 
238%--- determine  var names to read
239function varin=varname1(vel_type,FieldNames)
240%------------------------------------------------------------------------
241testder=0;
242C1='';
243C2='';
244for ilist=1:length(FieldNames)
245    if ~isempty(FieldNames{ilist})
246    switch FieldNames{ilist}
247        case 'ima_cor' %image correlation corresponding to a vel vector
248            C1='vec_C';
249            C2='vec2_C';
250        case 'error'
251            C1='vec_E';
252            C2='vec2_E';
253        case {'vort','div','strain'}
254            testder=1;
255    end
256    end
257end     
258switch vel_type
259    case 'civ1'
260        varin={'vec_X','vec_Y','vec_Z','vec_U','vec_V','vec_W',C1,'vec_F','vec_FixFlag'};
261    case 'interp1'
262        varin={'vec_patch_X','vec_patch_Y','','vec_patch0_U','vec_patch0_V','','','',''};
263    case 'filter1'
264        varin={'vec_patch_X','vec_patch_Y','','vec_patch_U','vec_patch_V','','','',''};
265    case 'civ2'
266        varin={'vec2_X','vec2_Y','vec2_Z','vec2_U','vec2_V','vec2_W',C2,'vec2_F','vec2_FixFlag'};
267    case 'interp2'
268        varin={'vec2_patch_X','vec2_patch_Y','vec2_patch_Z','vec2_patch0_U','vec2_patch0_V','vec2_patch0_W','','',''};
269    case 'filter2'
270        varin={'vec2_patch_X','vec2_patch_Y','vec2_patch_Z','vec2_patch_U','vec2_patch_V','vec2_patch0_W','','',''};
271end
272if testder
273     switch vel_type
274        case 'filter1'
275            varin=[varin {'vec_patch_DUDX','vec_patch_DVDX','vec_patch_DUDY','vec_patch_DVDY'}];
276        case 'filter2'
277            varin=[varin {'vec2_patch_DUDX','vec2_patch_DVDX','vec2_patch_DUDY','vec2_patch_DVDY'}];
278    end   
279end
Note: See TracBrowser for help on using the repository browser.