source: trunk/src/read_field.m @ 1042

Last change on this file since 1042 was 1041, checked in by sommeria, 6 years ago

psmn parameters updated + miscenaleous
File size: 15.1 KB
Line 
1%'read_field': read the fields from files in different formats (netcdf files, images, video)
2%--------------------------------------------------------------------------
3%  function [Field,ParamOut,errormsg] = read_field(FileName,FileType,ParamIn,num)
4%
5% OUTPUT:
6% Field: matlab structure representing the field
7% ParamOut: structure representing parameters:
8%        .FieldName; field name
9%        .VelType
10%        .CivStage: stage of civx processing (=0, not Civx, =1 (civ1), =2  (fix1)....     
11%        .Npx,.Npy: for images, nbre of pixels in x and y
12% errormsg: error message, ='' by default
13%
14%INPUT
15% FileName: name of the input file
16% FileType: type of file, as determined by the function get_file_info.m
17% ParamIn: movie object or Matlab structure of input parameters
18%     .FieldName: name (char string) of the input field (for Civx data)
19%     .VelType: char string giving the type of velocity data ('civ1', 'filter1', 'civ2'...)
20%     .ColorVar: variable used for vector color
21%     .Npx, .Npy: nbre of pixels along x and y (used for .vol input files)
22%     .TimeDimName: name of the dimension considered as 'time', selected index value then set by input 'num'   
23% num: frame number for movies
24%
25% see also read_image.m,read_civxdata.m,read_civdata.m,
26
27%=======================================================================
28% Copyright 2008-2018, LEGI UMR 5519 / CNRS UGA G-INP, Grenoble, France
29%   http://www.legi.grenoble-inp.fr
30%   Joel.Sommeria - Joel.Sommeria (A) legi.cnrs.fr
31%
32%     This file is part of the toolbox UVMAT.
33%
34%     UVMAT is free software; you can redistribute it and/or modify
35%     it under the terms of the GNU General Public License as published
36%     by the Free Software Foundation; either version 2 of the license,
37%     or (at your option) any later version.
38%
39%     UVMAT is distributed in the hope that it will be useful,
40%     but WITHOUT ANY WARRANTY; without even the implied warranty of
41%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
42%     GNU General Public License (see LICENSE.txt) for more details.
43%=======================================================================
44
45function [Field,ParamOut,errormsg] = read_field(FileName,FileType,ParamIn,num)
46%% default output and check input
47Field=[];
48if ~exist('num','var')
49    num=1;
50end
51if isempty(num)
52    num=1;
53end
54if ~exist('ParamIn','var')
55    ParamIn=[];
56end
57ParamOut=ParamIn;%default
58errormsg='';
59if isempty(regexp(FileName,'^http://'))&& ~exist(FileName,'file')
60    errormsg=['input file ' FileName ' does not exist'];
61    return
62end
63A=[];
64InputField={};
65check_colorvar=0;
66if isstruct(ParamIn)
67    if isfield(ParamIn,'FieldName')
68        if ischar(ParamIn.FieldName)
69            InputField={ParamIn.FieldName};
70        else
71            InputField= ParamIn.FieldName;
72        end
73    end
74    check_colorvar=zeros(size(InputField));
75    if isfield(ParamIn,'ColorVar')&&~isempty(ParamIn.ColorVar)
76        InputField=[ParamIn.FieldName {ParamIn.ColorVar}];
77        check_colorvar(numel(InputField))=1;
78    end
79end
80
81%% distingush different input file types
82switch FileType
83    case 'civdata'% new format for civ results
84        [Field,ParamOut.VelType,errormsg]=read_civdata(FileName,InputField,ParamIn.VelType);
85        if ~isempty(errormsg),errormsg=['read_civdata / ' errormsg];return,end
86        ParamOut.CivStage=Field.CivStage;
87    case 'pivdata_fluidimage'
88        [Field,ParamOut.VelType,errormsg]=read_pivdata_fluidimage(FileName,InputField,ParamIn.VelType);
89        ParamOut.CivStage=Field.CivStage;
90    case 'civx'% old (obsolete) format for civ results
91        ParamOut.FieldName='velocity';%Civx data found, set .FieldName='velocity' by default
92        [Field,ParamOut.VelType,errormsg]=read_civxdata(FileName,InputField,ParamIn.VelType);
93        if ~isempty(errormsg),errormsg=['read_civxdata / ' errormsg];return,end
94        ParamOut.CivStage=Field.CivStage;
95    case 'netcdf'% general netcdf file (not recognized as civ)
96        ListVar={};
97        Role={};
98        ProjModeRequest={};
99        ListInputField={};
100        ListOperator={};
101        checkU=0;
102        checkV=0;
103        % scan the list InputField
104        for ilist=1:numel(InputField)
105            % look for input variables to read
106            r=regexp(InputField{ilist},'(?<Operator>(^vec|^norm))\((?<UName>.+),(?<VName>.+)\)$','names');
107            if isempty(r)%  no operator used
108                if isempty(find(strcmp(InputField{ilist},ListVar),1))
109                    ListVar=[ListVar InputField(ilist)];%append the variable name if not already in the list
110                    ListInputField=[ListInputField InputField(ilist)];
111                    ListOperator=[ListOperator {''}];
112                end
113                if check_colorvar(ilist)
114                    if isempty(find(strcmp(InputField{ilist},ListVar),1))
115                    Role{numel(ListVar)}='ancillary';% not projected with interpolation
116                    ProjModeRequest{numel(ListVar)}='';
117                    end
118                else
119                    Role{numel(ListVar)}='scalar';
120                    ProjModeRequest{numel(ListVar)}='interp_lin';%scalar field (requires interpolation for plot)
121                end
122                if isfield(ParamIn,'Coord_y')
123                    if ~isempty(strcmp(InputField{ilist},ParamIn.Coord_y))
124                        Role{numel(ListVar)}='coord_y';
125                    end
126                end
127            else  % an operator 'vec' or 'norm' is used
128                if ~check_colorvar(ilist) && strcmp(r.Operator,'norm')
129                    ProjModeRequestVar='interp_lin';%scalar field (requires interpolation for plot)
130                else
131                    ProjModeRequestVar='';
132                end
133                ind_var_U=find(strcmp(r.UName,ListVar));%check previous listing of variable r.UName
134                ind_var_V=find(strcmp(r.VName,ListVar));%check previous listing of variable r.VName
135                if isempty(ind_var_U)
136                    ListVar=[ListVar {r.UName}]; % append the variable in the list if not previously listed
137                    Role=[Role {'vector_x'}];
138                    ProjModeRequest=[ProjModeRequest {ProjModeRequestVar}];
139                    ListInputField=[ListInputField InputField(ilist)];
140                else
141                    checkU=1;
142                end
143                if isempty(ind_var_V)
144                    ListVar=[ListVar {r.VName}];% append the variable in the list if not previously listed
145                    Role=[Role {'vector_y'}];
146                    ProjModeRequest=[ProjModeRequest {ProjModeRequestVar}];
147                    ListInputField=[ListInputField InputField(ilist)];               
148                else
149                    checkV=1;
150                end
151            end
152        end
153        if ~isfield(ParamIn,'Coord_z')
154            ParamIn.Coord_z=[];
155        end
156        NbCoord=~isempty(ParamIn.Coord_x)+~isempty(ParamIn.Coord_y)+~isempty(ParamIn.Coord_z);
157        if isfield(ParamIn,'TimeDimName')% case of reading of a single time index in a multidimensional array
158            [Field,var_detect,ichoice,errormsg]=nc2struct(FileName,'TimeDimName',ParamIn.TimeDimName,num,[ParamIn.Coord_x ParamIn.Coord_y ParamIn.Coord_z ListVar]);
159        elseif isfield(ParamIn,'TimeVarName')% case of reading of a single time  in a multidimensional array
160            [Field,var_detect,ichoice,errormsg]=nc2struct(FileName,'TimeVarName',ParamIn.TimeVarName,num,[ParamIn.Coord_x ParamIn.Coord_y ParamIn.Coord_z ListVar]);
161            if numel(num)~=1
162                NbCoord=NbCoord+1;% adds time coordinate, except if a single time has been selected
163            end
164        else
165            [Field,var_detect,ichoice,errormsg]=nc2struct(FileName,[ParamIn.Coord_x ParamIn.Coord_y ParamIn.Coord_z ListVar]);
166        end
167        if ~isempty(errormsg)
168            return
169        end
170        CheckStructured=1;
171        %scan all the variables beyond the two first NbCoord ones describing the coordinates.
172        for ilist=NbCoord+1:numel(Field.VarDimName)
173            if isequal(Field.VarDimName{1},Field.VarDimName{ilist}) % if a variable has the same dimension as the coordinate, it denotes a field with unstructured coordinates
174                Field.VarAttribute{1}.Role='coord_x';%unstructured coordinates
175                Field.VarAttribute{2}.Role='coord_y';
176                if NbCoord>=3
177                    Field.VarAttribute{3}.Role='coord_z';
178                end
179                CheckUnstructured=0;
180                break
181            end
182        end
183        if CheckStructured
184            for ilist=NbCoord+1:numel(Field.VarDimName)
185                if numel(Field.VarDimName{ilist})==NbCoord
186                    rank(1)=find(strcmp(ParamIn.Coord_x,Field.VarDimName{ilist}));
187                    rank(2)=find(strcmp(ParamIn.Coord_y,Field.VarDimName{ilist}));
188                    if NbCoord==3
189                        rank(3)=find(strcmp(ParamIn.Coord_z,Field.VarDimName{ilist}));
190                    end
191                    rank=rank(end:-1:1);
192                    VarName=Field.ListVarName{ilist};
193                    Field.(VarName)=permute(Field.(VarName),rank);
194                    Field.VarDimName{ilist}=Field.VarDimName{ilist}(rank);% permute the order of dimensions
195                end
196            end
197        end
198        NormName='';
199        UName='';
200        VName='';
201        if numel(Field.ListVarName)>NbCoord % if there are variables beyond coord (1 D plots)
202            for ilist=1:numel(ListVar)
203                Field.VarAttribute{ilist+NbCoord}.Role=Role{ilist};
204                Field.VarAttribute{ilist+NbCoord}.ProjModeRequest=ProjModeRequest{ilist};
205                if isfield(ParamIn,'FieldName')
206                    Field.VarAttribute{ilist+NbCoord}.FieldName=ListInputField{ilist};
207                end
208                r=regexp(ListInputField{ilist},'(?<Operator>(^vec|^norm))\((?<UName>.+),(?<VName>.+)\)$','names');
209                if ~isempty(r)&& strcmp(r.Operator,'norm')
210                    NormName='norm';
211                    if ~isempty(find(strcmp(ListVar,'norm')))
212                        NormName='norm_1';
213                    end
214                    Field.ListVarName=[Field.ListVarName {NormName}];
215                    ilistmax=numel(Field.ListVarName);
216                    Field.VarDimName{ilistmax}=Field.VarDimName{ilist+2};
217                    Field.VarAttribute{ilistmax}.Role='scalar';
218                    Field.(NormName)=Field.(r.UName).*Field.(r.UName)+Field.(r.VName).*Field.(r.VName);
219                    Field.(NormName)=sqrt(Field.(NormName));
220                    UName=r.UName;
221                    VName=r.VName;
222                end
223            end
224           
225            if ~isempty(NormName)% remove U and V if norm has been calculated and U and V are not needed as variables
226                ind_var_U=find(strcmp(UName,ListVar));%check previous listing of variable r.UName
227                ind_var_V=find(strcmp(VName,ListVar));%check previous listing of variable r.VName
228                if ~checkU && ~checkV
229                    Field.ListVarName([ind_var_U+2 ind_var_V+2])=[];
230                    Field.VarDimName([ind_var_U+2 ind_var_V+2])=[];
231                    Field.VarAttribute([ind_var_U+2 ind_var_V+2])=[];
232                elseif ~checkU
233                    Field.ListVarName(ind_var_U+2)=[];
234                    Field.VarDimName(ind_var_U+2)=[];
235                    Field.VarAttribute(ind_var_U+2 )=[];
236                elseif ~checkV
237                    Field.ListVarName(ind_var_V+2)=[];
238                    Field.VarDimName(ind_var_V+2)=[];
239                    Field.VarAttribute(ind_var_V+2 )=[];
240                end
241            end
242            % insert coordinates as indices in case of plots vs matrix index
243            if isfield(ParamIn,'CheckCoordIndex') && ParamIn.CheckCoordIndex
244                Field.ListVarName=[Field.ListDimName Field.ListVarName];
245                Field.VarDimName=[Field.ListDimName Field.VarDimName];
246                for idim=1:numel(Field.ListDimName)
247                    CoordName=Field.ListDimName{idim};
248                    Field.(CoordName)=1:Field.DimValue(idim);
249                end
250                Field.VarAttribute=[cell(1,numel(Field.ListDimName)) Field.VarAttribute];
251            end
252           
253        end
254    case 'video'
255        if strcmp(class(ParamIn),'VideoReader')
256            A=read(ParamIn,num);
257        else
258            ParamOut=VideoReader(FileName);
259            A=read(ParamOut,num);
260        end
261    case 'mmreader'
262        if strcmp(class(ParamIn),'mmreader')
263            A=read(ParamIn,num);
264        else
265            ParamOut=mmreader(FileName);
266            A=read(ParamOut,num);
267        end
268    case 'vol'
269        A=imread(FileName);
270        Npz=size(A,1)/ParamIn.Npy;
271        A=reshape(A',ParamIn.Npx,ParamIn.Npy,Npz);
272        A=permute(A,[3 2 1]);
273    case 'multimage'
274        %  warning 'off'
275        A=imread(FileName,num);
276    case 'image'
277        A=imread(FileName);
278    case 'rdvision'
279        [A,FileInfo,timestamps]=read_rdvision(FileName,num);
280    case 'image_DaVis'
281        Input=readimx(FileName);
282        if numel(Input.Frames)==1
283            num=1;
284        end
285        A=Input.Frames{num}.Components{1}.Planes{1}';
286        for ilist=1:numel(Input.Frames{1}.Attributes)
287            if strcmp(Input.Frames{1}.Attributes{ilist}.Name,'AcqTimeSeries')
288        timestamps=str2num(Input.Frames{1}.Attributes{ilist}.Value(1:end-3))/1000000;
289        break
290            end
291        end
292    case 'cine_phantom'
293        [A,FileInfo] = read_cine_phantom(FileName,num );
294    otherwise
295        errormsg=[ FileType ': invalid input file type for uvmat'];
296end
297
298%% case of image
299if ~isempty(A)
300    if strcmp(FileType,'rdvision')||strcmp(FileType,'image_DaVis')
301        Field.Time=timestamps;
302    end
303    if isstruct(ParamOut)
304        ParamOut.FieldName='image';
305    end
306    Npz=1;%default
307    npxy=size(A);
308    Field.NbDim=2;%default
309    Field.AName='image';
310    Field.ListVarName={'Coord_y','Coord_x','A'}; %
311    Field.VarAttribute{1}.Unit='pixel';
312    Field.VarAttribute{2}.Unit='pixel';
313    Field.VarAttribute{3}.Role='scalar';
314    if ndims(A)==3
315        if Npz==1;%color
316            Field.VarDimName={'Coord_y','Coord_x',{'Coord_y','Coord_x','rgb'}}; %
317            Field.Coord_y=[npxy(1)-0.5 0.5];
318            Field.Coord_x=[0.5 npxy(2)-0.5]; % coordinates of the first and last pixel centers
319            if isstruct(ParamOut)
320                ParamOut.Npx=npxy(2);% display image size on the interface
321                ParamOut.Npy=npxy(1);
322            end
323        else
324            Field.NbDim=3;
325            Field.ListVarName=['AZ' Field.ListVarName];
326            Field.VarDimName={'AZ','Coord_y','Coord_x',{'AZ','Coord_y','Coord_x'}};
327            Field.AZ=[npxy(1)-0.5 0.5];
328            Field.Coord_y=[npxy(2)-0.5 0.5];
329            Field.Coord_x=[0.5 npxy(3)-0.5]; % coordinates of the first and last pixel centers
330            if isstruct(ParamOut)
331                ParamOut.Npx=npxy(3);% display image size on the interface
332                ParamOut.Npy=npxy(2);
333            end
334        end
335    else
336        Field.VarDimName={'Coord_y','Coord_x',{'Coord_y','Coord_x'}}; %
337        Field.Coord_y=[npxy(1)-0.5 0.5];
338        Field.Coord_x=[0.5 npxy(2)-0.5]; % coordinates of the first and last pixel centers
339    end
340    Field.A=A;
341    Field.CoordUnit='pixel'; %used for mouse_motion
342end
343
344
345
Note: See TracBrowser for help on using the repository browser.