source: trunk/src/read_field.m @ 749

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

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