source: trunk/src/read_field.m @ 905

Last change on this file since 905 was 896, checked in by sommeria, 9 years ago

bug corrected in civ_series/patch

File size: 13.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-2014, LEGI UMR 5519 / CNRS UJF 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 ~exist('ParamIn','var')
52    ParamIn=[];
53end
54ParamOut=ParamIn;%default
55errormsg='';
56if ~exist(FileName,'file')
57    errormsg=['input file ' FileName ' does not exist'];
58    return
59end
60A=[];
61InputField={};
62check_colorvar=0;
63if isstruct(ParamIn)
64    if isfield(ParamIn,'FieldName')
65        if ischar(ParamIn.FieldName)
66            InputField={ParamIn.FieldName};
67        else
68            InputField= ParamIn.FieldName;
69        end
70    end
71    check_colorvar=zeros(size(InputField));
72    if isfield(ParamIn,'ColorVar')&&~isempty(ParamIn.ColorVar)
73        InputField=[ParamIn.FieldName {ParamIn.ColorVar}];
74        check_colorvar(numel(InputField))=1;
75    end
76end
77
78%% distingush different input file types
79switch FileType
80    case 'civdata'% new format for civ results
81        [Field,ParamOut.VelType,errormsg]=read_civdata(FileName,InputField,ParamIn.VelType);
82        if ~isempty(errormsg),errormsg=['read_civdata / ' errormsg];return,end     
83        ParamOut.CivStage=Field.CivStage;
84    case 'civx'% old (obsolete) format for civ results
85        ParamOut.FieldName='velocity';%Civx data found, set .FieldName='velocity' by default
86        [Field,ParamOut.VelType,errormsg]=read_civxdata(FileName,InputField,ParamIn.VelType);
87        if ~isempty(errormsg),errormsg=['read_civxdata / ' errormsg];return,end
88        ParamOut.CivStage=Field.CivStage;
89    case 'netcdf'% general netcdf file (not recognized as civ)
90        ListVar={};
91        Role={};
92        ProjModeRequest={};
93        ListInputField={};
94        ListOperator={};
95        checkU=0;
96        checkV=0;
97        for ilist=1:numel(InputField)
98            % look for input variables to read
99            r=regexp(InputField{ilist},'(?<Operator>(^vec|^norm))\((?<UName>.+),(?<VName>.+)\)$','names');
100            if isempty(r)%  no operator used
101                if isempty(find(strcmp(InputField{ilist},ListVar)))
102                    ListVar=[ListVar InputField(ilist)];%append the variable name if not already in the list
103                    ListInputField=[ListInputField InputField(ilist)];
104                    ListOperator=[ListOperator {''}];
105                end
106                if check_colorvar(ilist)
107                    Role{numel(ListVar)}='ancillary';% not projected with interpolation
108                    ProjModeRequest{numel(ListVar)}='';
109                else
110                    Role{numel(ListVar)}='scalar';
111                    ProjModeRequest{numel(ListVar)}='interp_lin';%scalar field (requires interpolation for plot)
112                end
113            else  % an operator 'vec' or 'norm' is used
114                if ~check_colorvar(ilist) && strcmp(r.Operator,'norm')
115                    ProjModeRequestVar='interp_lin';%scalar field (requires interpolation for plot)
116                else
117                    ProjModeRequestVar='';
118                end
119                ind_var_U=find(strcmp(r.UName,ListVar));%check previous listing of variable r.UName
120                ind_var_V=find(strcmp(r.VName,ListVar));%check previous listing of variable r.VName
121                if isempty(ind_var_U)
122                    ListVar=[ListVar {r.UName}]; % append the variable in the list if not previously listed
123                    Role=[Role {'vector_x'}];
124                    ProjModeRequest=[ProjModeRequest {ProjModeRequestVar}];
125                    ListInputField=[ListInputField InputField(ilist)];
126                else
127                    checkU=1;
128                end
129                if isempty(ind_var_V)
130                    ListVar=[ListVar {r.VName}];% append the variable in the list if not previously listed
131                    Role=[Role {'vector_y'}];
132                    ProjModeRequest=[ProjModeRequest {ProjModeRequestVar}];
133                    ListInputField=[ListInputField InputField(ilist)];
134                   
135                else
136                    checkV=1;
137                end
138            end
139        end
140        if ~isfield(ParamIn,'Coord_z')
141            ParamIn.Coord_z=[];
142        end
143        NbCoord=~isempty(ParamIn.Coord_x)+~isempty(ParamIn.Coord_y)+~isempty(ParamIn.Coord_z);
144        if isfield(ParamIn,'TimeDimName')% case of reading of a single time index in a multidimensional array
145            [Field,var_detect,ichoice,errormsg]=nc2struct(FileName,'TimeDimName',ParamIn.TimeDimName,num,[ParamIn.Coord_x ParamIn.Coord_y ParamIn.Coord_z ListVar]);
146        elseif isfield(ParamIn,'TimeVarName')% case of reading of a single time  in a multidimensional array
147            [Field,var_detect,ichoice,errormsg]=nc2struct(FileName,'TimeVarName',ParamIn.TimeVarName,num,[ParamIn.Coord_x ParamIn.Coord_y ParamIn.Coord_z ListVar]);
148            NbCoord=NbCoord+1;% adds time coordinate
149        else
150            [Field,var_detect,ichoice,errormsg]=nc2struct(FileName,[ParamIn.Coord_x ParamIn.Coord_y ParamIn.Coord_z ListVar]);
151        end
152        if ~isempty(errormsg)
153            return
154        end
155        %scan all the variables beyond the two first ones, ParamIn.Coord_x and ParamIn.Coord_y.
156        for ilist=NbCoord+1:numel(Field.VarDimName)
157            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
158                Field.VarAttribute{1}.Role='coord_x';%unstructured coordinates
159                Field.VarAttribute{2}.Role='coord_y';
160                if NbCoord>=3
161                    Field.VarAttribute{3}.Role='coord_z';
162                end
163                break
164            end
165        end
166        NormName='';
167        UName='';
168        VName='';
169        if numel(Field.ListVarName)>NbCoord % if there are variables beyond coord (1 D plots)
170        for ilist=1:numel(ListVar)
171            Field.VarAttribute{ilist+NbCoord}.Role=Role{ilist};
172            Field.VarAttribute{ilist+NbCoord}.ProjModeRequest=ProjModeRequest{ilist};
173            if isfield(ParamIn,'FieldName')
174                Field.VarAttribute{ilist+NbCoord}.FieldName=ListInputField{ilist};
175            end
176            r=regexp(ListInputField{ilist},'(?<Operator>(^vec|^norm))\((?<UName>.+),(?<VName>.+)\)$','names');
177            if ~isempty(r)&& strcmp(r.Operator,'norm')
178                NormName='norm';
179                if ~isempty(find(strcmp(ListVar,'norm')))
180                    NormName='norm_1';
181                end
182                Field.ListVarName=[Field.ListVarName {NormName}];
183                ilistmax=numel(Field.ListVarName);
184                Field.VarDimName{ilistmax}=Field.VarDimName{ilist+2};
185                Field.VarAttribute{ilistmax}.Role='scalar';
186                Field.(NormName)=Field.(r.UName).*Field.(r.UName)+Field.(r.VName).*Field.(r.VName);
187                Field.(NormName)=sqrt(Field.(NormName));
188                UName=r.UName;
189                VName=r.VName;
190            end
191        end
192
193        if ~isempty(NormName)% remove U and V if norm has been calculated and U and V are not needed as variables
194            ind_var_U=find(strcmp(UName,ListVar));%check previous listing of variable r.UName
195            ind_var_V=find(strcmp(VName,ListVar));%check previous listing of variable r.VName
196            if ~checkU && ~checkV
197                Field.ListVarName([ind_var_U+2 ind_var_V+2])=[];
198                Field.VarDimName([ind_var_U+2 ind_var_V+2])=[];
199                Field.VarAttribute([ind_var_U+2 ind_var_V+2])=[];
200            elseif ~checkU
201                Field.ListVarName(ind_var_U+2)=[];
202                Field.VarDimName(ind_var_U+2)=[];
203                Field.VarAttribute(ind_var_U+2 )=[];
204            elseif ~checkV
205                Field.ListVarName(ind_var_V+2)=[];
206                Field.VarDimName(ind_var_V+2)=[];
207                Field.VarAttribute(ind_var_V+2 )=[];
208            end
209        end
210        % insert coordinates as indices in case of plots vs matrix index
211        if isfield(ParamIn,'CheckCoordIndex') && ParamIn.CheckCoordIndex
212            Field.ListVarName=[Field.ListDimName Field.ListVarName];
213            Field.VarDimName=[Field.ListDimName Field.VarDimName];
214            for idim=1:numel(Field.ListDimName)
215                CoordName=Field.ListDimName{idim};
216                Field.(CoordName)=1:Field.DimValue(idim);
217            end
218            Field.VarAttribute=[cell(1,numel(Field.ListDimName)) Field.VarAttribute]
219        end
220                end
221    case 'video'
222        if strcmp(class(ParamIn),'VideoReader')
223            A=read(ParamIn,num);
224        else
225            ParamOut=VideoReader(FileName);
226            A=read(ParamOut,num);
227        end
228    case 'mmreader'
229        if strcmp(class(ParamIn),'mmreader')
230            A=read(ParamIn,num);
231        else
232            ParamOut=mmreader(FileName);
233            A=read(ParamOut,num);
234        end
235    case 'vol'
236        A=imread(FileName);
237        Npz=size(A,1)/ParamIn.Npy;
238        A=reshape(A',ParamIn.Npx,ParamIn.Npy,Npz);
239        A=permute(A,[3 2 1]);
240    case 'multimage'
241      %  warning 'off'
242        A=imread(FileName,num);
243    case 'image'
244        A=imread(FileName);
245    case 'rdvision'
246        [A,FileInfo,timestamps]=read_rdvision(FileName,num);
247    otherwise
248        errormsg=[ FileType ': invalid input file type for uvmat'];
249end
250
251%% case of image
252if ~isempty(A)
253    if strcmp(FileType,'rdvision')
254        Field.Time=timestamps;
255    end
256    if isstruct(ParamOut)
257        ParamOut.FieldName='image';
258    end
259    Npz=1;%default
260    npxy=size(A);
261    %     Rangx=[0.5 npxy(2)-0.5]; % coordinates of the first and last pixel centers
262    %     Rangy=[npxy(1)-0.5 0.5]; %
263    Field.NbDim=2;%default
264    Field.AName='image';
265    Field.ListVarName={'Coord_y','Coord_x','A'}; %
266    if ndims(A)==3
267        if Npz==1;%color
268            Field.VarDimName={'Coord_y','Coord_x',{'Coord_y','Coord_x','rgb'}}; %
269            Field.Coord_y=[npxy(1)-0.5 0.5];
270            Field.Coord_x=[0.5 npxy(2)-0.5]; % coordinates of the first and last pixel centers
271            if isstruct(ParamOut)
272                ParamOut.Npx=npxy(2);% display image size on the interface
273                ParamOut.Npy=npxy(1);
274            end
275%             Field.VarAttribute{3}.Mesh=1;
276        else
277            Field.NbDim=3;
278            Field.ListVarName=['AZ' Field.ListVarName];
279            Field.VarDimName={'AZ','Coord_y','Coord_x',{'AZ','Coord_y','Coord_x'}};
280            Field.AZ=[npxy(1)-0.5 0.5];
281            Field.Coord_y=[npxy(2)-0.5 0.5];
282            Field.Coord_x=[0.5 npxy(3)-0.5]; % coordinates of the first and last pixel centers
283            if isstruct(ParamOut)
284                ParamOut.Npx=npxy(3);% display image size on the interface
285                ParamOut.Npy=npxy(2);
286            end
287%             Field.VarAttribute{4}.Mesh=1;
288        end
289    else
290        Field.VarDimName={'Coord_y','Coord_x',{'Coord_y','Coord_x'}}; %
291        Field.Coord_y=[npxy(1)-0.5 0.5];
292        Field.Coord_x=[0.5 npxy(2)-0.5]; % coordinates of the first and last pixel centers
293%        ParamOut.Npx=npxy(2);% display image size on the interface
294 %       ParamOut.Npy=npxy(1);
295%         Field.VarAttribute{3}.Mesh=1;
296    end
297    Field.A=A;
298    Field.CoordUnit='pixel'; %used for mouse_motion
299       
300end
301
302
303
Note: See TracBrowser for help on using the repository browser.