source: trunk/src/read_field.m @ 723

Last change on this file since 723 was 693, checked in by sommeria, 11 years ago

bug repaired geometry_calib/Detect grid, + corrections in get_field

File size: 10.4 KB
RevLine 
[497]1%'read_field': read the fields from files in different formats (netcdf files, images, video)
[181]2%--------------------------------------------------------------------------
[450]3%  function [Field,ParamOut,errormsg] = read_field(FileName,FileType,ParamIn,num)
[181]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
[466]15% FileName: name of the input file
[497]16% FileType: type of file, as determined by the function get_file_type.m
[450]17% ParamIn: movie object or Matlab structure of input parameters
[445]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'...)
[181]20%     .ColorVar: variable used for vector color
21%     .Npx, .Npy: nbre of pixels along x and y (used for .vol input files)
[693]22%     .TimeDimName: name of the dimension considered as 'time', selected index value then set by input 'num'   
[466]23% num: frame number for movies
[497]24%
25% see also read_image.m,read_civxdata.m,read_civdata.m,
[466]26
[450]27function [Field,ParamOut,errormsg] = read_field(FileName,FileType,ParamIn,num)
[181]28Field=[];
[354]29if ~exist('num','var')
30    num=1;
31end
32if ~exist('ParamIn','var')
33    ParamIn=[];
34end
35ParamOut=ParamIn;%default
[181]36errormsg='';
[635]37if ~exist(FileName,'file')
38    erromsg=['input file ' FileName ' does not exist'];
39    return
40end
[334]41A=[];
[575]42InputField={};
43check_colorvar=0;
[526]44if isstruct(ParamIn)
[575]45    if isfield(ParamIn,'FieldName')
46        if ischar(ParamIn.FieldName)
[576]47            InputField={ParamIn.FieldName};
[575]48        else
49            InputField= ParamIn.FieldName;
50        end
[526]51    end
[654]52    check_colorvar=zeros(size(InputField));
[684]53    if isfield(ParamIn,'ColorVar')&&~isempty(ParamIn.ColorVar)
[526]54        InputField=[ParamIn.FieldName {ParamIn.ColorVar}];
[654]55        check_colorvar(numel(InputField))=1;
[526]56    end
[521]57end
[575]58
[334]59%% distingush different input file types
[527]60switch FileType
61    case 'civdata'
62        [Field,ParamOut.VelType,errormsg]=read_civdata(FileName,InputField,ParamIn.VelType);
63        if ~isempty(errormsg),errormsg=['read_civdata / ' errormsg];return,end
64        ParamOut.CivStage=Field.CivStage;
65    case 'civx'
66        ParamOut.FieldName='velocity';%Civx data found, set .FieldName='velocity' by default
67        [Field,ParamOut.VelType,errormsg]=read_civxdata(FileName,InputField,ParamIn.VelType);
68        if ~isempty(errormsg),errormsg=['read_civxdata / ' errormsg];return,end
69        ParamOut.CivStage=Field.CivStage;
70    case 'netcdf'
71        ListVar={};
[675]72        Role={};
73        ProjModeRequest={};
[684]74        ListInputField={};
75        ListOperator={};
[681]76        checkU=0;
77        checkV=0;
[527]78        for ilist=1:numel(InputField)
79            r=regexp(InputField{ilist},'(?<Operator>(^vec|^norm))\((?<UName>.+),(?<VName>.+)\)$','names');
[654]80            if isempty(r)%  no operator used
[681]81                if isempty(find(strcmp(InputField{ilist},ListVar)))
[684]82                    ListVar=[ListVar InputField(ilist)];%append the variable name if not already in the list
83                    ListInputField=[ListInputField InputField(ilist)];
84                    ListOperator=[ListOperator {''}];
[681]85                end
[654]86                if check_colorvar(ilist)
[667]87                    Role{numel(ListVar)}='ancillary';% not projected with interpolation
88                    ProjModeRequest{numel(ListVar)}='';
[654]89                else
[667]90                    Role{numel(ListVar)}='scalar';
91                    ProjModeRequest{numel(ListVar)}='interp_lin';%scalar field (requires interpolation for plot)
[654]92                end
93            else  % an operator 'vec' or 'norm' is used
[684]94                %Operator=r.Operator;
[667]95                if ~check_colorvar(ilist) && strcmp(r.Operator,'norm')
[675]96                    ProjModeRequestVar='interp_lin';%scalar field (requires interpolation for plot)
[667]97                else
[675]98                    ProjModeRequestVar='';
[667]99                end
[675]100                ind_var_U=find(strcmp(r.UName,ListVar));%check previous listing of variable r.UName
101                ind_var_V=find(strcmp(r.VName,ListVar));%check previous listing of variable r.VName
102                if isempty(ind_var_U)
[684]103                    ListVar=[ListVar {r.UName}]; % append the variable in the list if not previously listed
[675]104                    Role=[Role {'vector_x'}];
105                    ProjModeRequest=[ProjModeRequest {ProjModeRequestVar}];
[684]106                    ListInputField=[ListInputField InputField(ilist)];
107                    %ListOperator=[ListOperator {[r.Operator '_U']}];
[681]108                else
109                    checkU=1;
[675]110                end
111                if isempty(ind_var_V)
[684]112                    ListVar=[ListVar {r.VName}];% append the variable in the list if not previously listed
[675]113                    Role=[Role {'vector_y'}];
114                    ProjModeRequest=[ProjModeRequest {ProjModeRequestVar}];
[684]115                    ListInputField=[ListInputField {''}];
116                    %ListOperator=[ListOperator {[r.Operator '_V']}];
[681]117                else
118                    checkV=1;
[675]119                end
[493]120            end
[527]121        end
[648]122        if isfield(ParamIn,'TimeDimName')% case of reading of a single time index in a multidimensional array
123            [Field,var_detect,ichoice]=nc2struct(FileName,'TimeDimName',ParamIn.TimeDimName,num,[ParamIn.Coord_x (ParamIn.Coord_y) ListVar]);
124        else
[667]125            [Field,var_detect,ichoice]=nc2struct(FileName,[ParamIn.Coord_x (ParamIn.Coord_y) ListVar]);
[648]126        end
[534]127        if isfield(Field,'Txt')
128            errormsg=Field.Txt;
129            return
130        end
[675]131        for ilist=3:numel(Field.VarDimName)
132            if isequal(Field.VarDimName{1},Field.VarDimName{ilist})
133                Field.VarAttribute{1}.Role='coord_x';%unstructured coordinates
134                Field.VarAttribute{2}.Role='coord_y';
135                break
136            end
137        end
[684]138        NormName='';
139        UName='';
140        VName='';
141        for ilist=1:numel(ListVar)
142            Field.VarAttribute{ilist+2}.Role=Role{ilist};
143            Field.VarAttribute{ilist+2}.ProjModeRequest=ProjModeRequest{ilist};
[648]144            if isfield(ParamIn,'FieldName')
[684]145                Field.VarAttribute{ilist+2}.FieldName=ListInputField{ilist};
[648]146            end
[684]147            r=regexp(ListInputField{ilist},'(?<Operator>(^vec|^norm))\((?<UName>.+),(?<VName>.+)\)$','names');
148            if ~isempty(r)&& strcmp(r.Operator,'norm')
149                NormName='norm';
150                if ~isempty(find(strcmp(ListVar,'norm')))
151                    NormName='norm_1';
152                end
153                Field.ListVarName=[Field.ListVarName {NormName}];
154                ilistmax=numel(Field.ListVarName);
155                Field.VarDimName{ilistmax}=Field.VarDimName{ilist+2};
156                Field.VarAttribute{ilistmax}.Role='scalar';
157                Field.(NormName)=Field.(r.UName).*Field.(r.UName)+Field.(r.VName).*Field.(r.VName);
158                Field.(NormName)=sqrt(Field.(NormName));
159                UName=r.UName;
160                VName=r.VName;
[681]161            end
162        end
[684]163        if ~isempty(NormName)% remove U and V if norm has been calculated and U and V are not needed as variables
164            ind_var_U=find(strcmp(UName,ListVar));%check previous listing of variable r.UName
165            ind_var_V=find(strcmp(VName,ListVar));%check previous listing of variable r.VName
166             if ~checkU && ~checkV
167                    Field.ListVarName([ind_var_U+2 ind_var_V+2])=[];
168                    Field.VarDimName([ind_var_U+2 ind_var_V+2])=[];
169                    Field.VarAttribute([ind_var_U+2 ind_var_V+2])=[];
170                elseif ~checkU
171                    Field.ListVarName(ind_var_U+2)=[];
172                    Field.VarDimName(ind_var_U+2)=[];
173                    Field.VarAttribute(ind_var_U+2 )=[];
174                elseif ~checkV
175                    Field.ListVarName(ind_var_V+2)=[];
176                    Field.VarDimName(ind_var_V+2)=[];
177                    Field.VarAttribute(ind_var_V+2 )=[];
178             end
179        end
[527]180    case 'video'
181        if strcmp(class(ParamIn),'VideoReader')
182            A=read(ParamIn,num);
183        else
184            ParamOut=VideoReader(FileName);
185            A=read(ParamOut,num);
186        end
187    case 'mmreader'
188        if strcmp(class(ParamIn),'mmreader')
189            A=read(ParamIn,num);
190        else
191            ParamOut=mmreader(FileName);
192            A=read(ParamOut,num);
193        end
194    case 'vol'
195        A=imread(FileName);
196        Npz=size(A,1)/ParamIn.Npy;
197        A=reshape(A',ParamIn.Npx,ParamIn.Npy,Npz);
198        A=permute(A,[3 2 1]);
199    case 'multimage'
200        warning 'off'
201        A=imread(FileName,num);
202    case 'image'
203        A=imread(FileName);
204end
205if ~isempty(errormsg)
206    errormsg=[FileType ' input: ' errormsg];
[445]207    return
[334]208end
[445]209
[334]210%% case of image
211if ~isempty(A)
[452]212    if isstruct(ParamOut)
[580]213        ParamOut.FieldName='image';
[452]214    end
[182]215    Npz=1;%default
[181]216    npxy=size(A);
[580]217    %     Rangx=[0.5 npxy(2)-0.5]; % coordinates of the first and last pixel centers
218    %     Rangy=[npxy(1)-0.5 0.5]; %
[181]219    Field.NbDim=2;%default
220    Field.AName='image';
221    Field.ListVarName={'AY','AX','A'}; %
222    if ndims(A)==3
223        if Npz==1;%color
224            Field.VarDimName={'AY','AX',{'AY','AX','rgb'}}; %
225            Field.AY=[npxy(1)-0.5 0.5];
226            Field.AX=[0.5 npxy(2)-0.5]; % coordinates of the first and last pixel centers
[452]227            if isstruct(ParamOut)
[580]228                ParamOut.Npx=npxy(2);% display image size on the interface
229                ParamOut.Npy=npxy(1);
[452]230            end
[186]231            Field.VarAttribute{3}.Mesh=1;
[181]232        else
233            Field.NbDim=3;
234            Field.ListVarName=['AZ' Field.ListVarName];
235            Field.VarDimName={'AZ','AY','AX',{'AZ','AY','AX'}};
236            Field.AZ=[npxy(1)-0.5 0.5];
237            Field.AY=[npxy(2)-0.5 0.5];
238            Field.AX=[0.5 npxy(3)-0.5]; % coordinates of the first and last pixel centers
[452]239            if isstruct(ParamOut)
[580]240                ParamOut.Npx=npxy(3);% display image size on the interface
241                ParamOut.Npy=npxy(2);
242            end
[186]243            Field.VarAttribute{4}.Mesh=1;
[181]244        end
245    else
246        Field.VarDimName={'AY','AX',{'AY','AX'}}; %
247        Field.AY=[npxy(1)-0.5 0.5];
248        Field.AX=[0.5 npxy(2)-0.5]; % coordinates of the first and last pixel centers
[182]249        ParamOut.Npx=npxy(2);% display image size on the interface
250        ParamOut.Npy=npxy(1);
[186]251        Field.VarAttribute{3}.Mesh=1;
[181]252    end
253    Field.A=A;
254    Field.CoordUnit='pixel'; %used for mouse_motion
255end
256
257
[334]258
Note: See TracBrowser for help on using the repository browser.