Changeset 943 for trunk/src/series/merge_proj_polar.m
 Timestamp:
 May 2, 2016, 9:49:22 AM (8 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/src/series/merge_proj_polar.m
r941 r943 1 %'merge_proj': concatene several fields from series, can project them on a regular grid in phys coordinates1 %'merge_proj': concatene several fields from series, project on a polar grid 2 2 % 3 3 % function ParamOut=merge_proj(Param) … … 104 104 HeadData.azimuth=azimuth_arclength; 105 105 thresh2=16; % square of the interpolation range 106 106 107 %%%%%%%%%%%% STANDARD PART (DO NOT EDIT) %%%%%%%%%%%% 107 108 ParamOut=[]; %default output … … 160 161 161 162 %% prepare output file content 162 TimeData.ListGlobalAttribute={'Conventions','Project','CoordUnit','TimeUnit','ZPos'};163 TimeData.Conventions='uvmat';164 TimeData.Project='2016_Circumpolar';165 TimeData.CoordUnit='cm';166 TimeData.TimeUnit='s';167 TimeData.ZPos=0;168 TimeData.ListVarName={'time','radius','azimuth','U','V','curl','div'};169 TimeData.VarDimName={'time','radius','azimuth',{'time','radius','azimuth'},{'time','radius','azimuth'}...170 {'time','radius','azimuth'},{'time','radius','azimuth'}};171 TimeData.VarAttribute{1}.Role='';172 TimeData.VarAttribute{2}.Role='';173 TimeData.VarAttribute{3}.Role='';174 TimeData.VarAttribute{4}.Role='vector_x';175 TimeData.VarAttribute{5}.Role='vector_y';176 TimeData.VarAttribute{6}.Role='scalar';177 TimeData.VarAttribute{7}.Role='scalar';178 TimeData.time=nan(1,NbField);179 TimeData.radius=radius_shifted;180 TimeData.azimuth=azimuth_arclength;181 nby=numel(radius);182 nbx=numel(azimuth);183 TimeData.U=nan(NbField,nby,nbx);184 TimeData.V=nan(NbField,nby,nbx);185 TimeData.curl=nan(NbField,nby,nbx);186 TimeData.div=nan(NbField,nby,nbx);163 % TimeData.ListGlobalAttribute={'Conventions','Project','CoordUnit','TimeUnit','ZPos'}; 164 % TimeData.Conventions='uvmat'; 165 % TimeData.Project='2016_Circumpolar'; 166 % TimeData.CoordUnit='cm'; 167 % TimeData.TimeUnit='s'; 168 % TimeData.ZPos=0; 169 % TimeData.ListVarName={'time','radius','azimuth','U','V','curl','div'}; 170 % TimeData.VarDimName={'time','radius','azimuth',{'time','radius','azimuth'},{'time','radius','azimuth'}... 171 % {'time','radius','azimuth'},{'time','radius','azimuth'}}; 172 % TimeData.VarAttribute{1}.Role=''; 173 % TimeData.VarAttribute{2}.Role=''; 174 % TimeData.VarAttribute{3}.Role=''; 175 % TimeData.VarAttribute{4}.Role='vector_x'; 176 % TimeData.VarAttribute{5}.Role='vector_y'; 177 % TimeData.VarAttribute{6}.Role='scalar'; 178 % TimeData.VarAttribute{7}.Role='scalar'; 179 % TimeData.time=nan(1,NbField); 180 % TimeData.radius=radius_shifted; 181 % TimeData.azimuth=azimuth_arclength; 182 % nby=numel(radius); 183 % nbx=numel(azimuth); 184 % TimeData.U=nan(NbField,nby,nbx); 185 % TimeData.V=nan(NbField,nby,nbx); 186 % TimeData.curl=nan(NbField,nby,nbx); 187 % TimeData.div=nan(NbField,nby,nbx); 187 188 188 189 % if ~isempty(timeread) … … 254 255 for iview =2:numel(XmlData) 255 256 if ~(isfield(XmlData{iview},'GeometryCalib')&& isequal(XmlData{iview}.GeometryCalib.SliceCoord,ProjObjectCoord))... 256 disp('error: geometric calibration missing ')257 disp('error: geometric calibration missing or inconsistent plane positions') 257 258 return 258 259 end … … 285 286 286 287 % %% output file type 287 % FileExtOut='.nc'; %netcdf output 288 % if isempty(j1_series{1}) 289 % NomTypeOut='_1'; 290 % else 291 % NomTypeOut='_1_1'; 292 % end 293 % RootFileOut=RootFile{1}; 294 % for iview=2:NbView 295 % if ~strcmp(RootFile{iview},RootFile{1}) 296 % RootFileOut='mproj'; 297 % break 298 % end 299 % end 288 if isempty(j1_series{1}) 289 NomTypeOut='_1'; 290 else 291 NomTypeOut='_1_1'; 292 end 293 RootFileOut=RootFile{1}; 294 for iview=2:NbView 295 if ~strcmp(RootFile{iview},RootFile{1}) 296 RootFileOut='mproj'; 297 break 298 end 299 end 300 300 301 301 … … 313 313 for index=1:NbField 314 314 disp(['index=' num2str(index)]) 315 disp(['ellapsed time ' num2str(toc(tstart)/60, 2) ' minutes'])315 disp(['ellapsed time ' num2str(toc(tstart)/60,4) ' minutes']) 316 316 update_waitbar(WaitbarHandle,index/NbField) 317 317 if ~isempty(RUNHandle) && ~strcmp(get(RUNHandle,'BusyAction'),'queue') … … 319 319 return 320 320 end 321 322 %%% generating the name of the merged field323 %i1=i1_series{1}(index);324 %if ~isempty(i2_series{end})325 %i2=i2_series{end}(index);326 %else327 %i2=i1;328 %end329 %j1=1;330 %j2=1;331 %if ~isempty(j1_series{1})332 %j1=j1_series{1}(index);333 %if ~isempty(j2_series{end})334 %j2=j2_series{end}(index);335 %else336 %j2=j1;337 %end338 %end339 % % OutputFile=fullfile_uvmat(RootPath{1},OutputDir,RootFileOut,FileExtOut,NomTypeOut,i1,i2,j1,j2);340 %if ~CheckOverwrite && exist(OutputFile,'file')341 %disp(['existing output file ' OutputFile ' already exists, skip to next field'])342 %continue% skip iteration if the mode overwrite is desactivated and the result file already exists343 %end321 322 %% generating the name of the merged field 323 i1=i1_series{1}(index); 324 if ~isempty(i2_series{end}) 325 i2=i2_series{end}(index); 326 else 327 i2=i1; 328 end 329 j1=1; 330 j2=1; 331 if ~isempty(j1_series{1}) 332 j1=j1_series{1}(index); 333 if ~isempty(j2_series{end}) 334 j2=j2_series{end}(index); 335 else 336 j2=j1; 337 end 338 end 339 OutputFile=fullfile_uvmat(RootPath{1},OutputDir,RootFileOut,'.nc',NomTypeOut,i1,i2,j1,j2); 340 if ~CheckOverwrite && exist(OutputFile,'file') 341 disp(['existing output file ' OutputFile ' already exists, skip to next field']) 342 continue% skip iteration if the mode overwrite is desactivated and the result file already exists 343 end 344 344 345 345 %% z position 346 346 ZIndex=mod(i1_series{1}(index)1,NbSlice_calib{1})+1;%Zindex for phys transform 347 ZPos =ProjObjectCoord(ZIndex,3);347 ZPosNew=ProjObjectCoord(ZIndex,3); 348 348 if index==1 349 TimeData.ZPos=ZPos;349 ZPos=ZPosNew; 350 350 else 351 if ZPos ~=TimeData.ZPos351 if ZPosNew~=ZPos 352 352 disp('inconsistent z positions in the series') 353 353 return 354 354 end 355 355 end 356 % radius of the topography section at z pos 357 TopoRadius=0; 356 % radius of the topography section at z position 358 357 ind_mask=[]; 359 358 if ZPos<20 … … 372 371 return 373 372 end 374 disp([filecell{iview,index} ' read'])375 373 ListVar=Data{iview}.ListVarName; 376 374 for ilist=1:numel(ListVar) … … 438 436 end 439 437 438 440 439 %% time of the merged field: take the average of the different views 441 440 if ~isempty(time) … … 456 455 MergeData.div(ind_mask)=NaN; 457 456 end 458 TimeData.time(index)=timeread; 459 TimeData.U(index,:,:)=Unew; 460 TimeData.V(index,:,:)=Vnew; 461 TimeData.curl(index,:,:)=MergeData.curl; 462 TimeData.div(index,:,:)=MergeData.div; 463 464 %% recording the merged field 465 466 % MergeData.ListGlobalAttribute={'Conventions','Project','InputFile_1','InputFile_end','NbCoord','NbDim','CoordUnit','ZPos'}; 467 % MergeData.Conventions='uvmat'; 468 % MergeData.NbCoord=2; 469 % MergeData.NbDim=2; 470 % MergeData.CoordUnit=CoordUnit; 471 % MergeData.ZPos=ZPos; 472 % % time interval of PIV 473 % Dt=[]; 474 % if isfield(Data{1},'Dt')&& isnumeric(Data{1}.Dt) 475 % Dt=Data{1}.Dt; 476 % end 477 % for iview =2:numel(Data) 478 % if ~(isfield(Data{iview},'Dt')&& isequal(Data{iview}.Dt,Dt)) 479 % Dt=[];%dt not the same for all fields 480 % end 481 % end 482 % if ~isempty(timeread) 483 % MergeData.ListGlobalAttribute=[MergeData.ListGlobalAttribute {'Time'}]; 484 % MergeData.Time=timeread; 485 % end 486 % 487 % % time unit 488 % if isfield(Data{1},'TimeUnit') 489 % TimeUnit=Data{1}.TimeUnit; 490 % for iview =2:numel(Data) 491 % if ~(isfield(Data{iview},'TimeUnit')&& isequal(Data{iview}.TimeUnit,TimeUnit)) 492 % TimeUnit=[];%TimeUnit not the same for all fields 493 % end 494 % end 495 % if ~isempty(TimeUnit) 496 % MergeData.ListGlobalAttribute=[MergeData.ListGlobalAttribute {'TimeUnit'}]; 497 % MergeData.TimeUnit=TimeUnit; 498 % end 499 % end 500 501 end 502 error=struct2nc(OutputFile,TimeData);%save result file 503 if isempty(error) 504 disp(['output file ' OutputFile ' written']) 505 else 506 disp(error) 507 end 457 [npy,npx]=size(Unew); 458 459 %% create the output file for the first iteration of the loop 460 if index==1 461 TimeData.ListGlobalAttribute={'Conventions','Project','CoordUnit','TimeUnit','ZPos','Time'}; 462 TimeData.Conventions='uvmat'; 463 TimeData.Project='2016_Circumpolar'; 464 TimeData.CoordUnit='cm'; 465 TimeData.TimeUnit='s'; 466 TimeData.ZPos=ZPos; 467 TimeData.ListVarName={'radius','azimuth','U','V','curl','div'}; 468 TimeData.VarDimName={'radius','azimuth',{'radius','azimuth'},{'radius','azimuth'}... 469 {'radius','azimuth'},{'radius','azimuth'}}; 470 TimeData.VarAttribute{1}.Role=''; 471 TimeData.VarAttribute{2}.Role=''; 472 TimeData.VarAttribute{3}.Role='vector_x'; 473 TimeData.VarAttribute{4}.Role='vector_y'; 474 TimeData.VarAttribute{5}.Role='scalar'; 475 TimeData.VarAttribute{6}.Role='scalar'; 476 477 TimeData.radius=radius_shifted; 478 TimeData.azimuth=azimuth_arclength; 479 end 480 481 %% append data to the netcdf file for next iterations 482 TimeData.Time=timeread; 483 TimeData.U=Unew; 484 TimeData.V=Vnew; 485 TimeData.curl=MergeData.curl; 486 TimeData.div=MergeData.div; 487 488 [error,ncid]=struct2nc(OutputFile,TimeData);%save result file 489 if isempty(error) 490 disp(['output file ' OutputFile ' written']) 491 else 492 disp(error) 493 end 494 ellapsed_time=toc(tstart); 495 disp(['total ellapsed time ' num2str(ellapsed_time/60,2) ' minutes']) 496 end 497 508 498 ellapsed_time=toc(tstart); 509 499 disp(['total ellapsed time ' num2str(ellapsed_time/60,2) ' minutes'])
Note: See TracChangeset
for help on using the changeset viewer.