Changeset 606 for trunk/src/series/merge_proj.m
 Timestamp:
 Apr 7, 2013, 10:14:45 AM (8 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/src/series/merge_proj.m
r605 r606 41 41 function ParamOut=merge_proj(Param) 42 42 43 %% set the input elements needed on the GUI series when the function is selected in the menu ActionName 43 %% set the input elements needed on the GUI series when the function is selected in the menu ActionName or InputTable refreshed 44 44 if isstruct(Param) && isequal(Param.Action.RUN,0) 45 ParamOut.AllowInputSort='off';...% allow alphabetic sorting of the list of input file SubDir (options 'off'/'on', 'off' by default) 46 ParamOut.WholeIndexRange='on';...% prescribes the file index ranges from min to max (options 'off'/'on', 'off' by default) 47 ParamOut.NbSlice='off'; ...%nbre of slices ('off' by default) 48 ParamOut.VelType='one';...% menu for selecting the velocity type (options 'off'/'one'/'two', 'off' by default) 49 ParamOut.FieldName='one';...% menu for selecting the field (s) in the input file(options 'off'/'one'/'two', 'off' by default) 50 ParamOut.FieldTransform = 'on';...%can use a transform function 51 ParamOut.ProjObject='on';...%can use projection object(option 'off'/'on', 52 ParamOut.Mask='off';...%can use mask option (option 'off'/'on', 'off' by default) 45 ParamOut.AllowInputSort='off';% allow alphabetic sorting of the list of input file SubDir (options 'off'/'on', 'off' by default) 46 ParamOut.WholeIndexRange='on';% prescribes the file index ranges from min to max (options 'off'/'on', 'off' by default) 47 ParamOut.NbSlice='off'; %nbre of slices ('off' by default) 48 ParamOut.VelType='one';% menu for selecting the velocity type (options 'off'/'one'/'two', 'off' by default) 49 ParamOut.FieldName='one';% menu for selecting the field (s) in the input file(options 'off'/'one'/'two', 'off' by default) 50 ParamOut.FieldTransform = 'on';%can use a transform function 51 ParamOut.TransformPath=fullfile(fileparts(which('uvmat')),'transform_field');% path to transform functions (needed for compilation only) 52 ParamOut.ProjObject='on';%can use projection object(option 'off'/'on', 53 ParamOut.Mask='off';%can use mask option (option 'off'/'on', 'off' by default) 53 54 ParamOut.OutputDirExt='.mproj';%set the output dir extension 54 55 ParamOut.OutputFileMode='NbInput';% '=NbInput': 1 output file per input file index, '=NbInput_i': 1 file per input file index i, '=NbSlice': 1 file per slice … … 57 58 msgbox_uvmat('WARNING','the first input file does not exist') 58 59 elseif isequal(size(Param.InputTable,1),1) && ~isfield(Param,'ProjObject') 59 60 end 61 return60 msgbox_uvmat('WARNING','a projection object of type plane needs to be introduced for merge_proj') 61 end 62 return 62 63 end 63 64 64 65 %%%%%%%%%%%% STANDARD PART (DO NOT EDIT) %%%%%%%%%%%% 65 66 67 66 %% read input parameters from an xml file if input is a file name (batch mode) 68 67 checkrun=1; … … 72 71 end 73 72 74 ParamOut= Param;%default output73 ParamOut=[] %default output 75 74 if ~isfield(Param,'InputFields') 76 75 Param.InputFields.FieldName=''; 77 76 end 78 Output SubDir=[Param.OutputSubDir Param.OutputDirExt];% subdirectory for output files77 OutputDir=[Param.OutputSubDir Param.OutputDirExt];% subdirectory for output files 79 78 80 79 %% root input file type … … 93 92 % i1_series(iview,fileindex) expresses the same indices as a 1D array in file indices 94 93 %%%%%%%%%%%% 95 NbSlice=1;%default 96 if isfield(Param.IndexRange,'NbSlice')&&~isempty(Param.IndexRange.NbSlice) 97 NbSlice=Param.IndexRange.NbSlice; 98 end 99 nbview=numel(i1_series);%number of input file series (lines in InputTable) 100 nbfield_j=size(i1_series{1},1); %nb of fields for the j index (bursts or volume slices) 101 nbfield_i=size(i1_series{1},2); %nb of fields for the i index 102 nbfield=nbfield_j*nbfield_i; %total number of fields 103 nbfield_i=floor(nbfield/NbSlice);%total number of indexes in a slice (adjusted to an integer number of slices) 104 nbfield=nbfield_i*NbSlice; %total number of fields after adjustement 94 % NbSlice=1;%default 95 % if isfield(Param.IndexRange,'NbSlice')&&~isempty(Param.IndexRange.NbSlice) 96 % NbSlice=Param.IndexRange.NbSlice; 97 % end 98 NbView=numel(i1_series);%number of input file series (lines in InputTable) 99 NbField_j=size(i1_series{1},1); %nb of fields for the j index (bursts or volume slices) 100 NbField_i=size(i1_series{1},2); %nb of fields for the i index 101 NbField=NbField_j*NbField_i; %total number of fields 105 102 106 103 %determine the file type on each line from the first input file 107 104 ImageTypeOptions={'image','multimage','mmreader','video'}; 108 105 NcTypeOptions={'netcdf','civx','civdata'}; 109 for iview=1: nbview106 for iview=1:NbView 110 107 if ~exist(filecell{iview,1}','file') 111 108 displ_uvmat('ERROR',['the first input file ' filecell{iview,1} ' does not exist'],checkrun) … … 122 119 end 123 120 124 125 121 %% calibration data and timing: read the ImaDoc files 126 122 [XmlData,NbSlice_calib,time,errormsg]=read_multimadoc(RootPath,SubDir,RootFile,FileExt,i1_series,i2_series,j1_series,j2_series); … … 133 129 134 130 %% coordinate transform or other user defined transform 135 transform_fct='';%default 131 % transform_fct='';%default fct handle 132 % if isfield(Param,'FieldTransform')&&~isempty(Param.FieldTransform.TransformName) 133 % if isdeployed 134 % transform_fct=Param.FieldTransform.TransformName; 135 % dd=phys([]);%activate phys for compilation 136 % else 137 % currentdir=pwd; 138 % cd(Param.FieldTransform.TransformPath) 139 % transform_fct=str2func(Param.FieldTransform.TransformName); 140 % cd (currentdir) 141 % end 142 % end 143 transform_fct='';%default fct handle 136 144 if isfield(Param,'FieldTransform')&&~isempty(Param.FieldTransform.TransformName) 137 addpath(Param.FieldTransform.TransformPath)138 transform_fct=str2func(Param.FieldTransform.TransformName);139 rmpath(Param.FieldTransform.TransformPath)140 end 141 145 currentdir=pwd; 146 cd(Param.FieldTransform.TransformPath) 147 transform_fct=str2func(Param.FieldTransform.TransformName); 148 cd (currentdir) 149 end 142 150 %%%%%%%%%%%% END STANDARD PART %%%%%%%%%%%% 143 151 % EDIT FROM HERE … … 152 160 return 153 161 end 154 for iview=1: nbview162 for iview=1:NbView 155 163 if ~isequal(CheckImage{iview},CheckImage{1})~isequal(CheckNc{iview},CheckNc{1}) 156 164 displ_uvmat('ERROR','input set of input series: need either netcdf either image series',checkrun) … … 169 177 %% MAIN LOOP ON SLICES 170 178 %%%%%%%%%%%%% STANDARD PART (DO NOT EDIT) %%%%%%%%%%%% 171 for i_slice=1:NbSlice172 index_slice=i_slice:NbSlice:nbfield;% select file indices of the slice173 nbfiles=0;174 nbmissing=0;179 % for i_slice=1:NbSlice 180 % index_slice=i_slice:NbSlice:NbField;% select file indices of the slice 181 % NbFiles=0; 182 % nbmissing=0; 175 183 176 184 %%%%%%%%%%%%%%%% loop on field indices %%%%%%%%%%%%%%%% 177 for index=index_slice 178 if checkrun 179 stopstate=get(Param.RUNHandle,'BusyAction'); 180 update_waitbar(Param.WaitbarHandle,index/nbfield) 181 else 182 stopstate='queue'; 183 end 184 if ~isequal(stopstate,'queue')% enable STOP command 185 return 186 end 187 %%%%%%%%%%%%%%%% loop on views (input lines) %%%%%%%%%%%%%%%% 188 Data=cell(1,nbview);%initiate the set Data 189 nbtime=0; 190 for iview=1:nbview 191 %% reading input file(s) 192 [Data{iview},tild,errormsg] = read_field(filecell{iview,index},FileType{iview},Param.InputFields,frame_index{iview}(index)); 185 for index=1:NbField 186 if checkrun 187 stopstate=get(Param.RUNHandle,'BusyAction'); 188 update_waitbar(Param.WaitbarHandle,index/NbField) 189 else 190 stopstate='queue'; 191 end 192 if ~isequal(stopstate,'queue')% enable STOP command 193 return 194 end 195 %%%%%%%%%%%%%%%% loop on views (input lines) %%%%%%%%%%%%%%%% 196 Data=cell(1,NbView);%initiate the set Data 197 nbtime=0; 198 for iview=1:NbView 199 %% reading input file(s) 200 [Data{iview},tild,errormsg] = read_field(filecell{iview,index},FileType{iview},Param.InputFields,frame_index{iview}(index)); 201 if ~isempty(errormsg) 202 errormsg=['merge_proj/read_field/' errormsg]; 203 display(errormsg) 204 break 205 end 206 timeread(iview)=0; 207 if isfield(Data{iview},'Time') 208 timeread(iview)=Data{iview}.Time; 209 nbtime=nbtime+1; 210 end 211 if ~isempty(NbSlice_calib) 212 Data{iview}.ZIndex=mod(i1_series{iview}(index)1,NbSlice_calib{iview})+1;%Zindex for phys transform 213 end 214 215 %% transform the input field (e.g; phys) if requested 216 if ~isempty(transform_fct) 217 if nargin(transform_fct)>=2 218 Data{iview}=transform_fct(Data{iview},XmlData{iview}); 219 else 220 Data{iview}=transform_fct(Data{iview}); 221 end 222 end 223 % Data{iview}=phys(Data{iview},XmlData{iview}); 224 %% check whether tps is needed, then calculate tps coefficients if needed 225 check_tps=0; 226 if isfield(Param.InputFields,'FieldName') 227 if ischar(Param.InputFields.FieldName) 228 Param.InputFields.FieldName={Param.InputFields.FieldName}; 229 end 230 else 231 Param.InputFields.FieldName={}; 232 end 233 for ilist=1:numel(Param.InputFields.FieldName) 234 switch Param.InputFields.FieldName{ilist} 235 case {'vort','div','strain'} 236 check_tps=1; 237 end 238 end 239 240 %% calculate tps coeff if needed 241 check_proj_tps= ~isempty(Param.ProjObject)&& strcmp(Param.ProjObject.ProjMode,'filter')&&~isfield(Data{iview},'Coord_tps'); 242 Data{iview}=tps_coeff_field(Data{iview},check_proj_tps); 243 244 %% projection on object (gridded plane) 245 if Param.CheckObject 246 [Data{iview},errormsg]=proj_field(Data{iview},Param.ProjObject); 193 247 if ~isempty(errormsg) 194 errormsg=['merge_proj/read_field/' errormsg]; 195 display(errormsg) 196 break 197 end 198 timeread(iview)=0; 199 if isfield(Data{iview},'Time') 200 timeread(iview)=Data{iview}.Time; 201 nbtime=nbtime+1; 202 end 203 if ~isempty(NbSlice_calib) 204 Data{iview}.ZIndex=mod(i1_series{iview}(index)1,NbSlice_calib{iview})+1;%Zindex for phys transform 205 end 206 207 %% transform the input field (e.g; phys) if requested 208 if ~isempty(transform_fct) 209 if nargin(transform_fct)>=2 210 Data{iview}=transform_fct(Data{iview},XmlData{iview}); 211 else 212 Data{iview}=transform_fct(Data{iview}); 213 end 214 end 215 216 %% check whether tps is needed, then calculate tps coefficients if needed 217 check_tps=0; 218 if isfield(Param.InputFields,'FieldName') 219 if ischar(Param.InputFields.FieldName) 220 Param.InputFields.FieldName={Param.InputFields.FieldName}; 221 end 222 else 223 Param.InputFields.FieldName={}; 224 end 225 for ilist=1:numel(Param.InputFields.FieldName) 226 switch Param.InputFields.FieldName{ilist} 227 case {'vort','div','strain'} 228 check_tps=1; 229 end 230 end 231 232 %% calculate tps coeff if needed 233 check_proj_tps= ~isempty(Param.ProjObject)&& strcmp(Param.ProjObject.ProjMode,'filter')&&~isfield(Data{iview},'Coord_tps'); 234 Data{iview}=tps_coeff_field(Data{iview},check_proj_tps); 235 236 %% projection on object (gridded plane) 237 if Param.CheckObject 238 [Data{iview},errormsg]=proj_field(Data{iview},Param.ProjObject); 239 if ~isempty(errormsg) 240 displ_uvmat('ERROR',['error in merge_proge/proj_field: ' errormsg],checkrun) 241 return 242 end 243 end 244 end 245 %END LOOP ON VIEWS 246 247 %% merge the nbview fields 248 MergeData=merge_field(Data); 249 if isfield(MergeData,'Txt') 250 displ_uvmat('ERROR',MergeData.Txt,checkrun) 251 return 252 end 253 254 % time of the merged field: 255 if ~isempty(time)% time defined from ImaDoc 256 timeread=time(:,index); 257 end 258 timeread=mean(timeread); 259 260 % generating the name of the merged field 261 i1=i1_series{iview}(index); 262 if ~isempty(i2_series{iview}) 263 i2=i2_series{iview}(index); 264 else 265 i2=i1; 266 end 267 j1=1; 268 j2=1; 269 if ~isempty(j1_series{iview}) 270 j1=j1_series{iview}(index); 271 if ~isempty(j2_series{iview}) 272 j2=j2_series{iview}(index); 273 else 274 j2=j1; 275 end 276 end 277 OutputFile=fullfile_uvmat(RootPath{1},OutputSubDir,RootFile{1},FileExtOut,NomType{1},i1,i2,j1,j2); 278 279 % recording the merged field 280 if CheckImage{1} %in case of input images an image is produced 281 if isa(MergeData.A,'uint8') 282 bitdepth=8; 283 elseif isa(MergeData.A,'uint16') 284 bitdepth=16; 285 end 286 imwrite(MergeData.A,OutputFile,'BitDepth',bitdepth); 287 %write xml calibration file 288 siz=size(MergeData.A); 289 npy=siz(1); 290 npx=siz(2); 291 if isfield(MergeData,'VarAttribute')&&isfield(MergeData.VarAttribute{1},'Coord_2')&&isfield(MergeData.VarAttribute{1},'Coord_1') 292 Rangx=MergeData.VarAttribute{1}.Coord_2; 293 Rangy=MergeData.VarAttribute{1}.Coord_1; 294 elseif isfield(MergeData,'AX')&& isfield(MergeData,'AY') 295 Rangx=[MergeData.AX(1) MergeData.AX(end)]; 296 Rangy=[MergeData.AY(1) MergeData.AY(end)]; 297 else 298 Rangx=[0.5 npx0.5]; 299 Rangy=[npy0.5 0.5];%default 300 end 301 pxcmx=(npx1)/(Rangx(2)Rangx(1)); 302 pxcmy=(npy1)/(Rangy(1)Rangy(2)); 303 T_x=pxcmx*Rangx(1)+0.5; 304 T_y=pxcmy*Rangy(2)+0.5; 305 GeometryCal.focal=1; 306 GeometryCal.R=[pxcmx,0,0;0,pxcmy,0;0,0,1]; 307 GeometryCal.Tx_Ty_Tz=[T_x T_y 1]; 308 ImaDoc.GeometryCalib=GeometryCal; 309 % t=struct2xml(ImaDoc); 310 % t=set(t,1,'name','ImaDoc'); 311 % save(t,[filebase_merge '.xml']) 312 % display([filebase_merge '.xml saved']) 313 else 314 MergeData.ListGlobalAttribute={'Conventions','Project','InputFile_1','InputFile_end','nb_coord','nb_dim','dt','Time','civ'}; 315 MergeData.Conventions='uvmat'; 316 MergeData.nb_coord=2; 317 MergeData.nb_dim=2; 318 dt=[]; 319 if isfield(Data{1},'dt')&& isnumeric(Data{1}.dt) 320 dt=Data{1}.dt; 321 end 322 for iview =2:numel(Data) 323 if ~(isfield(Data{iview},'dt')&& isequal(Data{iview}.dt,dt)) 324 dt=[];%dt not the same for all fields 325 end 326 end 327 if isempty(dt) 328 MergeData.ListGlobalAttribute(6)=[]; 329 else 330 MergeData.dt=dt; 331 end 332 MergeData.Time=timeread; 333 error=struct2nc(OutputFile,MergeData);%save result file 334 if isempty(error) 335 display(['output file ' OutputFile ' written']) 336 else 337 display(error) 338 end 339 end 340 end 341 end 248 displ_uvmat('ERROR',['error in merge_proge/proj_field: ' errormsg],checkrun) 249 return 250 end 251 end 252 end 253 %END LOOP ON VIEWS 254 255 %% merge the NbView fields 256 MergeData=merge_field(Data); 257 if isfield(MergeData,'Txt') 258 displ_uvmat('ERROR',MergeData.Txt,checkrun) 259 return 260 end 261 262 % time of the merged field: 263 if ~isempty(time)% time defined from ImaDoc 264 timeread=time(:,index); 265 end 266 timeread=mean(timeread); 267 268 % generating the name of the merged field 269 i1=i1_series{iview}(index); 270 if ~isempty(i2_series{iview}) 271 i2=i2_series{iview}(index); 272 else 273 i2=i1; 274 end 275 j1=1; 276 j2=1; 277 if ~isempty(j1_series{iview}) 278 j1=j1_series{iview}(index); 279 if ~isempty(j2_series{iview}) 280 j2=j2_series{iview}(index); 281 else 282 j2=j1; 283 end 284 end 285 OutputFile=fullfile_uvmat(RootPath{1},OutputDir,RootFile{1},FileExtOut,NomType{1},i1,i2,j1,j2); 286 287 % recording the merged field 288 if CheckImage{1} %in case of input images an image is produced 289 if isa(MergeData.A,'uint8') 290 bitdepth=8; 291 elseif isa(MergeData.A,'uint16') 292 bitdepth=16; 293 end 294 imwrite(MergeData.A,OutputFile,'BitDepth',bitdepth); 295 %write xml calibration file 296 siz=size(MergeData.A); 297 npy=siz(1); 298 npx=siz(2); 299 if isfield(MergeData,'VarAttribute')&&isfield(MergeData.VarAttribute{1},'Coord_2')&&isfield(MergeData.VarAttribute{1},'Coord_1') 300 Rangx=MergeData.VarAttribute{1}.Coord_2; 301 Rangy=MergeData.VarAttribute{1}.Coord_1; 302 elseif isfield(MergeData,'AX')&& isfield(MergeData,'AY') 303 Rangx=[MergeData.AX(1) MergeData.AX(end)]; 304 Rangy=[MergeData.AY(1) MergeData.AY(end)]; 305 else 306 Rangx=[0.5 npx0.5]; 307 Rangy=[npy0.5 0.5];%default 308 end 309 pxcmx=(npx1)/(Rangx(2)Rangx(1)); 310 pxcmy=(npy1)/(Rangy(1)Rangy(2)); 311 T_x=pxcmx*Rangx(1)+0.5; 312 T_y=pxcmy*Rangy(2)+0.5; 313 GeometryCal.focal=1; 314 GeometryCal.R=[pxcmx,0,0;0,pxcmy,0;0,0,1]; 315 GeometryCal.Tx_Ty_Tz=[T_x T_y 1]; 316 ImaDoc.GeometryCalib=GeometryCal; 317 % t=struct2xml(ImaDoc); 318 % t=set(t,1,'name','ImaDoc'); 319 % save(t,[filebase_merge '.xml']) 320 % display([filebase_merge '.xml saved']) 321 else 322 MergeData.ListGlobalAttribute={'Conventions','Project','InputFile_1','InputFile_end','nb_coord','nb_dim','dt','Time','civ'}; 323 MergeData.Conventions='uvmat'; 324 MergeData.nb_coord=2; 325 MergeData.nb_dim=2; 326 dt=[]; 327 if isfield(Data{1},'dt')&& isnumeric(Data{1}.dt) 328 dt=Data{1}.dt; 329 end 330 for iview =2:numel(Data) 331 if ~(isfield(Data{iview},'dt')&& isequal(Data{iview}.dt,dt)) 332 dt=[];%dt not the same for all fields 333 end 334 end 335 if isempty(dt) 336 MergeData.ListGlobalAttribute(6)=[]; 337 else 338 MergeData.dt=dt; 339 end 340 MergeData.Time=timeread; 341 error=struct2nc(OutputFile,MergeData);%save result file 342 if isempty(error) 343 display(['output file ' OutputFile ' written']) 344 else 345 display(error) 346 end 347 end 348 end 349 342 350 343 351 %'merge_field': concatene fields … … 351 359 MergeData=Data{1};%default 352 360 error=0; 353 nbview=length(Data);354 if nbview==1361 NbView=length(Data); 362 if NbView==1 355 363 return 356 364 end … … 382 390 for ivar=VarIndex 383 391 VarName=MergeData.ListVarName{ivar}; 384 for iview=1: nbview392 for iview=1:NbView 385 393 MergeData.(VarName)=[MergeData.(VarName); Data{iview}.(VarName)]; 386 394 end … … 389 397 else 390 398 testFF=0; 391 for iview=2: nbview399 for iview=2:NbView 392 400 for ivar=VarIndex 393 401 VarName=MergeData.ListVarName{ivar}; … … 401 409 end 402 410 if testFF 403 nbaver= nbviewMergeData.FF;411 nbaver=NbViewMergeData.FF; 404 412 indgood=find(nbaver>0); 405 413 for ivar=VarIndex … … 410 418 for ivar=VarIndex 411 419 VarName=MergeData.ListVarName{ivar}; 412 MergeData.(VarName)=double(MergeData.(VarName))./ nbview;420 MergeData.(VarName)=double(MergeData.(VarName))./NbView; 413 421 end 414 422 end
Note: See TracChangeset
for help on using the changeset viewer.